Topics

resampling raster after clipping/masking


Eyal Saiet
 

Hello,
I have found references in rasterio.readthedocs.io both on how to crop/clip a raster and how to resample. I generally understand the steps and sought to combine the two steps into one:

#Clipping Raster
with rasterio.open(dem_raster_pd['rasters'][0]) as raster:
    #create a mask
    out_image,out_transform=rio.mask.mask(raster,shapes,crop=True)
    #adding meta of source
    out_meta=raster.meta


    #updading meta
    out_meta.update({"driver":"GTiff","height":out_image.shape[1],
                                     "width":out_image.shape[2],
                                    "transform":out_transform})

#Because I am resampling next I skipped writing crp_raster
   ##Rescaling to 1 m
    ras_xy=1 #(m)

   # figuring out the factor
    rs_factor=((ras_xy/out_image.shape[1],ras_xy/out_image.shape[2])
    rs_in_factor=(1/rs_factor[0],1/rs_factor[1])

    #But because I don't have a crp_raster raster object I cannot use read and apply yhr out_shape
     data_resample=crp_raster.read(out_shape=(crp_raster.count,int(crp_raster.height*rs_in_factor[0]),int(crp_raster.width*rs_in_factor[1])),resampling=Resampling.bilinear,masked=True)#,masked=True


    #Scaling image trnasform
    transform_resmp=out_transform * out_transform.scale(
       (out_image.shape[2]/data_resample.shape[-1]),
        (out_image.shape[1]/data_resample.shape[-2])
        )

    #updating meta
    out_meta_res=crp_raster.meta
    out_meta_res.update({"height":data_resample.shape[-2],"width":data_resample.shape[-1],"transform":transform_resmp})


    #New filename
    new_ras_nm=products_dir+os.sep+raster.name.split('.')[0]+"_crp_rsmp_1m.tif"

    #Saving
with rasterio.open(new_ras_nm,"w",**out_meta_res) as dest_crp_resamp:
    print('Resampled')
    dest_crp_resamp.write(data_resample)     

I am sure I am stuck over something obvious that currently I cannot see. If I could create a raster object,  I tried crp_raster=rasterio.open("w"..), then I could use the raster_objec.read(out_shape=).