David David - 2 months ago 69
Python Question

Writing to a Shapefile

I'm having trouble writing/reading a Shapefile in python. I have an array of points that I would like to write to a polygon using pyshp. The relevant parts of the code are:

dividedRects = [(7598325.0, 731579.0, 7698325.0, 631579.0), (7598325.0, 631579.0, 7698325.0, 611641.0), (7698325.0, 731579.0, 7728636.0, 631579.0), (7698325.0, 631579.0, 7728636.0, 611641.0)]

def createPolys(dividedRects):
w = shapefile.Writer(shapefile.POLYGON)
for i in range(0, len(dividedRects)):
print i
topLeft = [dividedRects[i][0],dividedRects[i][1]]
topRight = [dividedRects[i][2], dividedRects[i][1]]
bottomRight = [dividedRects[i][2], dividedRects[i][3]]
bottomLeft = [dividedRects[i][0], dividedRects[i][3]]
w.poly(parts=[[topLeft,topRight,bottomRight,bottomLeft]])
w.field("ID", "C", "40")
w.field("Events", "C", "40")
w.record(str(i), str(0))
w.save('cellFile')

createPolys(dividedRects)


This causes an error:

IndexError Traceback (most recent call last)
<ipython-input-36-503affbe838b> in <module>()
----> 1 createPolys(dividedRects)

<ipython-input-35-4c552ae29bc7> in createPolys(dividedRects)
10 w.field("ID", "C", "40")
11 w.field("Events", "C", "40")
---> 12 w.record(str(i), str(0))
13 w.save('cellFile')
14 # topLeft = [dividedRects[1][0],dividedRects[1][1]]

C:\Users\Me\Anaconda2\lib\site-packages\shapefile.pyc in record(self, *recordList, **recordDict)
967 if self.fields[0][0].startswith("Deletion"): fieldCount -= 1
968 if recordList:
--> 969 [record.append(recordList[i]) for i in range(fieldCount)]
970 elif recordDict:
971 for field in self.fields:

IndexError: tuple index out of range


If I remove the
field
and
records
lines from
createPolys
:

def createPolys(dividedRects):
w = shapefile.Writer(shapefile.POLYGON)
for i in range(0, len(dividedRects)):
print i
topLeft = [dividedRects[i][0],dividedRects[i][1]]
topRight = [dividedRects[i][2], dividedRects[i][1]]
bottomRight = [dividedRects[i][2], dividedRects[i][3]]
bottomLeft = [dividedRects[i][0], dividedRects[i][3]]
w.poly(parts=[[topLeft,topRight,bottomRight,bottomLeft]])
# w.field("ID", "C", "40")
# w.field("Events", "C", "40")
# w.record(str(i), str(0))
w.save('cellFile')


Then I get an assertionerror when reading the records from the file:

createPolys(dividedRects)

sf2 = shapefile.Reader("cellFile")
print sf2.records()
shapes = sf2.shapes()
bbox = shapes[1].bbox
#['%.3f' % coord for coord in bbox]
print bbox
points = shapes[1].points
print points

AssertionError Traceback (most recent call last)
<ipython-input-37-597af0b882ba> in <module>()
1 sf2 = shapefile.Reader("cellFile")
----> 2 print sf2.records()
3 shapes = sf2.shapes()
4 bbox = shapes[1].bbox
5 #['%.3f' % coord for coord in bbox]

C:\Users\Me\Anaconda2\lib\site-packages\shapefile.pyc in records(self)
528 """Returns all records in a dbf file."""
529 if not self.numRecords:
--> 530 self.__dbfHeader()
531 records = []
532 f = self.__getFileObj(self.dbf)

C:\Users\Me\Anaconda2\lib\site-packages\shapefile.pyc in __dbfHeader(self)
464 self.fields.append(fieldDesc)
465 terminator = dbf.read(1)
--> 466 assert terminator == b("\r")
467 self.fields.insert(0, ('DeletionFlag', 'C', 1, 0))
468

AssertionError:


When I removed the loop and wrote one record it seemed to work ok. What's going on?

Answer

I do not know about pyshp library, but I'll try to help anyways.

The two w.field() commands occur in the for loop. This may cause the two columns "ID" and "Events" to be defined multiple times. When you write only one record (polygon), it works fine (i.e. the w.record() command contains two values). After the first iteration, there would be 4, 6, etc columns. That would explain the behavior you describe.

Try moving the two w.field() lines before the for loop.

When you comment the w.record(), you obtain a shp (and shx) file with a number of records different from the corresponding dbffile. That explains the assertion error when reading.

Unrelated to your problem, you could also simplify your code with enumerate (built-in function).

w = shapefile.Writer(shapefile.POLYGON)
w.field("ID", "C", "40")
w.field("Events", "C", "40")    
for i,rect1 in enumerate(dividedRects):
    print i
    topLeft = [rect1[0],rect1[1]]
    topRight = [rect1[2], rect1[1]]
    bottomRight = [rect1[2], rect1[3]]
    bottomLeft = [rect1[0], rect1[3]]
    ....

(I cannot test since I do not have pyshp) Good luck!