I want to create a numpy array of size 2160 * 4320 with each grid cell containing the size of a 5' grid cell on earth.
Is there a python library or function to do this?
I am using GDAL and rasterio.
-
I suspect that there are many. Are you already using ArcPy or any other common Python libraries for GIS?– PolyGeo ♦Sep 7, 2016 at 4:30
-
@PolyGeo, I am not using ArcPy, can use gdal or rasterio though.– user1186Sep 7, 2016 at 4:32
-
2I'm not sure I totally understand what you are asking...do you need to load data into this numpy array or are you just trying to create it? Numpy by itself has no sense of spatial information, so creating a 2160x4320 array is simply that.– Evil GeniusSep 7, 2016 at 11:25
Add a comment
|
1 Answer
Numpy doesn't support any spatial info. You can use pyshp library for the same. https://github.com/GeospatialPython/pyshp
import shapefile as shp
import math
minx,maxx,miny,maxy = xxxxx, xxxxx, xxxxx, xxxxx (specify the extents)
dx = xx (size of grid)
dy = xx (size of grid)
nx = int(math.ceil(abs(maxx - minx)/dx))
ny = int(math.ceil(abs(maxy - miny)/dy))
w = shp.Writer(shp.POLYGON)
w.autoBalance = 1
w.field("ID")
id=0
for i in range(ny):
for j in range(nx):
id+=1
vertices = []
parts = []
vertices.append([min(minx+dx*j,maxx),max(maxy-dy*i,miny)])
vertices.append([min(minx+dx*(j+1),maxx),max(maxy-dy*i,miny)])
vertices.append([min(minx+dx*(j+1),maxx),max(maxy- dy*(i+1),miny)])
vertices.append([min(minx+dx*j,maxx),max(maxy-dy*(i+1),miny)])
parts.append(vertices)
w.poly(parts)
w.record(id)
w.save('polygon_grid')
