Como realizar cálculos GIS em bancos de dados não GIS

Introdução:

Ocasionalmente, você pode encontrar a necessidade de executar funções geoespaciais em seu aplicativo, como mapear localizações de usuários ou analisar dados geográficos. Existem inúmeras bibliotecas específicas de linguagem disponíveis para essas tarefas, como GDAL, Shapely e Geopandas para Python.

Alternativamente, a funcionalidade geoespacial pode ser implementada através de bancos de dados; por exemplo, você pode usar a extensão PostGIS com um banco de dados relacional como o PostgreSQL ou aproveitar o suporte nativo para tipos de dados espaciais em um banco de dados distribuído como o Azure CosmosDB.

No entanto, se o seu armazenamento de back-end, como Redis ou Google Spanner, não oferece suporte nativo a consultas espaciais e você precisa lidar com consultas geoespaciais em grande escala, este artigo foi feito sob medida para você.

Quais são minhas opções?

Você sempre pode criar outro microsserviço para lidar com dados espaciais, mas essa opção geralmente envolve a sobrecarga de manutenção de um aplicativo adicional. Outra abordagem é usar bibliotecas de indexação geográfica como S2 e H3. S2, desenvolvido pelo Google, é baseado na curva de Hilbert, enquanto H3, desenvolvido pela Uber, é baseado em um sistema de grade global geodésico discreto. S2 e H3 compartilham muitas semelhanças: ambos dividem uma determinada região em células e usam números inteiros de 64 bits para indexar essas células.

Porém, a principal diferença está no formato das células; S2 usa células quadradas, enquanto H3 usa células hexagonais. Para algumas aplicações, o H3 pode oferecer melhor desempenho. No entanto, no geral, qualquer uma das bibliotecas deve ser suficiente. Neste artigo usaremos S2, mas você pode executar funções semelhantes usando H3.

Conceitos básicos da biblioteca Google S2

• Células: S2 divide a esfera em células, cada uma com um identificador exclusivo de 64 bits.

• Níveis de células: A hierarquia permite diferentes níveis de detalhe, desde grandes regiões até pequenas áreas precisas. Cada nível representa uma resolução diferente:

  
  

  • Nível 0: As células maiores, cobrindo uma porção significativa da superfície da Terra.

    
    

  • Níveis mais altos: As células são progressivamente subdivididas em quadrantes menores. Por exemplo, cada célula de Nível 1 é dividida em quatro células de Nível 2 e assim por diante.

    
    

  • Resolução e Área: Níveis mais altos correspondem a resoluções mais finas e áreas de células menores. Essa hierarquia permite indexação e consulta precisas em vários níveis de detalhe.

    
    

A tabela abaixo mostra vários níveis de células juntamente com suas áreas correspondentes.

A partir da tabela fornecida, é evidente que você pode obter uma precisão de mapeamento de até 0,44 cm^2 usando S2. Dentro de cada quadrado de uma célula S2, existe uma célula filha que compartilha o mesmo pai, indicando uma estrutura hierárquica. O nível da célula pode ser um valor estático (mesmo nível aplicado a todas as células) ou pode ser dinâmico onde S2 decide qual resolução funciona melhor.

Calculando os vizinhos mais próximos

Vamos começar com um exemplo. Considere que estamos escrevendo um aplicativo que fornece recursos semelhantes a serviços de proximidade para a área de Seattle. Queremos retornar uma lista de cafeterias nas proximidades. Para realizar estas operações, dividiremos esta tarefa em 4 subtarefas:

• Carregando mapa de Seattle
• Visualize células S2 no mapa de Seattle
• Armazene alguns locais de cafeterias no banco de dados
• Consultar cafeterias mais próximas

Carregando mapa de Seattle

Para carregar um mapa do Google, usaríamos a biblioteca gmplot. Esta biblioteca requer uma chave de API do Google Maps para carregar. Para gerar a chave API, siga as instruções aqui .

 import gmplot import const # plot seattle with zoom level 13 gmap = gmplot.GoogleMapPlotter(47.6061, -122.3328, 13, apikey=const.API_KEY) # Draw the map to an HTML file: gmap.draw('map.html')

O código acima gera um arquivo map.html conforme mostrado abaixo:

Visualize células S2 no mapa de Seattle

Agora que temos o mapa, vamos desenhar algumas células S2 para mapas:

 from s2 import * import gmplot # plot seattle with zoom level 13 gmap = gmplot.GoogleMapPlotter(47.6061, -122.3328, 13, apikey=const.API_KEY) areatobeplotted = [ (47.64395531736767,-122.43597221319135), (47.51369277846956,-122.43597221319135), (47.51369277846956,-122.24156866779164), (47.64395531736767,-122.24156866779164), (47.64395531736767,-122.43597221319135) ] region_rect = S2LatLngRect( S2LatLng.FromDegrees(47.51369277846956,-122.43597221319135), S2LatLng.FromDegrees(47.64395531736767, -122.24156866779164)) coverer = S2RegionCoverer() coverer.set_min_level(8) coverer.set_max_level(15) covering = coverer.GetCovering(region_rect) geoms = 0 for cellid in covering: new_cell = S2Cell(cellid) vertices = [] for i in range(0, 4): vertex = new_cell.GetVertex(i) latlng = S2LatLng(vertex) vertices.append((latlng.lat().degrees(), latlng.lng().degrees())) gmap.polygon(*zip(*vertices), face_color='pink', edge_color='cornflowerblue', edge_width=5) geoms+=1 gmap.polygon(*zip(*areatobeplotted), face_color='red', edge_color='green', edge_width=5) print(f"Total Geometries: {geoms}") gmap.draw('/tmp/map.html')

 Output: Total Geometries: 273

No código acima, primeiro centralizamos o plotter do Google Map em torno da área de Seattle. Em S2RegionCoverer , inicializamos o cobridor de região para ter níveis dinâmicos entre um nível mínimo de 8 e um nível máximo de 15. Isso permite que S2 ajuste dinamicamente todas as células em tamanhos de células específicos para o melhor ajuste. O método GetCovering retorna a cobertura de um retângulo ao redor da área de Seattle.

Em seguida, iteramos sobre cada célula, calculando os vértices das células e plotando-os no mapa. Mantemos a contagem de células geradas em torno de 273. Finalmente, traçamos o retângulo de entrada em vermelho. Este código irá plotar as células S2 no mapa de Seattle em /tmp/map.html , conforme mostrado abaixo:

Armazene alguns locais de cafeterias no banco de dados

Vamos gerar um banco de dados de cafeterias junto com seus identificadores de células S2. Você pode armazenar essas células em qualquer banco de dados de sua escolha. Para este tutorial, usaremos um banco de dados SQLite. No exemplo de código abaixo, nos conectamos ao banco de dados SQLite para criar uma tabela CoffeeShops com 3 campos Id , name e cell_id .

Semelhante ao exemplo anterior, usamos S2RegionCoverer para calcular as células, mas desta vez usamos um nível fixo para plotar pontos. Finalmente, o ID calculado é convertido em uma string e armazenado no banco de dados.

 import sqlite3 from s2 import S2CellId,S2LatLng,S2RegionCoverer # Connect to SQLite database conn = sqlite3.connect('/tmp/sqlite_cells.db') cursor = conn.cursor() # Create a table to store cell IDs cursor.execute('''CREATE TABLE IF NOT EXISTS CoffeeShops ( id INTEGER PRIMARY KEY AUTOINCREMENT, name TEXT, cell_id TEXT )''') coverer = S2RegionCoverer() # Function to generate S2 cell ID for a given latitude and longitude def generate_cell_id(latitude, longitude, level=16): cell=S2CellId(S2LatLng.FromDegrees(latitude, longitude)) return str(cell.parent(level)) # Function to insert cell IDs into the database def insert_cell_ids(name,lat,lng): cell_id = generate_cell_id(lat, lng) cursor.execute("INSERT INTO CoffeeShops (name, cell_id) VALUES (?, ?)", (name, cell_id)) conn.commit() # Insert cell IDs into the database insert_cell_ids("Overcast Coffee", 47.616656277302155, -122.31156460382837) insert_cell_ids("Seattle Sunshine", 47.67366852914391, -122.29051997415843) insert_cell_ids("Sip House", 47.6682364706238, -122.31328618043693) insert_cell_ids("Victoria Coffee",47.624408595334536, -122.3117362652041) # Close connection conn.close()

Neste ponto, temos um banco de dados que armazena cafeterias junto com seus IDs de células, determinados pela resolução selecionada para o nível de célula.

Consulte as cafeterias mais próximas

Por fim, vamos consultar cafeterias na área do Distrito Universitário.

 import sqlite3 from s2 import S2RegionCoverer,S2LatLngRect, S2LatLng # Connect to SQLite database conn = sqlite3.connect('/tmp/sqlite_cells.db') cursor = conn.cursor() # Function to query database for cells intersecting with the given polygon def query_intersecting_cells(start_x,start_y,end_x,end_y): # Create S2RegionCoverer region_rect = S2LatLngRect( S2LatLng.FromDegrees(start_x,start_y), S2LatLng.FromDegrees(end_x,end_y)) coverer = S2RegionCoverer() coverer.set_min_level(8) coverer.set_max_level(15) covering = coverer.GetCovering(region_rect) # Query for intersecting cells intersecting_cells = set() for cell_id in covering: cursor.execute("SELECT name FROM CoffeeShops WHERE cell_id >= ? and cell_id<=?", (str(cell_id.range_min()),str(cell_id.range_max()),)) intersecting_cells.update(cursor.fetchall()) return intersecting_cells # Query for intersecting cells intersecting_cells = query_intersecting_cells(47.6527847,-122.3286438,47.6782181, -122.2797203) # Print intersecting cells print("Intersecting cells:") for cell_id in intersecting_cells: print(cell_id[0]) # Close connection conn.close()

 Output: Intersecting cells: Sip House Seattle Sunshine

Abaixo está uma representação visual das células. Para manter a brevidade, o código de visualização abaixo não foi adicionado.

Como todas as células filhas e pais compartilham um prefixo, podemos consultar intervalos de células entre mínimo e máximo para obter todas as células entre esses dois valores. Em nosso exemplo, usamos o mesmo princípio para consultar a cafeteria

Conclusão:

Neste artigo, demonstramos como usar a indexação geográfica para armazenar e consultar dados geoespaciais em bancos de dados que não oferecem suporte a consultas geoespaciais. Isso pode ser estendido a vários casos de uso, como cálculo de roteamento entre 2 pontos ou obtenção de vizinhos mais próximos.

Normalmente, para consultas de bancos de dados indexados geograficamente, você teria que realizar algum pós-processamento adicional nos dados. É necessária uma consideração cuidadosa da lógica de consulta e pós-processamento para garantir que não sobrecarregaremos o nó.

Referências:

• Plotador de mapas do Google - https://github.com/gmplot/gmplot/wiki/GoogleMapPlotter
• Guia do desenvolvedor S2 - http://s2geometry.io/devguide/
• SQLite - https://docs.python.org/3/library/sqlite3.html
• Geoespacial do Azure Cosmos DB - https://learn.microsoft.com/en-us/azure/cosmos-db/nosql/query/geospatial?tabs=javascript
• GDAL - https://gdal.org/index.html
• Shapley- https://shapely.readthedocs.io/en/stable/manual.html
• Geopandas - https://geopandas.org
• Posgis - https://postgis.net/
• Chave inglesa do Google - https://cloud.google.com/spanner?hl=en
• Redis - https://redis.io/