Search code examples
pythonjsongoogle-mapsgeojsonarcpy

Can polygonzo for google maps handle interior rings?


I've put together an python script with ArcGIS arcpy for creating polygonzo json polygons (http://code.google.com/p/polygonzo/). Here is my python script...

 import os, string, arcpy
 arcpy.env.overwriteOutput = True


 layer = "C:\\Other\\Shapefiles\\Geo500K_JSON\\GEOLOGY_500K_Project.shp"

 output = "C:\\Other\\Shapefiles\\Geo500K_JSON\\"

 outfile = output + "Geo500K.json"
 jsonFile = open(outfile,'w')
 jsonFile.write('var geo = {\n')
 jsonFile.write('\t"type": "FeatureCollection",\n')
 jsonFile.write('\t"features": [\n')

 idfield = "ORIG_LABEL"
 shape_field = arcpy.Describe(layer).shapeFieldName

 rows = arcpy.SearchCursor(layer,"","","",idfield + " A")
 row = rows.next()

 while row:

     geostring = '' #for each lat/lng pt
     geolist = [] # array for storing individual geostrings
     ringList = [] #array for storing geolist array with geostrings separated by commas
     partList = [] #array for storing partlist, final array used
     shapeString = ''

     jsonFile.write('\t\t{"type": "Feature", ')

     extent = row.Shape.extent
     ne = str(extent.XMax) + ',' + str(extent.YMax)
     sw = str(extent.XMin) + ',' + str(extent.YMin)
     jsonFile.write('"bbox": [' + sw + ', ' + ne + '],')

     jsonFile.write('"properties":{')

     geoLabel = str(row.getValue(idfield))
     jsonFile.write('"label": "' + geoLabel + '", ')    

     geoName = str(row.getValue("FM_NAME"))
     jsonFile.write('"name": "' + geoName + '", ')

     lithType = str(row.getValue("LithType"))
     jsonFile.write('"lithType": "' + lithType + '", ')

     rank = str(row.getValue("Rank"))
     jsonFile.write('"rank": "' + rank + '", ')

     lithName = str(row.getValue("LithName"))
     jsonFile.write('"lithName": "' + lithName + '", ')

     ageType = str(row.getValue("AgeType"))
     jsonFile.write('"ageType": "' + ageType + '", ')

     minAge = str(row.getValue("MinAge"))
     jsonFile.write('"minAge": "' + minAge + '", ')

     maxAge = str(row.getValue("MaxAge"))
     jsonFile.write('"maxAge": "' + maxAge + '", ')

     part = row.getValue(shape_field).centroid
     jsonFile.write('"center":[' + str(part.X) + ',' + str(part.Y) + '],')

     jsonFile.write('"centroid":[' + str(part.X) + ',' + str(part.Y) + ']},')


     jsonFile.write('"geometry":{"type":"MultiPolygon","coordinates":[[[')

     feat = row.shape
     for p in range(feat.partCount):
         pInt = p
         part = feat.getPart(p)
         pt = part.next()
         while pt:
             lat = str(round(pt.Y,6))
             lon = str(round(pt.X,6))

             geostring = '[' + lon + ',' + lat + ']'
             geolist.append(geostring)

             pt = part.next()

             #if now following point go to the next part which should be an interior ring.     
             if not pt:
                 ringList.append(',' .join(geolist))
                 geostring = ''
                 geolist = []
                 pt = part.next()
                 if pt:
                     print 'Interior Ring: ' + geoLabel

         partList.append(',' .join(ringList))
         ringList = []

     shapeString = ']], [[' .join(partList)
     jsonFile.write(shapeString)
     jsonFile.write(']]]}},\n')
     row = rows.next()

 #jsonFile.seek(-1, os.SEEK_END)
 #jsonFile.truncate()
 jsonFile.write('\t]\n')
 jsonFile.write('}')
 jsonFile.close()
 del row, rows

When the script encounters interior rings it only prints a warning. I don't know how to handle them. Unfortunately many of the polygons I work with have interior rings. I put together a test map using one polygon that has interior rings. Here is what it looks like... http://www.geology.ar.gov/test/test-polygonzo.html

Can polygonzo handle interior rings?

UPDATE: I really appreciate your response Mr. Michael Geary! However, I could not get your python script using the json module to work. There were a few bugs in it and I edited it above, but it spits out a blank document. Maybe I didn't try hard enough. After reviewing your example of what a multipolyon with interior rings should look like in json format, I went back to working on my python script (and yes it was slightly difficult getting the json valid withouth using the json module). I've added more comments, so maybe, if you have time, you could get your script working using the json module - I would like to see a working example. Here is my final python script....

import os, string, arcpy
arcpy.env.overwriteOutput = True


layer = "C:\\Other\\Shapefiles\\Geo500K_JSON\\GEOLOGY_500K_kn.shp"

output = "C:\\Other\\Shapefiles\\Geo500K_JSON\\"

outfile = output + "Geo500K_knTest.json"
jsonFile = open(outfile,'w')
jsonFile.write('var geo = {\n')
jsonFile.write('\t"type": "FeatureCollection",\n')
jsonFile.write('\t"features": [\n')

idfield = "ORIG_LABEL"
shape_field = arcpy.Describe(layer).shapeFieldName

rows = arcpy.SearchCursor(layer,"","","",idfield + " A")
row = rows.next()
#loop through the attribute table
while row:    

    jsonFile.write('\t\t{"type": "Feature", \n')

    extent = row.Shape.extent
    ne = str(extent.XMax) + ',' + str(extent.YMax)
    sw = str(extent.XMin) + ',' + str(extent.YMin)
    jsonFile.write('\t\t"bbox": [' + sw + ', ' + ne + '],\n')

    jsonFile.write('\t\t"properties":{\n')

    geoLabel = str(row.getValue(idfield))
    jsonFile.write('\t\t\t"label": "' + geoLabel + '", \n')    

    geoName = str(row.getValue("FM_NAME"))
    jsonFile.write('\t\t\t"name": "' + geoName + '", \n')

    lithType = str(row.getValue("LithType"))
    jsonFile.write('\t\t\t"lithType": "' + lithType + '", \n')

    rank = str(row.getValue("Rank"))
    jsonFile.write('\t\t\t"rank": "' + rank + '", \n')

    lithName = str(row.getValue("LithName"))
    jsonFile.write('\t\t\t"lithName": "' + lithName + '", \n')

    ageType = str(row.getValue("AgeType"))
    jsonFile.write('\t\t\t"ageType": "' + ageType + '", \n')

    minAge = str(row.getValue("MinAge"))
    jsonFile.write('\t\t\t"minAge": "' + minAge + '", \n')

    maxAge = str(row.getValue("MaxAge"))
    jsonFile.write('\t\t\t"maxAge": "' + maxAge + '", \n')

    centroid = row.getValue(shape_field).centroid
    jsonFile.write('\t\t\t"center":[' + str(centroid.X) + ',' + str(centroid.Y) + '], \n')
    jsonFile.write('\t\t\t"centroid":[' + str(centroid.X) + ',' + str(centroid.Y) + '] \n')

    jsonFile.write('\t\t\t}, \n') #end of properties

    jsonFile.write('\t\t"geometry":{\n\t\t\t"type":"MultiPolygon",\n\t\t\t"coordinates":[\n')

    feat = row.shape #get the shape/geography of the row in the attribute table
    partnum = 1

    #loop through the parts of the polygon (some may have more that one part)
    for p in range(feat.partCount):
        jsonFile.write('\t\t\t\t[\n\t\t\t\t\t[\n')
        jsonFile.write('\t\t\t\t\t\t//Part ' + str(partnum) + '\n')
        jsonFile.write('\t\t\t\t\t\t//Outer ring of Part ' + str(partnum) + '\n')

        part = feat.getPart(p) #return an array of point objects for particular part

        pt = part.next() #return specific pt object of array
        innerRingNum = 1

        #loop through each pt object/vertex of part
        while pt:
            lat = round(pt.Y,7) #get latitude of pt object and round to 7 decimal places
            lon = round(pt.X,7) #get longitude of pt object and round to 7 decimal places

            jsonFile.write('\t\t\t\t\t\t[' + str(lon) + ',' + str(lat) + '],\n') #assemble [lon,lat]

            pt = part.next() #go to next pt object to continue loop

            #if no following point go to the next part which should be an interior ring.
            if not pt:
                #we've got an interior ring so let's loop through the vertices of the ring
                pt = part.next()

                if pt:
                    jsonFile.seek(-3, os.SEEK_END)
                    jsonFile.truncate() #remove trailing comma
                    jsonFile.write('\n\t\t\t\t\t],\n')
                    jsonFile.write('\t\t\t\t\t[\n')
                    jsonFile.write('\t\t\t\t\t\t//Inner ring ' + str(innerRingNum) + ' of Part ' + str(partnum) + '\n')
                    print 'Interior Ring: ' + geoLabel
                    innerRingNum += 1



        partnum += 1
        jsonFile.seek(-3, os.SEEK_END)
        jsonFile.truncate() #remove trailing comma
        jsonFile.write('\n\t\t\t\t\t]\n\t\t\t\t],\n')

    jsonFile.seek(-3, os.SEEK_END)
    jsonFile.truncate() #remove trailing comma
    jsonFile.write('\n\t\t\t]\n\t\t\t}\n\t\t},\n')
    row = rows.next()

jsonFile.seek(-3, os.SEEK_END)
jsonFile.truncate() #remove trailing comma
jsonFile.write('\n\t]\n')
jsonFile.write('}')
jsonFile.close()
del row, rows

Let me also add, that I'm really impressed with polygonzo as well as your willingness to share it with others. However, the javascript and python that you provide really could use more comments for quicker understanding of it all.


Solution

  • PolyGonzo author here, sorry I didn't run across your question until just now.

    I don't know if this will still be relevant, but I looked at your test page.

    PolyGonzo does support interior rings, but there are no interior rings in your GeoJSON data.

    There's an example of an interior ring in the MultiPolygon example in the GeoJSON spec. Unfortunately it's formatted poorly, so here's an indented and commented version:

    {
        "type": "MultiPolygon",
        "coordinates": [
            // First polygon of the multipolygon
            [
                // Outer ring of the first polygon (there is no inner ring)
                [
                    [ 102.0, 2.0 ],
                    [ 103.0, 2.0 ],
                    [ 103.0, 3.0 ],
                    [ 102.0, 3.0 ],
                    [ 102.0, 2.0 ]
                ]
            ],
            // Second polygon of the multipolygon
            [
                // Outer ring of the second polygon
                [
                    [ 100.0, 0.0 ],
                    [ 101.0, 0.0 ],
                    [ 101.0, 1.0 ],
                    [ 100.0, 1.0 ],
                    [ 100.0, 0.0 ]
                ],
                // Inner ring of the second polygon
                [
                    [ 100.2, 0.2 ],
                    [ 100.8, 0.2 ],
                    [ 100.8, 0.8 ],
                    [ 100.2, 0.8 ],
                    [ 100.2, 0.2 ]
                ]
                // You could have additional inner rings here
            ]
        ]
    }
    

    Put another way, the coordinates property for a MultiPolygon is an array of polygons. Each polygon is in turn an array of rings. Each of those rings is an array of coordinate pairs (or triplets if you have height information, etc.).

    For a given polygon, the first ring is the outer ring, and any additional rings are inner rings.

    In your MultiPolygon, each of its polygons has only a single ring, so PolyGonzo interprets that as the outer ring.

    After looking at your data some more, I can see that this is what has happened: for each of the polygons within the MultiPolygon, you have all of the points for both the outer ring and all of the points in any inner rings all in one big array.

    The third polygon in your file is a good example. This is the largish area north of Hope. Looking through this part of the GeoJSON data, I found four inner rings. (I found them by zooming way in on the map to places that had spurious lines, and also by looking through the coordinates for large jumps.)

    I manually split the array in the GeoJSON so that these inner rings have their own arrays, and it fixes things up real nice.

    Here is a fiddle with the corrected data for the area north of Hope. The GeoJSON data is included inline in the JavaScript code for the fiddle, so you can see there what I changed. I also see a similar problem in the area around Arkadelphia but didn't correct it.

    Now about the way you're generating your JSON data...

    I strongly recommend that you not generate JSON by pasting together a bunch of bits and pieces of JSON text as the code is doing now.

    Instead, create an object (dict) in Python that represents your entire GeoJSON structure, and then call json.dump() or json.dumps() to convert the entire structure into JSON all at once.

    This makes things much easier. It will automatically guarantee that your JSON is valid - which it already is, but I'll bet you had to work at it to get there, right? ;-) json.dump() makes that trivial.

    It should also make it easier to avoid issues like this case where you meant to put out individual arrays for the outer ring and inner rings, but they accidentally got all jammed into a single array.

    Here's your code partly converted to use this technique. I didn't do the whole thing since I'm not familiar with arcpy, but this should give you the idea:

    import os, string, arcpy, json
    arcpy.env.overwriteOutput = True
    
    
    layer = "C:\\Other\\Shapefiles\\Geo500K_JSON\\GEOLOGY_500K_kn.shp"
    
    output = "C:\\Other\\Shapefiles\\Geo500K_JSON\\"
    
    outfile = output + "Geo500K_knTest_C2.json"
    
    idfield = "ORIG_LABEL"
    shape_field = arcpy.Describe(layer).shapeFieldName
    
    features = []
    geojson = {
        'type': 'FeatureCollection',
        'features': features
        }
    
    
    rows = arcpy.SearchCursor(layer,"","","",idfield + " A")
    for row in rows:
        geostring = '' #for each lat/lng pt
        geolist = [] # array for storing individual geostrings
        ringList = [] #array for storing geolist array with geostrings separated by commas
        partList = [] #array for storing partlist, final array used
        shapeString = ''
    
        extent = row.Shape.extent
        centroid = row.getValue(shape_field).centroid
    
        coordinates = []
    
        feature = {
            'type': 'Feature',
            'bbox': [ extent.XMin, extent.YMin, extent.XMax, extent.YMax ],
            'properties': {
                'label': str(row.getValue(idfield)),
                'name': str(row.getValue("FM_NAME")),
                'lithType': str(row.getValue("LithType")),
                'rank': str(row.getValue("Rank")),
                'lithName': str(row.getValue("LithName")),
                'ageType': str(row.getValue("AgeType")),
                'minAge': str(row.getValue("MinAge")),
                'maxAge': str(row.getValue("MaxAge")),
                'center': [ centroid.X, centroid.Y ],
                'centroid': [ centroid.X, centroid.Y ]
            },
            'geometry': {
                'type': 'MultiPolygon',
                'coordinates': coordinates
            }
        }
    
        feat = row.shape
        for p in range(feat.partCount):        
            part = feat.getPart(p)
            pt = part.next()
            while pt:
                lat = str(round(pt.Y,6))
                lon = str(round(pt.X,6))
    
                geostring = '[' + lon + ',' + lat + ']'
                geolist.append(geostring)
    
                pt = part.next()
    
                #if no following point go to the next part which should be an interior ring.     
                if not pt:
                    ringList.append(',' .join(geolist))
                    geostring = ''
                    geolist = []
                    pt = part.next()
                    if pt:
                        print 'Interior Ring: ' 
    
            partList.append('],[' .join(ringList))
            ringList = []        
    
    
        features.append( feature )
        shapeString = ']],[[' .join(partList)
        coordinates.append(shapeString)
    
    
    
    with open(outfile,'wb') as jsonFile:
        json.dump( geojson, jsonFile )