Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

I'm using the Python OGR API to add new features to a PostGIS layer. However, despite setting the geometry spatial reference, when the geometry is written out the spatial reference is missing. This is a problem as the destination table has check constraint forcing geom_srid to be a particular type.

For example:

ogr.UseExceptions()

# Create a feature in an existing layer (where layer is a PostGIS table)
feature = ogr.Feature(layer.GetLayerDefn())

# Create spatial reference
out_srs = ogr.osr.SpatialReference()
out_srs.ImportFromEPSG(27700)

# Assign to geometry
geom.AssignSpatialReference(out_srs)

# Set feature geometry
feature.SetGeometry(geom)

# Create the feature
layer.CreateFeature(feature)
feature.Destroy()

Returns:

Error: new row for relation "table" violates check constraint 
"enforce_srid_geom"
Command: INSERT INTO "table" ("geom", "attr1") VALUES (GeomFromEWKT('SRID=-1;POINT (5000, 1000)'::TEXT), 'row1')

At this point it is clear that OGR is converting the geometry to an EWKT string without pulling in the spatial reference attribute.

Specific questions:

1) Does anyone know how I can force the geometry output from SetGeometry() to include an SRID when it is converted to EWKT?

2) Has anyone found this with the CPP bindings?

3) Is this a bug in OGR? From OGR source the output string is being created in ogrpgtablelayer.cpp where nSRSID is being set to -1 if no valid SRS input is found.

Many thanks,

Tom

share|improve this question
1  
As an update it would appear that OGR assigns the SRS from the target layer and not from individual geometries. So as long as OGR can find the target tables SRS (i.e it's listed in geometry_columns) then the above code will work. However, note that this means no checking of individual geometry SRS is done on insert, meaning that geomgetry in the wrong coordinate system/format can be inserted even in enforce_geomtype_geom check is applied! –  Tomas Feb 3 '12 at 12:25

1 Answer 1

Your approach seems correct and I am not sure where the problem is. As an alternate, you can try using ImportFromProj4() and supply a proj4 string for EPSG 27700 and see if that works.

share|improve this answer
    
Thanks for your answer, but no luck I'm afraid. It's definitely not the assignment that's the problem (I can do geom.GetSpatialReference() to print the assigned SRS before feature.SetGeometry()). I've submitted a ticket to trac: trac.osgeo.org/gdal/ticket/4487 –  Tomas Feb 2 '12 at 22:10

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.