CDO

Hi Pao,

if you know the gridcell indices then you can use the selgridcell operator.

-Karin

Dear Karin,

thank you for the kind answer. Unfortunately I do not know the gridcell indices.

But I have the shapefile of the area from which I want to extract the netcdf data.

Here below the links for downloading both shapefile and netcdf file:

shapefile [[https://www.dropbox.com/s/vrhlzwk71ldpl0b/shape1.zip?dl=0]]
netcdf [[https://www.dropbox.com/s/03j1w4ncmcijkes/nc_file.nc?dl=0]]

I tried to follow this post [[https://code.mpimet.mpg.de/boards/2/topics/6693]] but I was unsuccessful. After adapting the NCL script and run it I got the following error: shapefile_mask_data: Error: not a valid rectilinear, curvilinear, or unstructured grid
Maybe the problem is that my netcdf variable ua has 4 dimensions: time, plev, lat, lon, and the NCL script only allows 3 dimensions in total. I do not need plev and I tried to remove it with NCO following these discussions [[https://sourceforge.net/p/nco/discussion/9830/thread/3d35b42d/]] but it did not work.

Is there any way I can extract form the netcdf file only the data which are contained in the shapefile boundaries?

Thanks for any suggestion.

You can reduce a single value dimension with

cdo --reduce_dim -copy nc_file.nc nc_file_3d.nc

Thanks, it does help. But my output is still not correct.

Here what I have done:

1) remove dimension

cdo --reduce_dim -copy nc_file.nc nc_file_3d.nc

2) updated the NCL script with input netcdf = nc_file_3d.nc and remove one dimension from ua variable, then create the netcdf mask:

ncl create_nc_mask.ncl

3) apply the mask to nc_file_3d.nc:

cdo div nc_file_3d.nc shape1_mask.nc nc_file_masked.nc

Attached below you can see the output plotted (ua grid-point mean), which is similar, but not exactly the same as the original shapefile.

Any suggestion on what it is wrong?
Thanks

The NCL script has to be modified:

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

shpname  = "./shape1/shape1.shp" 
maskname = "./shape1/shape1_mask.nc" 

print_shapefile_info(shpname)

;-- open data file to get the grid to be used
f   = addfile("nc_file.nc","r")
var = f->ua(0,0,:,:)
lat = f->lat
lon = f->lon

delete(lon@bounds)
delete(lat@bounds)

;-- 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   = "lat" 
mask_array!1   = "lon" 
mask_array&lon = lon
mask_array&lat = lat

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

fout->mask_array =  mask_array

Then you can do the following:

cdo --reduce_dim -copy nc_file.nc tmp1.nc 

cdo -remapbil,"./shape1/shape1_mask.nc" tmp1.nc tmp2.nc

cdo -mul tmp2.nc ./shape1/shape1_mask.nc outfile.nc

-Karin

Hi Karin,

thank you for your help and time.

I followed your instructions but I still get a wrong output file, with the final area different from the shapefile's one.

In simple words, my question is:

how can I extract CMIP6 data from climate models netcdf outputs that belong to SREX region 11 (North Europe, NEU)?
[[https://www.ipcc-data.org/guidelines/pages/ar5_regions.html]]

Thanks for any help.

The outfile.nc shows the correct data within the NEU region. But maybe you want to set the 0 values from the mask to missing value then do

cdo -mul tmp2.nc -setctomiss,0 ./shape1/shape1_mask.nc outfile.nc

I am sorry, but our outputs do not correspond to SREX region NEU (11).

Attached you can find the original NEU region boundaries. The pink area is the NEU area.

Thanks

But the shapefile has the above region.

ok. I made some progress. the outputs did not cover the British Isles because the netcdf file longitude coordinates needed to be converted to WGS84 -180,180.

Here the complete steps I followed:

1. remove plev dimension from netcdf file
   cdo --reduce_dim -copy nc_file.nc nc_file2.nc

1. convert longitude to -180, 180
   cdo sellonlatbox,-180,180,-90,90 nc_file2.nc nc_file3.nc

1. create mask (see NCL script updated [[https://www.dropbox.com/s/ivaazzimlj4fafc/create_nc_mask.ncl?dl=0]])
   ncl create_nc_mask_NEU.ncl

1. apply mask
   cdo div nc_file3.nc shape1_mask.nc nc_file4.nc #################

Attached the output. What is missing from the original SREX 11 area is its southern part (e.g. Denmark and southern England are included) and a correct representation/boundary of the southern diagonal line.