crear una característica geométrica puntual temporal a partir de coordenadas en OGR con Python

Question

Tengo dos archivos shape que contienen dos tipos de características geométricas: una contiene polígonos y la otra cadenas de líneas.

Porque me gustaría utilizar el Dentro del método para comprobar si el final (o el principio) de una cadena de líneas se encuentra dentro de un polígono situado en el otro archivo shape, necesito acceder a los puntos primero y último (como características geométricas). Por el momento itero a través de la característica linestring y puedo identificar los puntos utilizando geometry.GetPoint(x) que devuelve las coordenadas x e y de ese punto.

Pero cuando me gustaría usar within para comprobar si un punto determinado cae en el objeto polígono no puedo porque como explicó Paolo en esta pregunta Necesito la referencia geométrica (utilizando GetGeometryRef())de un punto que se puede recuperar como:

...
pointshp = driver.Open('test.shp')
pointlyr = pointshp.getLayer(0)
point = pointlyr.GetNextFeature()
point_geom = point.GetGeometryRef()
point_geom.Within(polygon_geom)

pero en realidad estoy abriendo un shp que tiene características linestring. Así que después de hojear Utah State University Diapositivas del SIG se me ocurrió esto:

def within(pointcoord):
driver = ogr.GetDriverByName('test.shp')
ds = driver.CreateDataSource('test.shp')
layer = ds.CreateLayer('test',geom_type=ogr.wkbPoint)
fieldDefn = ogr.FieldDefn('id',ogr.OFTInteger)
layer.CreateFiled(fieldDefn)
featureDefn = layer.GetLayerDefn()
feature = ogr.feature(featureDefn)

point = ogr.Geometry(ogr.wkbPoint)
point.AddPoint(pointcoord[0],pointcoord[1])
feature.SetGeometry(point)
feature.SetField('id',1)
layer.CreateFeature(feature)
pointlyr = layer.GetLayer(0)
point = pointlyr.GetNextFeature()
point_geom = point.GetGeometryRef()
for x in xrange(len(polygons))
        point_geom.Within(Network.polygons[x].GetGeometryRef())

Lo que me parece realmente ineficiente. Así que sólo por el bien de la utilización de un método(dentro) Tengo que crear no sólo una capa, pero un shapefile, cuando todo lo que realmente necesito es una simple característica de punto desde el que puedo obtener la referencia de geometría para compararlo con el de los polígonos.

La pregunta es, ¿cómo crear una característica de punto (temporal, para ser destruido) a partir de coordenadas x e y? ¿O hay una manera más sencilla de hacerlo que el ejemplo anterior?

Por otro lado, ¿hay alguna otra forma de obtener la referencia geométrica del punto inicial o final (como elemento geométrico) de un elemento de línea?

Answer 1

Utilizando las coordenadas de la característica de geometría linestring es posible crear simplemente un objeto de geometría como:

 point = ogr.Geometry(ogr.wkbPoint)
 point.AddPoint(g.GetX(j), g.GetY(j))#where j is the iteration in the linestring i.e. 0

que luego se puede utilizar con: point.Within(polygon)

Answer 2

La API de Python no está realmente bien documentada en este momento, pero es paralela a las API de C/C++. 1

Dicho esto, has enlazado una referencia de la API de Perl que no es de código directo (en tkk.fi), puede o no ser lo que estás buscando. Sin embargo, la API de C++ para LineString (en gdal.osgeo.org/ogr/classOGRLineString.html ) tiene métodos StartPoint y EndPoint, por lo que deberían reflejarse en Python, pero no es así debido a la forma en que se puede hacer referencia a objetos en Python.

Así que lo que tienes que hacer es averiguar cómo crear LineStrings en los bindings de python, utilizando tus datos. Mirando tus diapositivas, puedes ver que un ejemplo de LineString está definido en la diapositiva 4, pero vas a necesitar cambiarlo para importarlo desde tu archivo.

Lo bueno de python, es que usted puede encender el REPL, ver qué cosas son:

>>> import ogr
>>> dir(ogr.Geometry)
 ... # clipped
>>> dir(ogr.Geometry.Within)
 ... # clipped
>>> print ogr.Geometry.Within
<unbound method Geometry.Within>
>>> print ogr.Geometry.Within.__doc__
Within(self, Geometry other) -> bool

int OGR_G_Within(OGRGeometryH hThis,
OGRGeometryH hOther) 
>>> line = ogr.Geometry(ogr.wkbLineString) # the constant value in the constructor here determines the type of geometry it is
>>> print line
LINESTRING EMPTY
>>> line.AddPoint(10,10)
>>> line.AddPoint(20,20)
>>> line.AddPoint(30,30)
>>> line.GetPoint(1)
(20.0, 20.0, 0.0)
>>> line.GetPoint(1) # start point
(20.0, 20.0, 0.0)
>>> line.GetPoint(line.GetPointCount()-1) # end point
(30.0, 30.0, 0.0)

Extrapolando un poco, puedes llamar a Within en los puntos inicial y final de la línea para cada objeto poligonal.