[FWTools] setting the projection on a Shape file

[Frank Warmerdam]

Moreau, Jean-Sebastien wrote:
> Hi list!
> 
> I'm trying to create a shape file with a polygon in it and I want to set the
> coordinate system to NAD83. Here's my python script:
> ----------------------------------
> import ogr,osr
> 
> wkt = 'POLYGON ((%f %f, %f %f, %f %f, %f %f, %f %f))' %
> (llx,lly,ulx,uly,urx,ury,lrx,lry,llx,lly)
> 
> driver=ogr.GetDriverByName('ESRI Shapefile')
> datasource=driver.CreateDataSource(shapePath)
> 
> layer=datasource.CreateLayer('cadre',geom_type=ogr.wkbPolygon)
> geom=ogr.CreateGeometryFromWkt(wkt)
> 
> spatial=osr.SpatialReference()
> spatial.SetWellKnownGeogCS("NAD83")
> 
> #Assignation de la géoréference à l'objet géométrique
> 		
> geom.AssignSpatialReference(spatial)
> 		
> #Activation de la géométrie dans le feature puis dans le layer
> 		
> feat=ogr.Feature(layer.GetLayerDefn())
> feat.SetGeometry(geom)
> layer.CreateFeature(feat)
> 		
> #Fermeture du fichier
> 		
> datasource.Destroy()
> -----------------------------------------
> Assuming the llx,lly, etc are correct coordinates, I should get a shape file
> with the correct coordinate system. However when I open it in ArcMap and go
> to the properties, it says the coordinate system is NAD27.
> 
> Am I doing something wrong?

JS,

For shapefiles the spatial reference set on the geometry is ignored (and
this is true of essentially everything else in OGR), but you can pass
the spatial ref to the layer creation and it will be written out as a .prj
file.

eg.

  spatial=osr.SpatialReference()
  spatial.SetWellKnownGeogCS("NAD83")
  layer=datasource.CreateLayer('cadre',geom_type=ogr.wkbPolygon,
                               srs = spatial )

Good luck,
-- 
---------------------------------------+--------------------------------------
I set the clouds in motion - turn up   | Frank Warmerdam, warmerdam at pobox.com
light and sound - activate the windows | http://pobox.com/~warmerdam
and watch the world go round - Rush    | President OSGF, http://osgeo.org