Friday, November 6, 2015

Creating Data Driven Pages with a Maximum Scale

I often use Data Driven Pages to create sets of maps. Data Driven Pages iterates through a feature class and sets ArcMap's extent to match each feature of the feature class, and generates a PDF map for each. Sometimes the features in the feature class driving the Data Driven Pages can be small, meaning the scale of the map will be quite large. Large scale maps are sometimes not much use, as they don't have enough information to let the user know what the map represents. A while ago I creates a script that creates a new feature class from an input feature class. The script checks the scale of each feature in the feature class that will be generated with Data Driven Pages. If the scale is too large, it replaces the feature with a new one at the indicated maximum scale.

Inputs for the script are the feature class to check, the maximum scale desired in Layout View (for example, 1:2000 is input as 2000), the percent margin (corresponding to the percent margin in the Data Driven Pages setup window), the new feature class name and the new feature class location, and finally the map document you'll be using the new Data Driven Pages feature class with.

Here's an example of the output, with California's Bay Area divided into regions:


The input is the green features, and the output is the purple-bordered polygons. Smaller features of the input have been replaced by rectangular features at the maximum scale.

Using the output feature class to drive our Data Driven Pages, the maps created won't zoom in closer than desired, even on the smaller features:


The code:

#----------------------------#
#-------- User Inputs -------#
#----------------------------#
 
 
#Feature class used to generate DDP
DataDrivenPageFeature = r"C:\Users\e1b8\Desktop\E1B8\Blog\BlogData.gdb\Input"
#Maximum scale
MaxScale = 20000
#Percent margin map extent
MMEx = 125
 
#New feature class name
NewDDPFeat = "DDP_Ref_Output"
#New feature class GDB location
NewFeatGDB = r"C:\Users\e1b8\Desktop\E1B8\Blog\BlogData.gdb"
 
#Map document
MXDoc = r"C:\Users\e1b8\Desktop\E1B8\Blog\DdpMaxScale.mxd"
 
#----------------------------#
#---- End of User Inputs ----#
#----------------------------#
 
 
import arcpy
import os
 
arcpy.env.overwriteOutput = True
 
mxd = arcpy.mapping.MapDocument(MXDoc)
df = arcpy.mapping.ListDataFrames(mxd)[0]
mxd.activeView = "PAGE_LAYOUT"
 
#Data drive page feature type
fcShapeType = arcpy.Describe (DataDrivenPageFeature).shapeType
#Create temp feature class
DDPFeature = os.path.join ("in_memory", "tempFc")
if fcShapeType == "Polygon":
    print "Copying feature class to memory"
    arcpy.CopyFeatures_management (DataDrivenPageFeature,
                                   DDPFeature)
else:
    #Create minimum bounding rectangles
    print "Converting features to minimum bounding rectangles"
    arcpy.MinimumBoundingGeometry_management (DataDrivenPageFeature,
                                              DDPFeature, "ENVELOPE")
 
#Create layer for selection
arcpy.MakeFeatureLayer_management (DDPFeature, "DDPlyr")
 
#Get values for new features at maximum scale
df.scale = MaxScale / ((float(MMEx) / 100))
XAdjust = ((df.extent.XMax) - (df.extent.XMin)) / 2
YAdjust = ((df.extent.YMax) - (df.extent.YMin)) / 2
 
#Create new feature class
print
print "Creating new feature class"
print
 
NewFeat = os.path.join (NewFeatGDB, NewDDPFeat)
sr = arcpy.Describe(DDPFeature).SpatialReference
arcpy.CreateFeatureclass_management (NewFeatGDB,
                                     NewDDPFeat,
                                     "POLYGON",
                                     DDPFeature,
                                     "", "", sr)
 
#List fields
copyFeatureOIDList = []
 
shpFld = arcpy.Describe(DDPFeature).shapeFieldName
fieldNameList = [fld.name for fld in
                 arcpy.ListFields (DDPFeature)
                 if fld.editable and
                 fld.name != shpFld]
 
#Create cursor object
cursor = arcpy.da.SearchCursor(DDPFeature, ["SHAPE@",
                                            "SHAPE@TRUECENTROID",
                                            "OID@"]
                               + fieldNameList)
 
FeatureCount = 0
for row in cursor:
    df.extent = row[0].extent
    #Check if scale is larger than maximum
    if df.scale * (float(MMEx) / 100) >= MaxScale:
        #Add OID of feature to list if feature is large enough
        copyFeatureOIDList.append(str(row[2]))
        continue
     
    FeatureCount += 1
    print "Creating new feature", FeatureCount
    #Make x,y for four corners of extent
    array = arcpy.Array()
    FeatPnt1 = (row[1][0] + XAdjust, row[1][1] + YAdjust)
    FeatPnt2 = (row[1][0] + XAdjust, row[1][1] - YAdjust)
    FeatPnt3 = (row[1][0] - XAdjust, row[1][1] - YAdjust)
    FeatPnt4 = (row[1][0] - XAdjust, row[1][1] + YAdjust)
     
    pointList = [FeatPnt1, FeatPnt2, FeatPnt3, FeatPnt4]
    #Create each point and add to array
    for point in pointList:
        pnt = arcpy.Point(point[0], point[1])
        array.add(pnt)
    #Add first point to array to close polygon
    array.add(array.getObject(0))
    #Create a new polygon from array
    newPolygon = arcpy.Polygon(array)
    #Add polygon to new feature class
    arcpy.Append_management (newPolygon, NewFeat, "NO_TEST")
    #Update new feature with values
    oidFld = arcpy.Describe (NewFeat).OIDFieldName
    delimOidFld = arcpy.AddFieldDelimiters (NewFeat,
                                            oidFld)
    sql = "{0} = {1}".format (delimOidFld, FeatureCount)
    with arcpy.da.UpdateCursor(NewFeat,
                               ["OID@",
                                "SHAPE@",
                                "SHAPE@TRUECENTROID"]
                               + fieldNameList,
                               sql) as updateCursor:
        for udrow in updateCursor:
            #Populate fields
            if fieldNameList:
                udrow[3:] = row[3:]
            #Shift new rectangle from centroid
            #(center of gravity) to spatial center
            FeatmidpointX = (row[0].extent.XMin
                             + row[0].extent.XMax) / 2
            FeatmidpointY = (row[0].extent.YMin
                             + row[0].extent.YMax) / 2
            NewmidpointX = (udrow[1].extent.XMin
                            + udrow[1].extent.XMax) / 2
            NewmidpointY = (udrow[1].extent.YMin
                            + udrow[1].extent.YMax) / 2
            DiffX = NewmidpointX - FeatmidpointX
            DiffY = NewmidpointY - FeatmidpointY
            udrow[2] = udrow[2][0] - DiffX, udrow[2][1] - DiffY
            updateCursor.updateRow(udrow)
                 
    del updateCursor
     
del row
del cursor
 
#Select all features in list
if copyFeatureOIDList:
    if len(copyFeatureOIDList) == 1:
        print "Copying feature"
    else:
        print "Copying features"
    oidStr = ", ".join (map (str, copyFeatureOIDList))
    oidFld = arcpy.Describe("DDPlyr").OIDFieldName
    sql = "{0} IN ({1})".format (oidFld, oidStr)
    arcpy.SelectLayerByAttribute_management ("DDPlyr", "", sql)
     
#Copy selection
if arcpy.Describe("DDPlyr").FIDSet:
    arcpy.Append_management ("DDPlyr", NewFeat)
 
#delete temp feature
arcpy.Delete_management (DDPFeature)
 
print "Created:", NewFeat
print "Done"
profile for Emil Brundage at Geographic Information Systems Stack Exchange, Q&A for cartographers, geographers and GIS professionals

2 comments:

  1. Hi Emil.. Im a novice at coding - trying to get your code to work on my data - getting

    Runtime error
    Traceback (most recent call last):
    File "", line 36, in
    AttributeError: DescribeData: Method shapeType does not exist

    any chance you can tell me how to fix this?

    ReplyDelete