ROMS
Loading...
Searching...
No Matches
get_nudgcoef.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#ifndef ANA_NUDGCOEF
5!
6!git $Id$
7!================================================== Hernan G. Arango ===
8! Copyright (c) 2002-2025 The ROMS Group !
9! Licensed under a MIT/X style license !
10! See License_ROMS.md !
11!=======================================================================
12! !
13! This module reads grid time scales for nudging to climatology file !
14! using either the standard NetCDF library or the Parallel-IO (PIO) !
15! library. !
16! !
17!=======================================================================
18!
19 USE mod_param
20 USE mod_parallel
21 USE mod_clima
22 USE mod_grid
23 USE mod_iounits
24 USE mod_ncparam
25 USE mod_scalars
26!
28# ifdef SOLVE3D
30# endif
31# ifdef DISTRIBUTE
33# ifdef SOLVE3D
36# endif
37# endif
38 USE nf_fread2d_mod, ONLY : nf_fread2d
39# ifdef SOLVE3D
40 USE nf_fread3d_mod, ONLY : nf_fread3d
41# endif
43!
44 implicit none
45!
46 PUBLIC :: get_nudgcoef
47 PRIVATE :: get_nudgcoef_nf90
48# if defined PIO_LIB && defined DISTRIBUTE
49 PRIVATE :: get_nudgcoef_pio
50# endif
51!
52 CONTAINS
53!
54!***********************************************************************
55 SUBROUTINE get_nudgcoef (ng, tile, model)
56!***********************************************************************
57!
58! Imported variable declarations.
59!
60 integer, intent(in) :: ng, tile, model
61!
62! Local variable declarations.
63!
64 integer :: lbi, ubi, lbj, ubj
65!
66 character (len=*), parameter :: myfile = &
67 & __FILE__
68!
69!-----------------------------------------------------------------------
70! Read in nudging coefficients file according to IO type.
71!-----------------------------------------------------------------------
72!
73 lbi=bounds(ng)%LBi(tile)
74 ubi=bounds(ng)%UBi(tile)
75 lbj=bounds(ng)%LBj(tile)
76 ubj=bounds(ng)%UBj(tile)
77!
78 SELECT CASE (nud(ng)%IOtype)
79 CASE (io_nf90)
80 CALL get_nudgcoef_nf90 (ng, tile, model, &
81 & lbi, ubi, lbj, ubj)
82
83# if defined PIO_LIB && defined DISTRIBUTE
84 CASE (io_pio)
85 CALL get_nudgcoef_pio (ng, tile, model, &
86 & lbi, ubi, lbj, ubj)
87# endif
88 CASE DEFAULT
89 IF (master) WRITE (stdout,10) nud(ng)%IOtype
90 exit_flag=2
91 END SELECT
92 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
93!
94 10 FORMAT (' GET_NUDGCOEF - Illegal input file type, io_type = ',i0, &
95 & /,16x,'Check KeyWord ''INP_LIB'' in ''roms.in''.')
96!
97 RETURN
98 END SUBROUTINE get_nudgcoef
99!
100!***********************************************************************
101 SUBROUTINE get_nudgcoef_nf90 (ng, tile, model, &
102 & LBi, UBi, LBj, UBj)
103!***********************************************************************
104!
105 USE mod_netcdf
106!
107! Imported variable declarations.
108!
109 integer, intent(in) :: ng, tile, model
110 integer, intent(in) :: lbi, ubi, lbj, ubj
111!
112! Local variable declarations.
113!
114# ifdef SOLVE3D
115 logical :: got_generic
116 logical :: got_specific(nt(ng))
117!
118# endif
119 integer :: gtype, i, ic, itrc
120 integer :: nvatt, nvdim, status, vindex
121 integer :: vsize(4)
122# ifdef CHECKSUM
123 integer(i8b) :: fhash
124# endif
125!
126 real(dp) :: fscl
127 real(r8) :: fmax, fmin
128!
129 character (len=40 ) :: tunits
130 character (len=256) :: ncname
131
132 character (len=*), parameter :: myfile = &
133 & __FILE__//", get_nudgcoef_nf90"
134!
135 sourcefile=myfile
136!
137!-----------------------------------------------------------------------
138! Inquire about the contents of grid NetCDF file: Inquire about
139! the dimensions and variables. Check for consistency.
140!-----------------------------------------------------------------------
141!
142 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
143 ncname=nud(ng)%name
144!
145! Open grid NetCDF file for reading.
146!
147 IF (nud(ng)%ncid.eq.-1) THEN
148 CALL netcdf_open (ng, model, ncname, 0, nud(ng)%ncid)
149 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
150 WRITE (stdout,10) trim(ncname)
151 RETURN
152 END IF
153 END IF
154!
155! Check grid file dimensions for consitency
156!
157 CALL netcdf_check_dim (ng, model, ncname, &
158 & ncid = nud(ng)%ncid)
159 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
160!
161! Inquire about the variables.
162!
163 CALL netcdf_inq_var (ng, model, ncname, &
164 & ncid = nud(ng)%ncid)
165 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
166!
167!-----------------------------------------------------------------------
168! Check if required variables are available.
169!-----------------------------------------------------------------------
170!
171! Nudging coefficients for 2D momentum.
172!
173 IF (lnudgem2clm(ng)) THEN
174 IF (.not.find_string(var_name,n_var,vname(1,idm2nc), &
175 & nud(ng)%Vid(idm2nc))) THEN
176 IF (master) WRITE (stdout,20) trim(vname(1,idm2nc)), &
177 & trim(ncname)
178 exit_flag=2
179 RETURN
180 END IF
181 END IF
182
183# ifdef SOLVE3D
184!
185! Nudging coefficients for 3D momentum.
186!
187 IF (lnudgem3clm(ng)) THEN
188 IF (.not.find_string(var_name,n_var,vname(1,idm3nc), &
189 & nud(ng)%Vid(idm3nc))) THEN
190 IF (master) WRITE (stdout,20) trim(vname(1,idm3nc)), &
191 & trim(ncname)
192 exit_flag=2
193 RETURN
194 END IF
195 END IF
196!
197! Tracers nudging coefficients.
198!
199 IF (any(lnudgetclm(:,ng))) THEN
200 got_generic=find_string(var_name,n_var,vname(1,idgtnc), &
201 & nud(ng)%Vid(idgtnc))
202 DO itrc=1,nt(ng)
203 IF (lnudgetclm(itrc,ng)) THEN
204 vindex=idtnud(itrc)
205 got_specific(itrc)=find_string(var_name,n_var, &
206 & vname(1,vindex), &
207 & nud(ng)%Vid(vindex))
208 IF (.not.(got_generic.or.got_specific(itrc))) THEN
209 IF (master) WRITE (stdout,30) trim(vname(1,idgtnc)), &
210 & trim(vname(1,vindex)), &
211 & trim(ncname)
212 exit_flag=2
213 RETURN
214 END IF
215 END IF
216 END DO
217 END IF
218# endif
219!
220!-----------------------------------------------------------------------
221! Read in nudging coefficients.
222!-----------------------------------------------------------------------
223!
224! Set Vsize to zero to deativate interpolation of input data to model
225! grid in "nf_fread2d".
226!
227 DO i=1,4
228 vsize(i)=0
229 END DO
230!
231! If appropriate, read nudging coefficients for 2D momentum (inverse
232! time scales, 1/day).
233!
234 IF (master) WRITE (stdout,'(1x)')
235!
236 IF (lnudgem2clm(ng)) THEN
237 gtype=r2dvar
238 vindex=idm2nc
239 fscl=1.0_dp/day2sec ! default units: 1/day
240!
241 CALL netcdf_inq_var (ng, model, ncname, &
242 & ncid = nud(ng)%ncid, &
243 & myvarname = trim(vname(1,vindex)), &
244 & varid = nud(ng)%Vid(vindex), &
245 & nvardim = nvdim, &
246 & nvaratt = nvatt)
247 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
248 DO i=1,nvatt
249 IF (trim(var_aname(i)).eq.'units') THEN
250 tunits=trim(var_achar(i))
251 IF (tunits(1:3).eq.'day') THEN
252 fscl=1.0_dp/day2sec
253 ELSE IF (tunits(1:6).eq.'second') THEN
254 fscl=1.0_dp
255 END IF
256 END IF
257 END DO
258!
259 status=nf_fread2d(ng, model, ncname, nud(ng)%ncid, &
260 & vname(1,vindex), nud(ng)%Vid(vindex), &
261 & 0, gtype, vsize, &
262 & lbi, ubi, lbj, ubj, &
263 & fscl, fmin, fmax, &
264# ifdef MASKING
265 & grid(ng) % rmask, &
266# endif
267# ifdef CHECKSUM
268 & clima(ng) % M2nudgcof, &
269 & checksum = fhash)
270# else
271 & clima(ng) % M2nudgcof)
272# endif
273 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
274 exit_flag=2
275 ioerror=status
276 RETURN
277 ELSE
278 IF (master) THEN
279 WRITE (stdout,40) trim(vname(2,vindex))//': '// &
280 & trim(vname(1,vindex)), &
281 & ng, trim(ncname), fmin, fmax
282# ifdef CHECKSUM
283 WRITE (stdout,50) fhash
284# endif
285 END IF
286 END IF
287!
288 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
289 CALL exchange_r2d_tile (ng, tile, &
290 & lbi, ubi, lbj, ubj, &
291 & clima(ng) % M2nudgcof)
292 END IF
293# ifdef DISTRIBUTE
294 CALL mp_exchange2d (ng, tile, model, 1, &
295 & lbi, ubi, lbj, ubj, &
296 & nghostpoints, &
297 & ewperiodic(ng), nsperiodic(ng), &
298 & clima(ng) % M2nudgcof)
299# endif
300 END IF
301
302# ifdef SOLVE3D
303!
304! If appropriate, read nudging coefficients for 3D momentum (inverse
305! time scales, 1/day).
306!
307 IF (lnudgem3clm(ng)) THEN
308 gtype=r3dvar
309 vindex=idm3nc
310 fscl=1.0_dp/day2sec ! default units: 1/day
311!
312 CALL netcdf_inq_var (ng, model, ncname, &
313 & ncid = nud(ng)%ncid, &
314 & myvarname = trim(vname(1,vindex)), &
315 & varid = nud(ng)%Vid(vindex), &
316 & nvardim = nvdim, &
317 & nvaratt = nvatt)
318 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
319 DO i=1,nvatt
320 IF (trim(var_aname(i)).eq.'units') THEN
321 tunits=trim(var_achar(i))
322 IF (tunits(1:3).eq.'day') THEN
323 fscl=1.0_dp/day2sec
324 ELSE IF (tunits(1:6).eq.'second') THEN
325 fscl=1.0_dp
326 END IF
327 END IF
328 END DO
329!
330 status=nf_fread3d(ng, model, ncname, nud(ng)%ncid, &
331 & vname(1,vindex), nud(ng)%Vid(vindex), &
332 & 0, gtype, vsize, &
333 & lbi, ubi, lbj, ubj, 1, n(ng), &
334 & fscl, fmin, fmax, &
335# ifdef MASKING
336 & grid(ng) % rmask, &
337# endif
338# ifdef CHECKSUM
339 & clima(ng) % M3nudgcof, &
340 & checksum = fhash)
341# else
342 & clima(ng) % M3nudgcof)
343# endif
344 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
345 exit_flag=2
346 ioerror=status
347 RETURN
348 ELSE
349 IF (master) THEN
350 WRITE (stdout,40) trim(vname(2,vindex))//': '// &
351 & trim(vname(1,vindex)), &
352 & ng, trim(ncname), fmin, fmax
353# ifdef CHECKSUM
354 WRITE (stdout,50) fhash
355# endif
356 END IF
357 END IF
358!
359 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
360 CALL exchange_r3d_tile (ng, tile, &
361 & lbi, ubi, lbj, ubj, 1, n(ng), &
362 & clima(ng) % M3nudgcof)
363 END IF
364# ifdef DISTRIBUTE
365 CALL mp_exchange3d (ng, tile, model, 1, &
366 & lbi, ubi, lbj, ubj, 1, n(ng), &
367 & nghostpoints, &
368 & ewperiodic(ng), nsperiodic(ng), &
369 & clima(ng) % M3nudgcof)
370# endif
371 END IF
372!
373! If appropriate, read nudging coefficients for tracer variables
374! (inverse time scales, 1/day). If the nudging coefficients for a
375! specific tracer are available, it will read that NetCDF variable.
376! If NOT AND the generic coefficients are available, it will process
377! those values instead.
378!
379 ic=0
380 DO itrc=1,nt(ng)
381 IF (lnudgetclm(itrc,ng)) THEN
382 ic=ic+1
383 gtype=r3dvar
384 IF (got_specific(itrc)) THEN
385 vindex=idtnud(itrc) ! specific variable
386 fscl=1.0_dp/day2sec ! default units: 1/day
387 ELSE IF (got_generic) THEN
388 vindex=idgtnc ! generic variable
389 fscl=1.0_dp/day2sec ! default units: 1/day
390 END IF
391!
392 CALL netcdf_inq_var (ng, model, ncname, &
393 & ncid = nud(ng)%ncid, &
394 & myvarname = trim(vname(1,vindex)), &
395 & varid = nud(ng)%Vid(vindex), &
396 & nvardim = nvdim, &
397 & nvaratt = nvatt)
398 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
399 DO i=1,nvatt
400 IF (trim(var_aname(i)).eq.'units') THEN
401 tunits=trim(var_achar(i))
402 IF (tunits(1:3).eq.'day') THEN
403 fscl=1.0_dp/day2sec
404 ELSE IF (tunits(1:6).eq.'second') THEN
405 fscl=1.0_dp
406 END IF
407 END IF
408 END DO
409!
410 status=nf_fread3d(ng, model, ncname, nud(ng)%ncid, &
411 & vname(1,vindex), nud(ng)%Vid(vindex), &
412 & 0, gtype, vsize, &
413 & lbi, ubi, lbj, ubj, 1, n(ng), &
414 & fscl, fmin, fmax, &
415# ifdef MASKING
416 & grid(ng) % rmask, &
417# endif
418# ifdef CHECKSUM
419 & clima(ng) % Tnudgcof(:,:,:,ic), &
420 & checksum = fhash)
421# else
422 & clima(ng) % Tnudgcof(:,:,:,ic))
423# endif
424 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
425 exit_flag=2
426 ioerror=status
427 RETURN
428 ELSE
429 IF (master) THEN
430 WRITE (stdout,40) trim(vname(2,vindex))//': '// &
431 & trim(vname(1,vindex)), &
432 & ng, trim(ncname), fmin, fmax
433# ifdef CHECKSUM
434 WRITE (stdout,50) fhash
435# endif
436 END IF
437 END IF
438!
439 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
440 CALL exchange_r3d_tile (ng, tile, &
441 & lbi, ubi, lbj, ubj, 1, n(ng), &
442 & clima(ng) % Tnudgcof(:,:,:,ic))
443 END IF
444 END IF
445 END DO
446# ifdef DISTRIBUTE
447 IF (any(lnudgetclm(:,ng))) THEN
448 CALL mp_exchange4d (ng, tile, model, 1, &
449 & lbi, ubi, lbj, ubj, 1, n(ng), 1, ntclm(ng), &
450 & nghostpoints, &
451 & ewperiodic(ng), nsperiodic(ng), &
452 & clima(ng) % Tnudgcof)
453 END IF
454# endif
455# endif
456!
457 10 FORMAT (/,' GET_NUDGCOEF_NF90 - unable to open nudging NetCDF', &
458 & ' file: ',a)
459 20 FORMAT (/,' GET_NUDGCOEF_NF90 - unable to find nudging', &
460 & ' variable: ',a,/,21x,'in NetCDF file: ',a)
461 30 FORMAT (/,' GET_NUDGCOEF_NF90 - unable to find nudging', &
462 & ' variable: ',a, &
463 & /,21x,'or generic nudging variable: ',a, &
464 & /,21x,'in nudging NetCDF file: ',a)
465 40 FORMAT(2x,' GET_NUDGCOEF_NF90 - ',a,/,21x, &
466 & '(Grid = ',i2.2,', File: ',a,')',/,21x, &
467 & '(Min = ', 1p,e15.8,0p,' Max = ',1p,e15.8,0p,')')
468# ifdef CHECKSUM
469 50 FORMAT (19x,'(CheckSum = ',i0,')')
470# endif
471!
472 RETURN
473 END SUBROUTINE get_nudgcoef_nf90
474
475# if defined PIO_LIB && defined DISTRIBUTE
476!
477!***********************************************************************
478 SUBROUTINE get_nudgcoef_pio (ng, tile, model, &
479 & LBi, UBi, LBj, UBj)
480!***********************************************************************
481!
483!
484! Imported variable declarations.
485!
486 integer, intent(in) :: ng, tile, model
487 integer, intent(in) :: lbi, ubi, lbj, ubj
488!
489! Local variable declarations.
490!
491# ifdef SOLVE3D
492 logical :: got_generic
493 logical :: got_specific(nt(ng))
494!
495# endif
496 integer :: gtype, i, ic, ifield, itrc
497 integer :: nvatt, nvdim, status, vindex
498 integer :: vsize(4)
499# ifdef CHECKSUM
500 integer(i8b) :: fhash
501# endif
502!
503 real(dp) :: fscl
504 real(r8) :: fmax, fmin
505!
506 character (len=40 ) :: tunits
507 character (len=256) :: ncname
508
509 character (len=*), parameter :: myfile = &
510 & __FILE__//", get_nudgcoef_nf90"
511!
512 TYPE (io_desc_t), pointer :: iodesc
513!
514 sourcefile=myfile
515!
516!-----------------------------------------------------------------------
517! Inquire about the contents of grid NetCDF file: Inquire about
518! the dimensions and variables. Check for consistency.
519!-----------------------------------------------------------------------
520!
521 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
522 ncname=nud(ng)%name
523!
524! Open grid NetCDF file for reading.
525!
526 IF (nud(ng)%pioFile%fh.eq.-1) THEN
527 CALL pio_netcdf_open (ng, model, ncname, 0, nud(ng)%pioFile)
528 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
529 WRITE (stdout,10) trim(ncname)
530 RETURN
531 END IF
532 END IF
533!
534! Check grid file dimensions for consitency
535!
536 CALL pio_netcdf_check_dim (ng, model, ncname, &
537 & piofile = nud(ng)%pioFile)
538 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
539!
540! Inquire about the variables.
541!
542 CALL pio_netcdf_inq_var (ng, model, ncname, &
543 & piofile = nud(ng)%pioFile)
544 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
545!
546!-----------------------------------------------------------------------
547! Check if required variables are available.
548!-----------------------------------------------------------------------
549!
550! Nudging coefficients for 2D momentum.
551!
552 IF (lnudgem2clm(ng)) THEN
553 IF (.not.find_string(var_name,n_var,vname(1,idm2nc), &
554 & vindex)) THEN
555 IF (master) WRITE (stdout,20) trim(vname(1,idm2nc)), &
556 & trim(ncname)
557 exit_flag=2
558 RETURN
559 END IF
560 nud(ng)%pioVar(idm2nc)%vd=var_desc(vindex)
561 END IF
562
563# ifdef SOLVE3D
564!
565! Nudging coefficients for 3D momentum.
566!
567 IF (lnudgem3clm(ng)) THEN
568 IF (.not.find_string(var_name,n_var,vname(1,idm3nc), &
569 & vindex)) THEN
570 IF (master) WRITE (stdout,20) trim(vname(1,idm3nc)), &
571 & trim(ncname)
572 exit_flag=2
573 RETURN
574 END IF
575 nud(ng)%pioVar(idm3nc)%vd=var_desc(vindex)
576 END IF
577!
578! Tracers nudging coefficients.
579!
580 IF (any(lnudgetclm(:,ng))) THEN
581 IF (find_string(var_name,n_var,vname(1,idgtnc),vindex)) THEN
582 got_generic=.true.
583 nud(ng)%pioVar(idgtnc)%vd=var_desc(vindex)
584 ELSE
585 got_generic=.false.
586 END IF
587!
588 DO itrc=1,nt(ng)
589 IF (lnudgetclm(itrc,ng)) THEN
590 ifield=idtnud(itrc)
591 IF (find_string(var_name,n_var,vname(1,ifield),vindex)) THEN
592 got_specific(itrc)=.true.
593 nud(ng)%pioVar(ifield)%vd=var_desc(vindex)
594 ELSE
595 got_specific(itrc)=.false.
596 END IF
597!
598 IF (.not.(got_generic.or.got_specific(itrc))) THEN
599 IF (master) WRITE (stdout,30) trim(vname(1,idgtnc)), &
600 & trim(vname(1,ifield)), &
601 & trim(ncname)
602 exit_flag=2
603 RETURN
604 END IF
605 END IF
606 END DO
607 END IF
608# endif
609!
610!-----------------------------------------------------------------------
611! Read in nudging coefficients.
612!-----------------------------------------------------------------------
613!
614! Set Vsize to zero to deativate interpolation of input data to model
615! grid in "nf_fread2d".
616!
617 DO i=1,4
618 vsize(i)=0
619 END DO
620!
621! If appropriate, read nudging coefficients for 2D momentum (inverse
622! time scales, 1/day).
623!
624 IF (master) WRITE (stdout,'(1x)')
625!
626 IF (lnudgem2clm(ng)) THEN
627 ifield=idm2nc
628 fscl=1.0_dp/day2sec ! default units: 1/day
629 nud(ng)%pioVar(ifield)%gtype=r2dvar
630!
631 CALL pio_netcdf_inq_var (ng, model, ncname, &
632 & piofile = nud(ng)%pioFile, &
633 & myvarname = trim(vname(1,ifield)), &
634 & piovar = nud(ng)%pioVar(ifield)%vd, &
635 & nvardim = nvdim, &
636 & nvaratt = nvatt)
637 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
638!
639 DO i=1,nvatt
640 IF (trim(var_aname(i)).eq.'units') THEN
641 tunits=trim(var_achar(i))
642 IF (tunits(1:3).eq.'day') THEN
643 fscl=1.0_dp/day2sec
644 ELSE IF (tunits(1:6).eq.'second') THEN
645 fscl=1.0_dp
646 END IF
647 END IF
648 END DO
649!
650 IF (kind(clima(ng)%M2nudgcof).eq.8) THEN
651 nud(ng)%pioVar(ifield)%dkind=pio_double
652 iodesc => iodesc_dp_r2dvar(ng)
653 ELSE
654 nud(ng)%pioVar(ifield)%dkind=pio_real
655 iodesc => iodesc_sp_r2dvar(ng)
656 END IF
657
658 status=nf_fread2d(ng, model, ncname, nud(ng)%pioFile, &
659 & vname(1,ifield), nud(ng)%pioVar(ifield), &
660 & 0, iodesc, vsize, &
661 & lbi, ubi, lbj, ubj, &
662 & fscl, fmin, fmax, &
663# ifdef MASKING
664 & grid(ng) % rmask, &
665# endif
666# ifdef CHECKSUM
667 & clima(ng) % M2nudgcof, &
668 & checksum = fhash)
669# else
670 & clima(ng) % M2nudgcof)
671# endif
672 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
673 exit_flag=2
674 ioerror=status
675 RETURN
676 ELSE
677 IF (master) THEN
678 WRITE (stdout,40) trim(vname(2,ifield))//': '// &
679 & trim(vname(1,ifield)), &
680 & ng, trim(ncname), fmin, fmax
681# ifdef CHECKSUM
682 WRITE (stdout,50) fhash
683# endif
684 END IF
685 END IF
686!
687 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
688 CALL exchange_r2d_tile (ng, tile, &
689 & lbi, ubi, lbj, ubj, &
690 & clima(ng) % M2nudgcof)
691 END IF
692# ifdef DISTRIBUTE
693 CALL mp_exchange2d (ng, tile, model, 1, &
694 & lbi, ubi, lbj, ubj, &
695 & nghostpoints, &
696 & ewperiodic(ng), nsperiodic(ng), &
697 & clima(ng) % M2nudgcof)
698# endif
699 END IF
700
701# ifdef SOLVE3D
702!
703! If appropriate, read nudging coefficients for 3D momentum (inverse
704! time scales, 1/day).
705!
706 IF (lnudgem3clm(ng)) THEN
707 ifield=idm3nc
708 fscl=1.0_dp/day2sec ! default units: 1/day
709 nud(ng)%pioVar(ifield)%gtype=r3dvar
710!
711 CALL pio_netcdf_inq_var (ng, model, ncname, &
712 & piofile = nud(ng)%pioFile, &
713 & myvarname = trim(vname(1,ifield)), &
714 & piovar = nud(ng)%pioVar(ifield)%vd, &
715 & nvardim = nvdim, &
716 & nvaratt = nvatt)
717 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
718!
719 DO i=1,nvatt
720 IF (trim(var_aname(i)).eq.'units') THEN
721 tunits=trim(var_achar(i))
722 IF (tunits(1:3).eq.'day') THEN
723 fscl=1.0_dp/day2sec
724 ELSE IF (tunits(1:6).eq.'second') THEN
725 fscl=1.0_dp
726 END IF
727 END IF
728 END DO
729!
730 IF (kind(clima(ng)%M3nudgcof).eq.8) THEN
731 nud(ng)%pioVar(ifield)%dkind=pio_double
732 iodesc => iodesc_dp_r3dvar(ng)
733 ELSE
734 nud(ng)%pioVar(ifield)%dkind=pio_real
735 iodesc => iodesc_sp_r3dvar(ng)
736 END IF
737
738 status=nf_fread3d(ng, model, ncname, nud(ng)%pioFile, &
739 & vname(1,ifield), nud(ng)%pioVar(ifield), &
740 & 0, iodesc, vsize, &
741 & lbi, ubi, lbj, ubj, 1, n(ng), &
742 & fscl, fmin, fmax, &
743# ifdef MASKING
744 & grid(ng) % rmask, &
745# endif
746# ifdef CHECKSUM
747 & clima(ng) % M3nudgcof, &
748 & checksum = fhash)
749# else
750 & clima(ng) % M3nudgcof)
751# endif
752 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
753 exit_flag=2
754 ioerror=status
755 RETURN
756 ELSE
757 IF (master) THEN
758 WRITE (stdout,40) trim(vname(2,ifield))//': '// &
759 & trim(vname(1,ifield)), &
760 & ng, trim(ncname), fmin, fmax
761# ifdef CHECKSUM
762 WRITE (stdout,50) fhash
763# endif
764 END IF
765 END IF
766!
767 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
768 CALL exchange_r3d_tile (ng, tile, &
769 & lbi, ubi, lbj, ubj, 1, n(ng), &
770 & clima(ng) % M3nudgcof)
771 END IF
772# ifdef DISTRIBUTE
773 CALL mp_exchange3d (ng, tile, model, 1, &
774 & lbi, ubi, lbj, ubj, 1, n(ng), &
775 & nghostpoints, &
776 & ewperiodic(ng), nsperiodic(ng), &
777 & clima(ng) % M3nudgcof)
778# endif
779 END IF
780!
781! If appropriate, read nudging coefficients for tracer variables
782! (inverse time scales, 1/day). If the nudging coefficients for a
783! specific tracer are available, it will read that NetCDF variable.
784! If NOT AND the generic coefficients are available, it will process
785! those values instead.
786!
787 ic=0
788 DO itrc=1,nt(ng)
789 IF (lnudgetclm(itrc,ng)) THEN
790 ic=ic+1
791 IF (got_specific(itrc)) THEN
792 ifield=idtnud(itrc) ! specific variable
793 fscl=1.0_dp/day2sec ! default units: 1/day
794 ELSE IF (got_generic) THEN
795 ifield=idgtnc ! generic variable
796 fscl=1.0_dp/day2sec ! default units: 1/day
797 END IF
798 nud(ng)%pioVar(ifield)%gtype=r3dvar
799!
800 CALL pio_netcdf_inq_var (ng, model, ncname, &
801 & piofile = nud(ng)%pioFile, &
802 & myvarname = trim(vname(1,ifield)), &
803 & piovar = nud(ng)%pioVar(ifield)%vd, &
804 & nvardim = nvdim, &
805 & nvaratt = nvatt)
806 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
807!
808 DO i=1,nvatt
809 IF (trim(var_aname(i)).eq.'units') THEN
810 tunits=trim(var_achar(i))
811 IF (tunits(1:3).eq.'day') THEN
812 fscl=1.0_dp/day2sec
813 ELSE IF (tunits(1:6).eq.'second') THEN
814 fscl=1.0_dp
815 END IF
816 END IF
817 END DO
818!
819 IF (kind(clima(ng)%Tnudgcof).eq.8) THEN
820 nud(ng)%pioVar(ifield)%dkind=pio_double
821 iodesc => iodesc_dp_r3dvar(ng)
822 ELSE
823 nud(ng)%pioVar(ifield)%dkind=pio_real
824 iodesc => iodesc_sp_r3dvar(ng)
825 END IF
826
827 status=nf_fread3d(ng, model, ncname, nud(ng)%pioFile, &
828 & vname(1,ifield), nud(ng)%pioVar(ifield), &
829 & 0, iodesc, vsize, &
830 & lbi, ubi, lbj, ubj, 1, n(ng), &
831 & fscl, fmin, fmax, &
832# ifdef MASKING
833 & grid(ng) % rmask, &
834# endif
835# ifdef CHECKSUM
836 & clima(ng) % Tnudgcof(:,:,:,ic), &
837 & checksum = fhash)
838# else
839 & clima(ng) % Tnudgcof(:,:,:,ic))
840# endif
841 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
842 exit_flag=2
843 ioerror=status
844 RETURN
845 ELSE
846 IF (master) THEN
847 WRITE (stdout,40) trim(vname(2,ifield))//': '// &
848 & trim(vname(1,ifield)), &
849 & ng, trim(ncname), fmin, fmax
850# ifdef CHECKSUM
851 WRITE (stdout,50) fhash
852# endif
853 END IF
854 END IF
855!
856 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
857 CALL exchange_r3d_tile (ng, tile, &
858 & lbi, ubi, lbj, ubj, 1, n(ng), &
859 & clima(ng) % Tnudgcof(:,:,:,ic))
860 END IF
861 END IF
862 END DO
863# ifdef DISTRIBUTE
864 IF (any(lnudgetclm(:,ng))) THEN
865 CALL mp_exchange4d (ng, tile, model, 1, &
866 & lbi, ubi, lbj, ubj, 1, n(ng), 1, ntclm(ng), &
867 & nghostpoints, &
868 & ewperiodic(ng), nsperiodic(ng), &
869 & clima(ng) % Tnudgcof)
870 END IF
871# endif
872# endif
873!
874 10 FORMAT (/,' GET_NUDGCOEF_PIO - unable to open nudging NetCDF', &
875 & ' file: ',a)
876 20 FORMAT (/,' GET_NUDGCOEF_PIO - unable to find nudging', &
877 & ' variable: ',a,/,20x,'in NetCDF file: ',a)
878 30 FORMAT (/,' GET_NUDGCOEF_PIO - unable to find nudging', &
879 & ' variable: ',a, &
880 & /,20x,'or generic nudging variable: ',a, &
881 & /,20x,'in nudging NetCDF file: ',a)
882 40 FORMAT(2x,' GET_NUDGCOEF_PIO - ',a,/,21x, &
883 & '(Grid = ',i2.2,', File: ',a,')',/,21x, &
884 & '(Min = ', 1p,e15.8,0p,' Max = ',1p,e15.8,0p,')')
885# ifdef CHECKSUM
886 50 FORMAT (19x,'(CheckSum = ',i0,')')
887# endif
888!
889 RETURN
890 END SUBROUTINE get_nudgcoef_pio
891# endif
892#endif
893 END MODULE get_nudgcoef_mod
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine, private get_nudgcoef_pio(ng, tile, model, lbi, ubi, lbj, ubj)
subroutine, public get_nudgcoef(ng, tile, model)
subroutine, private get_nudgcoef_nf90(ng, tile, model, lbi, ubi, lbj, ubj)
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer ioerror
type(t_io), dimension(:), allocatable nud
integer stdout
character(len=256) sourcefile
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer idm3nc
integer, parameter io_pio
Definition mod_ncparam.F:96
integer idm2nc
integer idgtnc
integer, dimension(:), allocatable idtnud
character(len=maxlen), dimension(6, 0:nv) vname
character(len=1024), dimension(nvara) var_achar
Definition mod_netcdf.F:183
subroutine, public netcdf_check_dim(ng, model, ncname, ncid)
Definition mod_netcdf.F:538
subroutine, public netcdf_open(ng, model, ncname, omode, ncid)
character(len=100), dimension(mvars) var_name
Definition mod_netcdf.F:169
integer n_var
Definition mod_netcdf.F:152
character(len=100), dimension(nvara) var_aname
Definition mod_netcdf.F:181
subroutine, public netcdf_inq_var(ng, model, ncname, ncid, myvarname, searchvar, varid, nvardim, nvaratt)
logical master
integer, dimension(:), allocatable ntclm
Definition mod_param.F:494
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer, parameter r3dvar
Definition mod_param.F:721
integer nghostpoints
Definition mod_param.F:710
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer, parameter r2dvar
Definition mod_param.F:717
type(io_desc_t), dimension(:), pointer iodesc_dp_r3dvar
type(var_desc_t), dimension(:), pointer var_desc
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dvar
subroutine, public pio_netcdf_inq_var(ng, model, ncname, piofile, myvarname, searchvar, piovar, nvardim, nvaratt)
subroutine, public pio_netcdf_open(ng, model, ncname, omode, piofile)
subroutine, public pio_netcdf_check_dim(ng, model, ncname, piofile)
type(io_desc_t), dimension(:), pointer iodesc_sp_r3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dvar
real(dp), parameter day2sec
logical, dimension(:), allocatable lnudgem2clm
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
integer exit_flag
logical, dimension(:), allocatable lnudgem3clm
logical, dimension(:,:), allocatable lnudgetclm
integer noerror
subroutine mp_exchange4d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, lbt, ubt, nghost, ew_periodic, ns_periodic, a, b, c)
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
logical function, public find_string(a, asize, string, aindex)
Definition strings.F:417
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52