subdividing a rectangle using numpy

Seb spluque at gmail.com
Thu Sep 10 23:05:07 EDT 2015


Hello,

The code below is what I came up with to solve the problem:

1. We're given geographic coordinates for two opposite vertices of a
   rectangle.
2. If any of the sides of the rectangle is larger than some number of
   degrees, then subdivide the rectangle into squares/rectangles such
   that all sub-units have sides smaller than the specified amount.
3. Define each sub-unit by providing a list of all vertices, where the
   first and last vertices are identical so as to close the polygon.

I've ignored the "if" part of the problem to simplify.

The key to my solution was to use numpy's meshgrid to generate the
coordinates for defining the sub-units.  However, it seems awfully
complex and contrived, and am wondering if there's a simpler solution,
or perhaps some package offers this functionality.  I couldn't find any,
so any tips appreciated.

---<--------------------cut here---------------start------------------->---
# Let's say we have these coordinates
lons = [-96, -51.4]
lats = [60, 72]
# And we want to create 10x10 (max) degree polygons covering the rectangle
# defined by these coordinates
step = 10
# Calculate how many samples we need for linspace.  The ceiling is required
# to cover the last step, and then add two to accomodate for the inclusion
# of the end points... what a pain.
xn = np.ceil((lons[1] - lons[0]) / step) + 2
yn = np.ceil((lats[1] - lats[0]) / step) + 2
xgrd = np.linspace(lons[0], lons[1], xn)
ygrd = np.linspace(lats[0], lats[1], yn)
# Create grids of longitudes and latitudes with dimension (yn, xn).  The
# elements of the longitude grid are the longitude coordinates along the
# rows, where rows are identical.  The elements of the latitude grid are
# the latitude coordinates along the columns, where columns are identical.
longrd, latgrd = np.meshgrid(xgrd, ygrd, sparse=False)
for i in range(int(xn) - 1):
    for j in range(int(yn) - 1):
        print [(longrd[j, i], latgrd[j, i]), # lower left
               (longrd[j, i + 1], latgrd[j, i]), # lower right
               (longrd[j, i + 1], latgrd[j + 1, i]), # upper right
               (longrd[j, i], latgrd[j + 1, i]),     # upper left
               (longrd[j, i], latgrd[j, i])]         # close at lower left
---<--------------------cut here---------------end--------------------->---


-- 
Seb




More information about the Python-list mailing list