Random polygons

Question

Hello guys.

I need to create 20X30m polygons at random positions within a predefined boundary. And I also need the number of polygons created equal 25% of the area.

Any suggestion?

duvida.jpg

DanPatterson_Retired · Answer

Nice!  Now if you want Xander Bakker‌ try option 4 (spaced points ) in my Point Tools for ArcGIS Pro if you are looking to give 2.1 a test ride.  Point spacing should be guaranteed, then buffer by 1/2 the distance should also ensure no overlaps

XanderBakker · Answer

I thought I would have a go at it, but there are a number of things that would need to be resolved. * rectangles should not be reused * rectangles have common origin XY * minimum area of 25% should be covered This is the code I came up with, but it is far from fast and should not be used for a "large" set of data: import arcpy def main(): arcpy.env.overwriteOutput = True # settings fc_in = r'C:\GeoNet\RandRectangles\data.gdb\areas' fc_tmp = r'C:\GeoNet\RandRectangles\data.gdb ectangles_tmp' fc_sel = r'C:\GeoNet\RandRectangles\data.gdb ectangles_tmp_sel' fc_out = r'C:\GeoNet\RandRectangles\data.gdb ectangles_v01' rect_width = 20 rect_height = 30 min_area_frac = 0.25 # Some asumptions: # rectangles should not be reused # rectangles have common origin XY # minimum area of 25% should be covered sr = arcpy.Describe(fc_in).spatialReference ext = arcpy.Describe(fc_in).extent # create a tmp featureclass with the rectangles ext_fin = ReshapeExtent(ext, rect_width, rect_height) feats = [] for c in range(int(ext_fin.width / rect_width)): for r in range(int(ext_fin.height / rect_height)): polygon = CreateRectanglePolygon(ext_fin.XMin + c * rect_width, ext_fin.YMin + r * rect_height, ext_fin.XMin + (c+1) * rect_width, ext_fin.YMin + (r+1) * rect_height, sr) feats.append(polygon) arcpy.CopyFeatures_management(feats, fc_tmp) # now select only the ones that overlap any area polygon arcpy.MakeFeatureLayer_management(fc_in, 'area_lyr') arcpy.MakeFeatureLayer_management(fc_tmp, 'rect_lyr') arcpy.SelectLayerByLocation_management("rect_lyr", "INTERSECT", "area_lyr", None, "NEW_SELECTION") arcpy.CopyFeatures_management('rect_lyr', fc_sel) # create a dictionary from the selected reactangles dct_rect = {r[0]: r[1] for r in arcpy.da.SearchCursor(fc_sel, ('OID@', 'SHAPE@'))} # loop through polygons and determine which to use lst_oid_taken = [] lst_results = [] with arcpy.da.SearchCursor(fc_in, ('OID@', 'SHAPE@')) as curs: for row in curs: oid_area = row[0] area_pol = row[1] dct_overlap = GetDictOfRectanglesForArea(area_pol, dct_rect, lst_oid_taken) lst_selection, lst_oid_taken = GetRandomRectanglesForArea(area_pol, dct_overlap, lst_oid_taken, min_area_frac) lst_results.append([oid_area, lst_selection]) for result in lst_results: oid_area = result[0] lst_selection = result[1] for oid_rect in lst_selection: print "{0} {1}".format(oid_rect, oid_area) def GetRandomRectanglesForArea(area, dct_overlap, lst_oid_taken, min_area_frac): import random lst_selection = [] lst_oid = dct_overlap.keys() polygon_area = area.area tot_rect_area = 0 fraction = tot_rect_area / polygon_area cnt = 0 print " GetRandomRectanglesForArea" while fraction < min_area_frac: cnt += 1 print " - cnt:", cnt i = random.randrange(0, len(lst_oid) - 1) print " - i:", i oid = lst_oid[i] print" - oid:", oid if not oid in lst_oid_taken: rectangle = dct_overlap[oid] overlap_area = GetOverlapArea(area, rectangle) print " - overlap_area:", overlap_area tot_rect_area += overlap_area print " - tot_rect_area:", tot_rect_area fraction = tot_rect_area / polygon_area print " - fraction:", fraction lst_oid_taken.append(oid) lst_selection.append(oid) if cnt > 100: break return lst_selection, lst_oid_taken def GetOverlapArea(area, rectangle): if rectangle.within(area): return rectangle.area else: overlap = area.intersect(rectangle, 4) return overlap.area def GetDictOfRectanglesForArea(area, dct, lst_done): dct_overlap = {} for oid, rectangle in dct.items(): if not oid in lst_done: if rectangle.within(area): dct_overlap[oid] = rectangle elif rectangle.overlaps(area): dct_overlap[oid] = rectangle return dct_overlap def CreateRectanglePolygon(xmin, ymin, xmax, ymax, sr): points = [arcpy.Point(xmin, ymin), arcpy.Point(xmin, ymax), arcpy.Point(xmax, ymax), arcpy.Point(xmax, ymin), arcpy.Point(xmin, ymin)] return arcpy.Polygon(arcpy.Array(points), sr) def ReshapeExtent(ext, width, height): xmin_new = divmod(ext.XMin, width)[0] * width ymin_new = divmod(ext.YMin, height)[0] * height xmax_new = divmod(ext.XMax, width)[0] * width + width ymax_new = divmod(ext.YMax, height)[0] * height + height return arcpy.Extent(xmin_new, ymin_new, xmax_new, ymax_new) if __name__ == '__main__': main() ‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍ The result is list of rectangle OIDs with the corresponding area OIDs which I joined back to the rectangles (since I was too lazy to create a proper output featureclass. The steps are more or less the following: 1) a tmp featureclass with rectangles is created using lower left XY which is a multiple of the width and height and covers the entire area: 2) The rectangle features that overlap the areas are selected and read into a dictionary: 3) random rectangles are selected until at least 25% overlap is obtained (taking only the overlapping area of the rectangles into account): There are still some situations that are probably not desired. Such as rectangles that have a very small overlap or have a larger overlap with another polygon and that might not be what you want: So it will probably need to be optimized further.

DanPatterson_Retired · Answer

That might entail 'humanware' to provide some guidance.  You could however, determine the extent of the base polygon then use that extent to set the bounds of the fishnet.  You would just want to ensure that your fishnet has sufficient rows and columns that you could get a 25% sample of rectangles within the base polygon.  That is not a problem since 'select by spatial location' with a 'completely within' control would ensure that you select only polygons within the base polygon, then continue to select polygons perhaps incrementally, until the 25% areal coverage is met.