horiz_interp/horiz_interp.F90

   1 !***********************************************************************
   2 !*                   GNU Lesser General Public License
   3 !*
   4 !* This file is part of the GFDL Flexible Modeling System (FMS).
   5 !*
   6 !* FMS is free software: you can redistribute it and/or modify it under
   7 !* the terms of the GNU Lesser General Public License as published by
   8 !* the Free Software Foundation, either version 3 of the License, or (at
   9 !* your option) any later version.
  10 !*
  11 !* FMS is distributed in the hope that it will be useful, but WITHOUT
  12 !* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  13 !* FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  14 !* for more details.
  15 !*
  16 !* You should have received a copy of the GNU Lesser General Public
  17 !* License along with FMS.  If not, see <http://www.gnu.org/licenses/>.
  18 !***********************************************************************
  19 !> @defgroup horiz_interp_mod horiz_interp_mod
  20 !> @ingroup horiz_interp
  21 !> @brief Performs spatial interpolation between grids.
  22 !!
  23 !> @author Zhi Liang, Bruce Wyman
  24 !!
  25 !! This module can interpolate data from any logically rectangular grid
  26 !! to any logically rectangular grid. Four interpolation schems are used here:
  27 !! conservative, bilinear, bicubic and inverse of square distance weighted.
  28 !! The four interpolation schemes are implemented seperately in
  29 !! horiz_interp_conserver_mod, horiz_interp_blinear_mod, horiz_interp_bicubic_mod
  30 !! and horiz_interp_spherical_mod. bicubic interpolation requires the source grid
  31 !! is regular lon/lat grid. User can choose the interpolation method in the
  32 !! public interface horiz_interp_new through optional argument interp_method,
  33 !! with acceptable value "conservative", "bilinear", "bicubic" and "spherical".
  34 !! The default value is "conservative". There is an optional mask field for
  35 !! missing input data. An optional output mask field may be used in conjunction with
  36 !! the input mask to show where output data exists.
  37
  38 module horiz_interp_mod
  39
  40 !-----------------------------------------------------------------------
  41 !
  42 !        Performs spatial interpolation between grids.
  43 !
  44 !-----------------------------------------------------------------------
  45
  46 use fms_mod,                    only: write_version_number, fms_error_handler
  47 use fms_mod,                    only: check_nml_error
  48 use mpp_mod,                    only: mpp_error, FATAL, stdout, stdlog, mpp_min
  49 use mpp_mod,                    only: input_nml_file, WARNING, mpp_pe, mpp_root_pe
  50 use constants_mod,              only: pi
  51 use horiz_interp_type_mod,      only: horiz_interp_type, assignment(=)
  52 use horiz_interp_type_mod,      only: CONSERVE, BILINEAR, SPHERICA, BICUBIC
  53 use horiz_interp_conserve_mod,  only: horiz_interp_conserve_init, horiz_interp_conserve
  54 use horiz_interp_conserve_mod,  only: horiz_interp_conserve_new, horiz_interp_conserve_del
  55 use horiz_interp_bilinear_mod,  only: horiz_interp_bilinear_init, horiz_interp_bilinear
  56 use horiz_interp_bilinear_mod,  only: horiz_interp_bilinear_new, horiz_interp_bilinear_del
  57 use horiz_interp_bicubic_mod,   only: horiz_interp_bicubic_init, horiz_interp_bicubic
  58 use horiz_interp_bicubic_mod,   only: horiz_interp_bicubic_new, horiz_interp_bicubic_del
  59 use horiz_interp_spherical_mod, only: horiz_interp_spherical_init, horiz_interp_spherical
  60 use horiz_interp_spherical_mod, only: horiz_interp_spherical_new, horiz_interp_spherical_del
  61 use platform_mod,               only: r4_kind, r8_kind
  62
  63  implicit none
  64  private
  65
  66 !---- interfaces ----
  67
  68  public   horiz_interp_type, horiz_interp, horiz_interp_new, horiz_interp_del, &
  69           horiz_interp_init, horiz_interp_end, assignment(=)
  70
  71 !> Allocates space and initializes a derived-type variable
  72 !! that contains pre-computed interpolation indices and weights.
  73 !!
  74 !> Allocates space and initializes a derived-type variable
  75 !! that contains pre-computed interpolation indices and weights
  76 !! for improved performance of multiple interpolations between
  77 !! the same grids. This routine does not need to be called if you
  78 !! are doing a single grid-to-grid interpolation.
  79 !!
  80 !! @param lon_in
  81 !!      Longitude (in radians) for source data grid. You can pass 1-D lon_in to
  82 !!      represent the geographical longitude of regular lon/lat grid, or just
  83 !!      pass geographical longitude(lon_in is 2-D). The grid location may be
  84 !!      located at grid cell edge or center, decided by optional argument "grid_at_center".
  85 !!
  86 !! @param lat_in
  87 !!      Latitude (in radians) for source data grid. You can pass 1-D lat_in to
  88 !!      represent the geographical latitude of regular lon/lat grid, or just
  89 !!      pass geographical latitude(lat_in is 2-D). The grid location may be
  90 !!      located at grid cell edge or center, decided by optional argument "grid_at_center".
  91 !!
  92 !! @param lon_out
  93 !!      Longitude (in radians) for destination data grid. You can pass 1-D lon_out to
  94 !!      represent the geographical longitude of regular lon/lat grid, or just
  95 !!      pass geographical longitude(lon_out is 2-D). The grid location may be
  96 !!      located at grid cell edge or center, decided by optional argument "grid_at_center".
  97 !!
  98 !! @param lat_out
  99 !!      Latitude (in radians) for destination data grid. You can pass 1-D lat_out to
 100 !!      represent the geographical latitude of regular lon/lat grid, or just
 101 !!      pass geographical latitude(lat_out is 2-D). The grid location may be
 102 !!      located at grid cell edge or center, decided by optional argument "grid_at_center".
 103 !!
 104 !! @param verbose
 105 !!      Integer flag that controls the amount of printed output.
 106 !!      verbose = 0, no output; = 1, min,max,means; = 2, still more
 107 !!
 108 !! @param interp_method
 109 !!      interpolation method, = "conservative", using conservation scheme,
 110 !!      = "bilinear", using bilinear interpolation, = "spherical",using spherical regrid.
 111 !!      = "bicubic", using bicubic interpolation. The default value is "convervative".
 112 !!
 113 !! @param src_modulo
 114 !!      Indicate the source data grid is cyclic or not.
 115 !!
 116 !! @param grid_at_center
 117 !!      Indicate the data is on the center of grid box or the edge of grid box.
 118 !!      When true, the data is on the center of grid box. default vaule is false.
 119 !!      This option is only available when interp_method = "bilinear" or "bicubic".
 120 !!
 121 !! @param Interp
 122 !!      A derived-type variable containing indices and weights used for subsequent
 123 !!      interpolations. To reinitialize this variable for a different grid-to-grid
 124 !!      interpolation you must first use the "horiz_interp_del" interface.
 125  interface horiz_interp_new
 126     ! Source grid is 1d, destination grid is 1d
 127     module procedure horiz_interp_new_1d_r4
 128     module procedure horiz_interp_new_1d_r8
 129     ! Source grid is 1d, destination grid is 2d
 130     module procedure horiz_interp_new_1d_src_r4
 131     module procedure horiz_interp_new_1d_src_r8
 132     ! Source grid is 2d, destination grid is 2d
 133     module procedure horiz_interp_new_2d_r4
 134     module procedure horiz_interp_new_2d_r8
 135     ! Source grid is 2d, destination grid is 1d
 136     module procedure horiz_interp_new_1d_dst_r4
 137     module procedure horiz_interp_new_1d_dst_r8
 138  end interface
 139
 140
 141 !> Subroutine for performing the horizontal interpolation between two grids.
 142 !!
 143 !> Subroutine for performing the horizontal interpolation between
 144 !! two grids. There are two forms of this interface.
 145 !! Form A requires first calling horiz_interp_new, while Form B
 146 !! requires no initialization.
 147 !!
 148 !! @param Interp
 149 !!     Derived-type variable containing interpolation indices and weights.
 150 !!     Returned by a previous call to horiz_interp_new.
 151 !!
 152 !! @param data_in
 153 !!      Input data on source grid.
 154 !!
 155 !! @param verbose
 156 !!      flag for the amount of print output.
 157 !!               verbose = 0, no output; = 1, min,max,means; = 2, still more
 158 !!
 159 !! @param mask_in
 160 !!      Input mask, must be the same size as the input data. The real value of
 161 !!      mask_in must be in the range (0.,1.). Set mask_in=0.0 for data points
 162 !!      that should not be used or have missing data. It is Not needed for
 163 !!      spherical regrid.
 164 !!
 165 !! @param missing_value
 166 !!      Use the missing_value to indicate missing data.
 167 !!
 168 !! @param missing_permit
 169 !!      numbers of points allowed to miss for the bilinear interpolation. The value
 170 !!      should be between 0 and 3.
 171 !!
 172 !! @param lon_in, lat_in
 173 !!      longitude and latitude (in radians) of source grid. More explanation can
 174 !!      be found in the documentation of horiz_interp_new.
 175 !!
 176 !! @param lon_out, lat_out
 177 !!      longitude and latitude (in radians) of destination grid. More explanation can
 178 !!      be found in the documentation of horiz_interp_new.
 179 !!
 180 !! @param data_out
 181 !!      Output data on destination grid.
 182 !!
 183 !! @param mask_out
 184 !!      Output mask that specifies whether data was computed.
 185 !!
 186 !!
 187 !! @throws FATAL, size of input array incorrect
 188 !!      The input data array does not match the size of the input grid edges
 189 !!      specified. If you are using the initialization interface make sure you
 190 !!      have the correct grid size.
 191 !!
 192 !! @throws FATAL, size of output array incorrect
 193 !!      The output data array does not match the size of the input grid
 194 !!      edges specified. If you are using the initialization interface make
 195 !!      sure you have the correct grid size.
 196 !> @ingroup horiz_interp_mod
 197  interface horiz_interp
 198     module procedure horiz_interp_base_2d_r4
 199     module procedure horiz_interp_base_2d_r8
 200     module procedure horiz_interp_base_3d_r4
 201     module procedure horiz_interp_base_3d_r8
 202     module procedure horiz_interp_solo_1d_r4
 203     module procedure horiz_interp_solo_1d_r8
 204     module procedure horiz_interp_solo_1d_src_r4
 205     module procedure horiz_interp_solo_1d_src_r8
 206     module procedure horiz_interp_solo_2d_r4
 207     module procedure horiz_interp_solo_2d_r8
 208     module procedure horiz_interp_solo_1d_dst_r4
 209     module procedure horiz_interp_solo_1d_dst_r8
 210     module procedure horiz_interp_solo_old_r4
 211     module procedure horiz_interp_solo_old_r8
 212  end interface
 213
 214 !> Private helper routines
 215 interface is_lat_lon
 216     module procedure is_lat_lon_r4
 217     module procedure is_lat_lon_r8
 218 end interface
 219
 220 interface horiz_interp_solo_1d
 221   module procedure horiz_interp_solo_1d_r4
 222   module procedure horiz_interp_solo_1d_r8
 223 end interface
 224
 225
 226 !> @addtogroup horiz_interp_mod
 227 !> @{
 228
 229  logical :: reproduce_siena = .false. !< Set reproduce_siena = .true. to reproduce siena results.
 230                  !! Set reproduce_siena = .false. to decrease truncation error
 231                  !! in routine poly_area in file mosaic_util.c. The truncation error of
 232                  !! second order conservative remapping might be big for high resolution
 233                  !! grid.
 234
 235  namelist /horiz_interp_nml/ reproduce_siena
 236
 237 !-----------------------------------------------------------------------
 238 ! Include variable "version" to be written to log file.
 239 #include<file_version.h>
 240  logical            :: module_is_initialized = .FALSE.
 241 !-----------------------------------------------------------------------
 242
 243 contains
 244
 245 !#######################################################################
 246
 247   !> Initialize module and writes version number to logfile.out
 248   subroutine horiz_interp_init
 249   integer :: unit, ierr, io
 250
 251   if(module_is_initialized) return
 252   call write_version_number("HORIZ_INTERP_MOD", version)
 253
 254   read (input_nml_file, horiz_interp_nml, iostat=io)
 255   ierr = check_nml_error(io,'horiz_interp_nml')
 256   if (mpp_pe() == mpp_root_pe() ) then
 257      unit = stdlog()
 258      write (unit, nml=horiz_interp_nml)
 259   endif
 260
 261   if (reproduce_siena) then
 262     call mpp_error(FATAL, "horiz_interp_mod: You have overridden the default value of " // &
 263        "reproduce_siena and set it to .true. in horiz_interp_nml. This was a temporary workaround to " // &
 264        "allow for consistency in continuing experiments and is no longer supported. " // &
 265        "Please remove this namelist.")
 266   endif
 267
 268   call horiz_interp_conserve_init
 269   call horiz_interp_bilinear_init
 270   call horiz_interp_bicubic_init
 271   call horiz_interp_spherical_init
 272
 273   module_is_initialized = .true.
 274
 275   end subroutine horiz_interp_init
 276
 277 !> Deallocates memory used by "horiz_interp_type" variables.
 278 !! Must be called before reinitializing with horiz_interp_new.
 279  subroutine horiz_interp_del ( Interp )
 280
 281    type (horiz_interp_type), intent(inout) :: Interp !< A derived-type variable returned by previous
 282                                            !! call to horiz_interp_new. The input variable must have
 283                                            !! allocated arrays. The returned variable will contain
 284                                            !! deallocated arrays
 285
 286 !-----------------------------------------------------------------------
 287 !  releases space used by horiz_interp_type variables
 288 !  must be called before re-initializing the same variable
 289 !-----------------------------------------------------------------------
 290    select case(Interp % interp_method)
 291    case (CONSERVE)
 292       call horiz_interp_conserve_del(Interp )
 293    case (BILINEAR)
 294       call horiz_interp_bilinear_del(Interp )
 295    case (BICUBIC)
 296       call horiz_interp_bicubic_del(Interp )
 297    case (SPHERICA)
 298       call horiz_interp_spherical_del(Interp )
 299    end select
 300
 301    Interp%I_am_initialized = .false.
 302 !-----------------------------------------------------------------------
 303
 304  end subroutine horiz_interp_del
 305
 306  !#####################################################################
 307
 308  !> Dummy routine
 309  subroutine horiz_interp_end
 310  return
 311  end subroutine horiz_interp_end
 312
 313 #include "horiz_interp_r4.fh"
 314 #include "horiz_interp_r8.fh"
 315
 316 end module horiz_interp_mod
 317 !> @}
 318 ! close documentation grouping