7

I need to create polygons of a specific area based upon random points within bounding polygons, but haven't been able to find a tool able to do all these.

I've looked at growing Theissen polygons, but can't find a way to stop the growth to a specific area.

I've looked at buffering the points, but that doesn't provide a way to grow the polygon to the correct area within the bounding polygons.

I've attached a hand-drawn concept of what I need.

Point Buffer Mock-up

  • I Think you need Python for this. Iterate over each Point and buffer, clip with boundind polygon, check area, rebuffer with smaller/larger width, clip, check area etc. – BERA Mar 20 '18 at 14:09
  • 1
    Could you not do this in steps - run the buffer as a set distance, then run a clip using your bounding polygons and use the output? – MAJ742 Mar 20 '18 at 14:18
  • Ok, re-read the question - you need growth to a specific area, not just clipped to the bounding features so my comment above won't apply. Tricky one! – MAJ742 Mar 20 '18 at 14:55
  • Are all those blue shapes you've drawn supposed to be round (apart from where they run up against the green ones)? Some of them are a bit lumpy. If you're trying to create round buffers (blue) that are contained by the green polygons, BERA's and MAJ742's idea will work: buffer the point, clip the result by the green polygon, repeat as necessary to grow your polygons. – Dan C Mar 20 '18 at 16:22
3

I use Python and Buffer. Each Point gets buffered, clipped and area is checked. If the difference between clipped area and desired size is more than ok_diff, it will get buffered again with a slightly larger buffer distance and checked again.. When ok it gets appended to output fc. The code can be executed in the Python Window.

import arcpy, math
arcpy.env.overwriteOutput = True

#Change to match your data:
arcpy.env.workspace = r'C:\Default.gdb'
boundries = r'Boundries'
points = r'bufferpoints'
out_feature_class = r'GrowPoints'

#Change to match your desired size of polygon, ok difference and step in buffer radius increase
size = 7500
ok_diff = 50
increment = 1


arcpy.CreateFeatureclass_management(out_path=arcpy.env.workspace, out_name=out_feature_class, 
                                   geometry_type='POLYGON', 
                                   spatial_reference=arcpy.Describe(points).spatialReference)

arcpy.MakeFeatureLayer_management(in_features=boundries, out_layer='blyr')

with arcpy.da.SearchCursor(points,['OID@','SHAPE@']) as cursor:
    for row in cursor:
        bufferstart = math.sqrt((size/math.pi))
        sql = """{0} = {1}""".format(
            arcpy.AddFieldDelimiters(points,arcpy.Describe(points).OIDFieldName),row[0])
        arcpy.MakeFeatureLayer_management(in_features=points, out_layer='pointlyr',where_clause=sql)
        arcpy.SelectLayerByLocation_management(in_layer='blyr', overlap_type='INTERSECT', 
                                              select_features='pointlyr')
        if [i[0] for i in arcpy.da.SearchCursor('blyr','SHAPE@AREA')][0] > size:
            area = 1
            while abs(size-area)>ok_diff:
                print area
                arcpy.Buffer_analysis(in_features='pointlyr', out_feature_class=r'in_memory\point', 
                                     buffer_distance_or_field="{0} Meters".format(bufferstart))
                arcpy.Clip_analysis(in_features=r'in_memory\point', clip_features='blyr',out_feature_class=r'in_memory\clipbuffer')
                area = [i[0] for i in arcpy.da.SearchCursor(r'in_memory\clipbuffer','SHAPE@AREA')][0]
                bufferstart+=increment

            arcpy.Append_management(inputs=r'in_memory\clipbuffer', target=out_feature_class, schema_type='NO_TEST')
        else:
            print 'Impossible to fit buffer inside boundry for point number: ',row[0]

enter image description here

  • 1
    This looks like it will provide exactly what I need. Thanks for sharing the code! – Andrew Strassman Mar 20 '18 at 18:56
  • Set increment as small as possible (only downside is time taken) or it will possibly never leave the while loop if buffererd and clipped size gets as larger than size+ok_diff (as the code is currently written). – BERA Mar 21 '18 at 6:41
1

Rather than using a buffer I think you need this tool:

http://desktop.arcgis.com/en/arcmap/10.3/tools/production-mapping-toolbox/increase-polygon-area.htm

Crucially it specifies that...

You can optionally limit the features to be enlarged based on whether they intersect another feature class

I guess if you buffer by a really small amount to give a round polygon feature then use this tool it should do the trick.

Unfortunately it appears to only be available with the 'Production Mapping' extension/licence.

0

Script below is using bisection method to find solution:

import arcpy, math
from arcpy import env
env.overwriteOutput = True
pi=math.pi
points,polygons,buffers="POINTS","PGONS","BUFFERS"
tolerance=0.01
d=arcpy.Describe(polygons)
fld=d.OIDFieldName
outFC="in_memory/outFC"
# GET POLYGONS SHAPES
shapesDict={}
with arcpy.da.SearchCursor(polygons,(fld,"Shape@")) as cursor:
    for fid,shp in cursor:
        shapesDict[fid]=shp
curT = arcpy.da.InsertCursor(buffers,"Shape@")
# ITERATE THROUGH POINTS
with arcpy.da.SearchCursor(points,("Shape@","NEAR_FID","AREA",)) as cursor:
    for shp, pgonID, target in cursor:
        parent=shapesDict[pgonID]
        parentArea=parent.area
        if target>parentArea:
            insidePart=parent
        else:
            L = pow(parentArea/pi,0.5)/5
            H = parent.length
            while True:
                if (H-L)<tolerance:break
                M = (L+H)/2
                theBuffer=shp.buffer(M)
                cutter=theBuffer.boundary()
                insidePart=parent.cut(cutter)[1]
                if insidePart.area<target:
                    L=M
                else: H=M
        curT.insertRow((insidePart,))

Just name your layers in TOC of mxd accordingly. You'll need also to run near tool on points, so your points table looks like that:

enter image description here

Output: enter image description here

Unfortunately arcpy geometry method intersect doesn't do what it should so I was forsed to use cut method to find intersect between buffer and parent polygon.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.