Select to view content in your preferred language

Polygon Creation from text file (min xy, max xy) ISSUE

3106
2
02-20-2014 03:26 PM
RobertPotter
Deactivated User
We have been trying (for days) to create a polygon features from a CSV of min xy and max xy, held as GDA94 (lat,long) and it does not work for the majority of records.  Only one works.
I have added both the code and data.  We are creating the corresponding points to ensure the correct import and that works but when we Append to create the polygon it works but a feature of zero area is created.  There is no reported error in the script.

Python code

import arcpy, fileinput, os
from arcpy import env

env.overwriteOutput = True
infile = "C:/temp/M80/Metadata_RP.csv"

pt1 = arcpy.Point()
pt2 = arcpy.Point()
pt3 = arcpy.Point()
pt4 = arcpy.Point()
ptGeoms = []

point = arcpy.Point()
array = arcpy.Array()
featureList= []
curXY = arcpy.da.SearchCursor(infile, ["o_minx","o_miny","o_maxx","o_maxy"])

for row in curXY:
    message = "{0}, {1}, {2}, {3}".format(row[0], row[1], row[2], row[3])
    arcpy.AddMessage(message)
    o_minx= row[0]
    o_miny=row[1]
    o_maxx=row[2]
    o_maxy=row[3]
    message = "{0}, {1}, {2}, {3}".format(o_minx,o_miny,o_maxx,o_maxy)
   
    pt1.X = o_minx
    pt1.Y = o_miny
    ptGeoms.append(arcpy.PointGeometry(pt1))
    pt2.X = o_maxx
    pt2.Y = o_miny
    ptGeoms.append(arcpy.PointGeometry(pt2))
    pt3.X = o_maxx
    pt3.Y = o_maxy
    ptGeoms.append(arcpy.PointGeometry(pt3))
    pt4.X = o_minx
    pt4.Y = o_maxy
    ptGeoms.append(arcpy.PointGeometry(pt4))
    arcpy.AddMessage(message)
    array = arcpy.Array([pt1,pt2,pt3,pt4,pt1])

    polygon = arcpy.Polygon(array)
    featureList.append(polygon)
   
arcpy.CopyFeatures_management(featureList, r"C:/temp/M80/polys.shp")
arcpy.CopyFeatures_management(ptGeoms, r"C:/temp/M80/points.shp")

print arcpy.GetMessages()


CSV data
o_minx,o_miny,o_maxx,o_maxy
144.820900,-37.757700,144.827300,-37.755800
144.818700,-37.762200,144.823400,-37.759100
144.81900,-37.760800,144.824500,-37.758400
144.823800,-37.755800,144.829200,-37.753300
144.825500,-37.753400,144.832200,-37.75200
144.8209,-37.7577,144.8273,-37.7558
144.8187,-37.7622,144.8234,-37.7591
144.8268,-37.7524,144.8331,-37.7517
144.8292,-37.7506,144.8367,-37.7502
144.8203,-37.7621,144.8201,-37.761
144.8187,-37.7622,144.8234,-37.7591
144.8214,-37.7615,144.8224,-37.761
144.8192,-37.7622,144.8225,-37.7608
144.8209,-37.7577,144.8273,-37.7558
144.8187,-37.7622,144.8234,-37.7591
Tags (2)
0 Kudos
2 Replies
NeilAyres
MVP Alum
When creating features, it is important to specify the coordinate system.
The WKID of GDA_1994 is 4283.
So add this at the top of your script
SR = arcpy.SpatialReference(4283)

Then when you create the points and the polygon geometries, specify the SR as well.
ptGeoms.append(arcpy.PointGeometry(pt4, SR))
....

array = arcpy.Array([pt1,pt2,pt3,pt4,pt1])
polygon = arcpy.Polygon(array, SR)
featureList.append(polygon)

That fixed it for me.

Cheers,
Neil
0 Kudos
XanderBakker
Esri Esteemed Contributor
Hi Robert,

I tried the code and resulted with the same problems as you encountered; some polygons have a zero area and all that draw are overlapping exactly. I thought that this might be caused by a lacking coordinate system and even tried to project the polygon to a projected coordinate system. After this I notices that both seem to work. See code below:

import arcpy
arcpy.env.overwriteOutput = True

infile = r'C:\Project\_Forums\CSV2Polygon\csv\Metadata_RP.csv'
outshp_gcs = r'C:\Project\_Forums\CSV2Polygon\shp\poly_gcs.shp'
outshp_pcs = r'C:\Project\_Forums\CSV2Polygon\shp\poly_pcs.shp'

sr_in = arcpy.SpatialReference(4283) # GCS_GDA_1994
sr_out = arcpy.SpatialReference(28355) # GDA_1994_MGA_Zone_55

lstPolsGCS = []
lstPolsPCS = []
with arcpy.da.SearchCursor(infile, ["o_minx","o_miny","o_maxx","o_maxy"]) as curs:
    for row in curs:
        xvalues = [row[0], row[2]]
        yvalues = [row[1], row[3]]
        print "{0}, {1}, {2}, {3}".format(row[0], row[1], row[2], row[3])

        pnt1 = arcpy.Point(min(xvalues),min(yvalues))
        pnt2 = arcpy.Point(min(xvalues),max(yvalues))
        pnt3 = arcpy.Point(max(xvalues),max(yvalues))
        pnt4 = arcpy.Point(max(xvalues),min(yvalues))

        polygon_gcs = arcpy.Polygon(arcpy.Array([pnt1, pnt2, pnt3, pnt4, pnt1]), sr_in)
        lstPolsGCS.append(polygon_gcs)

        polygon_pcs = polygon_gcs.projectAs(sr_out)
        lstPolsPCS.append(polygon_pcs)

arcpy.CopyFeatures_management(lstPolsGCS, outshp_gcs)
arcpy.CopyFeatures_management(lstPolsPCS, outshp_pcs)



In my case I create two outputs; both polygon shapefiles, one for the geographic coordinates and one for the projected coordinates. Essential is providing the proper spatial reference when creating the geometry (see sr_in). I create a new polygon with the projected copy of the input polygon. Both are added to a list and written to their corresponding shapefiles.

BTW; the reason I first add the x and y coordinates of a row to lists, it to assure to get the minimum and maximum values and create polygons that are clockwise. There is one record (polygon # 10) where the xmin > xmax.

o_minx    o_miny    o_maxx    o_maxy
144,8203  -37,7621  144,8201  -37,761

Kind regards,

Xander
0 Kudos