import arcpy import os ################################################################################################### ### ### This tool moves the vertices on the 500-year polygons that are less than 2 feet from a 100-year ### polygon. It also removes polygons that are too small (based on input size). ### It creates a new SFHA layer with those vertices moved along with a file geodatabase that ### includes intermediary datasets that may be helpful for QA/QC purposes. ### ### Written by: Marianne Cardwell ### ################################################################################################### class Tool(object): def __init__(self): """Define the tool (tool name is the name of the class).""" self.label = 'Fix 500-year polygon slivers' self.description = 'Fix 500-year polygon slivers' self.canRunInBackground = False def getParameterInfo(self): """Define parameter definitions""" param0 = arcpy.Parameter( displayName="SFHA Feature Class", name="sfha_fc", datatype="DEFeatureClass", parameterType="Required", direction="Input" ) #param0.value = r'C:\Temp\self_intersecting_polys\Test_Data.gdb\SFHA_GCS\S_Fld_Haz_Ar_NewData' param1 = arcpy.Parameter( displayName="Working Geodatabase (created by tool)", name="working_gdb", datatype="DEWorkspace", parameterType="Required", direction="Output" ) #param1.value = r'C:\Temp\self_intersecting_polys\working' param2=arcpy.Parameter( displayName="Output Feature Class (created by tool)", name="output_fc", datatype="DEFeatureClass", parameterType="Required", direction="Output" ) #param2.value = r'C:\Temp\self_intersecting_polys\Test_Data.gdb\SFHA_GCS\S_Fld_Haz_Ar_NewData_New' param3=arcpy.Parameter( displayName="Minimum Polygon Size (Square Feet)", name="min_size", datatype="Long", parameterType="Required", direction="Input" ) param3.value = 1000 params = [param0,param1,param2,param3] return params def isLicensed(self): """Set whether tool is licensed to execute.""" return True def updateParameters(self, parameters): """Modify the values and properties of parameters before internal validation is performed. This method is called whenever a parameter has been changed.""" return def updateMessages(self, parameters): """Modify the messages created by internal validation for each tool parameter. This method is called after internal validation.""" return def execute(self, parameters, messages): # Settings arcpy.env.overwriteOutput = "True" spe = arcpy.SpatialReference(2965) # Get input parameter sfha_path = parameters[0].valueAsText working_gdb_path = parameters[1].valueAsText if working_gdb_path[-4:].upper() != '.GDB': working_gdb_path += '.gdb' working_gdb = arcpy.management.CreateFileGDB(os.path.dirname(working_gdb_path), os.path.basename(working_gdb_path)) #working_gdb_path = os.path.join(working_folder, 'working.gdb') output_path = parameters[2].valueAsText arcpy.AddMessage(f'Working geodatabase path: {working_gdb_path}') delete_size = parameters[3].value # Create feature class to track vertex changes sfha_sp = self.get_fc_spref(sfha_path) arcpy.management.CreateFeatureclass(working_gdb_path, 'vertex_changes', 'POLYLINE', has_m='DISABLED', has_z='DISABLED', spatial_reference=sfha_sp) vertex_changes_path = os.path.join(working_gdb_path, 'vertex_changes') vertex_cursor = arcpy.da.InsertCursor(vertex_changes_path, 'SHAPE@') # Go from multipart to single part geometry sfha_singlepart = os.path.join(working_gdb_path, 'single_part') arcpy.management.MultipartToSinglepart(sfha_path, sfha_singlepart) # Remove polygons smaller than delete_size sqft self.remove_small_polys(sfha_singlepart, "FLD_ZONE = 'X' And ZONE_SUBTY = '0.2 PCT ANNUAL CHANCE FLOOD HAZARD'", spe, delete_size) ########################## # Move 500-year vertices # ########################## # Extract 500-year vertices from the SFHA layer and measure distances to the ZONE A* features arcpy.AddMessage('Extracting 500-year polygons') sfha_500 = arcpy.conversion.FeatureClassToFeatureClass(sfha_singlepart, working_gdb, 'sfha_500', "FLD_ZONE = 'X' And ZONE_SUBTY = '0.2 PCT ANNUAL CHANCE FLOOD HAZARD'") arcpy.AddMessage('Extracting 100-year polygons') sfha_100 = arcpy.conversion.FeatureClassToFeatureClass(sfha_singlepart, working_gdb, 'sfha_100', "FLD_ZONE <> 'X'") arcpy.AddMessage('Extracting 500-year vertices') vertices_500 = arcpy.management.FeatureVerticesToPoints(sfha_500, os.path.join(working_gdb_path, 'sfha_500_vertices'), 'ALL') arcpy.AddMessage('Calculate distances between 500-year vertices and 100-year polygons') arcpy.analysis.Near(vertices_500, sfha_100, '2 Feet', 'LOCATION') # Find vertices on 500-year polygons that need to be moved move_500 = {} with arcpy.da.SearchCursor(vertices_500, ['SHAPE@X', 'SHAPE@Y','NEAR_X','NEAR_Y'], 'NEAR_DIST > 0') as cursor: for row in cursor: move_500[(row[0],row[1])] = arcpy.Point(row[2], row[3]) arcpy.AddMessage(move_500) # Move vertices from 500-year to 100-year arcpy.AddMessage('Moving the necessary 500-year flood vertices') lines = self.move_vertices(sfha_singlepart, "FLD_ZONE = 'X' And ZONE_SUBTY = '0.2 PCT ANNUAL CHANCE FLOOD HAZARD'", move_500) ################ # Create lines # ################ # Create feature class to track vertex changes sfha_sp = self.get_fc_spref(sfha_path) arcpy.management.CreateFeatureclass(working_gdb_path, 'vertex_changes', 'POLYLINE', has_m='DISABLED', has_z='DISABLED', spatial_reference=sfha_sp) vertex_changes_path = os.path.join(working_gdb_path, 'vertex_changes') vertex_cursor = arcpy.da.InsertCursor(vertex_changes_path, 'SHAPE@') # Populate lines feature class for line in lines: vertex_cursor.insertRow([line]) ############ # Clean up # ############ arcpy.AddMessage('Cleaning things up') arcpy.management.RepairGeometry(sfha_singlepart, 'DELETE_NULL', 'ESRI') arcpy.management.MultipartToSinglepart(sfha_singlepart, output_path) self.remove_small_polys(output_path, "FLD_ZONE = 'X' And ZONE_SUBTY = '0.2 PCT ANNUAL CHANCE FLOOD HAZARD'", spe, delete_size) # Add layers to map arcpy.AddMessage(os.path.join(working_gdb_path, 'sfha_500_vertices')) self.add_data_to_map(os.path.join(working_gdb_path, 'sfha_500_vertices')) arcpy.AddMessage(output_path) self.add_data_to_map(output_path) arcpy.AddMessage(vertex_changes_path) self.add_data_to_map(vertex_changes_path) def remove_small_polys(self, fc, filter, sp, area_in_sqft, area_unit='SQUARE_FEET_US'): arcpy.AddMessage(f'Deleting polygons with an area of less than {str(area_in_sqft)} {area_unit}') with arcpy.da.UpdateCursor(fc, ['OID@','SHAPE@AREA'], filter, spatial_reference=sp) as cursor: for row in cursor: if row[1] < area_in_sqft: arcpy.AddMessage(f'-- Deleted polygon with OID {str(row[0])} ({str(round(row[1],2))} {area_unit})') cursor.deleteRow() def add_data_to_map(self, path): aprx = arcpy.mp.ArcGISProject('CURRENT') mapx = aprx.activeMap mapx.addDataFromPath(path) def get_gdb_path(self, fc_path): desc = arcpy.Describe(fc_path) arcpy.AddMessage(desc.path) return desc.path def get_gdb_path(self, fc_path): '''Return the Geodatabase path from the input table or feature class. :param input_table: path to the input table or feature class From https://stackoverflow.com/questions/29191633/arcpy-get-database-path-of-feature-class-in-feature-dataset ''' workspace = os.path.dirname(fc_path) if [any(ext) for ext in ('.gdb', '.mdb', '.sde') if ext in os.path.splitext(workspace)]: return workspace else: return os.path.dirname(workspace) def get_fc_spref(self, fc_path): return arcpy.Describe(fc_path).spatialReference def move_vertices(self, fc_path, query, vertex_array): lines = [] # Start editing edit = arcpy.da.Editor(self.get_gdb_path(fc_path)) edit.startEditing(False, False) edit.startOperation() try: # Loop through the vertices of the 500-year polygons and move those that are in the vertex_array dictionary spref = self.get_fc_spref(fc_path) with arcpy.da.UpdateCursor(fc_path, ['SHAPE@', 'OID@'], query) as cursor: for row in cursor: # Loop through the polygon parts arcpy.AddMessage(f'-- Processing OID {str(row[1])}') part_num = 0 new_geom = [] geom = row[0] if geom.isMultipart: arcpy.AddWarning('Multi-Part polygon') # multi_part_array = arcpy.Array() # for part in geom: # array = arcpy.Array() # for pnt in part: # if pnt: # if (pnt.X, pnt.Y) in vertex_array: # new_pt = vertex_array[(pnt.X, pnt.Y)] # array.add(arcpy.Point(new_pt.X, new_pt.Y)) # arcpy.AddMessage(f'Moved a point from {pnt.X},{pnt.Y} to {new_pt.X},{new_pt.Y}') # # Add line to represent the moved vertex # insert_cursor.insertRow([[(pnt.X, pnt.Y),(new_pt.X, new_pt.Y)]]) # else: # array.add(arcpy.Point(pnt.X, pnt.Y)) # else: # arcpy.AddMessage('Got an island') # multi_part_array.add(array) # try: # polygon = arcpy.Polygon(multi_part_array, spref) # #arcpy.AddMessage(coords) # row[0] = polygon # cursor.updateRow(row) # except TypeError: # arcpy.AddMessage(f'Got an error on OID {str(row[1])}') else: coords = [] for part in geom: #arcpy.AddMessage('in part') for pnt in part: #arcpy.AddMessage('in point') if pnt: if (pnt.X, pnt.Y) in vertex_array: new_pt = vertex_array[(pnt.X, pnt.Y)] coords.append((new_pt.X, new_pt.Y)) arcpy.AddMessage(f'Moved a point from {pnt.X},{pnt.Y} to {new_pt.X},{new_pt.Y}') # Add line to represent the moved vertex lines.append([(pnt.X, pnt.Y),(new_pt.X, new_pt.Y)]) #insert_cursor.insertRow([[(pnt.X, pnt.Y),(new_pt.X, new_pt.Y)]]) else: coords.append((pnt.X, pnt.Y)) else: arcpy.AddMessage('Got an island') try: #arcpy.AddMessage(coords) row[0] = coords cursor.updateRow(row) except TypeError: arcpy.AddMessage(f'Got an error on OID {str(row[1])}') except: raise finally: # Stop editing edit.stopOperation() edit.stopEditing(True) return lines