FMS 2025.01-dev
Flexible Modeling System
Loading...
Searching...
No Matches
axis_utils2.inc
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 axis_utils2_mod axis_utils2_mod
20!> @ingroup axis_utils
21!> @brief A set of utilities for manipulating axes and extracting axis attributes.
22!! FMS2_IO equivalent version of @ref axis_utils_mod.
23!> @author M.J. Harrison
24
25!> @addtogroup axis_utils2_mod
26!> @{
27
28 !> get axis edge data from a given file
29 subroutine axis_edges_(fileobj, name, edge_data, reproduce_null_char_bug_flag)
30
31 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object to read from
32 character(len=*), intent(in) :: name !< Name of a given axis
33 real(FMS_AU_KIND_), dimension(:), intent(out) :: edge_data !< Returned edge data from given axis name
34 logical, intent(in), optional :: reproduce_null_char_bug_flag !< Flag indicating to reproduce
35 !! the mpp_io bug where the null characters were not removed
36 !! after reading a string attribute
37
38 integer :: ndims
39 character(len=128) :: buffer
40 integer, dimension(:), allocatable :: dim_sizes
41 real(kind=fms_au_kind_), dimension(:), allocatable :: r_var
42 real(kind=fms_au_kind_), dimension(:,:), allocatable :: r2d
43 integer :: i
44 integer :: n
45 logical :: reproduce_null_char_bug !< Local flag
46 !! indicating to reproduce the mpp_io bug where
47 !! the null characters were not removed after reading a string attribute
48 integer, parameter :: lkind = fms_au_kind_
49 integer :: edge_index(2) !< Index to use when reading the edges from the file
50 !! (/1, 2/) if the axis data is monotonically increasing
51 !! (/2, 1/) if the axis data is monotonically decreasing
52
53 ndims = get_variable_num_dimensions(fileobj, name)
54 allocate(dim_sizes(ndims))
55
56 call get_variable_size(fileobj, name, dim_sizes)
57
58 n = dim_sizes(1)
59 if (size(edge_data) .ne. n+1) then
60 call mpp_error(fatal, "axis_edge: incorrect size of edge_data array.")
61 endif
62 deallocate(dim_sizes)
63
64 reproduce_null_char_bug = .false.
65 if (present(reproduce_null_char_bug_flag)) reproduce_null_char_bug = reproduce_null_char_bug_flag
66
67 buffer = ""
68 if (variable_att_exists(fileobj, name, "edges")) then
69 !! If the reproduce_null_char_bug flag is turned on fms2io will not remove the null character
70 call get_variable_attribute(fileobj, name, "edges", buffer(1:128), &
71 reproduce_null_char_bug_flag=reproduce_null_char_bug)
72
73 !! Check for a null character here, if it exists *_bnds will be calculated instead of read in
74 if (reproduce_null_char_bug) then
75 i = 0
76 i = index(buffer, char(0))
77 if (i > 0) buffer = ""
78 endif
79 elseif (variable_att_exists(fileobj, name, "bounds")) then
80 !! If the reproduce_null_char_bug flag is turned on fms2io will not remove the null character
81 call get_variable_attribute(fileobj, name, "bounds", buffer(1:128), &
82 reproduce_null_char_bug_flag=reproduce_null_char_bug)
83
84 !! Check for a null character here, if it exists *_bnds will be calculated instead of read in
85 if (reproduce_null_char_bug) then
86 i = 0
87 i = index(buffer, char(0))
88 if (i > 0) buffer = ""
89 endif
90 endif
91 if (trim(buffer) .ne. "") then
92 ndims = get_variable_num_dimensions(fileobj, buffer)
93 allocate(dim_sizes(ndims))
94
95 call get_variable_size(fileobj, buffer, dim_sizes)
96
97 if (size(dim_sizes) .eq. 1) then
98 if (dim_sizes(1) .ne. n+1) then
99 call mpp_error(fatal, "axis_edges: incorrect size of edge data.")
100 endif
101
102 call read_data(fileobj, buffer, edge_data)
103
104 elseif (size(dim_sizes) .eq. 2) then
105 if (dim_sizes(1) .ne. 2) then
106 call mpp_error(fatal, "axis_edges: first dimension of edge must be of size 2")
107 endif
108 if (dim_sizes(2) .ne. n) then
109 call mpp_error(fatal, "axis_edges: incorrect size of edge data.")
110 endif
111
112 allocate(r2d(dim_sizes(1), dim_sizes(2)))
113 call read_data(fileobj, buffer, r2d)
114 edge_index = (/1, 2/)
115 if (r2d(1,1) .gt. r2d(1,2)) edge_index = (/2, 1 /)
116 edge_data(1:dim_sizes(2)) = r2d(edge_index(1),:)
117 edge_data(dim_sizes(2)+1) = r2d(edge_index(2),dim_sizes(2))
118 deallocate(r2d)
119 endif
120 deallocate(dim_sizes)
121 else
122 allocate(r_var(n))
123
124 call read_data(fileobj, name, r_var)
125
126 do i = 2, n
127 edge_data(i) = r_var(i-1) + 0.5_lkind*(r_var(i) - r_var(i-1))
128 enddo
129 edge_data(1) = r_var(1) - 0.5_lkind*(r_var(2) - r_var(1))
130 if (abs(edge_data(1)) .lt. 1.e-10_lkind) then
131 edge_data(1) = 0.0_lkind
132 endif
133 edge_data(n+1) = r_var(n) + 0.5_lkind*(r_var(n) - r_var(n-1))
134 deallocate(r_var)
135 endif
136 end subroutine axis_edges_
137
138 !> @brief Returns lon_strt <= longitude <= lon_strt+360
139 !! @return real lon_in_range */
140
141 function lon_in_range_(lon, l_strt)
142 real(kind=fms_au_kind_), intent(in) :: lon, l_strt
143 real(kind=fms_au_kind_) :: lon_in_range_
144 real(kind=fms_au_kind_) :: l_end
145 integer, parameter :: lkind = fms_au_kind_
146
147 lon_in_range_ = lon
148 l_end = l_strt + 360.0_lkind
149
150 if (abs(lon_in_range_ - l_strt) < 1.e-4_lkind) then
151 lon_in_range_ = l_strt
152 return
153 endif
154
155 if (abs(lon_in_range_ - l_end) < 1.e-4_lkind) then
156 lon_in_range_ = l_strt
157 return
158 endif
159
160 do
161 if (lon_in_range_ < l_strt) then
162 lon_in_range_ = real(lon_in_range_, fms_au_kind_) + real(f360, fms_au_kind_)
163 else if (lon_in_range_ > l_end) then
164 lon_in_range_ = real(lon_in_range_, fms_au_kind_) - real(f360, fms_au_kind_)
165 else
166 exit
167 end if
168 end do
169
170 end function lon_in_range_
171
172 !> @brief Returns monotonic array of longitudes s.t., lon_strt <= lon(:) < lon_strt+360.
173 !!
174 !! This may require that entries be moved from the beginning of the array to
175 !! the end. If no entries are moved (i.e., if lon(:) is already monotonic in
176 !! the range from lon_start to lon_start + 360), then istrt is set to 0. If
177 !! any entries are moved, then istrt is set to the original index of the entry
178 !! which becomes lon(1).
179 !!
180 !! e.g.,
181 !!
182 !! lon = 0 1 2 3 4 5 ... 358 359; lon_strt = 3
183 !! ==> lon = 3 4 5 6 7 8 ... 359 360 361 362; istrt = 4
184 !!
185 subroutine tranlon_(lon, lon_start, istrt)
186 real(kind=fms_au_kind_), intent(inout), dimension(:) :: lon
187 real(kind=fms_au_kind_), intent(in) :: lon_start
188 integer, intent(out) :: istrt
189 integer :: len, i
190 real(kind=fms_au_kind_) :: lon_strt, tmp(size(lon(:))-1)
191
192 len = size(lon(:))
193
194 do i = 1, len
195 lon(i) = lon_in_range(lon(i),lon_start)
196 enddo
197
198 istrt = 0
199 do i = 1,len-1
200 if (lon(i+1) < lon(i)) then
201 istrt = i+1
202 exit
203 endif
204 enddo
205
206 if (istrt>1) then ! grid is not monotonic
207 if (abs(lon(len)-lon(1)) < real(epsln, fms_au_kind_)) then
208 tmp = cshift(lon(1:len-1),istrt-1)
209 lon(1:len-1) = tmp
210 lon(len) = lon(1)
211 else
212 lon = cshift(lon,istrt-1)
213 endif
214
215 lon_strt = lon(1)
216 do i=2,len
217 lon(i) = lon_in_range(lon(i),lon_strt)
218 lon_strt = lon(i)
219 enddo
220 endif
221
222 return
223 end subroutine tranlon_
224
225
226 function frac_index_(rval, array)
227
228 integer :: ia, i, ii, iunit
229 real(kind=fms_au_kind_) :: rval !< arbitrary data...same units as elements in "array"
230 real(kind=fms_au_kind_) :: frac_index_
231 real(kind=fms_au_kind_), dimension(:) :: array !< array of data points (must be monotonically increasing)
232 logical :: keep_going
233 integer, parameter :: lkind = fms_au_kind_
234 ia = size(array(:))
235
236 do i = 2, ia
237 if (array(i) < array(i-1)) then
238 iunit = stdout()
239 write (iunit,*) '=> Error: "frac_index" array must be monotonically' &
240 & // 'increasing when searching for nearest value to ', rval
241 write (iunit,*) ' array(i) < array(i-1) for i=',i
242 write (iunit,*) ' array(i) for i=1..ia follows:'
243 do ii = 1, ia
244 write (iunit,*) 'i=',ii, ' array(i)=',array(ii)
245 enddo
246 call mpp_error(fatal,' "frac_index" array must be monotonically increasing.')
247 endif
248 enddo
249
250 if (rval < array(1) .or. rval > array(ia)) then
251 frac_index_ = -1.0_lkind
252 else
253 i = 1
254 keep_going = .true.
255 do while (i <= ia .and. keep_going)
256 i = i+1
257 if (rval <= array(i)) then
258 frac_index_ = real((i-1), lkind) + (rval-array(i-1)) / (array(i) - array(i-1))
259 keep_going = .false.
260 endif
261 enddo
262 endif
263 end function frac_index_
264
265 !> @brief Return index of nearest point along axis
266 !!
267 !> nearest_index = index of nearest data point within "array" corresponding to
268 !! "value".
269 !!
270 !! inputs:
271 !!
272 !! rval = arbitrary data...same units as elements in "array"
273 !! array = array of data points (must be monotonically)
274 !! ia = dimension of "array"
275 !!
276 !! output:
277 !!
278 !! nearest_index = index of nearest data point to "value"
279 !! if "value" is outside the domain of "array" then nearest_index = 1
280 !! or "ia" depending on whether array(1) or array(ia) is
281 !! closest to "value"
282 !!
283 !! note: if "array" is dimensioned array(0:ia) in the calling
284 !! program, then the returned index should be reduced
285 !! by one to account for the zero base.
286 !!
287 !! example:
288 !!
289 !! let model depths be defined by the following:
290 !! parameter (km=5)
291 !! dimension z(km)
292 !! data z /5.0, 10.0, 50.0, 100.0, 250.0/
293 !!
294 !! k1 = nearest_index (12.5, z, km)
295 !! k2 = nearest_index (0.0, z, km)
296 !!
297 !! k1 would be set to 2, and k2 would be set to 1 so that
298 !! z(k1) would be the nearest data point to 12.5 and z(k2) would
299 !! be the nearest data point to 0.0
300 !! @return integer nearest_index
301 function nearest_index_(rval, array)
302 real(kind=fms_au_kind_), intent(in) :: rval !< arbitrary data...same units as elements in "array"
303 real(kind=fms_au_kind_), intent(in), dimension(:) :: array !< array of data points (must be monotonic)
304
305 integer :: nearest_index_
306 integer :: ia !< dimension of "array"
307 integer :: i !< For looping through "array"
308
309 logical :: increasing !< .True. if the array is increasing
310
311 ia = SIZE(array(:))
312
313 ! check if array is increasing
314 increasing = .true.
315 DO i = 2, ia-1
316 IF( array(i) .lt. array(i-1)) then
317 increasing = .false.
318 exit
319 endif
320 END DO
321
322 if (.not. increasing) then
323 ! if not increasing, check that it is decreasing
324 DO i = 2, ia-1
325 IF( array(i) .gt. array(i-1)) &
326 call mpp_error(fatal, 'axis_utils2::nearest_index array is NOT monotonously ordered')
327 END DO
328 endif
329
330 array_is_increasing: if (increasing) then
331 !< Check if the rval is outside the range of the array
332 if (rval .le. array(1)) then
334 return
335 elseif (rval .ge. array(ia)) then
336 nearest_index_ = ia
337 return
338 endif
339
340 DO i = 2, ia
341 if (rval .le. array(i)) then
343 if (array(i) -rval .gt. rval - array(i-1)) nearest_index_ = i - 1
344 return
345 endif
346 END DO
347 else !array_is_decreasing
348 !< Check if the rval is outside the range of the array
349 if (rval .le. array(ia)) then
350 nearest_index_ = ia
351 return
352 elseif (rval .gt. array(1)) then
354 return
355 endif
356
357 DO i = 2, ia
358 if (rval .ge. array(i)) then
360 if (rval - array(i) .gt. array(i-1) -rval ) nearest_index_ = i - 1
361 return
362 endif
363 END DO
364 endif array_is_increasing
365 end function nearest_index_
366
367 !#############################################################################
368
369 subroutine interp_1d_linear_(grid1,grid2,data1,data2)
370
371 real(kind=fms_au_kind_), dimension(:), intent(in) :: grid1, data1, grid2
372 real(kind=fms_au_kind_), dimension(:), intent(inout) :: data2
373
374 integer :: n1, n2, i, n
375 real(kind=fms_au_kind_) :: w
376 integer, parameter :: lkind = fms_au_kind_
377
378 n1 = size(grid1(:))
379 n2 = size(grid2(:))
380
381
382 do i = 2, n1
383 if (grid1(i) <= grid1(i-1)) call mpp_error(fatal, 'grid1 not monotonic')
384 enddo
385
386 do i = 2, n2
387 if (grid2(i) <= grid2(i-1)) call mpp_error(fatal, 'grid2 not monotonic')
388 enddo
389
390 if (grid1(1) > grid2(1) ) call mpp_error(fatal, 'grid2 lies outside grid1')
391 if (grid1(n1) < grid2(n2) ) call mpp_error(fatal, 'grid2 lies outside grid1')
392
393 do i = 1, n2
394 n = nearest_index(grid2(i),grid1)
395
396 if (grid1(n) < grid2(i)) then
397 w = (grid2(i)-grid1(n))/(grid1(n+1)-grid1(n))
398 data2(i) = (1.0_lkind-w)*data1(n) + w*data1(n+1)
399 else
400 if(n==1) then
401 data2(i) = data1(n)
402 else
403 w = (grid2(i)-grid1(n-1))/(grid1(n)-grid1(n-1))
404 data2(i) = (1.0_lkind-w)*data1(n-1) + w*data1(n)
405 endif
406 endif
407 enddo
408
409
410 return
411
412 end subroutine interp_1d_linear_
413
414 !###################################################################
415 subroutine interp_1d_cubic_spline_(grid1, grid2, data1, data2, yp1, ypn)
416
417 real(kind=fms_au_kind_), dimension(:), intent(in) :: grid1, grid2, data1
418 real(kind=fms_au_kind_), dimension(:), intent(inout) :: data2
419 real(kind=fms_au_kind_), intent(in) :: yp1, ypn
420
421 real(kind=fms_au_kind_), dimension(size(grid1)) :: y2, u
422 real(kind=fms_au_kind_) :: sig, p, qn, un, h, a ,b
423 integer :: n, m, i, k, klo, khi
424 integer, parameter :: lkind = fms_au_kind_
425
426 n = size(grid1(:))
427 m = size(grid2(:))
428
429 do i = 2, n
430 if (grid1(i) <= grid1(i-1)) call mpp_error(fatal, 'grid1 not monotonic')
431 enddo
432
433 do i = 2, m
434 if (grid2(i) <= grid2(i-1)) call mpp_error(fatal, 'grid2 not monotonic')
435 enddo
436
437 if (grid1(1) > grid2(1) ) call mpp_error(fatal, 'grid2 lies outside grid1')
438 if (grid1(n) < grid2(m) ) call mpp_error(fatal, 'grid2 lies outside grid1')
439
440if (yp1>0.99e30_lkind) then
441 y2(1) = 0.0_lkind
442 u(1) = 0.0_lkind
443 else
444 y2(1) = -0.5_lkind
445 u(1) = (3.0_lkind)/(grid1(2)-grid1(1))*((data1(2)-data1(1))/(grid1(2)-grid1(1))-yp1)
446 endif
447
448 do i = 2, n-1
449 sig = (grid1(i)-grid1(i-1))/(grid1(i+1)-grid1(i-1))
450 p = sig*y2(i-1) + 2.0_lkind
451 y2(i) = (sig-1.0_lkind)/p
452 u(i) = (6.0_lkind*((data1(i+1)-data1(i))/(grid1(i+1)-grid1(i))-(data1(i)-data1(i-1)) &
453 /(grid1(i)-grid1(i-1)))/(grid1(i+1)-grid1(i-1))-sig*u(i-1))/p
454 enddo
455
456 if (ypn>0.99e30_lkind) then
457 qn = 0.0_lkind
458 un = 0.0_lkind
459 else
460 qn = 0.5_lkind
461 un = (3.0_lkind)/(grid1(n)-grid1(n-1))*(ypn-(data1(n)-data1(n-1))/ &
462 (grid1(n)-grid1(n-1)))
463 endif
464
465 y2(n) = (un-qn*u(n-1))/(qn*y2(n-1)+1.0_lkind)
466
467 do k = n-1,1,-1
468 y2(k) = y2(k)*y2(k+1)+u(k)
469 enddo
470
471 do k = 1, m
472 n = nearest_index(grid2(k),grid1)
473 if (grid1(n) < grid2(k)) then
474 klo = n
475 else
476 if(n==1) then
477 klo = n
478 else
479 klo = n -1
480 endif
481 endif
482
483 khi = klo+1
484 h = grid1(khi)-grid1(klo)
485 a = (grid1(khi) - grid2(k))/h
486 b = (grid2(k) - grid1(klo))/h
487 data2(k) = a*data1(klo) + b*data1(khi)+ ((a**3-a)*y2(klo) + (b**3-b)*y2(khi))*(h**2) &
488 /6.0_lkind
489 enddo
490
491 end subroutine interp_1d_cubic_spline_
492
493 !###################################################################
494
495 subroutine interp_1d_1d_(grid1,grid2,data1,data2, method, yp1, yp2)
496
497 real(kind=fms_au_kind_), dimension(:), intent(in) :: grid1, data1, grid2
498 real(kind=fms_au_kind_), dimension(:), intent(inout) :: data2
499 character(len=*), optional, intent(in) :: method
500 real(kind=fms_au_kind_), optional, intent(in) :: yp1, yp2
501
502 real(kind=fms_au_kind_) :: y1, y2
503 character(len=32) :: interp_method
504 integer :: k2, ks, ke
505 integer, parameter :: lkind = fms_au_kind_
506
507 k2 = size(grid2(:))
508
509 interp_method = "linear"
510 if(present(method)) interp_method = method
511 y1 = 1.0e30_lkind
512
513 if(present(yp1)) y1 = yp1
514 y2 = 1.0e30_lkind
515
516 if(present(yp2)) y2 = yp2
517 call find_index(grid1, grid2(1), grid2(k2), ks, ke)
518 select case(trim(interp_method))
519 case("linear")
520 call interp_1d_linear(grid1(ks:ke),grid2,data1(ks:ke),data2)
521 case("cubic_spline")
522 call interp_1d_cubic_spline(grid1(ks:ke),grid2,data1(ks:ke),data2, y1, y2)
523 case default
524 call mpp_error(fatal,"axis_utils: interp_method should be linear or cubic_spline")
525 end select
526
527 return
528
529 end subroutine interp_1d_1d_
530
531 !###################################################################
532
533
534 subroutine interp_1d_2d_(grid1,grid2,data1,data2)
535
536 real(kind=fms_au_kind_), dimension(:,:), intent(in) :: grid1, data1, grid2
537 real(kind=fms_au_kind_), dimension(:,:), intent(inout) :: data2
538
539 integer :: n1, n2, n, k2, ks, ke
540
541 n1 = size(grid1,1)
542 n2 = size(grid2,1)
543 k2 = size(grid2,2)
544
545 if (n1 /= n2) call mpp_error(fatal,'grid size mismatch')
546
547 do n = 1, n1
548 call find_index(grid1(n,:), grid2(n,1), grid2(n,k2), ks, ke)
549 call interp_1d_linear(grid1(n,ks:ke),grid2(n,:),data1(n,ks:ke),data2(n,:))
550 enddo
551
552 return
553
554 end subroutine interp_1d_2d_
555
556 !###################################################################
557
558 subroutine interp_1d_3d_(grid1,grid2,data1,data2, method, yp1, yp2)
559
560 real(fms_au_kind_), dimension(:,:,:), intent(in) :: grid1, data1, grid2
561 real(fms_au_kind_), dimension(:,:,:), intent(inout) :: data2
562 character(len=*), optional, intent(in) :: method
563 real(kind=fms_au_kind_), optional, intent(in) :: yp1, yp2
564
565 integer :: n1, n2, m1, m2, k2, n, m
566 real(kind=fms_au_kind_) :: y1, y2
567 character(len=32) :: interp_method
568 integer :: ks, ke
569 integer, parameter :: lkind = fms_au_kind_
570
571 n1 = size(grid1,1)
572 n2 = size(grid2,1)
573 m1 = size(grid1,2)
574 m2 = size(grid2,2)
575 k2 = size(grid2,3)
576
577 interp_method = "linear"
578 if(present(method)) interp_method = method
579 y1 = 1.0e30_lkind
580
581 if(present(yp1)) y1 = yp1
582 y2 = 1.0e30_lkind
583 if(present(yp2)) y2 = yp2
584
585 if (n1 /= n2 .or. m1 /= m2) call mpp_error(fatal,'grid size mismatch')
586
587 select case(trim(interp_method))
588 case("linear")
589 do m = 1, m1
590 do n = 1, n1
591 call find_index(grid1(n,m,:), grid2(n,m,1), grid2(n,m,k2), ks, ke)
592 call interp_1d_linear(grid1(n,m,ks:ke),grid2(n,m,:),data1(n,m,ks:ke),data2(n,m,:))
593 enddo
594 enddo
595
596 case("cubic_spline")
597 do m = 1, m1
598 do n = 1, n1
599 call find_index(grid1(n,m,:), grid2(n,m,1), grid2(n,m,k2), ks, ke)
600 call interp_1d_cubic_spline(grid1(n,m,ks:ke),grid2(n,m,:), data1(n,m,ks:ke),data2(n,m,:), y1, y2)
601 enddo
602 enddo
603
604 case default
605 call mpp_error(fatal,"axis_utils: interp_method should be linear or cubic_spline")
606 end select
607
608 return
609
610 end subroutine interp_1d_3d_
611
612
613 !#####################################################################
614 subroutine find_index_(grid1, xs, xe, ks, ke)
615 real(kind=fms_au_kind_), dimension(:), intent(in) :: grid1
616 real(kind=fms_au_kind_), intent(in) :: xs, xe
617 integer, intent(out) :: ks, ke
618
619 integer :: k, nk
620
621 nk = size(grid1(:))
622
623 ks = 0; ke = 0
624 do k = 1, nk-1
625 if(grid1(k) <= xs .and. grid1(k+1) > xs ) then
626 ks = k
627 exit
628 endif
629 enddo
630
631 do k = nk, 2, -1
632 if(grid1(k) >= xe .and. grid1(k-1) < xe ) then
633 ke = k
634 exit
635 endif
636 enddo
637
638 if(ks == 0 ) call mpp_error(fatal,' xs locate outside of grid1')
639 if(ke == 0 ) call mpp_error(fatal,' xe locate outside of grid1')
640
641 end subroutine find_index_
642 !> @}
643 ! close documentation grouping
subroutine axis_edges_(fileobj, name, edge_data, reproduce_null_char_bug_flag)
get axis edge data from a given file
real(kind=fms_au_kind_) function frac_index_(rval, array)
integer function nearest_index_(rval, array)
Return index of nearest point along axis.
subroutine tranlon_(lon, lon_start, istrt)
Returns monotonic array of longitudes s.t., lon_strt <= lon(:) < lon_strt+360.
real(kind=fms_au_kind_) function lon_in_range_(lon, l_strt)
Returns lon_strt <= longitude <= lon_strt+360.