Randomly moving, and rotating, polygons
5435
13
02-29-2012 12:33 PM
by
by
by
by
Hey all,
Thanks @Darren Wiens for this really cool code. I was helping my friend figuring out this code. The biggest issue to me is that the random polygon would not fall in the study area polygon. So I added a try/except and if/else statement to check and return until every polygon fall in the study area. Again, I have a simply polygons for test so it might not working if you have huge amount of polygon, and these polygons should be singlepart and no hollow inside. Also, I remove the "disjoint" function for checking if the result polygon is overlap or not. I hope this help 
import arcpy
from arcpy import da
import random
import math
#set-up environment
workspace = r"...:/.../..." #location of your files
pa = workspace + '\\' + "polygons.shp" # polygons to randomly rotate and move
sr = arcpy.Describe(pa).spatialReference # spatial reference
sa = workspace + '\\' + "study_area.shp" # study area
extent = arcpy.Describe(sa).extent # study area extent
sa_geom = [row[0] for row in arcpy.da.SearchCursor(sa,'SHAPE@',spatial_reference=sr)][0] # study area geometry
new_polys = [] # placeholder
# magic function
def rotate(poly,centroid):
#ang = random.random() * 2 * math.pi # random rotation angle in radians if you want to have unified angle to result polygons
new_array = arcpy.Array() # placeholder
rnd_x = random.uniform(extent.XMin,extent.XMax) # random x coordinate for new centroid
rnd_y = random.uniform(extent.YMin,extent.YMax) # random y coordinate for new centroid
rnd_centroid = arcpy.Point(rnd_x,rnd_y) # random centroid point
while True:
try:
for part in poly:#for each polygon part
ang = random.random() * 2 * math.pi # random rotation angle in radians # set as 0 if you don't need to rotate polygon
for pnt in part: # for each vertex
x_trans = pnt.X - centroid.X # normalize to zero
y_trans = pnt.Y - centroid.Y # normalize to zero
x_transprime = (math.cos(ang) * x_trans) - (math.sin(ang) * y_trans) # new x coord, from zero
y_transprime = (math.sin(ang) * x_trans) + (math.cos(ang) * y_trans) # new y coord, from zero
x_prime = x_transprime + rnd_centroid.X # move to new centroid x
y_prime = y_transprime + rnd_centroid.Y # move to new centroid y
new_array.add(arcpy.Point(x_prime, y_prime)) # add to array
temp_poly = arcpy.Polygon(new_array,sr) # make a temporarily polygon for test
if sa_geom.contains(temp_poly) is True:
print 'inside'
break
else:
raise Exception
except:
print 'outside sa, start over...'
return rotate(poly,centroid) #here you start over and run the loop again... until all the polygon fall in study area.
break
good_poly = arcpy.Polygon(new_array,sr) # make all good polygons (indside study area) from array
return (good_poly) # return the final good polygon
print "start processing..."
with arcpy.da.SearchCursor(pa,'SHAPE@',spatial_reference=sr) as cursor: # for each polygon
for row in cursor:
centroid = row[0].centroid # calculate centroid
new_polys.append(rotate(row[0],centroid)) # add returned good polygon to list
arcpy.CopyFeatures_management(new_polys, workspace + '\\' + "TheResult.shp") # write to disk
print "done"
by
