FMS  2025.04
Flexible Modeling System
horiz_interp.F90
1 !***********************************************************************
2 !* Apache License 2.0
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* Licensed under the Apache License, Version 2.0 (the "License");
7 !* you may not use this file except in compliance with the License.
8 !* You may obtain a copy of the License at
9 !*
10 !* http://www.apache.org/licenses/LICENSE-2.0
11 !*
12 !* FMS is distributed in the hope that it will be useful, but WITHOUT
13 !* WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied;
14 !* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
15 !* PARTICULAR PURPOSE. See the License for the specific language
16 !* governing permissions and limitations under the License.
17 !***********************************************************************
18 !> @defgroup horiz_interp_mod horiz_interp_mod
19 !> @ingroup horiz_interp
20 !> @brief Performs spatial interpolation between grids.
21 !!
22 !> @author Zhi Liang, Bruce Wyman
23 !!
24 !! This module can interpolate data from any logically rectangular grid
25 !! to any logically rectangular grid. Four interpolation schems are used here:
26 !! conservative, bilinear, bicubic and inverse of square distance weighted.
27 !! The four interpolation schemes are implemented separately in
28 !! horiz_interp_conserver_mod, horiz_interp_blinear_mod, horiz_interp_bicubic_mod
29 !! and horiz_interp_spherical_mod. bicubic interpolation requires the source grid
30 !! is regular lon/lat grid. User can choose the interpolation method in the
31 !! public interface horiz_interp_new through optional argument interp_method,
32 !! with acceptable value "conservative", "bilinear", "bicubic" and "spherical".
33 !! The default value is "conservative". There is an optional mask field for
34 !! missing input data. An optional output mask field may be used in conjunction with
35 !! the input mask to show where output data exists.
36 
37 module horiz_interp_mod
38 
39 !-----------------------------------------------------------------------
40 !
41 ! Performs spatial interpolation between grids.
42 !
43 !-----------------------------------------------------------------------
44 
45 use fms_mod, only: write_version_number, fms_error_handler
46 use fms_mod, only: check_nml_error
47 use mpp_mod, only: mpp_error, fatal, stdout, stdlog, mpp_min
48 use mpp_mod, only: input_nml_file, warning, mpp_pe, mpp_root_pe
49 use constants_mod, only: pi
50 use horiz_interp_type_mod, only: horiz_interp_type, assignment(=)
51 use horiz_interp_type_mod, only: conserve, bilinear, spherical, bicubic
52 use horiz_interp_conserve_mod, only: horiz_interp_conserve_init, horiz_interp_conserve
53 use horiz_interp_conserve_mod, only: horiz_interp_conserve_new, horiz_interp_conserve_del
54 use horiz_interp_bilinear_mod, only: horiz_interp_bilinear_init, horiz_interp_bilinear
55 use horiz_interp_bilinear_mod, only: horiz_interp_bilinear_new, horiz_interp_bilinear_del
56 use horiz_interp_bilinear_mod, only: horiz_interp_read_weights_bilinear
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 
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.
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  !> Subroutines for reading in weight files and using that to fill in the horiz_interp type instead
141  !! calculating it
143  module procedure horiz_interp_read_weights_r4
144  module procedure horiz_interp_read_weights_r8
145  end interface horiz_interp_read_weights
146 
147 !> Subroutine for performing the horizontal interpolation between two grids.
148 !!
149 !> Subroutine for performing the horizontal interpolation between
150 !! two grids. There are two forms of this interface.
151 !! Form A requires first calling horiz_interp_new, while Form B
152 !! requires no initialization.
153 !!
154 !! @param Interp
155 !! Derived-type variable containing interpolation indices and weights.
156 !! Returned by a previous call to horiz_interp_new.
157 !!
158 !! @param data_in
159 !! Input data on source grid.
160 !!
161 !! @param verbose
162 !! flag for the amount of print output.
163 !! verbose = 0, no output; = 1, min,max,means; = 2, still more
164 !!
165 !! @param mask_in
166 !! Input mask, must be the same size as the input data. The real value of
167 !! mask_in must be in the range (0.,1.). Set mask_in=0.0 for data points
168 !! that should not be used or have missing data. It is Not needed for
169 !! spherical regrid.
170 !!
171 !! @param missing_value
172 !! Use the missing_value to indicate missing data.
173 !!
174 !! @param missing_permit
175 !! numbers of points allowed to miss for the bilinear interpolation. The value
176 !! should be between 0 and 3.
177 !!
178 !! @param lon_in, lat_in
179 !! longitude and latitude (in radians) of source grid. More explanation can
180 !! be found in the documentation of horiz_interp_new.
181 !!
182 !! @param lon_out, lat_out
183 !! longitude and latitude (in radians) of destination grid. More explanation can
184 !! be found in the documentation of horiz_interp_new.
185 !!
186 !! @param data_out
187 !! Output data on destination grid.
188 !!
189 !! @param mask_out
190 !! Output mask that specifies whether data was computed.
191 !!
192 !!
193 !! @throws FATAL, size of input array incorrect
194 !! The input data array does not match the size of the input grid edges
195 !! specified. If you are using the initialization interface make sure you
196 !! have the correct grid size.
197 !!
198 !! @throws FATAL, size of output array incorrect
199 !! The output data array does not match the size of the input grid
200 !! edges specified. If you are using the initialization interface make
201 !! sure you have the correct grid size.
202 !> @ingroup horiz_interp_mod
203  interface horiz_interp
204  module procedure horiz_interp_base_2d_r4
205  module procedure horiz_interp_base_2d_r8
206  module procedure horiz_interp_base_3d_r4
207  module procedure horiz_interp_base_3d_r8
208  module procedure horiz_interp_solo_1d_r4
209  module procedure horiz_interp_solo_1d_r8
210  module procedure horiz_interp_solo_1d_src_r4
211  module procedure horiz_interp_solo_1d_src_r8
212  module procedure horiz_interp_solo_2d_r4
213  module procedure horiz_interp_solo_2d_r8
214  module procedure horiz_interp_solo_1d_dst_r4
215  module procedure horiz_interp_solo_1d_dst_r8
216  module procedure horiz_interp_solo_old_r4
217  module procedure horiz_interp_solo_old_r8
218  end interface
219 
220 !> Private helper routines
221 interface is_lat_lon
222  module procedure is_lat_lon_r4
223  module procedure is_lat_lon_r8
224 end interface
225 
227  module procedure horiz_interp_solo_1d_r4
228  module procedure horiz_interp_solo_1d_r8
229 end interface
230 
231 
232 !> @addtogroup horiz_interp_mod
233 !> @{
234 
235  logical :: reproduce_siena = .false. !< Set reproduce_siena = .true. to reproduce siena results.
236  !! Set reproduce_siena = .false. to decrease truncation error
237  !! in routine poly_area in file mosaic_util.c. The truncation error of
238  !! second order conservative remapping might be big for high resolution
239  !! grid.
240 
241  namelist /horiz_interp_nml/ reproduce_siena
242 
243 !-----------------------------------------------------------------------
244 ! Include variable "version" to be written to log file.
245 #include<file_version.h>
246  logical :: module_is_initialized = .false.
247 !-----------------------------------------------------------------------
248 
249 contains
250 
251 !#######################################################################
252 
253  !> Initialize module and writes version number to logfile.out
254  subroutine horiz_interp_init
255  integer :: iunit, ierr, io
256 
257  if(module_is_initialized) return
258  call write_version_number("HORIZ_INTERP_MOD", version)
259 
260  read (input_nml_file, horiz_interp_nml, iostat=io)
261  ierr = check_nml_error(io,'horiz_interp_nml')
262  if (mpp_pe() == mpp_root_pe() ) then
263  iunit = stdlog()
264  write (iunit, nml=horiz_interp_nml)
265  endif
266 
267  if (reproduce_siena) then
268  call mpp_error(fatal, "horiz_interp_mod: You have overridden the default value of " // &
269  "reproduce_siena and set it to .true. in horiz_interp_nml. This was a temporary workaround to " // &
270  "allow for consistency in continuing experiments and is no longer supported. " // &
271  "Please remove this namelist.")
272  endif
273 
274  call horiz_interp_conserve_init
275  call horiz_interp_bilinear_init
276  call horiz_interp_bicubic_init
277  call horiz_interp_spherical_init
278 
279  module_is_initialized = .true.
280 
281  end subroutine horiz_interp_init
282 
283 !> Deallocates memory used by "horiz_interp_type" variables.
284 !! Must be called before reinitializing with horiz_interp_new.
285  subroutine horiz_interp_del ( Interp )
286 
287  type (horiz_interp_type), intent(inout) :: interp !< A derived-type variable returned by previous
288  !! call to horiz_interp_new. The input variable must have
289  !! allocated arrays. The returned variable will contain
290  !! deallocated arrays
291 
292 !-----------------------------------------------------------------------
293 ! releases space used by horiz_interp_type variables
294 ! must be called before re-initializing the same variable
295 !-----------------------------------------------------------------------
296  select case(interp % interp_method)
297  case (conserve)
298  call horiz_interp_conserve_del(interp )
299  case (bilinear)
300  call horiz_interp_bilinear_del(interp )
301  case (bicubic)
302  call horiz_interp_bicubic_del(interp )
303  case (spherical)
304  call horiz_interp_spherical_del(interp )
305  end select
306 
307  interp%I_am_initialized = .false.
308 !-----------------------------------------------------------------------
309 
310  end subroutine horiz_interp_del
311 
312  !#####################################################################
313 
314  !> Dummy routine
315  subroutine horiz_interp_end
316  return
317  end subroutine horiz_interp_end
318 
319 #include "horiz_interp_r4.fh"
320 #include "horiz_interp_r8.fh"
321 
322 end module horiz_interp_mod
323 !> @}
324 ! close documentation grouping
integer function, public check_nml_error(IOSTAT, NML_NAME)
Checks the iostat argument that is returned after reading a namelist and determines if the error code...
Definition: fms.F90:523
subroutine, public write_version_number(version, tag, unit)
Prints to the log file (or a specified unit) the version id string and tag name.
Definition: fms.F90:701
logical function, public fms_error_handler(routine, message, err_msg)
Facilitates the control of fatal error conditions.
Definition: fms.F90:468
subroutine, public horiz_interp_bicubic_del(Interp)
Free memory from a horiz_interp_type used for bicubic interpolation (allocated via horiz_bicubic_new)
subroutine, public horiz_interp_bicubic_init
Initializes module and writes version number to logfile.out.
Creates a new horiz_interp_type for bicubic interpolation. Allocates space and initializes a derived-...
subroutine, public horiz_interp_bilinear_del(Interp)
Deallocates memory used by "horiz_interp_type" variables.
subroutine, public horiz_interp_bilinear_init
Initialize this module and writes version number to logfile.
Creates a horiz_interp_type for bilinear interpolation.
Subroutines for reading in weight files and using that to fill in the horiz_interp type instead calcu...
subroutine, public horiz_interp_conserve_init
Writes version number to logfile.
subroutine, public horiz_interp_conserve_del(Interp)
Deallocates memory used by "HI_KIND_TYPE" variables. Must be called before reinitializing with horiz_...
Allocates space and initializes a derived-type variable that contains pre-computed interpolation indi...
subroutine, public horiz_interp_del(Interp)
Deallocates memory used by "horiz_interp_type" variables. Must be called before reinitializing with h...
logical reproduce_siena
Set reproduce_siena = .true. to reproduce siena results. Set reproduce_siena = .false....
subroutine, public horiz_interp_init
Initialize module and writes version number to logfile.out.
subroutine, public horiz_interp_end
Dummy routine.
Subroutine for performing the horizontal interpolation between two grids.
subroutine, public horiz_interp_spherical_init
Initializes module and writes version number to logfile.out.
subroutine, public horiz_interp_spherical_del(Interp)
Deallocates memory used by "HI_KIND_TYPE" variables. Must be called before reinitializing with horiz_...
Holds data pointers and metadata for horizontal interpolations, passed between the horiz_interp modul...
integer function stdout()
This function returns the current standard fortran unit numbers for output.
Definition: mpp_util.inc:42
integer function stdlog()
This function returns the current standard fortran unit numbers for log messages. Log messages,...
Definition: mpp_util.inc:58
integer function mpp_pe()
Returns processor ID.
Definition: mpp_util.inc:406
Error handler.
Definition: mpp.F90:381
Reduction operations. Find the min of scalar a from the PEs in pelist result is also automatically br...
Definition: mpp.F90:559
Perform bicubic horizontal interpolation.
Allocates space and initializes a derived-type variable that contains pre-computed interpolation indi...
Subroutines for reading in weight files and using that to fill in the horiz_interp type instead calcu...
Private helper routines.