FMS 2025.01-dev
Flexible Modeling System
Loading...
Searching...
No Matches
column_diagnostics.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!> @file
20
21!> @brief initialize_diagnostic_columns returns the (i, j, lat, lon) coord-
22!! inates of any diagnostic columns that are located on the current
23!! processor.
24subroutine initialize_diagnostic_columns_ &
25 (module, num_diag_pts_latlon, num_diag_pts_ij, &
26 global_i , global_j , global_lat_latlon, &
27 global_lon_latlon, lonb_in, latb_in, &
28 do_column_diagnostics, &
29 diag_lon, diag_lat, diag_i, diag_j, diag_units)
30
31!---------------------------------------------------------------------
32! initialize_diagnostic_columns returns the (i, j, lat, lon) coord-
33! inates of any diagnostic columns that are located on the current
34! processor.
35!----------------------------------------------------------------------
36
37!---------------------------------------------------------------------
38character(len=*), intent(in) :: module !< module calling this subroutine
39integer, intent(in) :: num_diag_pts_latlon !< number of diagnostic columns specified
40 !! by lat-lon coordinates
41integer, intent(in) :: num_diag_pts_ij !< number of diagnostic columns specified
42 !! by global (i,j) coordinates
43integer, dimension(:), intent(in) :: global_i !< specified global i coordinates
44integer, dimension(:), intent(in) :: global_j !< specified global j coordinates
45real(FMS_CD_KIND_), dimension(:), intent(in) :: global_lat_latlon !< specified global lat coordinates
46real(FMS_CD_KIND_), dimension(:), intent(in) :: global_lon_latlon !< specified global lon coordinates
47real(FMS_CD_KIND_), dimension(:,:), intent(in) :: lonb_in, latb_in
48logical, dimension(:,:), intent(out) :: do_column_diagnostics !< is a diagnostic column in this jrow ?
49integer, dimension(:), intent(inout) :: diag_i !< processor i indices of diagnstic columns
50integer, dimension(:), intent(inout) :: diag_j !< processor j indices of diagnstic columns
51real(FMS_CD_KIND_), dimension(:), intent(out) :: diag_lat !< latitudes of diagnostic columns [ degrees ]
52real(FMS_CD_KIND_), dimension(:), intent(out) :: diag_lon !< longitudes of diagnostic columns [ degrees ]
53integer, dimension(:), intent(out) :: diag_units !< unit number for each diagnostic column
54!---------------------------------------------------------------------
55
56!---------------------------------------------------------------------
57! intent(in) variables:
58!
59! module module calling this subroutine
60! num_diag_pts_latlon number of diagnostic columns specified
61! by lat-lon coordinates
62! num_diag_pts_ij number of diagnostic columns specified
63! by global (i,j) coordinates
64! global_i specified global i coordinates
65! global_j specified global j coordinates
66! global_lat_latlon specified global lat coordinates
67! global_lon_latlon specified global lon coordinates
68!
69! intent(out) variables:
70!
71! do_column_diagnostics is a diagnostic column in this jrow ?
72! diag_i processor i indices of diagnstic columns
73! diag_j processor j indices of diagnstic columns
74! diag_lat latitudes of diagnostic columns
75! [ degrees ]
76! diag_lon longitudes of diagnostic columns
77! [ degrees ]
78! diag_units unit number for each diagnostic column
79!
80!---------------------------------------------------------------------
81
82!--------------------------------------------------------------------
83! local variables:
84
85 real(FMS_CD_KIND_), dimension(size(diag_i,1)) :: global_lat !< latitudes for all diagnostic columns [ degrees ]
86 real(FMS_CD_KIND_), dimension(size(diag_i,1)) :: global_lon !< longitudes for all diagnostic columns [ degrees ]
87 real(FMS_CD_KIND_), dimension(size(latb_in,1)-1, size(latb_in,2)-1) :: &
88 distance, distance_x, distance_y, &
89 distance_x2, distance2
90 real(FMS_CD_KIND_), dimension(size(latb_in,1), size(latb_in,2)) :: latb_deg
91 real(FMS_CD_KIND_), dimension(size(lonb_in,1), size(lonb_in,2)) :: lonb_deg
92 real(FMS_CD_KIND_) :: dellat, dellon
93 real(FMS_CD_KIND_) :: latb_max, latb_min, lonb_max, lonb_min
94
95 integer :: num_diag_pts !< total number of diagnostic columns
96 integer :: i !< do loop indices
97 integer :: j !< do loop indices
98 integer :: nn !< do loop indices
99 real(FMS_CD_KIND_) :: ref_lat
100 real(FMS_CD_KIND_) :: current_distance
101 character(len=8) :: char !< character string for diaganostic column index
102 character(len=FMS_FILE_LEN) :: filename !< filename for output file for diagnostic column
103 logical :: allow_ij_input
104 logical :: open_file
105 integer :: io
106
107 integer, parameter :: lkind=fms_cd_kind_
108 real(FMS_CD_KIND_) :: tmp
109!--------------------------------------------------------------------
110! local variables:
111!
112! global_lat latitudes for all diagnostic columns [ degrees ]
113! global_lon longitudes for all diagnostic columns
114! [ degrees ]
115! num_diag_pts total number of diagnostic columns
116! i, j, nn do loop indices
117! char character string for diaganostic column index
118! filename filename for output file for diagnostic column
119!
120!---------------------------------------------------------------------
121
122 if (.not. module_is_initialized) call column_diagnostics_init
123
124!--------------------------------------------------------------------
125! save the input lat and lon fields. define the delta of latitude
126! and longitude.
127!--------------------------------------------------------------------
128 latb_deg = real( real(latb_in,r8_kind)*radian, fms_cd_kind_) !< unit conversion in r8_kind
129 lonb_deg = real( real(lonb_in,r8_kind)*radian, fms_cd_kind_ ) !< unit conversion in r8_kind
130 dellat = latb_in(1,2) - latb_in(1,1)
131 dellon = lonb_in(2,1) - lonb_in(1,1)
132 latb_max = maxval(latb_deg(:,:))
133 latb_min = minval(latb_deg(:,:))
134 lonb_max = maxval(lonb_deg(:,:))
135 lonb_min = minval(lonb_deg(:,:))
136 if (lonb_min < 10.0_lkind .or. lonb_max > 350.0_lkind) then
137 lonb_min = 0.0_lkind
138 lonb_max = 360.0_lkind
139 endif
140
141 allow_ij_input = .true.
142 ref_lat = latb_in(1,1)
143 do i =2,size(latb_in,1)
144 if (latb_in(i,1) /= ref_lat) then
145 allow_ij_input = .false.
146 exit
147 endif
148 end do
149
150 if ( .not. allow_ij_input .and. num_diag_pts_ij /= 0) then
151 call error_mesg ('column_diagnostics_mod', &
152 'cannot specify column diagnostics column with (i,j) &
153 &coordinates when using cubed sphere -- must specify &
154 & lat/lon coordinates', fatal)
155 endif
156
157!----------------------------------------------------------------------
158
159!----------------------------------------------------------------------
160! initialize column_diagnostics flag and diag unit numbers. define
161! total number of diagnostic columns.
162!----------------------------------------------------------------------
163 do_column_diagnostics = .false.
164 diag_units(:) = -1
165 diag_i(:) = -99
166 diag_j(:) = -99
167 diag_lat(:) = -999.0_lkind
168 diag_lon(:) = -999.0_lkind
169 num_diag_pts = size(diag_i(:))
170
171!--------------------------------------------------------------------
172! define an array of lat-lon values for all diagnostic columns.
173!--------------------------------------------------------------------
174 do nn = 1, num_diag_pts_latlon
175 global_lat(nn) = global_lat_latlon(nn)
176 global_lon(nn) = global_lon_latlon(nn)
177 end do
178
179 do nn = 1, num_diag_pts_ij
180 tmp = (-0.5_lkind*acos(-1.0_lkind) + 0.5_lkind*dellat) + real(global_j(nn)-1,fms_cd_kind_)*dellat
181 global_lat(nn+num_diag_pts_latlon) = real( real(tmp,r8_kind)*radian, fms_cd_kind_ )
182 tmp = 0.5_lkind*dellon + real(global_i(nn)-1,fms_cd_kind_)*dellon
183 global_lon(nn+num_diag_pts_latlon) = real( real(tmp,r8_kind)*radian, fms_cd_kind_ )
184 end do
185
186!----------------------------------------------------------------------
187! loop over all diagnostic points to check for their presence on
188! this processor.
189!----------------------------------------------------------------------
190 do nn=1,num_diag_pts
191 open_file = .false.
192
193!----------------------------------------------------------------------
194! verify that the values of lat and lon are valid.
195!----------------------------------------------------------------------
196 if (global_lon(nn) >= 0.0_lkind .and. &
197 global_lon(nn) <= 360.0_lkind) then
198 else
199 call error_mesg ('column_diagnostics_mod', &
200 ' invalid longitude', fatal)
201 endif
202 if (global_lat(nn) >= -90.0_lkind .and. &
203 global_lat(nn) <= 90.0_lkind) then
204 else
205 call error_mesg ('column_diagnostics_mod', &
206 ' invalid latitude', fatal)
207 endif
208
209!--------------------------------------------------------------------
210! if the desired diagnostics column is within the current
211! processor's domain, define the total and coordinate distances from
212! each of the processor's grid points to the diagnostics point.
213!--------------------------------------------------------------------
214
215 if (global_lat(nn) .ge. latb_min .and. &
216 global_lat(nn) .le. latb_max) then
217 if (global_lon(nn) .ge. lonb_min .and.&
218 global_lon(nn) .le. lonb_max) then
219 do j=1,size(latb_deg,2) - 1
220 do i=1,size(lonb_deg,1) - 1
221 distance_y(i,j) = abs(global_lat(nn) - latb_deg(i,j))
222 distance_x(i,j) = abs(global_lon(nn) - lonb_deg(i,j))
223 distance_x2(i,j) = abs((global_lon(nn)-360.0_lkind) - lonb_deg(i,j))
224 distance(i,j) = (global_lat(nn)-latb_deg(i,j))**2 + (global_lon(nn)-lonb_deg(i,j))**2
225 distance2(i,j) = (global_lat(nn)-latb_deg(i,j))**2 + ((global_lon(nn)-360.0_lkind) - lonb_deg(i,j))**2
226 end do
227 end do
228
229!--------------------------------------------------------------------
230! find the grid point on the processor that is within the specified
231! critical distance and also closest to the requested diagnostics
232! column. save the (i,j) coordinates and (lon,lat) of this model
233! grid point. set a flag indicating that a disgnostics file should
234! be opened on this processor for this diagnostic point.
235!--------------------------------------------------------------------
236 current_distance = distance(1,1)
237 do j=1,size(latb_deg,2) - 1
238 do i=1,size(lonb_deg,1) - 1
239 if (distance_x(i,j) <= real(crit_xdistance,fms_cd_kind_) .and. &
240 distance_y(i,j) <= real(crit_ydistance,fms_cd_kind_)) then
241 if (distance(i,j) < current_distance) then
242 current_distance = distance(i,j)
243 do_column_diagnostics(i,j) = .true.
244 diag_j(nn) = j
245 diag_i(nn) = i
246 diag_lon(nn) = lonb_deg(i,j)
247 diag_lat(nn) = latb_deg(i,j)
248 open_file = .true.
249 endif
250 endif
251
252!---------------------------------------------------------------------
253! check needed because of the 0.0 / 360.0 longitude periodicity.
254!---------------------------------------------------------------------
255 if (distance_x2(i,j)<= real(crit_xdistance,fms_cd_kind_) .and. &
256 distance_y(i,j) <= real(crit_ydistance,fms_cd_kind_)) then
257 if (distance2(i,j) < current_distance) then
258 current_distance = distance2(i,j)
259 do_column_diagnostics(i,j) = .true.
260 diag_j(nn) = j
261 diag_i(nn) = i
262 diag_lon(nn) = lonb_deg(i,j)
263 diag_lat(nn) = latb_deg(i,j)
264 open_file = .true.
265 endif
266 endif
267 end do
268 end do
269
270!--------------------------------------------------------------------
271! if the point has been found on this processor, open a diagnostics
272! file.
273!--------------------------------------------------------------------
274 if (open_file) then
275 write (char, '(i2)') nn
276 filename = trim(module) // '_point' // &
277 trim(adjustl(char)) // '.out'
278 if(mpp_npes() > 10000) then
279 write( filename,'(a,i6.6)' )trim(filename)//'.', mpp_pe()-mpp_root_pe()
280 else
281 write( filename,'(a,i4.4)' )trim(filename)//'.', mpp_pe()-mpp_root_pe()
282 endif
283 open(newunit=diag_units(nn), file=trim(filename), action='WRITE', position='rewind', iostat=io)
284 if(io/=0) call error_mesg ('column_diagnostics_mod', 'Error in opening file '//trim(filename), fatal)
285 endif ! (open_file)
286 endif
287 endif
288 end do
289
290!---------------------------------------------------------------------
291
292
293end subroutine initialize_diagnostic_columns_
294
295
296
297
298!####################################################################
299!> @brief column_diagnostics_header writes out information concerning
300!! time and location of following data into the column_diagnostics
301!! output file.
302subroutine column_diagnostics_header_ &
303 (module, diag_unit, time, nn, diag_lon, &
304 diag_lat, diag_i, diag_j)
305
306!--------------------------------------------------------------------
307! column_diagnostics_header writes out information concerning
308! time and location of following data into the column_diagnostics
309! output file.
310!--------------------------------------------------------------------
311
312!--------------------------------------------------------------------
313character(len=*), intent(in) :: module !< module name calling this subroutine
314type(time_type), intent(in) :: Time !< current model time [ time_type ]
315integer, intent(in) :: diag_unit !< unit number for column_diagnostics output
316integer, intent(in) :: nn !< index of diagnostic column currently active
317real(FMS_CD_KIND_), dimension(:), intent(in) :: diag_lon !< longitude of current diagnostic column [ degrees ]
318real(FMS_CD_KIND_), dimension(:), intent(in) :: diag_lat !< latitude of current diagnostic column [ degrees ]
319integer, dimension(:), intent(in) :: diag_i !< i coordinate of current diagnostic column
320integer, dimension(:), intent(in) :: diag_j !< j coordinate of current diagnostic column
321
322!--------------------------------------------------------------------
323! intent(in) variables
324!
325! module module name calling this subroutine
326! Time current model time [ time_type ]
327! diag_unit unit number for column_diagnostics output
328! nn index of diagnostic column currently active
329! diag_lon longitude of current diagnostic column [ degrees ]
330! diag_lat latitude of current diagnostic column [ degrees ]
331! diag_i i coordinate of current diagnostic column
332! diag_j j coordinate of current diagnostic column
333!
334!---------------------------------------------------------------------
335
336
337!--------------------------------------------------------------------
338! local variables:
339
340 integer :: year !< integers defining the current time
341 integer :: month !< integers defining the current time
342 integer :: day !< integers defining the current time
343 integer :: hour !< integers defining the current time
344 integer :: minute !< integers defining the current time
345 integer :: second !< integers defining the current time
346 character(len=9) :: mon !< character string for the current month
347 character(len=64) :: header !< title for the output
348
349!--------------------------------------------------------------------
350! local variables:
351!
352! year, month, day, hour, minute, seconds
353! integers defining the current time
354! mon character string for the current month
355! header title for the output
356!
357!--------------------------------------------------------------------
358
359 if (.not. module_is_initialized) call column_diagnostics_init
360
361!--------------------------------------------------------------------
362! convert the time type to a date and time for printing. convert
363! month to a character string.
364!--------------------------------------------------------------------
365 call get_date (time, year, month, day, hour, minute, second)
366 mon = month_name(month)
367
368!---------------------------------------------------------------------
369! write timestamp and column location information to the diagnostic
370! columns output unit.
371!---------------------------------------------------------------------
372 write (diag_unit,'(a)') ' '
373 write (diag_unit,'(a)') ' '
374 write (diag_unit,'(a)') &
375 '======================================================'
376 write (diag_unit,'(a)') ' '
377 header = ' PRINTING ' // module // ' DIAGNOSTICS'
378 write (diag_unit,'(a)') header
379 write (diag_unit,'(a)') ' '
380 write (diag_unit,'(a, i6,2x, a,i4,i4,i4,i4)') ' time stamp:', &
381 year, trim(mon), day, &
382 hour, minute, second
383 write (diag_unit,'(a, i4)') &
384 ' DIAGNOSTIC POINT COORDINATES, point #', nn
385 write (diag_unit,'(a)') ' '
386 write (diag_unit,'(a,f8.3,a,f8.3)') ' longitude = ', &
387 diag_lon(nn), ' latitude = ', diag_lat(nn)
388 write (diag_unit,'(a, i6, a,i6,a,i6)') &
389 ' on processor # ', mpp_pe(), &
390 ' : processor i =', diag_i(nn), &
391 ' , processor j =', diag_j(nn)
392 write (diag_unit,'(a)') ' '
393
394!---------------------------------------------------------------------
395
396end subroutine column_diagnostics_header_
397!@}