0

Obviously it's difficult to reproduce a netCDF here on Stackexchange, especially a curvilinear/rotated grid type.

So let me start with what I already know. I know how to mask a regular netCDF with rectangular/square grid with a shapefile.

For this I use NCL programming and it's routine "shapefile_utils.ncl". This is my code:

;-- copy shapefile_utils.ncl from
;-- https://www.ncl.ucar.edu/Applications/Scripts/shapefile_utils.ncl
load "shapefile_utils.ncl"

shpname  = "gadm36_GER_0.shp"
maskname = "germany_mask.nc"

print_shapefile_info(shpname)

;-- open data file to get the grid to be used
f = addfile("Europe.nc","r")

;-- read variable
var                    =  f->tp(0,:,:)
var@lat2d              =  f->latitude
var@lon2d              =  f->longitude

;-- shapefile mask resources
opt             =  True
opt@return_mask =  True    ;-- this forces the return of a 0s and 1s mask array

;-- create the mask based on the given shapefile
mask_array             =  shapefile_mask_data(var, shpname, opt)
mask_array!0           = "y"
mask_array!1           = "x"
mask_array@coordinates = "latitude longitude"

;-- create new netCDF file and write mask array
system("rm -f " + maskname)
fout = addfile(maskname,"c")

fout->mask_array =  mask_array
fout->latitude       =  f->latitude
fout->longitude      =  f->longitude

So in this specific example of NCL code, gadm36_GER_0.shp is the shapefile of Germany, then the mask nc file I want to create is germany_mask.nc and the nc file I am going to mask is Europe.nc.

Then finally in I use CDO in Linux to mask the netCDF.

cdo div inputfile maskfile outputfile

Or in this example it the code would be:

cdo div Europe.nc germany_mask.nc Germany_masked.nc

My question would be: how to do this with a netCDF with a curvilinear/rotated grid.

All I can really do is show you some output of curvilinear/rotated grid nc file in R.

     5 variables (excluding dimension variables):
        double time_bnds[bnds,time]   (Chunking: [2,1])  
        double lon[rlon,rlat]   (Contiguous storage)  
            standard_name: longitude
            long_name: longitude
            units: degrees_east
            _CoordinateAxisType: Lon
        double lat[rlon,rlat]   (Contiguous storage)  
            standard_name: latitude
            long_name: latitude
            units: degrees_north
            _CoordinateAxisType: Lat
        char rotated_pole[]   (Contiguous storage)  
            grid_mapping_name: rotated_latitude_longitude
            grid_north_pole_latitude: 39.25
            grid_north_pole_longitude: -162
            long_name: coordinates of the rotated North Pole
        float pr[rlon,rlat,time]   (Chunking: [415,423,1])  
            standard_name: precipitation_flux
            long_name: Precipitation
            units: kg m-2 s-1
            grid_mapping: rotated_pole
            coordinates: lat lon
            _FillValue: 1.00000002004088e+20
            missing_value: 1.00000002004088e+20
            cell_methods: time: mean

     4 dimensions:
        time  Size:1   *** is unlimited *** 
            standard_name: time
            long_name: time
            bounds: time_bnds
            units: days since 1969-12-01T00:00:00Z
            calendar: proleptic_gregorian
            axis: T
        bnds  Size:2 (no dimvar)
        rlon  Size:415 
            standard_name: projection_x_coordinate
            long_name: longitude in rotated pole grid
            units: degrees
            axis: X
        rlat  Size:423 
            standard_name: projection_y_coordinate
            long_name: latitude in rotated pole grid
            units: degrees
            axis: Y

When I directly apply the NCL code shown above on the this specific curvilinear/rotated grid netCDF, it doesn't work. It doesn't throw an error, but the entire masked nc file has NA values. It basically masks the entire matrix and not the region I want to mask.

My question: Does anybody have any ideas on how to mask a curvilinear/rotated grid with a shapefile?

0

Your Answer

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