[FWTools] Get projection from shp using GDAL and Python

David Fawcett david.fawcett at gmail.com
Mon Oct 25 11:10:40 EST 2010


My guess is that you are trying to write out the spatial reference
object, not the text representation of that object.

Instead of this line:

 wksht.row(row).write(8, layer.GetSpatialRef())

Try something like:

srsObj = layer.GetSpatialRef()
srsText = srsObj.ExportToWkt()
wksht.row(row).write(8, layer.srsText)


This is a little verbose, but should give you the idea.

David.

On Mon, Oct 25, 2010 at 9:38 AM, Susana Iraiis Delgado Rodriguez
<susana.delgado_s at utzmg.edu.mx> wrote:
> Hello list!
>
> I designed a python module to get all the information I need from a
> shapefile. I've been using GDAL to work with my script and I already
> collected all the information I need. My application saves the values in an
> excel sheet, but I haven't gotten the projection info, the code is:
>
> import os, time, socket
> from xlwt import Workbook
> from osgeo import ogr,gdal,osr
> from dbf import *
> gdal.AllRegister()
> file_list = []
> folders = None
>  # look in this (root) folder for files with specified extension
> for root, folders, files in os.walk( "C:\\" ):
>  file_list.extend(os.path.join(root,fi) for fi in files if
> fi.endswith(".shp"))
> wrkbk = Workbook()
> wksht = wrkbk.add_sheet('shp')
> wksht.row(0).write(0,'ruta')
> wksht.row(0).write(1,'archivo')
> wksht.row(0).write(2,'x_min')
> wksht.row(0).write(3,'x_max')
> wksht.row(0).write(4,'y_min')
> wksht.row(0).write(5,'y_max')
> wksht.row(0).write(6,'geometria')
> wksht.row(0).write(7,'num_elem')
> wksht.row(0).write(8,'proyeccion')
> wksht.row(0).write(9,'prj')
> wksht.row(0).write(10,'estructura bd')
> wksht.row(0).write(11,'fecha_modificacion')
> wksht.row(0).write(12,'maquina_host')
> wksht.row(0).write(13,'usuario')
> for row, filepath in enumerate(file_list, start=1):
>  #Get path and filename
>  wksht.row(row).write(0, filepath)
>  #Get filename
>  (ruta, filename) = os.path.split(filepath)
>  wksht.row(row).write(1, filename)
>  # Get extent
>  shapeData = ogr.Open(filepath)
>  layer = shapeData.GetLayer()
>  feature = layer.GetNextFeature()
>  x_min, x_max, y_min, y_max = layer.GetExtent()
>  x_y = layer.GetExtent()
>  wksht.row(row).write(2, x_y[0])
>  wksht.row(row).write(3, x_y[1])
>  wksht.row(row).write(4, x_y[2])
>  wksht.row(row).write(5, x_y[3])
>  #Get geometry type;1=Point, 2=LineString, 3=Polygon
>  defn = layer.GetLayerDefn()
>  wksht.row(row).write(6, defn.GetGeomType())
>  #Get featurecount()
>  wksht.row(row).write(7, layer.GetFeatureCount())
>  #Verify prj from shp
>  n = os.path.splitext(filepath)
>  p = n[0]+'.prj'
>  if os.path.lexists(p):
>   wksht.row(row).write(9, 1)
>  else:
>   wksht.row(row).write(9, 0)
>  #Get host from the file
>  wksht.row(row).write(12, socket.gethostname())
>
>  #Extarct database
>  t = n[0]+'_bd.txt'
>  d = n[0]+'.dbf'
>  if os.path.lexists(d):
>   a = open (t,"w+")
>   dbf = Dbf(d,new=False)
>
>   for fldName in dbf.fieldDefs:
>    a.write(fldName.name)
>    a.write(" || ")
>    a.write(fldName.typeCode)
>    a.write("\n")
>   dbf.close()
>   a.close()
>   wksht.row(row).write(10, t)
>  else:
>   print "El archivo " +n[0]+".shp" " no tiene dbf"
>   wksht.row(row).write(10, "Sin bd")
>
> #Get projection, this is the part that fails
>  print layer.GetSpatialRef()
>  wksht.row(row).write(8, layer.GetSpatialRef())
>
>  #Get file last modification
>  t = time.strftime("%m/%d/%Y %I:%M:%S
> %p",time.localtime(os.path.getmtime(filepath)))
>  wksht.row(row).write(11, t)
>  #Get username
>  wksht.row(row).write(13,os.environ.get("USERNAME"))
>
> wrkbk.save('shp.xls')
>
> SEARCH_PATH = os.getcwd()
> TARGET_FILE = os.path.realpath('shp.xls')
> print "Buscando en", SEARCH_PATH, "and writing to", TARGET_FILE
> print "Encontrando archivos..."
> print "Escribiendo archivo de Excel..."
> print "Listo."
>
> I got the next output:
> Python 2.6.6 (r266:84297, Aug 24 2010, 18:46:32) [MSC v.1500 32 bit (Intel)]
> on
> win32
> Type "help", "copyright", "credits" or "license" for more information.
>>>> import crawler_shp
> None
> None
> PROJCS["WGS_1984_UTM_Zone_13N",
>     GEOGCS["GCS_WGS_1984",
>         DATUM["WGS_1984",
>             SPHEROID["WGS_1984",6378137,298.257223563]],
>         PRIMEM["Greenwich",0],
>         UNIT["Degree",0.017453292519943295]],
>     PROJECTION["Transverse_Mercator"],
>     PARAMETER["latitude_of_origin",0],
>     PARAMETER["central_meridian",-105],
>     PARAMETER["scale_factor",0.9996],
>     PARAMETER["false_easting",500000],
>     PARAMETER["false_northing",0],
>     UNIT["Meter",1]]
> Traceback (most recent call last):
>   File "<stdin>", line 1, in <module>
>   File "crawler_shp.py", line 85, in <module>
>     wksht.row(row).write(8, layer.GetSpatialRef())
>   File "C:\Python26\lib\site-packages\xlwt\Row.py", line 248, in write
>     raise Exception("Unexpected data type %r" % type(label))
> Exception: Unexpected data type <class 'osgeo.osr.SpatialReference'>
>>>>
> Any idea?
>
> _______________________________________________
> FWTools mailing list
> FWTools at lists.maptools.org
> http://lists.maptools.org/mailman/listinfo/fwtools
> http://fwtools.maptools.org/
>
>


More information about the FWTools mailing list