2

I have a point shapefile or feature class and I want to create bounding boxes around each point. The size of the bounding box will depend on user input, for example the user might want a bounding box of 50 x 50 meters.

The shapefiles are in WGS 84. This boils down to getting the x/y coordinates of a point x meters away from the existing point to the South East and to the North West (bearing 135 and 315 degrees).

How can I use arcpy to get the x/y coordinates of a point n meters away to the NW and the SE?

5
  • 2
    Why not just buffer the points by the desired amount and then use the Minimum Bounding Geometry tool to get the extent? Aug 19, 2015 at 13:58
  • A "box" in that coordinate reference shouldn't look rectangular (unless the point is on the equator), so you'll need at least four points (adding bearings 45 and 225). Large shapes (10+KM) should have even more vertices. This question is among the most common here at GIS SE
    – Vince
    Aug 19, 2015 at 13:59
  • I am creating the box so I can query a wms service that requires a bounding box parameter. I don't expect the size to be more than 5 km in diameter, but most should be much smaller (50 to 500 m).
    – bdwyer
    Aug 19, 2015 at 14:15
  • Based off Evil Genius' comment, I think I can run a buffer on the points and then get the arcpy geometry extent and use the XMin,XMax,YMin,YMax values to get my SE and NW points.
    – bdwyer
    Aug 19, 2015 at 14:26
  • 1
    As long as the distances aren't massive and you don't need high levels of accuracy, you could just use trigonometry. You could even project the points on the fly if you wanted.
    – Paul
    Aug 19, 2015 at 15:55

1 Answer 1

1

As @Paul suggested, you can use some simple trig to do this:

import arcpy

def getPoly(pt, w, h, wkid):
    """calculates polygon corners

    pt -- tuple of coordinates (x, y)
    w -- width of polygon (x length)
    h -- height of polygon (y length)
    sr -- spatial reference object
    """
    x, y = pt
    ll = (x - (w * 0.5), y - (h * 0.5))
    ul = (x - (w * 0.5), y + (h * 0.5))
    ur = (x + (w * 0.5), y + (h * 0.5))
    lr = (x + (w * 0.5), y - (h * 0.5))
    return arcpy.Polygon(arcpy.Array([arcpy.Point(*coords) for coords in [ll,ul,ur,lr,ll]]), sr)

coords =(-10377711.530988107, 5617827.075599414)
sr = arcpy.SpatialReference(102100)
poly = getPoly(coords, 50, 50, sr) # make a bounding box 50x50 meters
e = poly.extent
bbox = ','.join(map(str, [e.XMin, e.YMin, e.XMax, e.YMax])) #here's the bbox

The bounding box came out as:

'-10377736.531,5617802.0756,-10377686.531,5617852.0756'

2
  • crmackey, thanks for the post. Will the code work with a shapefile in WGS84 (4326)?
    – bdwyer
    Aug 20, 2015 at 15:12
  • Yep, it should work in any projection, just make sure you pass that WKID into the arcpy.SpatialReference() class.
    – crmackey
    Aug 20, 2015 at 15:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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