Set spatial extent to arcpy.da.SerachCursor and UpdateCursor

[yosukekimura]

Hi,

How can I iterate through feature layer or feature class by arcpy.da.SerachCursor and arcpy.da.UpdateCursor within a certain geographic extent?  Right now I am making feature class copying the content.  Is this most efficient method?  Basically I have giant feature class, and I want to visit geographic subdivision of the feature class and then do my processing there.

# I have a feature class.  I made it into feature layer
arcpy.MakeFeatureLayer_management(poly_fc, tmp_lyr)

# I then set extent to the part of data that I want to process
arcpy.env.extent = extent1

# if I count the feature from the layer, I got the count of feature in the extent. 
#  if I copy it into feature class I got the featueres inside the extent
n0 = int( arcpy.GetCount_management(tmp_lyr).getOutput(0)) # I got corrent n0
arcpy.CopyFeatures_management(tmp_lyr, selected_fc)  # I got fc of only the feature within extent1

# But iterating the tmp_lyr visits all featuers in tmp_lyr, ignoring the arcpy.env.extent
lst =  [ _[0] for _ in arcpy.da.SearchCursor(tmp_lyr, ['SHAPE@'])]   
#     lst above contains all features, not just the one in extent

# if I access through the copied fc, I get what I want, i.e., features within the extent
lst_small =  [ _[0] for _ in arcpy.da.SearchCursor(selected_fc, ['SHAPE@'])]   
#     lst_small above contains what I want

# so I just process on the copy....
with arcpy.da.SearchCursor(selected_fc, ['SHAPE@']) as cur:
  for row in cur:
    # my processing here...
    pass

[DanPatterson_Retired]

you can redirect your copy to in_memory if you want, if not, reuse the same filename during your process and delete when all done.

[yosukekimura]

Thanks, all my feature/layer are in_momory now.

[JoshuaBixby]

Use arcpy.SelectLayerByLocation_management first to create a selection on your feature layer, then create the search cursor.  If you are using ArcGIS 10.3 and already have an extent object created, you can pass the extent object's polygon property to arcpy.SelectLayerByLocation_management.  If you are using ArcGIS 10.2.2 or earlier, you can manually create a polygon using the extent object's bounds.

[yosukekimura]

This is what I wanted, seems like working the way I wanted.  Thanks!

[yosukekimura]

Actually, I am having trouble with Select Layer By Location approach.  I am getting different results with two method.  I dumped what's selected into arcmap and found that setting extent and copy method selects fearues I expected, which is any features which has shared point with the extent rectangle.  Somehow method 1 failed to choose some of feature.  My questions are

- Is "INTERSECT" correct method to mimic what extent does?

EDIT

Looking into further, this occurred for one small polygon near the edge of the extent (but completely contained by the extent rectangle), not being selected by Select Layer By Location method, whereas setting extent and then copying features catches this one feature.  I heve everything in same sapatial reference.  The extent polygon has dimension of tens of degree (WGS84), and the feature having trouble is an rectangle of  half minute in each side, off from the extent boundary by about half minutes.

I am using ArcGIS 10.2.2, by the way.

# prepare extent from precalculated corner coordinates
sr = arcpy.Describe(feat_lyr).spatialReference
myextent = arcpy.Extent( ((x0,y0), (x0,y1), (x1, y1), (x1, y0)), sr )

# Method 1, using Select By Location
arcpy.env.extent = None
geom_ext = arcpy.Polygon(arcpy.Array(
    [getattr( myextent, _ ) for _ in 
    ('lowerLeft', 'upperLeft', 'upperRight', 'lowerRight') ] )
arcpy.SelectLayerByLocation(feat_lyr, 'INTERSECT', geom_ext )
geoms1 = [row[0] for row in arcpy.SearchCursor(feat_lyr, ['SHAPE@'])]

# Method 2, copy feature
arcpy.env.extent = myextent
arcpy.CopyFeatures_management(feat_lyr, feat_fc)
geoms2 = [row[0] for row in arcpy.SearchCursor(feat_fc, ['SHAPE@'])]

[JoshuaBixby]

Regarding spatial reference, yes, it is always a good idea to specify a spatial reference when creating ArcPy geometry objects.  There are numerous, known precision issues that crop up with ArcPy geometry objects when spatial references aren't specified.  Is lack of a spatial reference the issue here?  We will know soon enough when you test it.

When you are talking about a "shared point with the extent rectangle," I am a bit uncertain whether that is what you mean.  "Shared point" could be taken to mean something very specific, like a polygon in question must share a vertex with the extent rectangle.  Try adding the spatial reference, and then report back.  Also, are the differences you are finding only along the periphery of the extent object?  Are any polygons more in the interior being handled different.

There are 15+ different types of overlaps one can specify.  I picked Intersect in this case because it is the default and most common.  It could be the CopyFeatures_management tool uses a different overlap type, or it could be a bug in one tool or the other.

[yosukekimura]

Thank you for help!

So my extent is one giant box.  my feature having trouble in particular this time is a small rectangle sitting inside the extent rectangle.  I understand that INTERSECT, CONTAINS etc has slight difference but I understand that they both select the blue feature in the pic below.  pink is the extent polygon

The extent is a lot larger than the feature.  Below show the extent, and tiny dot near top left corner is the feature causing trouble.

[JoshuaBixby]

If I understand you correctly, it is one polygon on or near the periphery that is not being handled the way you think.  A few questions.

• Is the missing polygon a single-part or multi-part polygon?
• You appear to be working at continental map scale, i.e., 1:60,000,000.  What spatial reference are you using?
• Which part specifically isn't working?  Does the arcpy.SelectLayerByLocation_management tool miss the polygon or does that step work and it is the search cursor that is missing it?
• Regarding cursors, you should really consider using the newer ArcPy Data Access (arcpy.da) cursors instead of the older/original ArcPy (arcpy) cursors.

[yosukekimura]

Thank you so much again for the help!

1.  It is a single part.  Most of features are, and this particular one definitely is.

2. I am using WGS84, and everything is made to WGS84.  I figured how to set reference to arcpy.Polygon, so the extent polygon knows the spatial reference now.

3. I actually tested what I am doing in my ArcMap.  Specifically, I run the script but I also exported the extent polygon, the feature class, feature layers, while I set the extent to be global first.  The image I had was showing those exported feature (shape file).  I then manually tried to select by geometry in two way.  One is from the ArcMap's tool (Selection ->Select By Location ...).  I also run "Select Layer By Location" tool from the tool box.  Neither of them did not select this blue feature, even though blue is included in the pink (extent).  Extent's top is 56.565833 dd, and the feature's top is 56.557500 dd , so off by 0.0833 degree (~ 0.5 minutes).

I also have a new observation.  I made the problematic feature into a feature class lone feature.  Then this select by location feature works, it picks up the feature.  But if I apply the tool to the original feature which has ~1600 features spreading across the continent, that particular feature is not selected.  The features are wide spread across the area, and each polygon is tiny in comparison to the extent.  So I am wondering ArcGIS somehow sets criteria to lower for the spatial accuracy and then this polygon got dropped somehow...?

I made three shape files, and posted here

https://utexas.box.com/s/a6q45925jrpxzknxdp357fremtjhprhz

It has three features, ext1, feat1 and feat1_sel.  ext1 is the extent polygon, feat1 the feature class of ~1600 features.  feat1_sel is the one having trouble. 

4. sorry for the type for cursor.  I always use arcpy.da.cursor.  They should have just overwritten old  cursor.

yosuke