Version 4 (modified by nicolasmartin, 3 years ago) (diff)

Interpolation with CDO

You can find here an example of bash script to interpolate forcing core2 files.

To interpolate files you need to do 3 steps:

  1. have good grid description file (if not you can create it, see instructions 1.) below)
  1. prepare files to be interpolated ((see 2.) below)
  1. remapping on ORCA2 grid ((see 3.) below)

NOTE: This script works with cdo version 1.4.0.1

1) Prepare grid description file for the ORCA curvilinear grid

You need to prepare a netCDF grid description file for the ORCA curvilinear grid.
This file must contain:

  • a longitude variable with exactly 2 dimensions and the attribute: units = "degrees_east"
  • a latitude variable with exactly 2 dimensions and the attribute: units = "degrees_north"
  • a dummy variable with the attribute: coordinates = "xxx yyy" where xxx and yyy are the name of the longitude and latitude variables mentioned just above

There is an example of a netCDF grid description file for the T grid of ORCA2:

netcdf grid_ORCA2_T {
dimensions:
        y = 149 ;
        x = 182 ;
variables:
        byte dummy(y, x) ;
                dummy:coordinates = "glamt gphit" ;
        double glamt(y, x) ;
                glamt:units = "degrees_east" ;
        double gphit(y, x) ;
                gphit:units = "degrees_north" ;
}

There is an example of nco commands to create this file from coordinates.nc (input file used by NEMO for its grid description)

# select glamt and gphit variables from the coordinates file
ncks -O -C -a -v glamt,gphit coordinates_ORCA_R2.nc grid_ORCA2_T.nc

make sure that coordinates variables contains only 2 dimension. Use ncwa -a to remove degenerated dimensions (with a size of 1) for example:

# remove degenerated dimention time (if existing)
ncwa -O -a time grid_ORCA2_T.nc  grid_ORCA2_T.nc
# remove degenerated dimention z (if existing)
ncwa -O -a z grid_ORCA2_T.nc  grid_ORCA2_T.nc
# add a dummy variable
ncap2 -O -s 'dummy[y,x]=1b' grid_ORCA2_T.nc grid_ORCA2_T.nc
# add needed attributes
ncatted -a coordinates,dummy,c,c,'glamt gphit' \
        -a units,glamt,c,c,'degrees_east'      \
        -a units,gphit,c,c,'degrees_north'  grid_ORCA2_T.nc

2) Prepare file to be interpolated

for a curvilinear grid:

  • you have to follow above spefication

For a regular lon/lat grid:

  • The coordinate variables must have the same name as the dimensions. You can chose any name. If you use the name x for the dimension then you have to use also x for the coordinate variable.
  • The coordinate variable must also have the attribute units = "degrees_east" or "degrees_north" to identify the lon/lat variables:
netcdf etopo5 {
 dimensions:
         x = 4320 ;
         y = 2161 ;
 variables:
         float x(x) ;
                  x:units = "degree_east" ;
         float y(y) ;
                 y:units = "degree_north" ;
         float bath(y, x) ;
 }

3) Remapping on ORCA2

Generate weights with bilinear interpolation, and then remapping for a scalar field:

    cdo genbil,grid_ORCA2_T.nc $file bil_orca2_weights.nc
    cdo remap,grid_ORCA2_T.nc,bil_orca2_weights.nc $file_fill.nc $file_orca2.nc

where

  • $file_fill.nc is a file whose missing values is filled (see below interpolation_loop.sh)
  • $file_orca2.nc is the output
  • and bil_orca2_weights.nc is a weight file generated (see below)

Generate weights with bilinear interpolation, and then remapping for a vectorial field:

    cdo genbic,grid_ORCA2_T.nc $file bic_orca2_weights.nc
    cdo remap,grid_ORCA2_T.nc,bic_orca2_weights.nc $file_fill.nc $file_orca2.nc

where

  • $file_fill.nc is a file whose missing values is filled (see below interpolation_loop.sh)
  • $file_orca2.nc is the output
  • and bic_orca2_weights.nc is a weight file generated (see below)

Example of "interpolation_loop.sh" script

This script needs to have input directory in your execution's directory to copy *.nc files and it calls 4 others sub-scripts: "create_orca2.sh", "prepare_mask.sh", "bilweights.sh" and "bicweights.sh"

NOTE: Before running "interpolation_loop.sh", you nedd to checks different things:

  1. control latitude of your mask file. If it goes from north to south and latitude in your file goes from south to north you have to invert it:
cdo invertlat ${file_mask_land}.nc ${file_mask_land}.SN.nc
  1. control format of land variable, if it is short and your file has float format you will have problems with cdo created missing values, and filling missing values
#convert in float (from short) land variable
ncap -O -s "land=float(land)" ${file_mask_land}.SN.nc ${file_mask_land}.SN.float.nc
  1. control values of your land variable, if it is 1 on earth and 0 on ocean is not good if you are interpolating files for ocean forcing
#convert land=0 on earth, 1 on ocean
ncap -O -s 'land=1-land' mask.nc mask.nc
#rename "land" in "mask"
ncrename -v land,mask mask.nc
  1. control values of calendar's attribute in your file, if it is "NOLEAP" it is better to convert it in a recognized format by nccated (ex. 365_day)
ncatted -a calendar,${variable},m,c,365_day ${file}.nc ${file}_365.nc

main programm: interpolation_loop.sh

#!/bin/bash
#script to interpolate files on ORCA2 grids

set -u
LC_ALL=C; export LC_ALL

type cdo
status_cdo=${?}
if [ ${status_cdo} != 0 ] ; then
        echo "ERROR: cdo not found in your PATH"
        exit 1
fi

file_orca2="grid_ORCA2_T.nc"
#build ORCA2_grid if necessary
if [ ! -f ${file_orca2} ] ; then
 /bin/sh create_orca2.sh
fi

#preparing mask file if necessary
if [ ! -f mask.nc ] ; then
 /bin/sh prepare_mask.sh
fi

dir=/PATH/FILES/TO/BE/INTERPOLATED

`cp ${dir}/*.nc input/.`
cd input

for file in $(echo `ls *.nc | sed -e "s/\.nc//"`) ; do
 
 # create missing values over land with e.g. 'ifthen':
 echo "create missing values ifthen mask"
 cdo ifthen ../mask.nc ${file}_365.nc ${file}_mask.nc
 status_ifthen=${?}
 if [ ${status_ifthen} != 0 ] ; then
         echo "ERROR ifthen mask"
         exit 1
 fi
 
 
 # fill missing value (here over land):
 echo "fill missing values"
 cdo fillmiss ${file}_mask.nc ${file}_fill.nc
 status_fillmiss=${?}
 if [ ${status_fillmiss} != 0 ] ; then
  echo "ERRORE cdo fillmiss"
         exit 1
 fi
 
 
 #remapping on ORCA2 grid
 echo "remapping "
 cd ..
# looking for first part of filename, to see if it is vectorial field, for interannual files
 name=$(echo $file | awk -F \_ {'print $1'} )
 if [ ${name} = "u" -o ${name} = "v" ] ; then
    echo "vectorial field"
    /bin/sh bicweights.sh input/${file}_fill.nc
    cdo remap,grid_ORCA2_T.nc,bic_orca2_weights.nc input/${file}_fill.nc input/${file}_orca2.nc
 else 
    echo "scalar field"
    /bin/sh bilweights.sh input/${file}_fill.nc
    cdo remap,grid_ORCA2_T.nc,bil_orca2_weights.nc input/${file}_fill.nc input/${file}_orca2.nc
 fi
 status_remap=${?}
 if [ ${status_remap} != 0 ] ; then
         echo "ERROR cdo remap"
         exit 1
 fi
 echo "status remap: ${status_remap}"
 cd input
done

first sub-programm: file create_orca2.sh

#!/bin/bash
#script to create ORCA2 grid

set -u
LC_ALL=C ; export LC_ALL

echo "create orca2 grid"
ncks -O -C -a -v glamt,gphit coordinates.nc grid_ORCA2_T.nc 

#remove not necessar axes:
ncwa -O -a time,z grid_ORCA2_T.nc grid_ORCA2_T.nc

#add dummy variable:
ncap -O -s 'dummy[y,x]=1b' grid_ORCA2_T.nc grid_ORCA2_T.nc

#add needed attributes
ncatted -a coordinates,dummy,c,c,'glamt gphit'  \
        -a units,glamt,c,c,'degrees_east'       \
        -a units,gphit,c,c,'degrees_north' grid_ORCA2_T.nc

Second sub-programm: prepare_mask.sh

#!/bin/bash
#preparing mask file

set -u
LC_ALL=C ; export LC_ALL

file_mask_land="land.sfc.gauss"
# invert latitude (from SOUTH to NORTH)
cdo invertlat ${file_mask_land}.nc ${file_mask_land}.SN.nc

#convert in float (from short) land variable
ncap -O -s "land=float(land)" ${file_mask_land}.SN.nc ${file_mask_land}.SN.float.nc 

# remove all attribute to land
ncatted -a ,land,d,, ${file_mask_land}.SN.float.nc ${file_mask_land}.SN.float.miss.nc

# extract land
#ncks -O -v land land.sfc.gauss.SN.float.nc mask.nc 
ncks -O -v land ${file_mask_land}.SN.float.miss.nc mask.nc 
ncwa -O -a time mask.nc mask.nc

#convert land=0 on earth, 1 on ocean
ncap -O -s 'land=1-land' mask.nc mask.nc

#rename "land" in "mask"
ncrename -v land,mask mask.nc
chmod 755 mask.nc

Third sub-programm: bicweights.sh

#!/bin/bash

#remapping on ORCA2 grid
echo "gen bicubic WEIGHTS"
export REMAP_EXTRAPOLATE=on
cdo genbic,grid_ORCA2_T.nc $1 bic_orca2_weights.nc

Forth sub-programm: bicweights.sh

#!/bin/bash

#remapping on ORCA2 grid
echo "gen bilinear WEIGHTS"
export REMAP_EXTRAPOLATE=on
cdo genbil,grid_ORCA2_T.nc $1 bil_orca2_weights.nc