Expand extent for purposes of clipping

5757
14
Jump to solution
12-09-2015 11:27 AM
RebeccaStrauch__GISP
MVP Emeritus

I'm wanting to clip a hydro layer (with a reselect) based on the extent of my study area extent....but expanded a bit (% to be determined). I have code that works, but it's a bit clunky. I grab the extent, create a poly, buffer, grab-new-extent, create the poly and then do my clip/query.  It's the getting the extent->poly->buffer->extent->poly that seems a bit overkill.  I'm wondering if I am missing an obvious arcpy option for expanding the extent for the clip.  Trying to keep it to arcpy and not get deep into the geometry of this one...but all suggestions will be considered.

Thanks.

This works:

import arcpy
from arcpy import env
import os
arcpy.env.overwriteOutput = "True"

inStudy =  r'C:\Prep.gdb\study_boundary'
srDesc = arcpy.Describe(inStudy).spatialReference
extent = arcpy.Describe(inStudy).extent
arcpy.env.outputCoordinateSystem = srDesc.PCSCode

# ---> Ugly way to do this, but it works and fairly fast.  Ok for now.
# get extent, a buff distance (1% of X), and create a temp polygon
extBuffDist = ((int(abs(extent.lowerLeft.X - extent.lowerRight.X))) * .01)

# Array to hold points for the bounding box for initial extent
origExtentPts = arcpy.Array()
origExtentPts.add(extent.lowerLeft)
origExtentPts.add(extent.lowerRight)
origExtentPts.add(extent.upperRight)
origExtentPts.add(extent.upperLeft)
origExtentPts.add(extent.lowerLeft)  # ensures polygon closes
# Create and buffer the polygon object
polygonTmp1 = arcpy.Polygon(origExtentPts)
arcpy.Buffer_analysis(polygonTmp1, "polyBuff", extBuffDist, "OUTSIDE_ONLY")
origExtentPts.removeAll()

#get extent of buffered poly
extent = arcpy.Describe("polyBuff").extent
# Array to hold points for the buffer bounding box
newExtentPts = arcpy.Array()
newExtentPts.add(extent.lowerLeft)
newExtentPts.add(extent.lowerRight)
newExtentPts.add(extent.upperRight)
newExtentPts.add(extent.upperLeft)
newExtentPts.add(extent.lowerLeft) # ensures polygon closes
polygonTmp2 = arcpy.Polygon(newExtentPts)

# save to disk
#extentPoly = "extentPoly" 
arcpy.CopyFeatures_management(polygonTmp2, "extentPoly")
del polygonTmp1, polygonTmp2

hydroNHDSource = r"\\myServer\NHDH_AK.gdb\Hydrography\NHDFlowline"
hydroNHDQuery = "(GNIS_Name >= '1' AND ( FType = 566 OR FType = 460 OR FType = 558))"

arcpy.Clip_analysis(hydroNHDSource, "extentPoly", "hydroExtent")   
arcpy.MakeFeatureLayer_management("hydroExtent", "HydroReselect", hydroNHDQuery)
arcpy.Buffer_analysis("HydroReselect", "hydroBuff", 500, "FULL", "ROUND", "ALL" )
hydroBuffStudy = arcpy.Clip_analysis("hydroBuff", inStudy)
0 Kudos
1 Solution

Accepted Solutions
DanPatterson_Retired
MVP Emeritus

Darren, just grab the points from the extent to create the polygon, unfortunately, one would assume that they were a variant of a list but you have to go through the process Rebecca of getting them by name to create the polygon.  This hasn't even changed in Pro Extent—ArcPy Classes | ArcGIS for Desktop but at least your buffer cuts out a step.

View solution in original post

14 Replies
DarrenWiens2
MVP Honored Contributor

You can make use of the buffer method of the polygon object to make a new buffered polygon geometry:

polygonTmp1 = arcpy.Polygon(origExtentPts)
buffPoly = polygonTmp1.buffer(extBuffDist)
newExtent = buffPoly.extent # this returns an 'extent' object. May have to convert to 'polygon'
DanPatterson_Retired
MVP Emeritus

Darren, just grab the points from the extent to create the polygon, unfortunately, one would assume that they were a variant of a list but you have to go through the process Rebecca of getting them by name to create the polygon.  This hasn't even changed in Pro Extent—ArcPy Classes | ArcGIS for Desktop but at least your buffer cuts out a step.

RebeccaStrauch__GISP
MVP Emeritus

Darren's answer was helpful, but as you mention Dan, I think I would still have to go thru the process of grabbing the LL, LR, etc corners into an array to make the poly.      ??    I guess it's the extent->poly that I was hoping to simplify, i.e. without the array.  Maybe not possible?

As a side note, the code Darren provided did allow me to remove two lines...

arcpy.Buffer_analysis(polygonTmp1, "polyBuff", extBuffDist, "OUTSIDE_ONLY")
#...
extent = arcpy.Describe("polyBuff").extent
0 Kudos
DanPatterson_Retired
MVP Emeritus

Yes you still have to get the values since they aren't variants of lists, but so you can't slice them but have to get them by name.

for row in arcpy.da.SearchCursor(in_FC, ["SHAPE@"]):
    extent = row[0].extent
    L,B,R,T = [extent.XMin, extent.YMin, extent.XMax, extent.YMax]
    frmt = "arcpy extent object {}\nScraped values for ... L {} B {} R {} T {}"
    print(frmt.format(extent, L, B, R, T))

but then at least you add/subtract from L,B,R,T to produce new extent objects and formulate a new point if you wanted to. or just use that after buffering the r[0] object above

0 Kudos
DarrenWiens2
MVP Honored Contributor

I suppose since you're expanding the extent by some number in all directions, you could just add/subtract that number to/from your original extent, rather than running a buffer process at all (I'm not sure if that's what Dan was getting at. If so, carry on).

0 Kudos
RebeccaStrauch__GISP
MVP Emeritus

I took bits from both you and Dan and was able to get it to something that works (and give the same results, as I needed).  I didn't need a cursor as Dan pointed out (overkill for this project) but reassigning the extent.lowerLeft etc so shorter variables and then creating the array with a for loop worked for me. I wanted to use a buffer since that was the easiest/cleanest to move the extent out equal, in all directions.

I could use the output of the temp polygon object created from the extent of the buffer in my clip, so all works fine.  Final snippet of code..

import arcpy
from arcpy import env
import os
arcpy.env.overwriteOutput = "True"

inStudy =  r'C:\Prep.gdb\study_boundary' # r"c:\temp\USA.gdb\States"
srDesc = arcpy.Describe(inStudy).spatialReference
arcpy.env.outputCoordinateSystem = srDesc.PCSCode

# get extent, a buff distance (1% of X), and create a temp polygon
extent = arcpy.Describe(inStudy).extent
ll, lr, ur, ul = [extent.lowerLeft, extent.lowerRight, extent.upperRight, extent.upperLeft]
extBuffDist = ((int(abs(ll.X - lr.X))) * .01)
origExtentPts = arcpy.Array()
# Array to hold points for the bounding box for initial extent
for coords in [ll, lr, ur, ul, ll]:
    origExtentPts.add(coords)

polygonTmp1 = arcpy.Polygon(origExtentPts)
# buffer the temporary poly by 1% of width of extent as caluculated above
buffPoly = polygonTmp1.buffer(extBuffDist)
newExtent = buffPoly.extent

# extent points for polygonTmp2
ll2, lr2, ur2, ul2 = [newExtent.lowerLeft, newExtent.lowerRight, newExtent.upperRight, newExtent.upperLeft]
newExtentPts = arcpy.Array()
for coords in [ll2, lr2, ur2, ul2, ll2]:
    newExtentPts.add(coords)

polygonTmp2 = arcpy.Polygon(newExtentPts)

hydroNHDSource = r"\\myServer\NHDH_AK.gdb\Hydrography\NHDFlowline"
hydroNHDQuery = "(GNIS_Name >= '1' AND ( FType = 566 OR FType = 460 OR FType = 558))"

arcpy.Clip_analysis(hydroNHDSource, polygonTmp2, "hydroExtent") 
arcpy.MakeFeatureLayer_management("hydroExtent", "HydroReselect", hydroNHDQuery)
arcpy.Buffer_analysis("HydroReselect", "hydroBuff", 500, "FULL", "ROUND", "ALL" )
hydroBuffStudy = arcpy.Clip_analysis("hydroBuff", inStudy)

At least the array is easier to read.  So, I think Dan get's the correct answer for pointing out that I couldn't really get around the array for the coordinates.

For anyone that cares, this is a part of a much larger project, but I needed to buffer a subset of the rivers within a study area, but since the NHD hydro in our area is still a bit fuzzy, and many of our study areas use the topo maps to create the boundaries (many times based on hydro), I needed to buffer the rivers before I clipped.

Temporary layers..

and the final result (study area and buffers)

Thanks Darren and Dan for you thoughts on this.  It helped a lot.

0 Kudos
DarrenWiens2
MVP Honored Contributor

I see you've got this figured out, but for this, why fight with the extents at all? Just buffer the study area, clip those streams in buffered area, buffer streams, clip back to study area. Perhaps it takes more time to buffer the study area than to clip a bunch of extra streams. Maybe it's six one way, half dozen the other, but I'd guess dealing with the fewest streams possible would be a good idea.

RebeccaStrauch__GISP
MVP Emeritus

Ahh, yes good thoughts, and in many cases that would be better, but I'm using the same clip boundary for a bunch of other layers I'm clipping in the same script (topo maps, DEM, etc) that I want to be that same extent+, so this was just a piece of what I needed.  I just simplified the code to keep that part out, but worth noting for others.

DanPatterson_Retired
MVP Emeritus

Rebecca Strauch, GISP​ alternate ways to get extent parameters to simplify code

>>> import numpy as np
>>> ext = extent.__str__()
>>> L, B, R, T, *ZZMM = np.fromstring(ext,dtype=float,sep= " ")
>>> l, b, r, t, *zzmm = [float(i) for i in ext.split(" ")]

Where lines 2 and 4 are the key lines are for pure python implementations and 2 and 3 are if numpy is used.... of course the letters represent their obvious values and *zzmm variants are the extra z, m min, max values if they are needed.

0 Kudos