Line data Source code
1 : !-------------------------------------------------------------------------------
2 : ! $Id$
3 : !===============================================================================
4 : module output_grads
5 :
6 :
7 : ! Description:
8 : ! This module contains structure and subroutine definitions to
9 : ! create GrADS output data files for one dimensional arrays.
10 : !
11 : ! The structure type (stat_file) contains all necessay information
12 : ! to generate a GrADS file and a list of variables to be output
13 : ! in the data file.
14 : !
15 : ! References:
16 : ! None
17 : !
18 : ! Original Author:
19 : ! Chris Golaz, updated 2/18/2003
20 : !-------------------------------------------------------------------------------
21 :
22 : use clubb_precision, only: &
23 : core_rknd ! Variable(s)
24 :
25 : implicit none
26 :
27 : public :: open_grads, write_grads
28 :
29 : private :: format_date, check_grads, &
30 : determine_time_inc
31 :
32 : ! Undefined value
33 : real( kind = core_rknd ), private, parameter :: undef = -9.99e33_core_rknd
34 :
35 : private ! Default scope
36 :
37 : contains
38 :
39 : !-------------------------------------------------------------------------------
40 0 : subroutine open_grads( iunit, fdir, fname, &
41 0 : ia, iz, nlat, nlon, z, &
42 0 : day, month, year, lat_vals, lon_vals, &
43 : time, &
44 : stats_metadata, nvar, &
45 : grads_file )
46 : ! Description:
47 : ! Opens and initialize variable components for derived type 'grads_file'
48 : ! If the GrADS file already exists, open_grads will overwrite it.
49 :
50 : ! References:
51 : ! None
52 : !-------------------------------------------------------------------------------
53 : use constants_clubb, only: &
54 : fstderr, & ! Constant(s)
55 : sec_per_min
56 :
57 : use stat_file_module, only: &
58 : stat_file ! Type
59 :
60 : use clubb_precision, only: &
61 : time_precision ! Variable
62 :
63 : use stats_variables, only: &
64 : stats_metadata_type
65 :
66 : implicit none
67 :
68 : ! Input Variables
69 :
70 : integer, intent(in) :: iunit ! File unit being written to [-]
71 :
72 : character(len=*), intent(in) :: &
73 : fdir, & ! Directory where file is stored [-]
74 : fname ! Name of file [-]
75 :
76 : integer, intent(in) :: &
77 : ia, & ! Lower Bound of z (altitude) [-]
78 : iz, & ! Upper Bound of z (altitude) [-]
79 : nlat, & ! Number of points in the y direction (latitude) [-]
80 : nlon ! Number of points in the x direction (longitude) [-]
81 :
82 : real( kind = core_rknd ), dimension(:), intent(in) :: &
83 : z ! Vertical levels [m]
84 :
85 : integer, intent(in) :: &
86 : day, & ! Day of Month at Model Start [dd]
87 : month, & ! Month of Year at Model start [mm]
88 : year ! Year at Model Start [yyyy]
89 :
90 : real( kind = core_rknd ), dimension(nlat), intent(in) :: &
91 : lat_vals ! Latitude [Degrees E]
92 :
93 : real( kind = core_rknd ), dimension(nlon), intent(in) :: &
94 : lon_vals ! Longitude [Degrees N]
95 :
96 : real( kind = time_precision ), intent(in) :: &
97 : time ! Time since Model start [s]
98 :
99 : type (stats_metadata_type), intent(in) :: &
100 : stats_metadata
101 :
102 : ! Number of GrADS variables to store [#]
103 : integer, intent(in) :: nvar
104 :
105 : ! Input/Output Variables
106 : type (stat_file), intent(inout) :: &
107 : grads_file ! File data [-]
108 :
109 : ! Local Variables
110 :
111 : integer :: k
112 : logical :: l_ctl, l_dat, l_error
113 :
114 : ! ---- Begin Code ----
115 :
116 : ! Define parameters for the GrADS ctl and dat files
117 :
118 0 : grads_file%iounit = iunit
119 0 : grads_file%fdir = fdir
120 0 : grads_file%fname = fname
121 0 : grads_file%ia = ia
122 0 : grads_file%iz = iz
123 :
124 : ! Determine if the altitudes are ascending or descending and setup the
125 : ! variable z accordingly.
126 0 : if ( ia <= iz ) then
127 0 : do k=1,iz-ia+1
128 0 : grads_file%z(k) = z(ia+k-1)
129 : end do
130 : else
131 0 : do k=1,ia-iz+1
132 0 : grads_file%z(k) = z(ia-k+1)
133 : end do
134 : end if
135 :
136 0 : grads_file%day = day
137 0 : grads_file%month = month
138 0 : grads_file%year = year
139 :
140 0 : grads_file%nlat = nlat
141 0 : grads_file%nlon = nlon
142 :
143 0 : allocate( grads_file%lat_vals(nlat), grads_file%lon_vals(nlon) )
144 :
145 0 : grads_file%lat_vals = lat_vals
146 0 : grads_file%lon_vals = lon_vals
147 :
148 0 : grads_file%dtwrite = stats_metadata%stats_tout
149 :
150 0 : grads_file%nvar = nvar
151 :
152 : ! Check to make sure the timestep is appropriate. GrADS does not support an
153 : ! output timestep less than 1 minute.
154 0 : if ( stats_metadata%stats_tout < sec_per_min ) then
155 : write(fstderr,*) "Warning: GrADS requires an output timestep of at least &
156 : &one minute, but the requested output timestep &
157 0 : &(stats_tout) is less than one minute."
158 0 : if (.not. stats_metadata%l_allow_small_stats_tout) then
159 : write(fstderr,*) "To override this warning, set l_allow_small_stats_tout = &
160 : &.true. in the stats_setting namelist in the &
161 0 : &appropriate *_model.in file."
162 0 : error stop "Fatal error in open_grads"
163 : end if
164 : end if
165 :
166 : ! Check whether GrADS files already exists
167 :
168 : ! We don't use this feature for the single-column model. The
169 : ! clubb_standalone program will simply overwrite existing data files if they
170 : ! exist. The restart function will create a new GrADS file starting from
171 : ! the restart time in the output directory.
172 :
173 : ! inquire( file=trim(fdir)//trim(fname)//'.ctl', exist=l_ctl )
174 : ! inquire( file=trim(fdir)//trim(fname)//'.dat', exist=l_dat )
175 :
176 0 : l_ctl = .false.
177 0 : l_dat = .false.
178 :
179 : ! If none of the files exist, set ntimes and nrecord and
180 : ! to initial values and return
181 :
182 : if ( .not.l_ctl .and. .not.l_dat ) then
183 :
184 0 : grads_file%time = time
185 0 : grads_file%ntimes = 0
186 0 : grads_file%nrecord = 1
187 : return
188 :
189 : ! If both files exists, attempt to append to existing files
190 :
191 : else if ( l_ctl .and. l_dat ) then
192 :
193 : ! Check existing ctl file
194 :
195 : call check_grads( iunit, fdir, fname, & ! intent(in)
196 : ia, iz, & ! intent(in)
197 : day, month, year, time, stats_metadata%stats_tout, & ! intent(in)
198 : nvar, & ! intent(in)
199 : l_error, grads_file%ntimes, grads_file%nrecord, & ! intent(out)
200 : grads_file%time ) ! intnet(out)
201 :
202 : if ( l_error ) then
203 : write(unit=fstderr,fmt=*) "Error in open_grads:"
204 : write(unit=fstderr,fmt=*) &
205 : "Attempt to append to existing files failed"
206 : ! call stopcode('open_grads')
207 : error stop 'open_grads'
208 : end if
209 :
210 : return
211 :
212 : ! If one file exists, but not the other, give up
213 :
214 : else
215 : write(unit=fstderr,fmt=*) 'Error in open_grads:'
216 : write(unit=fstderr,fmt=*) &
217 : "Attempt to append to existing files failed,"// &
218 : " because only one of the two GrADS files was found."
219 : error stop "open_grads"
220 :
221 : end if
222 :
223 : return
224 0 : end subroutine open_grads
225 :
226 : !-------------------------------------------------------------------------------
227 : subroutine check_grads( iunit, fdir, fname, &
228 : ia, iz, &
229 : day, month, year, time, dtwrite, &
230 : nvar, &
231 : l_error, ntimes, nrecord, time_grads )
232 : ! Description:
233 : ! Given a GrADS file that already exists, this subroutine will attempt
234 : ! to determine whether data can be safely appended to existing file.
235 : ! References:
236 : ! None
237 : !-------------------------------------------------------------------------------
238 : use stat_file_module, only: &
239 : grid_avg_variable ! Type
240 :
241 : use clubb_precision, only: &
242 : time_precision ! Variable
243 :
244 : use constants_clubb, only: &
245 : fstderr, & ! Variable
246 : sec_per_hr, &
247 : sec_per_min
248 :
249 : implicit none
250 :
251 : ! Input Variables
252 :
253 : integer, intent(in) :: &
254 : iunit, & ! Fortran file unit
255 : ia, iz, & ! First and last level
256 : day, month, year, & ! Day, month and year numbers
257 : nvar ! Number of variables in the file
258 :
259 : character(len=*), intent(in) :: &
260 : fdir, fname ! File directory and name
261 :
262 : real( kind = time_precision ), intent(in) :: &
263 : time ! Current model time [s]
264 :
265 : real( kind = core_rknd ), intent(in) :: &
266 : dtwrite ! Time interval between writes to the file [s]
267 :
268 : ! Output Variables
269 : logical, intent(out) :: &
270 : l_error
271 :
272 : integer, intent(out) :: &
273 : ntimes, nrecord
274 :
275 : real(kind=time_precision), intent(out) :: time_grads
276 :
277 : ! Local Variables
278 : logical :: l_done
279 : integer :: ierr
280 : character(len = 256) :: line, tmp, date, dt
281 :
282 : integer :: &
283 : i, nx, ny, nzmax, &
284 : ihour, imin, &
285 : ia_in, iz_in, ntimes_in, nvar_in, &
286 : day_in, month_in, year_in
287 :
288 : real( kind = core_rknd ) :: dtwrite_in
289 :
290 : real( kind = core_rknd ), dimension(:), allocatable :: z_in
291 :
292 : type (grid_avg_variable), dimension(:), allocatable :: var_in
293 :
294 : !-------------------------------------------------------------------------------
295 :
296 : ! ---- Begin Code ----
297 :
298 : ! Initialize logical variables
299 : l_error = .false.
300 : l_done = .false.
301 :
302 : ! Open control file
303 : open( unit = iunit, &
304 : file = trim( fdir )//trim( fname )//'.ctl', &
305 : status = 'old', iostat = ierr )
306 : if ( ierr < 0 ) l_done = .true.
307 :
308 : ! Read and process it
309 :
310 : read(unit=iunit,iostat=ierr,fmt='(a256)') line
311 : if ( ierr < 0 ) l_done = .true.
312 :
313 : do while ( .not. l_done )
314 :
315 : if ( index(line,'XDEF') > 0 ) then
316 :
317 : read(unit=line,fmt=*) tmp, nx
318 : if ( nx /= 1 ) then
319 : write(unit=fstderr,fmt=*) 'Error: XDEF can only be 1'
320 : l_error = .true.
321 : end if
322 :
323 : else if ( index(line,'YDEF') > 0 ) then
324 :
325 : read(unit=line,fmt=*) tmp, ny
326 : if ( ny /= 1 ) then
327 : write(unit=fstderr,fmt=*) "Error: YDEF can only be 1"
328 : l_error = .true.
329 : end if
330 :
331 : else if ( index(line,'ZDEF') > 0 ) then
332 :
333 : read(unit=line,fmt=*) tmp, iz_in
334 :
335 : if ( index(line,'LEVELS') > 0 ) then
336 : ia_in = 1
337 : allocate( z_in(ia_in:iz_in) )
338 : read(unit=iunit,fmt=*) (z_in(i),i=ia_in,iz_in)
339 : end if
340 :
341 : else if ( index(line,'TDEF') > 0 ) then
342 :
343 : read(unit=line,fmt=*) tmp, ntimes_in, tmp, date, dt
344 : read(unit=date(1:2),fmt=*) ihour
345 : read(unit=date(4:5),fmt=*) imin
346 : time_grads = real( ihour, kind=time_precision) * real(sec_per_hr,kind=time_precision) &
347 : + real( imin, kind=time_precision ) * real(sec_per_min,kind=time_precision)
348 : read(unit=date(7:8),fmt=*) day_in
349 : read(unit=date(12:15),fmt=*) year_in
350 :
351 : select case( date(9:11) )
352 : case( 'JAN' )
353 : month_in = 1
354 : case( 'FEB' )
355 : month_in = 2
356 : case( 'MAR' )
357 : month_in = 3
358 : case( 'APR' )
359 : month_in = 4
360 : case( 'MAY' )
361 : month_in = 5
362 : case( 'JUN' )
363 : month_in = 6
364 : case( 'JUL' )
365 : month_in = 7
366 : case( 'AUG' )
367 : month_in = 8
368 : case( 'SEP' )
369 : month_in = 9
370 : case( 'OCT' )
371 : month_in = 10
372 : case( 'NOV' )
373 : month_in = 11
374 : case( 'DEC' )
375 : month_in = 12
376 : case default
377 : write(unit=fstderr,fmt=*) "Unknown month: "//date(9:11)
378 : l_error = .true.
379 : end select
380 :
381 : read(unit=dt(1:len_trim(dt)-2),fmt=*) dtwrite_in
382 : dtwrite_in = dtwrite_in * sec_per_min
383 :
384 : else if ( index(line,'ENDVARS') > 0 ) then
385 :
386 : l_done = .true.
387 :
388 : else if ( index(line,'VARS') > 0 ) then
389 :
390 : read(line,*) tmp, nvar_in
391 : allocate( var_in(nvar_in) )
392 : do i=1, nvar_in
393 : read(unit=iunit,iostat=ierr,fmt='(a256)') line
394 : read(unit=line,fmt=*) var_in(i)%name, nzmax
395 : if ( nzmax /= iz_in ) then
396 : write(unit=fstderr,fmt=*) &
397 : "Error reading ", trim( var_in(i)%name )
398 : l_error = .true.
399 : end if ! nzmax /= iz_in
400 : end do ! 1..nvar_in
401 : end if
402 :
403 : read(unit=iunit,iostat=ierr,fmt='(a256)') line
404 : if ( ierr < 0 ) l_done = .true.
405 :
406 : end do ! while ( .not. l_done )
407 :
408 : close( unit=iunit )
409 :
410 : ! Perform some error check
411 :
412 : if ( abs(ia_in - iz_in) /= abs(ia - iz) ) then
413 : write(unit=fstderr,fmt=*) "check_grads: size mismatch"
414 : l_error = .true.
415 : end if
416 :
417 : if ( day_in /= day ) then
418 : write(unit=fstderr,fmt=*) "check_grads: day mismatch"
419 : l_error = .true.
420 : end if
421 :
422 : if ( month_in /= month ) then
423 : write(unit=fstderr,fmt=*) "check_grads: month mismatch"
424 : l_error = .true.
425 : end if
426 :
427 : if ( year_in /= year ) then
428 : write(unit=fstderr,fmt=*) "check_grads: year mismatch"
429 : l_error = .true.
430 : end if
431 :
432 : if ( int( time_grads ) + ntimes_in*int( dtwrite_in ) &
433 : /= int( time ) ) then
434 : write(unit=fstderr,fmt=*) "check_grads: time mismatch"
435 : l_error = .true.
436 : end if
437 :
438 : if ( int( dtwrite_in ) /= int( dtwrite) ) then
439 : write(unit=fstderr,fmt=*) 'check_grads: dtwrite mismatch'
440 : l_error = .true.
441 : end if
442 :
443 : if ( nvar_in /= nvar ) then
444 : write(unit=fstderr,fmt=*) 'check_grads: nvar mismatch'
445 : l_error = .true.
446 : end if
447 :
448 : if ( l_error ) then
449 : write(unit=fstderr,fmt=*) "check_grads diagnostic"
450 : write(unit=fstderr,fmt=*) "ia = ", ia_in, ia
451 : write(unit=fstderr,fmt=*) "iz = ", iz_in, iz
452 : write(unit=fstderr,fmt=*) "day = ", day_in, day
453 : write(unit=fstderr,fmt=*) "month = ", month_in, month
454 : write(unit=fstderr,fmt=*) "year = ", year_in, year
455 : write(unit=fstderr,fmt=*) "time_grads / time = ", time_grads, time
456 : write(unit=fstderr,fmt=*) "dtwrite = ", dtwrite_in, dtwrite
457 : write(unit=fstderr,fmt=*) "nvar = ", nvar_in, nvar
458 : end if
459 :
460 : ! Set ntimes and nrecord to append to existing files
461 :
462 : ntimes = ntimes_in
463 : nrecord = ntimes_in * nvar_in * iz_in + 1
464 :
465 : deallocate( z_in )
466 :
467 : ! The purpose of this statement is to avoid a compiler warning
468 : ! for tmp
469 : if (tmp =="") then
470 : end if
471 : ! Joshua Fasching June 2008
472 :
473 : return
474 : end subroutine check_grads
475 :
476 : !-------------------------------------------------------------------------------
477 0 : subroutine write_grads( grads_file )
478 :
479 : ! Description:
480 : ! Write part of a GrADS file to data (.dat) file update control file (.ctl.
481 : ! Can be called as many times as necessary
482 : ! References:
483 : ! None
484 : !-------------------------------------------------------------------------------
485 :
486 : use constants_clubb, only: &
487 : fstderr ! Variable(s)
488 :
489 : use model_flags, only: &
490 : l_byteswap_io ! Variable
491 :
492 : use endian, only: &
493 : big_endian, & ! Variable
494 : little_endian
495 :
496 : use stat_file_module, only: &
497 : stat_file ! Type
498 :
499 : ! use stat_file_module, only: &
500 : ! clubb_i, clubb_j ! Variable(s)
501 :
502 : implicit none
503 :
504 : ! External
505 : intrinsic :: selected_real_kind
506 :
507 : ! Constant parameters
508 : integer, parameter :: &
509 : r4 = selected_real_kind( p=5 ) ! Specify 5 decimal digits of precision
510 :
511 : ! Input Variables
512 : type (stat_file), intent(inout) :: &
513 : grads_file ! Contains all information on the files to be written to
514 :
515 : ! Local Variables
516 : integer :: &
517 : ivar, & ! Loop indices
518 : ios ! I/O status indicator
519 :
520 : character(len=15) :: date
521 :
522 : integer :: dtwrite_ctl ! Time increment for the ctl file
523 : character(len=2) :: dtwrite_units ! Units on dtwrite_ctl
524 :
525 : ! ---- Begin Code ----
526 : ! Check number of variables and write nothing if less than 1
527 :
528 0 : if ( grads_file%nvar < 1 ) return
529 :
530 : #include "recl.inc"
531 :
532 : ! Output data to file
533 : open( unit=grads_file%iounit, &
534 : file=trim( grads_file%fdir )//trim( grads_file%fname )//'.dat', &
535 : form='unformatted', access='direct', &
536 : recl=F_RECL*abs( grads_file%iz-grads_file%ia+1 )*grads_file%nlon*grads_file%nlat, &
537 0 : status='unknown', iostat=ios )
538 0 : if ( ios /= 0 ) then
539 : write(unit=fstderr,fmt=*) &
540 0 : "write_grads: error opening binary file"
541 0 : write(unit=fstderr,fmt=*) "iostat = ", ios
542 0 : error stop
543 : end if
544 :
545 0 : if ( grads_file%ia <= grads_file%iz ) then
546 0 : do ivar=1,grads_file%nvar
547 : write(grads_file%iounit,rec=grads_file%nrecord) &
548 0 : real( grads_file%grid_avg_var(ivar)%ptr(1:grads_file%nlon, &
549 0 : 1:grads_file%nlat,grads_file%ia:grads_file%iz), kind=r4)
550 0 : grads_file%nrecord = grads_file%nrecord + 1
551 : end do
552 :
553 : else
554 0 : do ivar=1, grads_file%nvar
555 : write(grads_file%iounit,rec=grads_file%nrecord) &
556 0 : real( grads_file%grid_avg_var(ivar)%ptr(1:grads_file%nlon, &
557 0 : 1:grads_file%nlat,grads_file%ia:grads_file%iz:-1), kind=r4)
558 0 : grads_file%nrecord = grads_file%nrecord + 1
559 : end do
560 :
561 : end if ! grads_file%ia <= grads_file%iz
562 :
563 0 : close( unit=grads_file%iounit, iostat = ios )
564 :
565 0 : if ( ios /= 0 ) then
566 : write(unit=fstderr,fmt=*) &
567 0 : "write_grads: error closing binary file"
568 0 : write(unit=fstderr,fmt=*) "iostat = ", ios
569 0 : error stop
570 : end if
571 :
572 0 : grads_file%ntimes = grads_file%ntimes + 1
573 :
574 : ! Write control file
575 :
576 : open(unit=grads_file%iounit, &
577 : file=trim( grads_file%fdir )//trim( grads_file%fname )//'.ctl', &
578 0 : status='unknown', iostat=ios)
579 0 : if ( ios > 0 ) then
580 : write(unit=fstderr,fmt=*) &
581 0 : "write_grads: error opening control file"
582 0 : write(unit=fstderr,fmt=*) "iostat = ", ios
583 0 : error stop
584 : end if
585 :
586 : ! Write file header
587 : if ( ( big_endian .and. .not. l_byteswap_io ) &
588 : .or. ( little_endian .and. l_byteswap_io ) ) then
589 : write(unit=grads_file%iounit,fmt='(a)') 'OPTIONS BIG_ENDIAN'
590 :
591 : else
592 0 : write(unit=grads_file%iounit,fmt='(a)') 'OPTIONS LITTLE_ENDIAN'
593 :
594 : end if
595 :
596 0 : write(unit=grads_file%iounit,fmt='(a)') 'DSET ^'//trim( grads_file%fname )//'.dat'
597 0 : write(unit=grads_file%iounit,fmt='(a,e12.5)') 'UNDEF ',undef
598 :
599 0 : if ( grads_file%nlon == 1 ) then ! Use linear for a singleton X dimesion
600 0 : write(unit=grads_file%iounit,fmt='(a,f8.3,a)') 'XDEF 1 LINEAR ', grads_file%lon_vals, ' 1.'
601 : else
602 0 : write(unit=grads_file%iounit,fmt='(a,i5,a)') 'XDEF', grads_file%nlon,' LEVELS '
603 0 : write(unit=grads_file%iounit,fmt='(6f13.4)') grads_file%lon_vals
604 : end if
605 :
606 0 : if ( grads_file%nlat == 1 ) then ! Use linear for a singleton Y dimension
607 0 : write(unit=grads_file%iounit,fmt='(a,f8.3,a)') 'YDEF 1 LINEAR ', grads_file%lat_vals, ' 1.'
608 : else
609 0 : write(unit=grads_file%iounit,fmt='(a,i5,a)') 'YDEF', grads_file%nlat,' LEVELS '
610 0 : write(unit=grads_file%iounit,fmt='(6f13.4)') grads_file%lat_vals
611 : end if
612 :
613 0 : if ( grads_file%ia == grads_file%iz ) then ! If ia == iz, then Z is also singleton
614 0 : write(unit=grads_file%iounit,fmt='(a)') 'ZDEF 1 LEVELS 0.'
615 0 : else if ( grads_file%ia < grads_file%iz ) then
616 : write(unit=grads_file%iounit,fmt='(a,i5,a)') &
617 0 : 'ZDEF', abs(grads_file%iz-grads_file%ia)+1,' LEVELS '
618 : write(unit=grads_file%iounit,fmt='(6f13.4)') &
619 0 : (grads_file%z(ivar-grads_file%ia+1),ivar=grads_file%ia,grads_file%iz)
620 : else
621 : write(unit=grads_file%iounit,fmt='(a,i5,a)') &
622 0 : 'ZDEF',abs(grads_file%iz-grads_file%ia)+1,' LEVELS '
623 0 : write(grads_file%iounit,'(6f13.4)') (grads_file%z(grads_file%ia-ivar+1), &
624 0 : ivar=grads_file%ia,grads_file%iz,-1)
625 : end if
626 :
627 : call format_date( grads_file%day, grads_file%month, grads_file%year, grads_file%time, & ! In
628 0 : date ) ! Out
629 :
630 : call determine_time_inc( grads_file%dtwrite, & ! In
631 0 : dtwrite_ctl, dtwrite_units ) ! Out
632 :
633 0 : write(unit=grads_file%iounit,fmt='(a,i6,a,a,i5,a)') 'TDEF ', &
634 0 : grads_file%ntimes, ' LINEAR ', date, dtwrite_ctl, dtwrite_units
635 :
636 : ! Variables description
637 0 : write(unit=grads_file%iounit,fmt='(a,i5)') 'VARS', grads_file%nvar
638 :
639 0 : do ivar=1, grads_file%nvar, 1
640 0 : write(unit=grads_file%iounit,fmt='(a,i5,a,a)') &
641 0 : grads_file%grid_avg_var(ivar)%name(1:len_trim(grads_file%grid_avg_var(ivar)%name)), &
642 0 : abs(grads_file%iz-grads_file%ia)+1,' 99 ', &
643 : grads_file%grid_avg_var(ivar)%description( &
644 0 : 1:len_trim(grads_file%grid_avg_var(ivar)%description))
645 : end do
646 :
647 0 : write(unit=grads_file%iounit,fmt='(a)') 'ENDVARS'
648 :
649 0 : close( unit=grads_file%iounit, iostat=ios )
650 0 : if ( ios > 0 ) then
651 : write(unit=fstderr,fmt=*) &
652 0 : "write_grads: error closing control file"
653 0 : write(unit=fstderr,fmt=*) "iostat = ",ios
654 0 : error stop
655 : end if
656 :
657 : return
658 : end subroutine write_grads
659 :
660 : !---------------------------------------------------------
661 0 : subroutine format_date( day_in, month_in, year_in, time_in, &
662 0 : date )
663 : !
664 : ! Description:
665 : ! This subroutine formats the current time of the model (given in seconds
666 : ! since the start time) to a date format usable as GrADS output.
667 : ! References:
668 : ! None
669 : !---------------------------------------------------------
670 : use clubb_precision, only: &
671 : time_precision ! Variable(s)
672 :
673 : use calendar, only: &
674 : compute_current_date ! Procedure(s)
675 :
676 : use calendar, only: &
677 : month_names ! Variable(s)
678 :
679 : use constants_clubb, only: &
680 : sec_per_hr, & ! Variable(s)
681 : min_per_hr
682 :
683 : implicit none
684 :
685 : ! Input Variables
686 : integer, intent(in) :: &
687 : day_in, & ! Day of the Month at Model Start [dd]
688 : month_in, & ! Month of the Year at Model Start [mm]
689 : year_in ! Year at Model Start [yyyy]
690 :
691 : real(kind=time_precision), intent(in) :: &
692 : time_in ! Time since Model Start [s]
693 :
694 : ! Output Variables
695 : character(len=15), intent(out) :: &
696 : date ! Current Date in format 'hh:mmZddmmmyyyy'
697 :
698 : ! Local Variables
699 : integer :: iday, imonth, iyear ! Day, month, year
700 : real(kind=time_precision) :: time ! time [s]
701 :
702 : ! ---- Begin Code ----
703 :
704 : ! Copy input arguments into local variables
705 :
706 0 : iday = day_in
707 0 : imonth = month_in
708 0 : iyear = year_in
709 0 : time = time_in
710 :
711 : call compute_current_date( day_in, month_in, & ! In
712 : year_in, & ! In
713 : time_in, & ! In
714 : iday, imonth, & ! Out
715 : iyear, & ! Out
716 0 : time ) ! Out
717 :
718 0 : date = 'hh:mmZddmmmyyyy'
719 0 : write(unit=date(7:8),fmt='(i2.2)') iday
720 0 : write(unit=date(9:11),fmt='(a3)') month_names(imonth)
721 0 : write(unit=date(12:15),fmt='(i4.4)') iyear
722 0 : write(unit=date(1:2),fmt='(i2.2)') floor(time/real(sec_per_hr,kind=time_precision ))
723 : write(unit=date(4:5),fmt='(i2.2)') &
724 0 : int( mod( nint( time ), nint(sec_per_hr) ) / nint(min_per_hr) )
725 :
726 0 : return
727 : end subroutine format_date
728 :
729 : !-------------------------------------------------------------------------------
730 0 : subroutine determine_time_inc( dtwrite_sec, &
731 0 : dtwrite_ctl, units )
732 : ! Description:
733 : ! Determine the units on the time increment, since GrADS only allows a 2 digit
734 : ! time increment.
735 : ! References:
736 : ! None
737 : !-------------------------------------------------------------------------------
738 : use constants_clubb, only: &
739 : sec_per_day, & ! Constants
740 : sec_per_hr, &
741 : sec_per_min
742 :
743 :
744 : implicit none
745 :
746 : ! External
747 : intrinsic :: max, floor
748 :
749 : ! Input Variables
750 : real(kind=core_rknd), intent(in) :: &
751 : dtwrite_sec ! Time increment in GrADS [s]
752 :
753 : ! Output Variables
754 : integer, intent(out) :: &
755 : dtwrite_ctl ! Time increment in GrADS [units vary]
756 :
757 : character(len=2), intent(out) :: units ! Units on dtwrite_ctl
758 :
759 : ! Local variables
760 : real(kind=core_rknd) :: &
761 : dtwrite_min, & ! Time increment [minutes]
762 : dtwrite_hrs, & ! Time increment [hours]
763 : dtwrite_days ! Time increment [days]
764 :
765 : ! ---- Begin Code ----
766 :
767 : ! Since GrADs can't handle a time increment of less than a minute we assume
768 : ! 1 minute output for an output frequency of less than a minute.
769 0 : dtwrite_min = real( floor( dtwrite_sec/sec_per_min ), kind=core_rknd )
770 0 : dtwrite_min = max( 1._core_rknd, dtwrite_min )
771 :
772 0 : if ( dtwrite_min <= 99._core_rknd ) then
773 0 : dtwrite_ctl = int( dtwrite_min )
774 0 : units = 'mn'
775 : else
776 0 : dtwrite_hrs = dtwrite_sec / sec_per_hr
777 0 : if ( dtwrite_hrs <= 99._core_rknd ) then
778 0 : dtwrite_ctl = int( dtwrite_hrs )
779 0 : units = 'hr'
780 : else
781 0 : dtwrite_days = dtwrite_sec / sec_per_day
782 0 : if ( dtwrite_days <= 99._core_rknd ) then
783 0 : dtwrite_ctl = int( dtwrite_days )
784 0 : units = 'dy'
785 : else
786 0 : error stop "Fatal error in determine_time_inc"
787 : end if ! dwrite_days <= 99.
788 : end if ! dtwrite_hrs <= 99.
789 : end if ! dtwrite_min <= 99.
790 :
791 0 : return
792 : end subroutine determine_time_inc
793 :
794 : end module output_grads
795 : !-------------------------------------------------------------------------------
|