FMS  2025.04
Flexible Modeling System
horiz_interp_conserve.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_conserve_mod horiz_interp_conserve_mod
19 !> @ingroup horiz_interp
20 !> @brief Performs spatial interpolation between grids using conservative interpolation
21 !!
22 !> @author Bruce Wyman, Zhi Liang
23 !!
24 !> This module can conservatively interpolate data from any logically rectangular grid
25 !! to any rectangular grid. The interpolation scheme is area-averaging
26 !! conservative scheme. There is an optional mask field for missing input data in both
27 !! horiz_interp__conserveinit and horiz_interp_conserve. For efficiency purpose, mask should only be
28 !! kept in horiz_interp_init (will remove the mask in horiz_interp in the future).
29 !! There are 1-D and 2-D version of horiz_interp_conserve_init for 1-D and 2-D grid.
30 !! There is a optional argument mask in horiz_interp_conserve_init_2d and no mask should
31 !! to passed into horiz_interp_conserv. optional argument mask will not be passed into
32 !! horiz_interp_conserve_init_1d and optional argument mask may be passed into
33 !! horiz_interp_conserve (For the purpose of reproduce Memphis??? results).
34 !! An optional output mask field may be used in conjunction with the input mask to show
35 !! where output data exists.
36 
37 module horiz_interp_conserve_mod
38 
39  use platform_mod, only: r4_kind, r8_kind
40  use mpp_mod, only: mpp_send, mpp_recv, mpp_pe, mpp_root_pe, mpp_npes
41  use mpp_mod, only: mpp_error, fatal, mpp_sync_self
42  use mpp_mod, only: comm_tag_1, comm_tag_2
43  use fms_mod, only: write_version_number
44  use grid2_mod, only: get_great_circle_algorithm
45  use constants_mod, only: pi
46  use horiz_interp_type_mod, only: horiz_interp_type, conserve
47 
48 
49  implicit none
50  private
51 
52  ! public interface
53 
54 
55  !> @brief Allocates space and initializes a derived-type variable
56  !! that contains pre-computed interpolation indices and weights.
57  !!
58  !> Allocates space and initializes a derived-type variable
59  !! that contains pre-computed interpolation indices and weights
60  !! for improved performance of multiple interpolations between
61  !! the same grids.
62  !! @param lon_in
63  !! Longitude (in radians) for source data grid.
64  !!
65  !! @param lat_in
66  !! Latitude (in radians) for source data grid.
67  !!
68  !! @param lon_out
69  !! Longitude (in radians) for destination data grid.
70  !!
71  !! @param lat_out
72  !! Latitude (in radians) for destination data grid.
73  !!
74  !! @param verbose
75  !! flag for the amount of print output.
76  !!
77  !! @param mask_in
78  !! Input mask. must be the size (size(lon_in)-1, size(lon. The real value of
79  !! mask_in must be in the range (0.,1.). Set mask_in=0.0 for data points
80  !! that should not be used or have missing data.
81  !!
82  !! @param mask_out
83  !! Output mask that specifies whether data was computed.
84  !!
85  !! @param Interp
86  !! A derived-type variable containing indices and weights used for subsequent
87  !! interpolations. To reinitialize this variable for a different grid-to-grid
88  !! interpolation you must first use the "horiz_interp_del" interface.
89  !!
90  !> @ingroup horiz_interp_conserve_mod
92  module procedure horiz_interp_conserve_new_1dx1d_r4
93  module procedure horiz_interp_conserve_new_1dx2d_r4
94  module procedure horiz_interp_conserve_new_2dx1d_r4
95  module procedure horiz_interp_conserve_new_2dx2d_r4
96  module procedure horiz_interp_conserve_new_1dx1d_r8
97  module procedure horiz_interp_conserve_new_1dx2d_r8
98  module procedure horiz_interp_conserve_new_2dx1d_r8
99  module procedure horiz_interp_conserve_new_2dx2d_r8
100  end interface
101 
103  module procedure horiz_interp_conserve_r4
104  module procedure horiz_interp_conserve_r8
105  end interface
106 
107 !> private helper routines
108  interface data_sum
109  module procedure data_sum_r4
110  module procedure data_sum_r8
111  end interface
112 
113  interface stats
114  module procedure stats_r4
115  module procedure stats_r8
116  end interface
117 
119  module procedure horiz_interp_conserve_version1_r8
120  module procedure horiz_interp_conserve_version1_r4
121  end interface
122 
124  module procedure horiz_interp_conserve_version2_r8
125  module procedure horiz_interp_conserve_version2_r4
126  end interface
127 
128 
129 
130  !> @addtogroup horiz_interp_conserve_mod
131  !> @{
134 
135  integer :: pe, root_pe
136  !-----------------------------------------------------------------------
137  ! Include variable "version" to be written to log file.
138 #include<file_version.h>
139  logical :: module_is_initialized = .false.
140 
141  logical :: great_circle_algorithm = .false.
142 
143 contains
144 
145  !> Writes version number to logfile.
147 
148  if(module_is_initialized) return
149  call write_version_number("HORIZ_INTERP_CONSERVE_MOD", version)
150 
151  great_circle_algorithm = get_great_circle_algorithm()
152 
153  module_is_initialized = .true.
154 
155  end subroutine horiz_interp_conserve_init
156 
157  !> Deallocates memory used by "HI_KIND_TYPE" variables.
158  !! Must be called before reinitializing with horiz_interp_new.
159  subroutine horiz_interp_conserve_del ( Interp )
160 
161  type (horiz_interp_type), intent(inout) :: interp !< A derived-type variable returned by
162  !! previous call to horiz_interp_new. The input variable must have
163  !! allocated arrays. The returned variable will contain deallocated arrays.
164 
165  select case(interp%version)
166  case (1)
167  if( interp%horizInterpReals8_type%is_allocated) then
168  if(allocated(interp%horizInterpReals8_type%area_src)) deallocate(interp%horizInterpReals8_type%area_src)
169  if(allocated(interp%horizInterpReals8_type%area_dst)) deallocate(interp%horizInterpReals8_type%area_dst)
170  if(allocated(interp%horizInterpReals8_type%facj)) deallocate(interp%horizInterpReals8_type%facj)
171  if(allocated(interp%jlat)) deallocate(interp%jlat)
172  if(allocated(interp%horizInterpReals8_type%faci)) deallocate(interp%horizInterpReals8_type%faci)
173  if(allocated(interp%ilon)) deallocate(interp%ilon)
174  else if( interp%horizInterpReals4_type%is_allocated) then
175  if(allocated(interp%horizInterpReals4_type%area_src)) deallocate(interp%horizInterpReals4_type%area_src)
176  if(allocated(interp%horizInterpReals4_type%area_dst)) deallocate(interp%horizInterpReals4_type%area_dst)
177  if(allocated(interp%horizInterpReals4_type%facj)) deallocate(interp%horizInterpReals4_type%facj)
178  if(allocated(interp%jlat)) deallocate(interp%jlat)
179  if(allocated(interp%horizInterpReals4_type%faci)) deallocate(interp%horizInterpReals4_type%faci)
180  if(allocated(interp%ilon)) deallocate(interp%ilon)
181  endif
182  case (2)
183  if( interp%horizInterpReals8_type%is_allocated) then
184  if(allocated(interp%i_src)) deallocate(interp%i_src)
185  if(allocated(interp%j_src)) deallocate(interp%j_src)
186  if(allocated(interp%i_dst)) deallocate(interp%i_dst)
187  if(allocated(interp%j_dst)) deallocate(interp%j_dst)
188  if(allocated(interp%horizInterpReals8_type%area_frac_dst)) &
189  deallocate(interp%horizInterpReals8_type%area_frac_dst)
190  else if( interp%horizInterpReals4_type%is_allocated ) then
191  if(allocated(interp%i_src)) deallocate(interp%i_src)
192  if(allocated(interp%j_src)) deallocate(interp%j_src)
193  if(allocated(interp%i_dst)) deallocate(interp%i_dst)
194  if(allocated(interp%j_dst)) deallocate(interp%j_dst)
195  if(allocated(interp%horizInterpReals4_type%area_frac_dst)) &
196  deallocate(interp%horizInterpReals4_type%area_frac_dst)
197  endif
198  end select
199  interp%horizInterpReals4_type%is_allocated = .false.
200  interp%horizInterpReals8_type%is_allocated = .false.
201 
202  end subroutine horiz_interp_conserve_del
203 
204 #include "horiz_interp_conserve_r4.fh"
205 #include "horiz_interp_conserve_r8.fh"
206 
207 end module horiz_interp_conserve_mod
208 !> @}
209 ! close documentation grouping
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 get_great_circle_algorithm()
Determine if we are using the great circle algorithm.
Definition: grid2.F90:188
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...
Holds data pointers and metadata for horizontal interpolations, passed between the horiz_interp modul...
subroutine mpp_sync_self(pelist, check, request, msg_size, msg_type)
This is to check if current PE's outstanding puts are complete but we can't use shmem_fence because w...
integer function mpp_npes()
Returns processor count for current pelist.
Definition: mpp_util.inc:420
integer function mpp_pe()
Returns processor ID.
Definition: mpp_util.inc:406
Error handler.
Definition: mpp.F90:381
Receive data from another PE.
Definition: mpp.F90:950
Send data to a receiving PE.
Definition: mpp.F90:1017