The interpolation

The following transformations, controlled by interp.f, are available in OASIS3.

• BLASOLD:

  BLASOLD (routine blasold.f) performs a linear combination of the current coupling field with other coupling fields or with a constant before the interpolation per se.

  This transformation requires at least one configuring line with two parameters:

   # BLASOLD operation 
       $XMULT     $NBFIELDS
  

  where $XMULT is the multiplicative coefficient of the current field, and $NBFIELDS the number of additional fields to be combined with the current field. For each additional field, an additional input line is required:

  # nbfields lines
       $CNAME  $AMULT
  

  where $CNAME and $AMULT are the symbolic name and the multiplicative coefficient for the additional field. To add a constant value to the original field, put $XMULT = 1, $NBFIELDS = 1, $CNAME = CONSTANT, $AMULT = value to add.

• SCRIPR:

  SCRIPR gathers the interpolation techniques offered by Los Alamos National Laboratory SCRIP 1.4 library¹⁵(1). SCRIPR routines are in oasis3/lib/scrip. See the SCRIP 1.4 documentation in oasis3/doc/SCRIPusers.pdf for more details on the interpolation algorithms.

  The following types of interpolations are available:

  • DISTWGT performs a distance weighted nearest-neighbour interpolation (N neighbours). All types of grids are supported.
    • Masked target grid points: the zero value is associated to masked target grid points.

    • Non-masked target grid points having some of the N source nearest neighbours masked: a nearest neighbour algorithm using the remaining non masked source nearest neighbours is applied.

    • Non-masked target grid points having all of the N source nearest neighbours masked: by default, the nearest non-masked source neighbour is used. Warning: this default behaviour has changed since release oasis3_prism_2_5¹⁶.

    The configuring line is:

     # SCRIPR/DISWGT 
      $CMETH $CGRS $CFTYP $REST $NBIN $NV $ASSCMP $PROJCART
    

    where:
    • $CMETH = DISTWGT.
    • $CGRS is the source grid type (LR, D or U)- see annexe A.

    • $CFTYP is the field type: SCALAR if the field is a scalar one, or VECTOR_I or VECTOR_J whether the field represents respectively the first or the second component of a vector field (see paragraph Support of vector fields below). The option VECTOR, which is fact leads to a scalar treatment of the field (as in the previous versions), is still accepted.

    • $REST is the search restriction type: LATLON or LATITUDE (see SCRIP 1.4 documentation SCRIPusers.pdf). Note that for D or U grid, the restriction may influence sligthly the result near the borders of the restriction bins. (XXX to be checked)
    • $NBIN the number of restriction bins (see SCRIP 1.4 documentation SCRIPusers.pdf).
    • $NV is the number of neighbours used.
    • $ASSCMP: optional, for VECTOR_I or VECTOR_J vector fields only; the source symbolic name of the associated vector component.
    • $PROJCART: optional, for vector fields only; should be PROJCART if the user wants the vector components to be projected in a Cartesian coordinate system before interpolation (see paragraph Support of vector fields below).

  • GAUSWGT performs a N nearest-neighbour interpolation weighted by their distance and a gaussian function. All grid types are supported.
    • Masked target grid points: the zero value is associated to masked target grid points.

    • Non-masked target grid points having some of the N source nearest neighbours masked: a nearest neighbour algorithm using the remaining non masked source nearest neighbours is applied.

    • Non-masked target grid points having their N nearest neighbours all masked: the zero value will be associated to these target points. The value 1.0E+20 will however be assigned to these non-masked target grid points if routines scriprmp.f or vector.F90 (for vector interpolation) in oasis3/lib/scrip/src/ are compiled with
      ll_weightot=.true.. Furthermore, if oasis3/lib/scrip/src/remap_gauswgt.f is compiled with ll_nnei=.true., the non-masked nearest neighbour value will be given to these target grid points.

    The configuring line is:

     # SCRIPR/GAUSWGT
    $CMETH  $CGRS  $CFTYP  $REST  $NBIN  $NV $VAR $ASSCMP $PROJCART
    

    where all entries are as for DISTWGT, except that:
    • $CMETH = GAUSWGT
    • $VAR, which must be given as a REAL value (e.g 2.0 and not 2), defines the weight given to a neighbour source grid point as proportional to $exp(-1/2 \cdot d^2/\sigma^2)$ where $d$ is the distance between the source and target grid points, and $\sigma^2 = \$VAR \cdot \overline{d}^2$ where $\overline{d}^2$ is the average distance between two source grid points (calculated automatically by OASIS3).

  • BILINEAR performs an interpolation based on a local bilinear approximation (see details in chapter 4 of SCRIP 1.4 documentation SCRIPusers.pdf)

  • BICUBIC performs an interpolation based on a local bicubic approximation (see details in chapter 5 of SCRIP 1.4 documentation SCRIPusers.pdf)

    For BILINEAR and BICUBIC, Logically-Rectangular (LR) and Reduced (D) source grid types are supported.

    • Masked target grid points: the zero value is associated to masked target grid points.

    • Non-masked target grid points having some of the source points normally used in the bilinear or bicubic interpolation masked: a N nearest neighbour algorithm using the remaining non masked source points is applied.

    • Non-masked target grid points having their N nearest neighbours all masked: the zero value will be associated to these target points. The value 1.0E+20 will however be assigned to these non-masked target grid points if routines scriprmp.f or vector.F90 (for vector interpolation) in oasis3/lib/scrip/src/ are compiled with
      ll_weightot=.true.. Furthermore, if oasis3/lib/scrip/src/remap_bicubic.f or remap_bilinear.f is compiled with ll_nnei=.true., the non-masked nearest neighbour value will be given to these target grid points.

    The configuring line is:

     # SCRIPR/BILINEAR or SCRIPR/BICUBIC
         $CMETH  $CGRS  $CFTYP  $REST  $NBIN $ASSCMP $PROJCART
    

    where:
    • $CMETH = BILINEAR or BICUBIC
    • $CGRS is the source grid type (LR or D)
    • $CFTYP, $NBIN, $ASSCMP $PROJCART are as for DISTWGT.
    • $REST is as for DISTWGT, except that only LATITUDE is possible for a Reduced (D) source grid.

  • CONSERV performs 1st or 2nd order conservative remapping, which means that the weight of a source cell is proportional to area intersected by the target cell.

    The configuring line is:

     # SCRIPR/CONSERV
    $CMETH  $CGRS  $CFTYP  $REST  $NBIN  $NORM  $ORDER $ASSCMP $PROJCART
    

    where:
    • $CMETH = CONSERV
    • $CGRS is the source grid type: LR, D and U are supported for first-order remapping if the grid corners are given by the user in the grid data file grids.nc ; only LR is supported if the grid corners are not available in grids.nc and therefore have to be calculated automatically by OASIS3. For second-order remapping, only LR is supported because the gradient of the coupling field used in the transformation has to be calculated automatically by OASIS3.
    • $CFTYP, $REST, $NBIN, $ASSCMP,and $PROJCART are as for DISTWGT.
    • $NORM is the NORMalization option:
      • FRACAREA: The sum of the non-masked source cell intersected areas is used to NORMalise each target cell field value: the flux is not locally conserved, but the flux value itself is reasonable.
      • DESTAREA: The total target cell area is used to NORMalise each target cell field value even if it only partly intersects non-masked source grid cells: local flux conservation is ensured, but unreasonable flux values may result.
      • FRACNNEI: as FRACAREA, except that at least the source nearest unmasked neighbour is used for unmasked target cells that intersect only masked source cells. Note that a zero value will be assigned to a target cell that does not intersect any source cells (masked or unmasked), even with FRACNNEI option.
    • $ORDER: FIRST or SECOND for first or second order remapping respectively (see SCRIP 1.4 documentation). Note that CONSERV/SECOND is not positive definite and has not been fully validated yet.

  Precautions related to the use of the SCRIPR/CONSERV remapping in particular

  • For the 1st order conservative remapping: the weight of a source cell is proportional to area of the source cell intersected by target cell. Using the divergence theorem, the SCRIP library evaluates this area with the line integral along the cell borders enclosing the area. As the real shape of the borders is not known (only the location of the 4 corners of each cell is known), the library assumes that the borders are linear in latitude and longitude between two corners. In general, this assumption becomes less valid closer to the pole and for latitudes above the north_thresh or below the south_thresh values specified in oasis3/lib/scrip/remap_conserv.F, the library evaluates the intersection between two border segments using a Lambert equivalent azimuthal projection. Problems have been observed in some cases for the grid cell located around this north_thresh or south_thresh latitude.

  • Another limitation of the SCRIP 1st order conservative remapping algorithm is that is also supposes, for line integral calculation, that $sin(latitude)$ is linear in longitude on the cell borders which again is in general not valid close to the pole.

  • For a proper consevative remapping, the corners of a cell have to coincide with the corners of its neighbour cell.

  • If two cells of a grid overlay, at least the one with the greater numerical index must be masked (they also can be both masked) for a proper conservative remapping. For example, if the grid line with i=1 overlaps the grid line with i=imax, it is the latter that must be masked. When this is not the case with the mask defined in masks.nc, OASIS must be compiled with the CPP key TREAT_OVERLAY which will ensure that these rules are respected. This CPP key was introduced in oasis3_3 version.

  • A target grid cell intersecting no source cell (either masked or non masked) at all i.e. falling in a ``hole'' of the source grid will in all cases get a zero value.

  • If a target grid cell intersects only masked source cells, it will still get a zero value unless:

    - the FRACNNEI normalisation option is used, in which case it will get the nearest non masked neighbour value, or

    - the routines oasis3/lib/scrip/src/scriprmp.f or vector.F90 - for vector interpolation - are compiled with ll_weightot=.true. in which case, the value 1.0E+20 will be assigned to these target grid cell intersecting only masked source cells (for easier identification).

  Precautions related to the use of the SCRIPR remappings in general

  • For using SCRIPR interpolations, linking with the NetCDF library is mandatory and the grid data files (see section 7.2) must be NetCDF files (binary files are not supported).

  • When the SCRIP library performs a remapping, it first checks if the file containing the corresponding remapping weights and addresses exists. If it exists, it reads them from the file; if not, it calculates them and store them in a file. The file is created in the working directory and is called rmp_srcg_to_tgtg_INTTYPE_NORMAOPT.nc, where srcg and tgtg are the acronyms of respetively the source and the target grids, INTTYPE is the interpolation type (i.e. DISTWGT, GAUSWGT, BILINEA, BICUBIC, or CONSERV) and NORMAOPT is the normalization option (i.e. DESTAREA, FRACAREA or FRACNNEI for CONSERV only). The problem comes from the fact that the weights and addresses will also differ whether or not the MASK and EXTRAP transformations are first activated during the pre-processing phase (see section 6.3) and this option is not stored in the remapping file name. Therefore, the remapping file used will be the one created for the first field having the same source grid, target grid, and interpolation type (and the same normalization option for CONSERV), even if the MASK and EXTRAP transformations are used or not for that field. (This inconsistency is however usually not a problem as the MASK and EXTRAP transformations are usually used for all fields having the same source grid, target grid, and interpolation type, or not at all.)

  Support of vector fields with the SCRIPR remappings

  SCRIPR supports 2D vector interpolation. Note however that this functionality has been validated only in a reduced number of test cases. The two vector components have to be identified by assigning VECTOR_I or VECTOR_J to $CFTYP and have to be associated by giving, for each component field, the source symbolic name of the associated vector component to $ASSCMP (see above). The grids of the two vector components can be different but have to have the same number of points, the same overlap, the same mask; the same interpolation must be used for the two components. A proper example of vector interpolation is given in the interpolator-only mode example testinterp (see section 8.4.1). The details of the vector treatment, performed by the routines
  scriprmp_vector.F90 and rotations.F90 in oasis3/lib/scrip/src are the following:

  • If the angles of the source grid local coordinate system are defined in the grids.nc data file (see section 7.2), an automatic rotation from the local to the geographic spherical coordinate system is performed.
  • If the two source vector components are not defined on the same source grid, one component is automatically interpolated on the grid of the other component. Note that if the components are not given in a Cartesian coordinate system, this operation is not exact as the coordinate system is not fixed spatially (for grids not covering the North or South poles, the error is however small at the scale of the grid mesh).
  • If the user puts the PROJCART keyword at the end of the SCRIPR configuring line (see above), projection of the two vector components in a Cartesian coordinate system, interpolation of the resulting 3 Cartesian components, and projection back in the spherical coordinate system are performed. In debug mode (compilation with __DEBUG pre-compiling key), the resulting vertical component in the spherical coordinate system after interpolation is written to a file projection.nc; as the source vector is horizontal, this component should be very close to 0.
  • If the user did not put the PROJCART keyword at the end of the SCRIPR configuring line, the two spherical coordinate system components are interpolated. Note that this operation is not exact as the coordinate system is not fixed spatially (in most cases, the error is however small at the scale of the grid mesh).
  • If the angles of the target grid local coordinate system are defined in the grids.nc data file (see section 7.2), an automatic rotation from the geographic spherical to the local coordinate system is performed.
  • The first and second components of the interpolated vector field are then present in the target fields associated respectively to the first and second source vector component. The target grids for the two vector components can be different.

  
  
  
  

• INTERP:

  INTERP gathers different techniques of interpolation controlled by routine fiasco.f. The following interpolations are available:

  • BILINEAR performs a bilinear interpolation using 4 neighbours.

  • BICUBIC performs a bicubic interpolation.

  • NNEIBOR performs a nearest-neighbour interpolation.

    These three interpolations are performed by routines in /oasis3/lib/fscint and support only A, B, G, L, Y, or Z grids (see appendix A). All sources grid points, masked or not, are used in the calculation. To avoid the `contamination' by masked source grid points, transformations MASK and EXTRAP should be used. Values are calculated for all target grid points, masked or not.

    The configuring line is as follows:

     # BILINEAR or BICUBIC or NNEIBOR interpolation 
        $CMETH $CGRS $CFTYP
    

    where
    • $CMETH = BILINEAR, BICUBIC or NNEIBOR
    • $CGRS is the source grid type (A, B, G, L, Y, or Z, see appendix A)
    • $CFTYP the field type (SCALAR or VECTOR). VECTOR has an effect for target grid points located near the pole: the sign of a source value located on the other side of the pole will be reversed.

  • SURFMESH (routines in /oasis3/lib/anaism) is a first-order conservative remapping from a fine to a coarse grid (the source grid must be finer over the whole domain) and supports only Lat-Lon grids (see appendix A). For a target grid cell, all the underlying not masked source grid cells are found and the target grid field value is the sum of the source grid field values weighted by the overlapped surfaces. No value is assigned to masked cells. Note that it is not recommended to use this interpolation anymore, as the more general SCRIPR/CONSERV remapping is now available. The configuring line is as follows:

     # SURFMESH remapping
         $CMETH $CGRS $CFTYP $NID $NV $NIO
    

    where
    • $CMETH = SURFMESH
    • $CGRS and $CFTYP are as for BILINEAR
    • $NID is the identificator for the weight-address dataset in all the different INTERP/SURFMESH datasets in the present coupling. This dataset will be calculated by OASIS3 if $NIO= 1, or only read if $NIO= 0.
    • $NV is the maximum number of source grid meshes used in the remapping.

  • GAUSSIAN (routines in /oasis3/lib/anaisg) is a gaussian weighted nearest-neighbour interpolation technique. The user can choose the variance of the function and the number of neighbours considered. The $NV non masked nearest source grid points are automatically used and no value are calculated for masked target grid points.

    The configuring line is:

        # GAUSSIAN interpolation
        $CMETH $CGRS $CFTYP $NID $NV $VAR  $NIO
    

    where
    • $CMETH = GAUSSIAN
    • $CGRS is the source grid type; all grids are supported¹⁷
    • $CFTYP is the field type SCALAR or VECTOR.
    • $NID is the identificator for the weight-address dataset in all the different INTERP/GAUSSIAN datasets in the present coupling. This weight-address dataset will be calculated by OASIS3 if $NIO= 1, or only read if $NIO= 0.
    • $NV is the number of neighbours used in the interpolation.
    • $VAR is as for SCRIPR/GAUSWGT (see above).

• MOZAIC:

  MOZAIC performs the mapping of a field from a source to a target grid. The grid-mapping dataset, i.e. the weights and addresses of the source grid points used to calculate the value of each target grid point are defined by the user in a file (see section 7.5). The configuring line is:

  # MOZAIC operation
       $CFILE  $NUMLU  $NID  $NV
  

  where
  • $CFILE and $NUMLU are the grid-mapping file name and associated logical unit on which the grid-mapping dataset is going to be read),
  • $NID the identificator for this grid-mapping dataset in all MOZAIC grid-mapping datasets in the present coupling
  • $NV is the maximum number of target grid points use in the mapping.

• NOINTERP:

  NOINTERP is the analysis that has to be chosen when no other transformation from the interpolation class is chosen. There is no configuring line.

• FILLING:

  FILLING (routine oasis3/src/filling.f) performs the blending of a regional data set with a climatological global one for a Sea Surface Temperature (SST) or a Sea Ice Extent field. This occurs when coupling a non-global ocean model with a global atmospheric model. FILLING can only handle fields on Logically Rectangular grid (LR, but also A, B, G, L, Y, and Z grids, see section A.

  The global data set has to be a set of SST given in Celsius degrees (for the filling of a Sea Ice Extent field, the presence or absence of ice is deduced from the value of the SST). The frequency of the global set can be interannual monthly, climatological monthly or yearly.

  The blending can be smooth or abrupt. If the blending is abrupt, only model values are used within the model domain, and only the global data set values are used outside. If the blending is smooth, a linear interpolation is performed between the two fields within the model domain over narrow bands along the boundaries. The linear interpolation can also be performed giving a different weight to the regional or and global fields.

  The smoothing is defined by parameters in oasis3/src/mod_smooth.F90. The lower smoothing band in the global model second dimension is defined by nsltb (outermost point) and nslte (innermost point); the upper smoothing band in the global model second dimension is defined by nnltb (outermost point) and nnlte (innermost point). The parameter qalfa controls the weights given to the regional and to the global fields in the linear interpolation. qalfa has to be $1/(nslte-nsltb)$ or $1/(nnltb-nnlte)$. For the outermost points (nsltb or nnltb) in the smoothing band, the weight given to the regional and global fields will respectively be 0 and 1; for the innermost points (nslte or nnlte) in the smoothing band, the weight given to the regional and global fields will respectively be 1 and 0; within the smoothing band, the weights will be a linear interpolation of the outermost and innermost weights.

  The smoothing band in the global model first dimension will be a band of nliss points following the coastline. To calculate this band, OASIS3 needs nwlgmx, the greater first dimension index of the lower coastline and nelgmx, the smaller first dimension index on the upper coastline. The parameter qbeta controls the weights given to the regional and to the global fields in the linear interpolation. qbeta has to be $1/(nliss-1)$. The weights given to the regional and global fields in the global model first dimension smoothing bands will be calculated as for the second dimension.

  The user must provide the climatological data file with a specific format described in 7.5. When one uses FILLING for SST with smooth blending, thermodynamics consistency requires to modify the heat fluxes over the blending regions. The correction term is proportional to the difference between the blended SST and the original SST interpolated on the atmospheric grid and can be written out on a file to be read later, for analysis CORRECT for example. In that case, the symbolic name of the flux correction term read through the input file namcouple must correspond in FILLING and CORRECT analyses.

  In case the regional ocean model includes a coastal part or islands, a sea-land mask mismatch may occur and a coastal point correction can be performed if the field has been previously interpolated with INTER/SURFMESH. In fact, the mismatch could result in the atmosphere undesirably ``seeing'' climatological SST's directly adjacent to ocean model SST's. Where this situation arises, the coastal correction consists in identifying the suitable ocean model grid points that can be used to extrapolate the field, excluding climatological grid points.

  This analysis requires one configuring line with 3, 4 or 6 arguments.

  1. If FILLING performs the blending of a regional data set with a global one for the Sea Ice Extent, the 3-argument input line is:

     # Sea Ice Extent FILLING operation
          $CFILE  $NUMLU  $CMETH
     

     where $CFILE is the file name for the global data set, $NUMLU the associated logical unit. $CMETH, the FILLING technique, is a CHARACTER*8 variable: the first 3 characters are either SMO, smooth filling, or RAW, no smoothing ; the next three characters must be SIE for a Sea Ice Extent filling operation; the last two define the time characteristics of the global data file, respectively MO, SE and AN for interannual monthly, climatological monthly and yearly. Note that in all cases, the global data file has to be a Sea Surface Temperature field in Celsius degrees.

  2. If FILLING performs the blending of a regional data set with a global one for the Sea Surface Temperature without any smoothing, the 4-argument input line is:

     #Sea Surface Temperature FILLING operation without smoothing
          $CFILE  $NUMLU  $CMETH  $NFCOAST
     

     where $CFILE, $NUMLU are as for the SIE filling. In this case however, $CMETH(1:3) = RAW, $CMETH(4:6) = SST, and the last two characters define the time characteristics of the global data file, as for the SIE filling. $NFCOAST is the flag for the calculation of the coastal correction ( 0 no, 1 yes).

  3. If FILLING performs the blending of a regional data set with a global one for the Sea Surface Temperature with smoothing, the 6-argument input line is:

     #Sea Surface Temperature FILLING operation with smoothing
          $CFILE  $NUMLU  $CMETH  $NFCOAST  $CNAME  $NUNIT
     

     where $CFILE, $NUMLU and $NFCOAST are as for the SST filling without smoothing. In this case, $CMETH(1:3) = SMO, $CMETH(4:6) = SST, and the last two characters define the time characteristics of the global data file, as for the SIE filling. $CNAME is the symbolic name for the correction term that is calculated by OASIS3 and $NUNIT the logical unit on which it is going to be written.

Footnotes

... library¹⁵

See the copyright statement in appendix 3.2.2.

...oasis3_prism_2_5¹⁶

To reproduce the previous default behaviour, one has to compile with CPP key NOT_NNEIGHBOUR. In this case, the zero value is associated to the target points having all of the N source nearest neighbours masked. However, the value 1.0E+20 will be assigned to these non-masked target grid points if routines scriprmp.f or vector.F90 (for vector interpolation) in oasis3/lib/scrip/src/ are compiled with ll_weightot=.true..

... supported¹⁷

Note that for ``U" grids, the average distance between two source grid point $\overline{d}^2$, used in the calculation of the Gaussian weighting function, is not exact but this has only a limited impact as this value is in all cases multiplied by the $VAR value defined by the user (see SCRIPR/GAUSWGT above).