Hi Peter,
On Mon, Feb 14, 2011 at 4:03 PM, Peter Tittmann <[email protected]> wrote:
> Hi all,
> I am in need of making maps from postgis data and am sort of stymied by the
> range of options and difficulty in finding an example. My apologies if its
> out there and I've just not been able to find it. My goal is to be able to
> plot multiple point, polygon, and possibly raster layers on a matplotlib
> figure.
> Iv'e tried Jose Gomez-Dans example
> (http://sites.google.com/site/jgomezdans/ogr,pythonymatplotlib) but he seems
> to be just using postgis to get the SRID of the display extent argument
> transformed.
>
> I've tried the example from
> here: http://trac.gispython.org/lab/browser/Shapely/trunk/examples/world.py?rev=749
> substituting an ExecuteSQL for GetLayersByName like this:
> import ogr
> import pylab
> from numpy import asarray
> from shapely.wkb import loads
> source = ogr.Open("PG:host=localhost dbname=fb user=*** password= ***")
> borders = source.ExecuteSQL("select std_id as stand, boundary from
> vaneck.vaneck_units")
>
> fig = pylab.figure(1, figsize=(4,2), dpi=300)
> while 1:
> feature = borders.GetNextFeature()
> if not feature:
> break
>
> geom = loads(feature.GetGeometryRef().ExportToWkb())
> a = asarray(geom)
> pylab.plot(a[:,0], a[:,1])
> pylab.show()
> but asarray cannot process the syntax:
> Traceback (most recent call last):
> File "<stdin>", line 20, in <module>
> IndexError: invalid index
>
> At this point I'm inclined to generate the wdb using a query via Psycopg
> (instead of the above OGR method) but am stuck trying to figure out how to
> plot the multi polygon returned to shapely using matplotlib.
> I'd be most grateful for any suggestions.
> Peter
The asarray function only works on geometric objects that have a
coordinate sequence attribute: line strings and rings. If your objects
are multipolygons or polygons, you'll need to access their parts and
rings:
for polygon in geom:
a = asarray(polygon.exterior)
pylab.plot(a[:,0], a[:,1])
for ring in polygon.interiors:
b = asarray(ring)
pylab.plot(b[:,0], b[:,1])
--
Sean
_______________________________________________
Community mailing list
[email protected]
http://lists.gispython.org/mailman/listinfo/community
