Como encontrar o polígono dentro do qual um ponto se encontra?

8

Eu tenho uma camada com recursos poligonais. Cada recurso possui atributos e valores. Também tenho uma lista de coordenadas e gostaria de saber em qual característica (ou polígono) as coordenadas se encontram.

Alguém poderia me orientar sobre como fazer isso? Existe uma função na API que pode me ajudar a atingir meu objetivo ou devo usar algum algoritmo de geometria computacional para fazer isso sozinho? Eu sei como fazer o último, mas me pouparia algum tempo se já houvesse uma função incorporada.

Obrigado.

Shubham Goyal
fonte

Respostas:

6

Poucos recursos

O que você provavelmente quer fazer é:

  1. Crie uma lista de QgsPoint a partir de suas coordenadas
  2. Iterar todos os recursos da sua camada (polígonos)
  3. Passa sobre a lista de pontos (dentro da iteração)
  4. Ligue para feature.geometry (). Contains (point) para verificar se você tem uma correspondência

É claro que agora você pode melhorar o desempenho se, por exemplo, souber que um ponto só pode ser contido por um polígono; você pode remover um ponto da lista, depois que o polígono apropriado for encontrado.

Muitos recursos (Usando SpatialIndex)

Como apontado nos comentários, um índice espacial pode ser usado para acelerar significativamente o processo.

Os passos aqui seriam

  1. Crie uma lista de QgsPoint a partir de suas coordenadas
  2. Crie um QgsSpatialIndex
  3. Iterar todos os recursos da sua camada (polígonos)
  4. Adicione os recursos ao seu índice com insertFeature
  5. Iterar todos os seus pontos
  6. Ligue para index.intersects (QgsRectangle (ponto, ponto)) para verificar se há uma correspondência

Há também um exemplo de código de NathanW

Matthias Kuhn
fonte
Ah, eu não sabia sobre a chamada feature.geometry.contains (point). Eu usei mathplotlib. pastebin.com/61LSeMWv Por favor, perdoe a bagunça do meu código. Estou com pressa e vou limpá-lo mais tarde.
Shubham Goyal
Eu não implementei sua solução e, portanto, não posso dizer com certeza se ela funciona ou não. Mas eu acredito que deveria e por isso estou marcando este como a resposta correta :)
Shubham Goyal
2
Isso poderia ser acelerado usando um QgsSpatialIndex?
Snorfalorpagus
1
Eu recomendo o uso de um QgsSpatailIndex. Existe um exemplo em nathanw.net/2013/01/04/…
Nathan W
6

Primeiro, você precisa importar a lista de coordenadas para o seu projeto. Este tutorial explica bem como fazer isso: http://qgis.spatialthoughts.com/2012/01/importing-spreadsheets-or-csv-files-to.html

Quando você tiver duas camadas (polígonos e pontos) em seu projeto, vá para vetor> ferramentas de gerenciamento de dados> associar atributos por local

insira a descrição da imagem aqui

Você obtém uma janela onde é possível definir quais camadas você deseja combinar:

insira a descrição da imagem aqui

  • Defina sua camada de pontos como a 'camada de vetor de destino'.
  • Defina sua camada de polígono como 'unir camada de vetor'.
  • Defina seu Shapefile de Saída (será criado um novo. Portanto, se você perdeu, os dados originais são preservados).
  • Você pode optar por manter todos os dados no novo arquivo de forma, mesmo que não haja correspondência com um polígono: marque 'manter todos os registros (incluindo registros de destino não correspondentes)'

Clique OK'. O novo shapefile é criado e você será perguntado 'Deseja adicionar a nova camada ao sumário?' Clique novamente em OK.

Abra a atribuição de atributo do novo shapefile adicionado e você verá que todos os recursos do polígono correspondente são adicionados ao ponto estabelecido nesse polígono.

PieterB
fonte
2

Uma maneira mais simples de fazer isso usando o PyQGIS. Imaginei que você pode construir um QgsRectangleobjeto com um único ponto e usá-lo QgsFeatureRequestpara filtrar recursos da camada que o intercepta.

point = QgsPoint(10, 10)
# Construct a QgsRectange with the same point
rect = QgsRectangle(point, point)
req = QgsFeatureRequest()
req.setFilterRect(rect)
# You get the feature that intersects the point
f = layer.getFeatures(req).next()
pensamentos espaciais
fonte
0

No QuantumGIS, você pode adicionar a lista de coordenadas com a função 'adicionar camada de texto delimitada' (se for um arquivo csv). Adicione também os polígonos. Então você pode fazer uma 'interseção' ou 'pontos no polígono'.

Stefan
fonte