1

I am trying to clip a raster using a geoJSON without writing to my disk.

The following code compiles and produces the correct result.

import tempfile
import fiona
import rasterio
from rasterio.io import MemoryFile
import rasterio.mask
from scipy import stats

def getSummStats(raster, shapefile): #Open Polygon with fiona.open(shapefile, "r") as shapefile: shapes = [feature["geometry"] for feature in shapefile]

#Clip the raster
out_image, out_transform = rasterio.mask.mask(src, shapes, crop=True)
out_meta = src.meta

#Setup to save and write the clipped raster in memory
with MemoryFile() as memfile:
    with memfile.open(driver = "GTiff", height = out_image.shape[1], width = out_image.shape[2], count=1, 
                      dtype=rasterio.float64, transform = out_transform) as dataset:
        #write the clipped raster to memory
        dataset.write(out_image)
        #read the clipped raster from memory
        ds = dataset.read()
        return stats.describe(ds, axis = None)    

with rasterio.open("B:/dd.weather.gc.ca/model_gem_global/15km/grib2/lat_lon/00/000/CMC_glb_TMP_TGL_2_latlon.15x.15_2020052500_P000_Copy.grib2") as src: print(getSummStats(src, "C:/Users/Tom/Downloads/canada_provinces/canada_provinces/canada_provinces.geojson"))

However, I would like to be able to call this script from the command line and input the raster and geoJSON when the script is called. I made some changes to accomplish this.

import tempfile
import fiona
import rasterio
from rasterio.io import MemoryFile
import rasterio.mask
from scipy import stats
import os, sys

def getSummStats(raster, geoJSON_path): #Open Polygon with fiona.open(geoJSON_path, "r") as shapefile: shapes = [feature["geometry"] for feature in shapefile]

#Clip the raster
out_image, out_transform = rasterio.mask.mask(src, shapes, crop=True)
out_meta = src.meta

#Setup to save and write the clipped raster in memory
with MemoryFile() as memfile:
    with memfile.open(driver = "GTiff", height = out_image.shape[1], width = out_image.shape[2], count=1, 
                      dtype=rasterio.float64, transform = out_transform) as dataset:
        #write the clipped raster to memory
        dataset.write(out_image)
        #read the clipped raster from memory
        ds = dataset.read()
        return stats.describe(ds, axis = None)    

def main(raster_path, geoJSON_path): with rasterio.open(raster_path) as src: print(getSummStats(src, geoJSON_path))

if name == 'main': # # example run : $ Polygon_CMDVersion.py /<full-path>/<raster-name>.grib2 /<full-path>/<geoJSON-name>.geoJSON # if len( sys.argv ) < 2: print ("[ ERROR ] you must two args. 1) the full raster path and 2) the full geoJSON path") sys.exit( 1 )

main(sys.argv[1], sys.argv[2])

However, this second method gives me the following error

File "C:\Users\Tom\Downloads\Polygon_CMDVersion.py" line 21, in getSummStats
out_image, out_transform = rasterio.mask.mask(src, shapes, crop=True)
NameError: name 'src' is not defined

If this helps, here is what I wrote in the command line

Polygon_CMDVersion.py B:/dd.weather.gc.ca/model_gem_global/15km/grib2/lat_lon/00/000/CMC_glb_TMP_TGL_2_latlon.15x.15_2020052500_P000_Copy.grib2 C:/Users/Tom/Downloads/canada_provinces/canada_provinces/canada_provinces.geoJSON

How can I resolve this issue?

notmenotme
  • 23
  • 5

1 Answers1

0

The comment made by @user2856 was the correct answer:

You haven't defined src in your getSummStats function, nor is it passed in as an argument. You have a raster arg but do not use it.

Aaron
  • 51,658
  • 28
  • 154
  • 317
notmenotme
  • 23
  • 5