[Gdal-dev] Append two shapefiles

Alex Mandel tech_dev at wildintellect.com
Mon Dec 5 20:49:30 EST 2005


Thanks Dylan -
I was putting full paths for the layers and using the -f, and for this 
time I also switched my slashes \ to / and it works great now.

The next questions is how can I get just the ogr2ogr component into my 
program. I need to distribute a very thin and lite python program and 
can't expect my users to have FWtools.

Ironically I think I'm one error away from writing a python module with 
gdal to do what I want, which would save me a ton of MB and give me 
control to copy only certain points.

There's something wrong with the geometry line though. boa says that 
feat2.SetGeometry(pnt) has invalid syntax
Anyone see my mistake right away, I figure it's something obvious that I 
just haven't learned yet. Keep in mind I'm only dealing with points, and 
will never need to handle any other geometry type.

Thanks Everyone - Alex
---------------------------------------

import gdal
import gdal.ogr as ogr

#Set ESRI Shapefile driver just in case
driver = ogr.GetDriverByName('ESRI Shapefile')
#Bring in 1st Shapefile, src
ds1 = driver.Open("C:/Boxes/sandbox/one/one.shp")
pLayer1 = ds1.GetLayer(0)
#initialize first feature read
feat = pLayer1.GetNextFeature()

#Bring in 2nd Shapefile, dest
ds2 = driver.Open("C:/Boxes/sandbox/two/two.shp")
pLayer2 = ds2.GetLayer(0)

#Count the numnber of fields
layerDef2 = pLayer2.GetLayerDefn()
field_count = layerDef.GetFieldCount()

#Create index of fields to copy
fld_list = {}

for i in range(layerDef2.GetFieldCount()):
     field = layerDef.GetFieldDefn(i)
     fld_list[field.GetName()] = i

while feat is not None:

     multi_geom = feat.GetGeometryRef()

     for i in range(multi_geom.GetGeometryCount()):
         pnt = multi_geom.GetGeometryRef(i)

         feat2 = ogr.Feature(feature_def=layerDef2)
#Copy Feature table data
         for fld_index in range(field_count):
             fld_name = layerDef.GetFieldDefn(fld_index)
             #if fld_list.has_key(fld_name.GetName()):
             feat2.SetField(fld_index, \
feat.GetField(fld_list.get(fld_name))

#Paste into destination
         feat2.SetGeometry(pnt)
         pLayer2.CreateFeature(feat2)
         feat2.Destroy()

#Clean up and iterate 	
	feat.Destroy()
	feat = pLayer1.GetNextFeature()


----------------------------------------

Dylan Beaudette wrote:
 > hi Alex!
 >
 > didn't know you were lurking here!
 >
 > welcome the GDAL land --
 >
 > lots of good stuff here.
 >
 >
 > here is how i do it:
 >
 > ogr2ogr -update -append $merge_file $this_file_to_append -nln 
$merge_layer
 > $this_file_to_append_layer
 >
 > where the 2nd set of 'layers' are missing the '.shp' suffix.
 >
 > cheers,
 >
 > Dylan
 >





More information about the Gdal-dev mailing list