FMS  2025.04
Flexible Modeling System
horiz_interp_spherical.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_spherical_mod horiz_interp_spherical_mod
19 !> @ingroup horiz_interp
20 !> @brief Performs spatial interpolation between grids using inverse-distance-weighted scheme.
21 !> This module can interpolate data from rectangular/tripolar grid
22 !! to rectangular/tripolar grid. The interpolation scheme is inverse-distance-weighted
23 !! scheme. There is an optional mask field for missing input data.
24 !! An optional output mask field may be used in conjunction with
25 !! the input mask to show where output data exists.
26 
27 !> @addtogroup horiz_interp_spherical_mod
28 !> @{
29 module horiz_interp_spherical_mod
30 
31  use platform_mod, only : r4_kind, r8_kind
32  use mpp_mod, only : mpp_error, fatal, warning, stdout
33  use mpp_mod, only : mpp_root_pe, mpp_pe
34  use mpp_mod, only : input_nml_file
35  use fms_mod, only : write_version_number
36  use fms_mod, only : check_nml_error
37  use constants_mod, only : pi
38  use horiz_interp_type_mod, only : horiz_interp_type, stats, spherical
39 
40  implicit none
41  private
42 
44  module procedure horiz_interp_spherical_r4
45  module procedure horiz_interp_spherical_r8
46  end interface
47 
49  module procedure horiz_interp_spherical_new_r4
50  module procedure horiz_interp_spherical_new_r8
51  end interface
52 
54  module procedure horiz_interp_spherical_wght_r4
55  module procedure horiz_interp_spherical_wght_r8
56  end interface
57 
60 
61  !> private helper routines
62  interface full_search
63  module procedure full_search_r4
64  module procedure full_search_r8
65  end interface
66 
67  interface radial_search
68  module procedure radial_search_r4
69  module procedure radial_search_r8
70  end interface
71 
73  module procedure spherical_distance_r4
74  module procedure spherical_distance_r8
75  end interface
76 
77  integer, parameter :: max_neighbors = 400
78  real(R8_KIND), parameter :: max_dist_default = 0.1_r8_kind ! radians
79  integer, parameter :: num_nbrs_default = 4
80  real(R8_KIND), parameter :: large=1.e20_r8_kind
81  real(R8_KIND), parameter :: epsln=1.e-10_r8_kind
82 
83  integer :: pe, root_pe
84 
85 
86  character(len=32) :: search_method = "radial_search" !< Namelist variable to indicate the searching
87  !! method to find the
88  !! nearest neighbor points. Its value can be "radial_search" and "full_search",
89  !! with default value "radial_search". when search_method is "radial_search",
90  !! the search may be not quite accurate for some cases. Normally the search will
91  !! be ok if you chose suitable max_dist. When search_method is "full_search",
92  !! it will be always accurate, but will be slower comparing to "radial_search".
93  !! Normally these two search algorithm will produce same results other than
94  !! order of operation. "radial_search" are recommended to use. The purpose to
95  !! add "full_search" is in case you think you interpolation results is
96  !! not right, you have other option to verify.
97 
98 !or "full_search"
99  namelist /horiz_interp_spherical_nml/ search_method
100 
101  !-----------------------------------------------------------------------
102  ! Include variable "version" to be written to log file.
103 #include<file_version.h>
104  logical :: module_is_initialized = .false.
105 
106 contains
107 
108  !#######################################################################
109 
110  !> Initializes module and writes version number to logfile.out
112  integer :: ierr, io
113 
114 
115  if(module_is_initialized) return
116  call write_version_number("horiz_interp_spherical_mod", version)
117  read (input_nml_file, horiz_interp_spherical_nml, iostat=io)
118  ierr = check_nml_error(io,'horiz_interp_spherical_nml')
119 
120  module_is_initialized = .true.
121 
122  end subroutine horiz_interp_spherical_init
123 
124  !#######################################################################
125 
126  !> Deallocates memory used by "HI_KIND_TYPE" variables.
127  !! Must be called before reinitializing with horiz_interp_spherical_new.
128  subroutine horiz_interp_spherical_del( Interp )
129 
130  type (horiz_interp_type), intent(inout) :: interp !< A derived-type variable returned by previous
131  !! call to horiz_interp_spherical_new. The input variable
132  !! must have allocated arrays. The returned variable will
133  !! contain deallocated arrays.
134 
135  if(interp%horizInterpReals4_type%is_allocated) then
136  if(allocated(interp%horizInterpReals4_type%src_dist)) deallocate(interp%horizInterpReals4_type%src_dist)
137  else if (interp%horizInterpReals8_type%is_allocated) then
138  if(allocated(interp%horizInterpReals8_type%src_dist)) deallocate(interp%horizInterpReals8_type%src_dist)
139  endif
140  if(allocated(interp%num_found)) deallocate(interp%num_found)
141  if(allocated(interp%i_lon)) deallocate(interp%i_lon)
142  if(allocated(interp%j_lat)) deallocate(interp%j_lat)
143 
144  interp%horizInterpReals4_type%is_allocated = .false.
145  interp%horizInterpReals8_type%is_allocated = .false.
146 
147  end subroutine horiz_interp_spherical_del
148 
149  !#######################################################################
150 
151 #include "horiz_interp_spherical_r4.fh"
152 #include "horiz_interp_spherical_r8.fh"
153 
154 end module horiz_interp_spherical_mod
155 !> @}
156 ! 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
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_...
character(len=32) search_method
Namelist variable to indicate the searching method to find the nearest neighbor points....
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 mpp_pe()
Returns processor ID.
Definition: mpp_util.inc:406
Error handler.
Definition: mpp.F90:381