ROMS
Loading...
Searching...
No Matches
get_data.F File Reference
#include "cppdefs.h"
Include dependency graph for get_data.F:

Go to the source code of this file.

Functions/Subroutines

subroutine get_data (ng)
 

Function/Subroutine Documentation

◆ get_data()

subroutine get_data ( integer, intent(in) ng)

Definition at line 2 of file get_data.F.

3!
4!git $Id$
5!================================================== Hernan G. Arango ===
6! Copyright (c) 2002-2025 The ROMS Group !
7! Licensed under a MIT/X style license !
8! See License_ROMS.md !
9!=======================================================================
10! !
11! This routine reads in forcing, climatology and other data from !
12! NetCDF files. If there is more than one time-record, data is !
13! loaded into global two-time record arrays. The interpolation !
14! is carried elsewhere. !
15! !
16! Currently, this routine is only executed in serial mode by the !
17! main thread. !
18! !
19!=======================================================================
20!
21 USE mod_param
22#if defined HYPOXIA_SRM || defined RED_TIDE
23 USE mod_biology
24#endif
25 USE mod_boundary
26 USE mod_clima
27 USE mod_forces
28 USE mod_grid
29 USE mod_iounits
30 USE mod_ncparam
31#if defined HYPOXIA_SRM || \
32 defined nlm_outer || \
33 defined rbl4dvar || \
34 defined rbl4dvar_ana_sensitivity || \
35 defined rbl4dvar_fct_sensitivity || \
36 defined red_tide || \
37 defined sp4dvar || \
38 defined tlm_check || \
39 defined tl_rbl4dvar
40 USE mod_ocean
41#endif
42 USE mod_scalars
43 USE mod_sources
44 USE mod_stepping
45!
46 USE strings_mod, ONLY : founderror
47!
48 implicit none
49!
50! Imported variable declarations.
51!
52 integer, intent(in) :: ng
53!
54! Local variable declarations.
55!
56 logical :: Lprocess
57
58 logical, save :: recordless = .false.
59
60 logical, dimension(3) :: update = (/ .false., .false., .false. /)
61!
62 integer :: ILB, IUB, JLB, JUB
63 integer :: LBi, UBi, LBj, UBj
64 integer :: i, ic, my_tile
65!
66 character (len=*), parameter :: MyFile = &
67 & __FILE__
68!
69! Lower and upper bounds for nontiled (global values) boundary arrays.
70!
71 my_tile=-1 ! for global values
72 ilb=bounds(ng)%LBi(my_tile)
73 iub=bounds(ng)%UBi(my_tile)
74 jlb=bounds(ng)%LBj(my_tile)
75 jub=bounds(ng)%UBj(my_tile)
76!
77! Lower and upper bounds for tiled arrays.
78!
79 lbi=lbound(grid(ng)%h,dim=1)
80 ubi=ubound(grid(ng)%h,dim=1)
81 lbj=lbound(grid(ng)%h,dim=2)
82 ubj=ubound(grid(ng)%h,dim=2)
83
84#ifdef PROFILE
85!
86!-----------------------------------------------------------------------
87! Turn on input data time wall clock.
88!-----------------------------------------------------------------------
89!
90 CALL wclock_on (ng, inlm, 3, __line__, myfile)
91#endif
92
93#ifndef ANA_PSOURCE
94!
95!=======================================================================
96! Point Sources/Sinks time dependent data.
97!=======================================================================
98!
99! Point Source/Sink vertically integrated mass transport.
100!
101 IF (luvsrc(ng).or.lwsrc(ng)) THEN
102 CALL get_ngfld (ng, inlm, idrtra, ssf(ng)%ncid, &
103# if defined PIO_LIB && defined DISTRIBUTE
104 & ssf(ng)%pioFile, &
105# endif
106 & 1, ssf(ng), recordless, update(1), &
107 & 1, nsrc(ng), 1, 2, 1, nsrc(ng), 1, &
108 & sources(ng) % QbarG)
109 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
110 END IF
111
112# ifdef SOLVE3D
113!
114! Tracer Sources/Sinks.
115!
116 DO i=1,nt(ng)
117 IF (ltracersrc(i,ng)) THEN
118 CALL get_ngfld (ng, inlm, idrtrc(i), ssf(ng)%ncid, &
119# if defined PIO_LIB && defined DISTRIBUTE
120 & ssf(ng)%pioFile, &
121# endif
122 & 1, ssf(ng), recordless, update(1), &
123 & 1, nsrc(ng), n(ng), 2, 1, nsrc(ng), n(ng), &
124 & sources(ng) % TsrcG(:,:,:,i))
125 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
126 END IF
127 END DO
128# endif
129#endif
130!
131!=======================================================================
132! Read in forcing data from FORCING NetCDF file(s).
133!=======================================================================
134!
135! Set switch to process surface atmospheric fields.
136!
137#if defined FOUR_DVAR && \
138 defined bulk_fluxes && defined prior_bulk_fluxes
139! In 4D-Var data assimilation applications, the user have the option
140! to fix the prior (background phase) surface fluxes in the successive
141! outer loops (Nouter>1) or the final analysis phase. In such case, the
142! fluxes are read from the background trajector
143!
144 IF (nrun.eq.1) THEN
145 lprocess=.true.
146 ELSE
147 lprocess=.false.
148 END IF
149#else
150 lprocess=.true.
151#endif
152
153#if !(defined ANA_WINDS || defined FRC_COUPLING) && \
154 (defined bulk_fluxes || defined ecosim)
155!
156! Surface wind components.
157!
158 IF (lprocess) THEN
159 CALL get_2dfld (ng, inlm, iduair, frcncid(iduair,ng), &
160# if defined PIO_LIB && defined DISTRIBUTE
161 & frcpiofile(iduair,ng), &
162# endif
163 & nffiles(ng), frc(1,ng), update(1), &
164 & lbi, ubi, lbj, ubj, 2, 1, &
165# ifdef MASKING
166 & grid(ng) % rmask, &
167# endif
168 & forces(ng) % UwindG)
169 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
170!
171 CALL get_2dfld (ng , inlm, idvair, frcncid(idvair,ng), &
172# if defined PIO_LIB && defined DISTRIBUTE
173 & frcpiofile(idvair,ng), &
174# endif
175 & nffiles(ng), frc(1,ng), update(1), &
176 & lbi, ubi, lbj, ubj, 2, 1, &
177# ifdef MASKING
178 & grid(ng) % rmask, &
179# endif
180 & forces(ng) % VwindG)
181 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
182 END IF
183#endif
184
185#if !(defined ANA_SMFLUX || defined FRC_COUPLING || defined BULK_FLUXES)
186!
187! Surface wind stress components.
188!
189 IF (lprocess) THEN
190 CALL get_2dfld (ng, inlm, idusms, frcncid(idusms,ng), &
191# if defined PIO_LIB && defined DISTRIBUTE
192 & frcpiofile(idusms,ng), &
193# endif
194 & nffiles(ng), frc(1,ng), update(1), &
195 & lbi, ubi, lbj, ubj, 2, 1, &
196# ifdef MASKING
197 & grid(ng) % umask, &
198# endif
199 & forces(ng) % sustrG)
200 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
201!
202 CALL get_2dfld (ng, inlm, idvsms, frcncid(idvsms,ng), &
203# if defined PIO_LIB && defined DISTRIBUTE
204 & frcpiofile(idvsms,ng), &
205# endif
206 & nffiles(ng), frc(1,ng), update(1), &
207 & lbi, ubi, lbj, ubj, 2, 1, &
208# ifdef MASKING
209 & grid(ng) % vmask, &
210# endif
211 & forces(ng) % svstrG)
212 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
213 END IF
214#endif
215
216#if !(defined ANA_PAIR || defined FRC_COUPLING) && \
217 (defined bulk_fluxes || defined ecosim || defined atm_press)
218!
219! Surface air pressure.
220!
221 IF (lprocess) THEN
222 CALL get_2dfld (ng, inlm, idpair, frcncid(idpair,ng), &
223# if defined PIO_LIB && defined DISTRIBUTE
224 & frcpiofile(idpair,ng), &
225# endif
226 & nffiles(ng), frc(1,ng), update(1), &
227 & lbi, ubi, lbj, ubj, 2, 1, &
228# ifdef MASKING
229 & grid(ng) % rmask, &
230# endif
231 & forces(ng) % PairG)
232 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
233 END IF
234#endif
235
236#if !defined ANA_WWAVE && defined WAVE_DATA
237!
238! Surface wind induced wave properties.
239!
240# ifdef WAVES_DIR
241 CALL get_2dfld (ng, inlm, idwdir, frcncid(idwdir,ng), &
242# if defined PIO_LIB && defined DISTRIBUTE
243 & frcpiofile(idwdir,ng), &
244# endif
245 & nffiles(ng), frc(1,ng), update(1), &
246 & lbi, ubi, lbj, ubj, 2, 1, &
247# ifdef MASKING
248 & grid(ng) % rmask, &
249# endif
250 & forces(ng) % DwaveG)
251 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
252# endif
253
254# ifdef WAVES_DIRP
255!
256 CALL get_2dfld (ng, inlm, idwdip, frcncid(idwdip,ng), &
257# if defined PIO_LIB && defined DISTRIBUTE
258 & frcpiofile(idwdip,ng), &
259# endif
260 & nffiles(ng), frc(1,ng), update(1), &
261 & lbi, ubi, lbj, ubj, 2, 1, &
262# ifdef MASKING
263 & grid(ng) % rmask, &
264# endif
265 & forces(ng) % DwavepG)
266 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
267# endif
268
269# ifdef WAVES_HEIGHT
270!
271 CALL get_2dfld (ng, inlm, idwamp, frcncid(idwamp,ng), &
272# if defined PIO_LIB && defined DISTRIBUTE
273 & frcpiofile(idwamp,ng), &
274# endif
275 & nffiles(ng), frc(1,ng), update(1), &
276 & lbi, ubi, lbj, ubj, 2, 1, &
277# ifdef MASKING
278 & grid(ng) % rmask, &
279# endif
280 & forces(ng) % HwaveG)
281 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
282# endif
283
284# ifdef WAVES_LENGTH
285!
286 CALL get_2dfld (ng, inlm, idwlen, frcncid(idwlen,ng), &
287# if defined PIO_LIB && defined DISTRIBUTE
288 & frcpiofile(idwlen,ng), &
289# endif
290 & nffiles(ng), frc(1,ng), update(1), &
291 & lbi, ubi, lbj, ubj, 2, 1, &
292# ifdef MASKING
293 & grid(ng) % rmask, &
294# endif
295 & forces(ng) % LwaveG)
296 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
297# endif
298
299# ifdef WAVES_LENGTHP
300!
301 CALL get_2dfld (ng, inlm, idwlep, frcncid(idwlep,ng), &
302# if defined PIO_LIB && defined DISTRIBUTE
303 & frcpiofile(idwlep,ng), &
304# endif
305 & nffiles(ng), frc(1,ng), update(1), &
306 & lbi, ubi, lbj, ubj, 2, 1, &
307# ifdef MASKING
308 & grid(ng) % rmask, &
309# endif
310 & forces(ng) % LwavepG)
311 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
312# endif
313
314# ifdef WAVES_TOP_PERIOD
315!
316 CALL get_2dfld (ng, inlm, idwptp, frcncid(idwptp,ng), &
317# if defined PIO_LIB && defined DISTRIBUTE
318 & frcpiofile(idwptp,ng), &
319# endif
320 & nffiles(ng), frc(1,ng), update(1), &
321 & lbi, ubi, lbj, ubj, 2, 1, &
322# ifdef MASKING
323 & grid(ng) % rmask, &
324# endif
325 & forces(ng) % Pwave_topG)
326 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
327# endif
328
329# ifdef WAVES_BOT_PERIOD
330!
331 CALL get_2dfld (ng, inlm, idwpbt, frcncid(idwpbt,ng), &
332# if defined PIO_LIB && defined DISTRIBUTE
333 & frcpiofile(idwpbt,ng), &
334# endif
335 & nffiles(ng), frc(1,ng), update(1), &
336 & lbi, ubi, lbj, ubj, 2, 1, &
337# ifdef MASKING
338 & grid(ng) % rmask, &
339# endif
340 & forces(ng) % Pwave_botG(:,:,1))
341 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
342# endif
343
344# ifdef WAVES_UB
345!
346 CALL get_2dfld (ng, inlm, idworb, frcncid(idworb,ng), &
347# if defined PIO_LIB && defined DISTRIBUTE
348 & frcpiofile(idworb,ng), &
349# endif
350 & nffiles(ng), frc(1,ng), update(1), &
351 & lbi, ubi, lbj, ubj, 2, 1, &
352# ifdef MASKING
353 & grid(ng) % rmask, &
354# endif
355 & forces(ng) % Uwave_rmsG)
356 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
357# endif
358
359# ifdef WAVES_DISS
360!
361 CALL get_2dfld (ng, inlm, idwdib, frcncid(idwdib,ng), &
362# if defined PIO_LIB && defined DISTRIBUTE
363 & frcpiofile(idwdib,ng), &
364# endif
365 & nffiles(ng), frc(1,ng), update(1), &
366 & lbi, ubi, lbj, ubj, 2, 1, &
367# ifdef MASKING
368 & grid(ng) % rmask, &
369# endif
370 & forces(ng) % Dissip_breakG)
371 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
372!
373 CALL get_2dfld (ng, inlm, idwdiw, frcncid(idwdiw,ng), &
374# if defined PIO_LIB && defined DISTRIBUTE
375 & frcpiofile(idwdiw,ng), &
376# endif
377 & nffiles(ng), frc(1,ng), update(1), &
378 & lbi, ubi, lbj, ubj, 2, 1, &
379# ifdef MASKING
380 & grid(ng) % rmask, &
381# endif
382 & forces(ng) % Dissip_wcapG)
383 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
384# endif
385
386# ifdef ROLLER_SVENDSEN
387!
388 CALL get_2dfld (ng, inlm, idwbrk, frcncid(idwbrk,ng), &
389# if defined PIO_LIB && defined DISTRIBUTE
390 & frcpiofile(idwbrk,ng), &
391# endif
392 & nffiles(ng), frc(1,ng), update(1), &
393 & lbi, ubi, lbj, ubj, 2, 1, &
394# ifdef MASKING
395 & grid(ng) % rmask, &
396# endif
397 & forces(ng) % Wave_breakG)
398 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
399# endif
400#endif
401
402#ifdef SOLVE3D
403
404# if !(defined ANA_CLOUD || defined FRC_COUPLING) && defined CLOUDS
405!
406! Cloud fraction.
407!
408 CALL get_2dfld (ng, inlm, idcfra, frcncid(idcfra,ng), &
409# if defined PIO_LIB && defined DISTRIBUTE
410 & frcpiofile(idcfra,ng), &
411# endif
412 & nffiles(ng), frc(1,ng), update(1), &
413 & lbi, ubi, lbj, ubj, 2, 1, &
414# ifdef MASKING
415 & grid(ng) % rmask, &
416# endif
417 & forces(ng) % cloudG)
418 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
419# endif
420
421# if !(defined ANA_SRFLUX || defined FRC_COUPLING) && defined SHORTWAVE
422!
423! Surface solar shortwave radiation.
424!
425 IF (lprocess) THEN
426 CALL get_2dfld (ng, inlm, idsrad, frcncid(idsrad,ng), &
427# if defined PIO_LIB && defined DISTRIBUTE
428 & frcpiofile(idsrad,ng), &
429# endif
430 & nffiles(ng), frc(1,ng), update(1), &
431 & lbi, ubi, lbj, ubj, 2, 1, &
432# ifdef MASKING
433 & grid(ng) % rmask, &
434# endif
435 & forces(ng) % srflxG)
436 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
437 END IF
438# endif
439
440# if defined RED_TIDE && defined DAILY_SHORTWAVE
441!
442! Daily-averaged Surface solar shortwave radiation.
443!
444 CALL get_2dfld (ng, inlm, idasrf, frcncid(idasrf,ng), &
445# if defined PIO_LIB && defined DISTRIBUTE
446 & frcpiofile(idasrf,ng), &
447# endif
448 & nffiles(ng), frc(1,ng), update(1), &
449 & lbi, ubi, lbj, ubj, 2, 1, &
450# ifdef MASKING
451 & grid(ng) % rmask, &
452# endif
453 & forces(ng) % srflxG_avg)
454 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
455# endif
456
457# if defined BULK_FLUXES && \
458 !(defined FRC_COUPLING || defined LONGWAVE || defined LONGWAVE_OUT)
459!
460! Surface net longwave radiation.
461!
462 IF (lprocess) THEN
463 CALL get_2dfld (ng, inlm, idlrad, frcncid(idlrad,ng), &
464# if defined PIO_LIB && defined DISTRIBUTE
465 & frcpiofile(idlrad,ng), &
466# endif
467 & nffiles(ng), frc(1,ng), update(1), &
468 & lbi, ubi, lbj, ubj, 2, 1, &
469# ifdef MASKING
470 & grid(ng) % rmask, &
471# endif
472 & forces(ng) % lrflxG)
473 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
474 END IF
475# endif
476
477# if defined BULK_FLUXES && defined LONGWAVE_OUT && \
478 !defined FRC_COUPLING
479!
480! Surface downwelling longwave radiation.
481!
482 IF (lprocess) THEN
483 CALL get_2dfld (ng, inlm, idldwn, frcncid(idldwn,ng), &
484# if defined PIO_LIB && defined DISTRIBUTE
485 & frcpiofile(idldwn,ng), &
486# endif
487 & nffiles(ng), frc(1,ng), update(1), &
488 & lbi, ubi, lbj, ubj, 2, 1, &
489# ifdef MASKING
490 & grid(ng) % rmask, &
491# endif
492 & forces(ng) % lrflxG)
493 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
494 END IF
495# endif
496
497# if !(defined ANA_TAIR || defined FRC_COUPLING) && \
498 ( defined bulk_fluxes || defined ecosim || \
499 (defined shortwave && defined ana_srflux && defined albedo) )
500!
501! Surface air temperature.
502!
503 IF (lprocess) THEN
504 CALL get_2dfld (ng, inlm, idtair, frcncid(idtair,ng), &
505# if defined PIO_LIB && defined DISTRIBUTE
506 & frcpiofile(idtair,ng), &
507# endif
508 & nffiles(ng), frc(1,ng), update(1), &
509 & lbi, ubi, lbj, ubj, 2, 1, &
510# ifdef MASKING
511 & grid(ng) % rmask, &
512# endif
513 & forces(ng) % TairG)
514 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
515 END IF
516# endif
517
518# if !(defined ANA_HUMIDITY || defined FRC_COUPLING) && \
519 (defined bulk_fluxes || defined ecosim)
520!
521! Surface air humidity.
522!
523 IF (lprocess) THEN
524 CALL get_2dfld (ng, inlm, idqair, frcncid(idqair,ng), &
525# if defined PIO_LIB && defined DISTRIBUTE
526 & frcpiofile(idqair,ng), &
527# endif
528 & nffiles(ng), frc(1,ng), update(1), &
529 & lbi, ubi, lbj, ubj, 2, 1, &
530# ifdef MASKING
531 & grid(ng) % rmask, &
532# endif
533 & forces(ng) % HairG)
534 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
535 END IF
536# endif
537
538# if !(defined ANA_RAIN || defined FRC_COUPLING) && defined BULK_FLUXES
539!
540! Rain fall rate.
541!
542 IF (lprocess) THEN
543 CALL get_2dfld (ng, inlm, idrain, frcncid(idrain,ng), &
544# if defined PIO_LIB && defined DISTRIBUTE
545 & frcpiofile(idrain,ng), &
546# endif
547 & nffiles(ng), frc(1,ng), update(1), &
548 & lbi, ubi, lbj, ubj, 2, 1, &
549# ifdef MASKING
550 & grid(ng) % rmask, &
551# endif
552 & forces(ng) % rainG)
553 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
554 END IF
555# endif
556
557# if !(defined ANA_STFLUX || defined FRC_COUPLING || \
558 defined bulk_fluxes)
559!
560! Surface net heat flux.
561!
562 IF (lprocess) THEN
563 CALL get_2dfld (ng, inlm, idtsur(itemp), &
564 & frcncid(idtsur(itemp),ng), &
565# if defined PIO_LIB && defined DISTRIBUTE
566 & frcpiofile(idtsur(itemp),ng), &
567# endif
568 & nffiles(ng), frc(1,ng), update(1), &
569 & lbi, ubi, lbj, ubj, 2, 1, &
570# ifdef MASKING
571 & grid(ng) % rmask, &
572# endif
573 & forces(ng) % stfluxG(:,:,:,itemp))
574 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
575 END IF
576# endif
577
578# if !defined ANA_SST && defined QCORRECTION
579!
580! Surface net heat flux correction fields: sea surface temperature
581! (SST).
582!
583 CALL get_2dfld (ng, inlm, idsstc, frcncid(idsstc,ng), &
584# if defined PIO_LIB && defined DISTRIBUTE
585 & frcpiofile(idsstc,ng), &
586# endif
587 & nffiles(ng), frc(1,ng), update(1), &
588 & lbi, ubi, lbj, ubj, 2, 1, &
589# ifdef MASKING
590 & grid(ng) % rmask, &
591# endif
592 & forces(ng) % sstG)
593 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
594# endif
595
596# if !defined ANA_DQDSST && defined QCORRECTION
597!
598! Surface net heat flux correction fields: heat flux sensitivity to
599! SST (dQdSST).
600!
601 CALL get_2dfld (ng, inlm, iddqdt, frcncid(iddqdt,ng), &
602# if defined PIO_LIB && defined DISTRIBUTE
603 & frcpiofile(iddqdt,ng), &
604# endif
605 & nffiles(ng), frc(1,ng), update(1), &
606 & lbi, ubi, lbj, ubj, 2, 1, &
607# ifdef MASKING
608 & grid(ng) % rmask, &
609# endif
610 & forces(ng) % dqdtG)
611 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
612# endif
613
614# ifndef ANA_BTFLUX
615!
616! Bottom net heat flux.
617!
618 CALL get_2dfld (ng, inlm, idtbot(itemp), &
619 & frcncid(idtbot(itemp),ng), &
620# if defined PIO_LIB && defined DISTRIBUTE
621 & frcpiofile(idtbot(itemp),ng), &
622# endif
623 & nffiles(ng), frc(1,ng), update(1), &
624 & lbi, ubi, lbj, ubj, 2, 1, &
625# ifdef MASKING
626 & grid(ng) % rmask, &
627# endif
628 & forces(ng) % btfluxG(:,:,:,itemp))
629 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
630# endif
631
632# if !defined ANA_SSFLUX && defined SALINITY
633# if !(defined BULK_FLUXES || defined EMINUSP || \
634 defined frc_coupling)
635!
636! Surface net freshwater flux: E-P from NetCDF variable "swflux".
637!
638 IF (lprocess) THEN
639 CALL get_2dfld (ng, inlm, idsfwf, frcncid(idsfwf,ng), &
640# if defined PIO_LIB && defined DISTRIBUTE
641 & frcpiofile(idsfwf,ng), &
642# endif
643 & nffiles(ng), frc(1,ng), update(1), &
644 & lbi, ubi, lbj, ubj, 2, 1, &
645# ifdef MASKING
646 & grid(ng) % rmask, &
647# endif
648 & forces(ng) % stfluxG(:,:,:,isalt))
649 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
650 END IF
651
652# elif defined BULK_FLUXES && !defined EMINUSP
653!
654! Surface net freshwater flux (E-P) from NetCDF variable "EminusP".
655!
656 IF (lprocess) THEN
657 CALL get_2dfld (ng, inlm, idempf, frcncid(idempf,ng), &
658# if defined PIO_LIB && defined DISTRIBUTE
659 & frcpiofile(idempf,ng), &
660# endif
661 & nffiles(ng), frc(1,ng), update(1), &
662 & lbi, ubi, lbj, ubj, 2, 1, &
663# ifdef MASKING
664 & grid(ng) % rmask, &
665# endif
666 & forces(ng) % stfluxG(:,:,:,isalt))
667 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
668 END IF
669# endif
670
671# if !defined ANA_SSS && (defined SCORRECTION || defined SRELAXATION)
672!
673! Surface net freshwater flux correction field: sea surface salinity.
674!
675 CALL get_2dfld (ng, inlm, idsssc, frcncid(idsssc,ng), &
676# if defined PIO_LIB && defined DISTRIBUTE
677 & frcpiofile(idsssc,ng), &
678# endif
679 & nffiles(ng), frc(1,ng), update(1), &
680 & lbi, ubi, lbj, ubj, 2, 1, &
681# ifdef MASKING
682 & grid(ng) % rmask, &
683# endif
684 & forces(ng) % sssG)
685 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
686# endif
687
688# ifndef ANA_BSFLUX
689!
690! Bottom net freshwater flux.
691!
692 CALL get_2dfld (ng, inlm, idtbot(isalt), &
693 & frcncid(idtbot(isalt),ng), &
694# if defined PIO_LIB && defined DISTRIBUTE
695 & frcpiofile(idtbot(isalt),ng), &
696# endif
697 & nffiles(ng), frc(1,ng), update(1), &
698 & lbi, ubi, lbj, ubj, 2, 1, &
699# ifdef MASKING
700 & grid(ng) % rmask, &
701# endif
702 & forces(ng) % btfluxG(:,:,:,isalt))
703 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
704# endif
705# endif
706
707# if defined BIOLOGY || defined SEDIMENT || defined T_PASSIVE
708# ifndef ANA_SPFLUX
709!
710! Passive tracers surface fluxes.
711!
712 DO i=nat+1,nt(ng)
713 CALL get_2dfld (ng, inlm, idtsur(i), frcncid(idtsur(i),ng), &
714# if defined PIO_LIB && defined DISTRIBUTE
715 & frcpiofile(idtsur(i),ng), &
716# endif
717 & nffiles(ng), frc(1,ng), update(1), &
718 & lbi, ubi, lbj, ubj, 2, 1, &
719# ifdef MASKING
720 & grid(ng) % rmask, &
721# endif
722 & forces(ng) % stfluxG(:,:,:,i))
723 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
724 END DO
725# endif
726
727# ifndef ANA_BPFLUX
728!
729! Passive tracers bottom fluxes.
730!
731 DO i=nat+1,nt(ng)
732 CALL get_2dfld (ng, inlm, idtbot(i), frcncid(idtbot(i),ng), &
733# if defined PIO_LIB && defined DISTRIBUTE
734 & frcpiofile(idtbot(i),ng), &
735# endif
736 & nffiles(ng), frc(1,ng), update(1), &
737 & lbi, ubi, lbj, ubj, 2, 1, &
738# ifdef MASKING
739 & grid(ng) % rmask, &
740# endif
741 & forces(ng) % btfluxG(:,:,:,i))
742 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
743 END DO
744# endif
745# endif
746#endif
747
748#if !defined ANA_RESPIRATION && defined HYPOXIA_SRM
749!
750! Total respiration rate for hypoxia.
751!
752 CALL get_3dfld (ng, inlm, idresr, frcncid(idresr,ng), &
753# if defined PIO_LIB && defined DISTRIBUTE
754 & frcpiofile(idresr,ng), &
755# endif
756 & nffiles(ng), frc(1,ng), update(1), &
757 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
758# ifdef MASKING
759 & grid(ng) % rmask, &
760# endif
761 & ocean(ng) % respirationG)
762 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
763#endif
764
765#ifdef RED_TIDE
766!
767! Red tide Observed Dissolved Inorganic Nutrient.
768!
769 CALL get_3dfld (ng, inlm, idodin, frcncid(idodin,ng), &
770# if defined PIO_LIB && defined DISTRIBUTE
771 & frcpiofile(idodin,ng), &
772# endif
773 & nffiles(ng), frc(1,ng), update(1), &
774 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
775# ifdef MASKING
776 & grid(ng) % rmask, &
777# endif
778 & ocean(ng) % DIN_obsG)
779 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
780#endif
781
782#if defined FOUR_DVAR && \
783 defined bulk_fluxes && defined prior_bulk_fluxes
784!
785!=======================================================================
786! In 4D-Var data assimilation algorithms, the user has the option to
787! impose the prior (background phase) surface fluxes in the successive
788! outer loops (Nouter>1) or the final analysis phase. Such fluxes were
789! computed by routine "bulk_fluxes" and stored in the NLM background
790! trajectory, BLK structure. Notice that "bulk_fluxes" is called in
791! "main3d" only during the prior trajectory computation.
792!
793! The strategy is to save the 2D surface fluxes in the quicksave
794! NetCDF file to allow frequent data records to resolve the daily
795! cycle while maintaining its size small.
796!=======================================================================
797!
798 IF (.not.lprocess) THEN
799!
800! Get prior surface wind stress components.
801!
802 CALL get_2dfld (ng, inlm, idusms, blk(ng)%ncid, &
803# if defined PIO_LIB && defined DISTRIBUTE
804 & blk(ng)%pioFile, &
805# endif
806 & 1, blk(ng), update(1), &
807 & lbi, ubi, lbj, ubj, 2, 1, &
808# ifdef MASKING
809 & grid(ng) % umask, &
810# endif
811 & forces(ng) % sustrG)
812 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
813!
814 CALL get_2dfld (ng, inlm, idvsms, blk(ng)%ncid, &
815# if defined PIO_LIB && defined DISTRIBUTE
816 & blk(ng)%pioFile, &
817# endif
818 & 1, blk(ng), update(1), &
819 & lbi, ubi, lbj, ubj, 2, 1, &
820# ifdef MASKING
821 & grid(ng) % vmask, &
822# endif
823 & forces(ng) % svstrG)
824 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
825
826# ifdef ATM_PRESS
827!
828! Get prior surface air pressure.
829!
830 CALL get_2dfld (ng, inlm, idpair, blk(ng)%ncid, &
831# if defined PIO_LIB && defined DISTRIBUTE
832 & blk(ng)%pioFile, &
833# endif
834 & 1, blk(ng), update(1), &
835 & lbi, ubi, lbj, ubj, 2, 1, &
836# ifdef MASKING
837 & grid(ng) % rmask, &
838# endif
839 & forces(ng) % PairG)
840 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
841# endif
842
843# ifdef SOLVE3D
844# ifdef SHORTWAVE
845!
846! Get prior shortwave radiation flux. For consistency, we need to
847! process the same values using in the computation of the net heat
848! flux since there is a time interpolation from snapshots.
849!
850 CALL get_2dfld (ng, inlm, idsrad, blk(ng)%ncid, &
851# if defined PIO_LIB && defined DISTRIBUTE
852 & blk(ng)%pioFile, &
853# endif
854 & 1, blk(ng), update(1), &
855 & lbi, ubi, lbj, ubj, 2, 1, &
856# ifdef MASKING
857 & grid(ng) % rmask, &
858# endif
859 & forces(ng) % srflxG)
860 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
861# endif
862!
863! Get prior surface net heat flux.
864!
865 CALL get_2dfld (ng, inlm, idtsur(itemp), blk(ng)%ncid, &
866# if defined PIO_LIB && defined DISTRIBUTE
867 & blk(ng)%pioFile, &
868# endif
869 & 1, blk(ng), update(1), &
870 & lbi, ubi, lbj, ubj, 2, 1, &
871# ifdef MASKING
872 & grid(ng) % rmask, &
873# endif
874 & forces(ng) % stfluxG(:,:,:,itemp))
875 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
876
877# ifdef SALINITY
878!
879! Get prior net freshwater flux (E-P). The kinematic net freshwater
880! will be computed in "set_vbc" by multiplying with surface salinity.
881!
882 CALL get_2dfld (ng, inlm, idempf, blk(ng)%ncid, &
883# if defined PIO_LIB && defined DISTRIBUTE
884 & blk(ng)%pioFile, &
885# endif
886 & 1, blk(ng), update(1), &
887 & lbi, ubi, lbj, ubj, 2, 1, &
888# ifdef MASKING
889 & grid(ng) % rmask, &
890# endif
891 & forces(ng) % stfluxG(:,:,:,isalt))
892 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
893# endif
894# endif
895 END IF
896#endif
897!
898!=======================================================================
899! Read in open boundary conditions from BOUNDARY NetCDF file.
900!=======================================================================
901
902#ifndef ANA_FSOBC
903!
904! Free-surface open boundary conditions.
905!
906 IF (lprocessobc(ng)) THEN
907 IF (lbc(iwest,isfsur,ng)%acquire) THEN
908 CALL get_ngfld (ng, inlm, idzbry(iwest), &
909 & bryncid(idzbry(iwest),ng), &
910# if defined PIO_LIB && defined DISTRIBUTE
911 & brypiofile(idzbry(iwest),ng), &
912# endif
913 & nbcfiles(ng), bry(1,ng), &
914 & recordless, update(1), &
915 & jlb, jub, 1, 2, 0, mm(ng)+1, 1, &
916 & boundary(ng) % zetaG_west)
917 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
918 END IF
919!
920 IF (lbc(ieast,isfsur,ng)%acquire) THEN
921 CALL get_ngfld (ng, inlm, idzbry(ieast), &
922 & bryncid(idzbry(ieast),ng), &
923# if defined PIO_LIB && defined DISTRIBUTE
924 & brypiofile(idzbry(ieast),ng), &
925# endif
926 & nbcfiles(ng), bry(1,ng), &
927 & recordless, update(1), &
928 & jlb, jub, 1, 2, 0, mm(ng)+1, 1, &
929 & boundary(ng) % zetaG_east)
930 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
931 END IF
932!
933 IF (lbc(isouth,isfsur,ng)%acquire) THEN
934 CALL get_ngfld (ng, inlm, idzbry(isouth), &
935 & bryncid(idzbry(isouth),ng), &
936# if defined PIO_LIB && defined DISTRIBUTE
937 & brypiofile(idzbry(isouth),ng), &
938# endif
939 & nbcfiles(ng), bry(1,ng), &
940 & recordless, update(1), &
941 & ilb, iub, 1, 2, 0, lm(ng)+1, 1, &
942 & boundary(ng) % zetaG_south)
943 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
944 END IF
945!
946 IF (lbc(inorth,isfsur,ng)%acquire) THEN
947 CALL get_ngfld (ng, inlm, idzbry(inorth), &
948 & bryncid(idzbry(inorth),ng), &
949# if defined PIO_LIB && defined DISTRIBUTE
950 & brypiofile(idzbry(inorth),ng), &
951# endif
952 & nbcfiles(ng), bry(1,ng), &
953 & recordless, update(1), &
954 & ilb, iub, 1, 2, 0, lm(ng)+1, 1, &
955 & boundary(ng) % zetaG_north)
956 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
957 END IF
958 END IF
959#endif
960
961#ifndef ANA_M2OBC
962!
963! 2D momentum components open boundary conditions.
964!
965 IF (lprocessobc(ng)) THEN
966 IF (lbc(iwest,isubar,ng)%acquire) THEN
967 CALL get_ngfld (ng, inlm, idu2bc(iwest), &
968 & bryncid(idu2bc(iwest),ng), &
969# if defined PIO_LIB && defined DISTRIBUTE
970 & brypiofile(idu2bc(iwest),ng), &
971# endif
972 & nbcfiles(ng), bry(1,ng), &
973 & recordless, update(1), &
974 & jlb, jub, 1, 2, 0, mm(ng)+1, 1, &
975 & boundary(ng) % ubarG_west)
976 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
977 END IF
978!
979 IF (lbc(iwest,isvbar,ng)%acquire) THEN
980 CALL get_ngfld (ng, inlm, idv2bc(iwest), &
981 & bryncid(idv2bc(iwest),ng), &
982# if defined PIO_LIB && defined DISTRIBUTE
983 & brypiofile(idv2bc(iwest),ng), &
984# endif
985 & nbcfiles(ng), bry(1,ng), &
986 & recordless, update(1), &
987 & jlb, jub, 1, 2, 1, mm(ng)+1, 1, &
988 & boundary(ng) % vbarG_west)
989 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
990 END IF
991!
992 IF (lbc(ieast,isubar,ng)%acquire) THEN
993 CALL get_ngfld (ng, inlm, idu2bc(ieast), &
994 & bryncid(idu2bc(ieast),ng), &
995# if defined PIO_LIB && defined DISTRIBUTE
996 & brypiofile(idu2bc(ieast),ng), &
997# endif
998 & nbcfiles(ng), bry(1,ng), &
999 & recordless, update(1), &
1000 & jlb, jub, 1, 2, 0, mm(ng)+1, 1, &
1001 & boundary(ng) % ubarG_east)
1002 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1003 END IF
1004!
1005 IF (lbc(ieast,isvbar,ng)%acquire) THEN
1006 CALL get_ngfld (ng, inlm, idv2bc(ieast), &
1007 & bryncid(idv2bc(ieast),ng), &
1008# if defined PIO_LIB && defined DISTRIBUTE
1009 & brypiofile(idv2bc(ieast),ng), &
1010# endif
1011 & nbcfiles(ng), bry(1,ng), &
1012 & recordless, update(1), &
1013 & jlb, jub, 1, 2, 1, mm(ng)+1, 1, &
1014 & boundary(ng) % vbarG_east)
1015 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1016 END IF
1017!
1018 IF (lbc(isouth,isubar,ng)%acquire) THEN
1019 CALL get_ngfld (ng, inlm, idu2bc(isouth), &
1020 & bryncid(idu2bc(isouth),ng), &
1021# if defined PIO_LIB && defined DISTRIBUTE
1022 & brypiofile(idu2bc(isouth),ng), &
1023# endif
1024 & nbcfiles(ng), bry(1,ng), &
1025 & recordless, update(1), &
1026 & ilb, iub, 1, 2, 1, lm(ng)+1, 1, &
1027 & boundary(ng) % ubarG_south)
1028 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1029 END IF
1030!
1031 IF (lbc(isouth,isvbar,ng)%acquire) THEN
1032 CALL get_ngfld (ng, inlm, idv2bc(isouth), &
1033 & bryncid(idv2bc(isouth),ng), &
1034# if defined PIO_LIB && defined DISTRIBUTE
1035 & brypiofile(idv2bc(isouth),ng), &
1036# endif
1037 & nbcfiles(ng), bry(1,ng), &
1038 & recordless, update(1), &
1039 & ilb, iub, 1, 2, 0, lm(ng)+1, 1, &
1040 & boundary(ng) % vbarG_south)
1041 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1042 END IF
1043!
1044 IF (lbc(inorth,isubar,ng)%acquire) THEN
1045 CALL get_ngfld (ng, inlm, idu2bc(inorth), &
1046 & bryncid(idu2bc(inorth),ng), &
1047# if defined PIO_LIB && defined DISTRIBUTE
1048 & brypiofile(idu2bc(inorth),ng), &
1049# endif
1050 & nbcfiles(ng), bry(1,ng), &
1051 & recordless, update(1), &
1052 & ilb, iub, 1, 2, 1, lm(ng)+1, 1, &
1053 & boundary(ng) % ubarG_north)
1054 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1055 END IF
1056!
1057 IF (lbc(inorth,isvbar,ng)%acquire) THEN
1058 CALL get_ngfld (ng, inlm, idv2bc(inorth), &
1059 & bryncid(idv2bc(inorth),ng), &
1060# if defined PIO_LIB && defined DISTRIBUTE
1061 & brypiofile(idv2bc(inorth),ng), &
1062# endif
1063 & nbcfiles(ng), bry(1,ng), &
1064 & recordless, update(1), &
1065 & ilb, iub, 1, 2, 0, lm(ng)+1, 1, &
1066 & boundary(ng) % vbarG_north)
1067 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1068 END IF
1069 END IF
1070#endif
1071
1072#ifdef SOLVE3D
1073# ifndef ANA_M3OBC
1074!
1075! 3D momentum components open boundary conditions.
1076!
1077 IF (lprocessobc(ng)) THEN
1078 IF (lbc(iwest,isuvel,ng)%acquire) THEN
1079 CALL get_ngfld (ng, inlm, idu3bc(iwest), &
1080 & bryncid(idu3bc(iwest),ng), &
1081# if defined PIO_LIB && defined DISTRIBUTE
1082 & brypiofile(idu3bc(iwest),ng), &
1083# endif
1084 & nbcfiles(ng), bry(1,ng), &
1085 & recordless, update(1), &
1086 & jlb, jub, n(ng), 2, 0, mm(ng)+1, n(ng), &
1087 & boundary(ng) % uG_west)
1088 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1089 END IF
1090!
1091 IF (lbc(iwest,isvvel,ng)%acquire) THEN
1092 CALL get_ngfld (ng, inlm, idv3bc(iwest), &
1093 & bryncid(idv3bc(iwest),ng), &
1094# if defined PIO_LIB && defined DISTRIBUTE
1095 & brypiofile(idv3bc(iwest),ng), &
1096# endif
1097 & nbcfiles(ng), bry(1,ng), &
1098 & recordless, update(1), &
1099 & jlb, jub, n(ng), 2, 1, mm(ng)+1, n(ng), &
1100 & boundary(ng) % vG_west)
1101 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1102 END IF
1103!
1104 IF (lbc(ieast,isuvel,ng)%acquire) THEN
1105 CALL get_ngfld (ng, inlm, idu3bc(ieast), &
1106 & bryncid(idu3bc(ieast),ng), &
1107# if defined PIO_LIB && defined DISTRIBUTE
1108 & brypiofile(idu3bc(ieast),ng), &
1109# endif
1110 & nbcfiles(ng), bry(1,ng), &
1111 & recordless, update(1), &
1112 & jlb, jub, n(ng), 2, 0, mm(ng)+1, n(ng), &
1113 & boundary(ng) % uG_east)
1114 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1115 END IF
1116!
1117 IF (lbc(ieast,isvvel,ng)%acquire) THEN
1118 CALL get_ngfld (ng, inlm, idv3bc(ieast), &
1119 & bryncid(idv3bc(ieast),ng), &
1120# if defined PIO_LIB && defined DISTRIBUTE
1121 & brypiofile(idv3bc(ieast),ng), &
1122# endif
1123 & nbcfiles(ng), bry(1,ng), &
1124 & recordless, update(1), &
1125 & jlb, jub, n(ng), 2, 1, mm(ng)+1, n(ng), &
1126 & boundary(ng) % vG_east)
1127 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1128 END IF
1129!
1130 IF (lbc(isouth,isuvel,ng)%acquire) THEN
1131 CALL get_ngfld (ng, inlm, idu3bc(isouth), &
1132 & bryncid(idu3bc(isouth),ng), &
1133# if defined PIO_LIB && defined DISTRIBUTE
1134 & brypiofile(idu3bc(isouth),ng), &
1135# endif
1136 & nbcfiles(ng), bry(1,ng), &
1137 & recordless, update(1), &
1138 & ilb, iub, n(ng), 2, 1, lm(ng)+1, n(ng), &
1139 & boundary(ng) % uG_south)
1140 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1141 END IF
1142!
1143 IF (lbc(isouth,isvvel,ng)%acquire) THEN
1144 CALL get_ngfld (ng, inlm, idv3bc(isouth), &
1145 & bryncid(idv3bc(isouth),ng), &
1146# if defined PIO_LIB && defined DISTRIBUTE
1147 & brypiofile(idv3bc(isouth),ng), &
1148# endif
1149 & nbcfiles(ng), bry(1,ng), &
1150 & recordless, update(1), &
1151 & ilb, iub, n(ng), 2, 0, lm(ng)+1, n(ng), &
1152 & boundary(ng) % vG_south)
1153 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1154 END IF
1155!
1156 IF (lbc(inorth,isuvel,ng)%acquire) THEN
1157 CALL get_ngfld (ng, inlm, idu3bc(inorth), &
1158 & bryncid(idu3bc(inorth),ng), &
1159# if defined PIO_LIB && defined DISTRIBUTE
1160 & brypiofile(idu3bc(inorth),ng), &
1161# endif
1162 & nbcfiles(ng), bry(1,ng), &
1163 & recordless, update(1), &
1164 & ilb, iub, n(ng), 2, 1, lm(ng)+1, n(ng), &
1165 & boundary(ng) % uG_north)
1166 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1167 END IF
1168!
1169 IF (lbc(inorth,isvvel,ng)%acquire) THEN
1170 CALL get_ngfld (ng, inlm, idv3bc(inorth), &
1171 & bryncid(idv3bc(inorth),ng), &
1172# if defined PIO_LIB && defined DISTRIBUTE
1173 & brypiofile(idv3bc(inorth),ng), &
1174# endif
1175 & nbcfiles(ng), bry(1,ng), &
1176 & recordless, update(1), &
1177 & ilb, iub, n(ng), 2, 0, lm(ng)+1, n(ng), &
1178 & boundary(ng) % vG_north)
1179 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1180 END IF
1181 END IF
1182# endif
1183
1184# ifndef ANA_TOBC
1185!
1186! Tracer variables open boundary conditions.
1187!
1188 IF (lprocessobc(ng)) THEN
1189 DO i=1,nt(ng)
1190 IF (lbc(iwest,istvar(i),ng)%acquire) THEN
1191 CALL get_ngfld (ng, inlm, idtbry(iwest,i), &
1192 & bryncid(idtbry(iwest,i),ng), &
1193# if defined PIO_LIB && defined DISTRIBUTE
1194 & brypiofile(idtbry(iwest,i),ng), &
1195# endif
1196 & nbcfiles(ng), bry(1,ng), &
1197 & recordless, update(1), &
1198 & jlb, jub, n(ng), 2, 0, mm(ng)+1, n(ng), &
1199 & boundary(ng) % tG_west(:,:,:,i))
1200 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1201 END IF
1202 END DO
1203!
1204 DO i=1,nt(ng)
1205 IF (lbc(ieast,istvar(i),ng)%acquire) THEN
1206 CALL get_ngfld (ng, inlm, idtbry(ieast,i), &
1207 & bryncid(idtbry(ieast,i),ng), &
1208# if defined PIO_LIB && defined DISTRIBUTE
1209 & brypiofile(idtbry(ieast,i),ng), &
1210# endif
1211 & nbcfiles(ng), bry(1,ng), &
1212 & recordless, update(1), &
1213 & jlb, jub, n(ng), 2, 0, mm(ng)+1, n(ng), &
1214 & boundary(ng) % tG_east(:,:,:,i))
1215 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1216 END IF
1217 END DO
1218!
1219 DO i=1,nt(ng)
1220 IF (lbc(isouth,istvar(i),ng)%acquire) THEN
1221 CALL get_ngfld (ng, inlm, idtbry(isouth,i), &
1222 & bryncid(idtbry(isouth,i),ng), &
1223# if defined PIO_LIB && defined DISTRIBUTE
1224 & brypiofile(idtbry(isouth,i),ng), &
1225# endif
1226 & nbcfiles(ng), bry(1,ng), &
1227 & recordless, update(1), &
1228 & ilb, iub, n(ng), 2, 0, lm(ng)+1, n(ng), &
1229 & boundary(ng) % tG_south(:,:,:,i))
1230 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1231 END IF
1232 END DO
1233!
1234 DO i=1,nt(ng)
1235 IF (lbc(inorth,istvar(i),ng)%acquire) THEN
1236 CALL get_ngfld (ng, inlm, idtbry(inorth,i), &
1237 & bryncid(idtbry(inorth,i),ng), &
1238# if defined PIO_LIB && defined DISTRIBUTE
1239 & brypiofile(idtbry(inorth,i),ng), &
1240# endif
1241 & nbcfiles(ng), bry(1,ng), &
1242 & recordless, update(1), &
1243 & ilb, iub, n(ng), 2, 0, lm(ng)+1, n(ng), &
1244 & boundary(ng) % tG_north(:,:,:,i))
1245 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1246 END IF
1247 END DO
1248 END IF
1249# endif
1250#endif
1251!
1252!=======================================================================
1253! Read in data from Climatology NetCDF file.
1254!=======================================================================
1255
1256#ifndef ANA_SSH
1257!
1258! Free-surface climatology.
1259!
1260 IF (lsshclm(ng)) THEN
1261 CALL get_2dfld (ng, inlm, idsshc, clmncid(idsshc,ng), &
1262# if defined PIO_LIB && defined DISTRIBUTE
1263 & clmpiofile(idsshc,ng), &
1264# endif
1265 & nclmfiles(ng), clm(1,ng), update(1), &
1266 & lbi, ubi, lbj, ubj, 2, 1, &
1267# ifdef MASKING
1268 & grid(ng) % rmask, &
1269# endif
1270 & clima(ng) % sshG)
1271 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1272 END IF
1273#endif
1274#ifndef ANA_M2CLIMA
1275!
1276! 2D momentum components climatology.
1277!
1278 IF (lm2clm(ng)) THEN
1279 CALL get_2dfld (ng, inlm, idubcl, clmncid(idubcl,ng), &
1280# if defined PIO_LIB && defined DISTRIBUTE
1281 & clmpiofile(idubcl,ng), &
1282# endif
1283 & nclmfiles(ng), clm(1,ng), update(1), &
1284 & lbi, ubi, lbj, ubj, 2, 1, &
1285# ifdef MASKING
1286 & grid(ng) % umask, &
1287# endif
1288 & clima(ng) % ubarclmG)
1289 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1290!
1291 CALL get_2dfld (ng, inlm, idvbcl, clmncid(idvbcl,ng), &
1292# if defined PIO_LIB && defined DISTRIBUTE
1293 & clmpiofile(idvbcl,ng), &
1294# endif
1295 & nclmfiles(ng), clm(1,ng), update(1), &
1296 & lbi, ubi, lbj, ubj, 2, 1, &
1297# ifdef MASKING
1298 & grid(ng) % vmask, &
1299# endif
1300 & clima(ng) % vbarclmG)
1301 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1302 END IF
1303#endif
1304#ifdef SOLVE3D
1305# ifndef ANA_M3CLIMA
1306!
1307! 3D momentum components climatology.
1308!
1309 IF (lm3clm(ng)) THEN
1310 CALL get_3dfld (ng, inlm, iduclm, clmncid(iduclm,ng), &
1311# if defined PIO_LIB && defined DISTRIBUTE
1312 & clmpiofile(iduclm,ng), &
1313# endif
1314 & nclmfiles(ng), clm(1,ng), update(1), &
1315 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
1316# ifdef MASKING
1317 & grid(ng) % umask, &
1318# endif
1319 & clima(ng) % uclmG)
1320 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1321!
1322 CALL get_3dfld (ng, inlm, idvclm, clmncid(idvclm,ng), &
1323# if defined PIO_LIB && defined DISTRIBUTE
1324 & clmpiofile(idvclm,ng), &
1325# endif
1326 & nclmfiles(ng), clm(1,ng), update(1), &
1327 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
1328# ifdef MASKING
1329 & grid(ng) % vmask, &
1330# endif
1331 & clima(ng) % vclmG)
1332 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1333 END IF
1334# endif
1335# ifndef ANA_TCLIMA
1336!
1337! Tracers variables climatology.
1338!
1339 ic=0
1340 DO i=1,nt(ng)
1341 IF (ltracerclm(i,ng)) THEN
1342 ic=ic+1
1343 CALL get_3dfld (ng, inlm, idtclm(i), &
1344 & clmncid(idtclm(i),ng), &
1345# if defined PIO_LIB && defined DISTRIBUTE
1346 & clmpiofile(idtclm(i),ng), &
1347# endif
1348 & nclmfiles(ng), clm(1,ng), update(1), &
1349 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
1350# ifdef MASKING
1351 & grid(ng) % rmask, &
1352# endif
1353 & clima(ng) % tclmG(:,:,:,:,ic))
1354 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1355 END IF
1356 END DO
1357# endif
1358#endif
1359
1360#ifdef TLM_CHECK
1361!
1362!=======================================================================
1363! If tangent linear model check, read in nonlinear forward solution
1364! to compute dot product with perturbated nonlinear solution. Time
1365! interpolation between snapshot is not required (see subroutine
1366! "nl_dotproduct").
1367!=======================================================================
1368!
1369 IF (outer.ge.1) THEN
1370!
1371! Read in free-surface.
1372!
1373 CALL get_2dfld (ng, inlm, idfsur, fwd(ng)%ncid, &
1374# if defined PIO_LIB && defined DISTRIBUTE
1375 & fwd(ng)%pioFile, &
1376# endif
1377 & 1, fwd(ng), update(1), &
1378 & lbi, ubi, lbj, ubj, 2, 1, &
1379# ifdef MASKING
1380 & grid(ng) % rmask, &
1381# endif
1382 & ocean(ng) % zetaG)
1383 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1384!
1385! Read 2D momentum.
1386!
1387 CALL get_2dfld (ng, inlm, idubar, fwd(ng)%ncid, &
1388# if defined PIO_LIB && defined DISTRIBUTE
1389 & fwd(ng)%pioFile, &
1390# endif
1391 & 1, fwd(ng), update(1), &
1392 & lbi, ubi, lbj, ubj, 2, 1, &
1393# ifdef MASKING
1394 & grid(ng) % umask, &
1395# endif
1396 & ocean(ng) % ubarG)
1397 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1398
1399 CALL get_2dfld (ng, inlm, idvbar, fwd(ng)%ncid, &
1400# if defined PIO_LIB && defined DISTRIBUTE
1401 & fwd(ng)%pioFile, &
1402# endif
1403 & 1, fwd(ng), update(1), &
1404 & lbi, ubi, lbj, ubj, 2, 1, &
1405# ifdef MASKING
1406 & grid(ng) % vmask, &
1407# endif
1408 & ocean(ng) % vbarG)
1409 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1410# ifdef SOLVE3D
1411!
1412! Read in 3D momentum.
1413!
1414 CALL get_3dfld (ng, inlm, iduvel, fwd(ng)%ncid, &
1415# if defined PIO_LIB && defined DISTRIBUTE
1416 & fwd(ng)%pioFile, &
1417# endif
1418 & 1, fwd(ng), update(1), &
1419 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
1420# ifdef MASKING
1421 & grid(ng) % umask, &
1422# endif
1423 & ocean(ng) % uG)
1424 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1425
1426 CALL get_3dfld (ng, inlm, idvvel, fwd(ng)%ncid, &
1427# if defined PIO_LIB && defined DISTRIBUTE
1428 & fwd(ng)%pioFile, &
1429# endif
1430 & 1, fwd(ng), update(1), &
1431 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
1432# ifdef MASKING
1433 & grid(ng) % vmask, &
1434# endif
1435 & ocean(ng) % vG)
1436 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1437!
1438! Read in 3D tracers.
1439!
1440 DO i=1,nt(ng)
1441 CALL get_3dfld (ng, inlm, idtvar(i), fwd(ng)%ncid, &
1442# if defined PIO_LIB && defined DISTRIBUTE
1443 & fwd(ng)%pioFile, &
1444# endif
1445 & 1, fwd(ng), update(1), &
1446 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
1447# ifdef MASKING
1448 & grid(ng) % rmask, &
1449# endif
1450 & ocean(ng) % tG(:,:,:,:,i))
1451 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1452 END DO
1453# endif
1454 END IF
1455#endif
1456
1457#if defined NLM_OUTER || \
1458 defined rbl4dvar || \
1459 defined rbl4dvar_ana_sensitivity || \
1460 defined rbl4dvar_fct_sensitivity || \
1461 defined sp4dvar || \
1462 defined tl_rbl4dvar
1463!
1464!=======================================================================
1465! Read weak contraint forcing snapshots. Notice that the forward
1466! basic state snapshops arrays are reused here.
1467!=======================================================================
1468!
1469 IF (frequentimpulse(ng)) THEN
1470!
1471! Read in free-surface.
1472!
1473 CALL get_2dfld (ng, inlm, idfsur, tlf(ng)%ncid, &
1474# if defined PIO_LIB && defined DISTRIBUTE
1475 & tlf(ng)%pioFile, &
1476# endif
1477 & 1, tlf(ng), update(1), &
1478 & lbi, ubi, lbj, ubj, 2, 1, &
1479# ifdef MASKING
1480 & grid(ng) % rmask, &
1481# endif
1482 & ocean(ng) % zetaG)
1483 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1484
1485# ifndef SOLVE3D
1486!
1487! Read 2D momentum.
1488!
1489 CALL get_2dfld (ng, inlm, idubar, tlf(ng)%ncid, &
1490# if defined PIO_LIB && defined DISTRIBUTE
1491 & tlf(ng)%pioFile, &
1492# endif
1493 & 1, tlf(ng), update(1), &
1494 & lbi, ubi, lbj, ubj, 2, 1, &
1495# ifdef MASKING
1496 & grid(ng) % umask, &
1497# endif
1498 & ocean(ng) % ubarG)
1499 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1500
1501 CALL get_2dfld (ng, inlm, idvbar, tlf(ng)%ncid, &
1502# if defined PIO_LIB && defined DISTRIBUTE
1503 & tlf(ng)%pioFile, &
1504# endif
1505 & 1, tlf(ng), update(1), &
1506 & lbi, ubi, lbj, ubj, 2, 1, &
1507# ifdef MASKING
1508 & grid(ng) % vmask, &
1509# endif
1510 & ocean(ng) % vbarG)
1511 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1512# else
1513!
1514! Read in 3D momentum.
1515!
1516 CALL get_3dfld (ng, inlm, iduvel, tlf(ng)%ncid, &
1517# if defined PIO_LIB && defined DISTRIBUTE
1518 & tlf(ng)%pioFile, &
1519# endif
1520 & 1, tlf(ng), update(1), &
1521 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
1522# ifdef MASKING
1523 & grid(ng) % umask, &
1524# endif
1525 & ocean(ng) % uG)
1526 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1527
1528 CALL get_3dfld (ng, inlm, idvvel, tlf(ng)%ncid, &
1529# if defined PIO_LIB && defined DISTRIBUTE
1530 & tlf(ng)%pioFile, &
1531# endif
1532 & 1, tlf(ng), update(1), &
1533 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
1534# ifdef MASKING
1535 & grid(ng) % vmask, &
1536# endif
1537 & ocean(ng) % vG)
1538 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1539!
1540! Read in 3D tracers.
1541!
1542 DO i=1,nt(ng)
1543 CALL get_3dfld (ng, inlm, idtvar(i), tlf(ng)%ncid, &
1544# if defined PIO_LIB && defined DISTRIBUTE
1545 & tlf(ng)%pioFile, &
1546# endif
1547 & 1, tlf(ng), update(1), &
1548 & lbi, ubi, lbj, ubj, 1, n(ng), 2, 1, &
1549# ifdef MASKING
1550 & grid(ng) % rmask, &
1551# endif
1552 & ocean(ng) % tG(:,:,:,:,i))
1553 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1554 END DO
1555# endif
1556 END IF
1557#endif
1558
1559#ifdef PROFILE
1560!
1561!-----------------------------------------------------------------------
1562! Turn off input data time wall clock.
1563!-----------------------------------------------------------------------
1564!
1565 CALL wclock_off (ng, inlm, 3, __line__, myfile)
1566#endif
1567!
1568 RETURN
subroutine get_2dfld(ng, model, ifield, ncid, piofile, nfiles, s, update, lbi, ubi, lbj, ubj, iout, irec, fmask, fout)
Definition get_2dfld.F:12
subroutine get_3dfld(ng, model, ifield, ncid, piofile, nfiles, s, update, lbi, ubi, lbj, ubj, lbk, ubk, iout, irec, fmask, fout)
Definition get_3dfld.F:14
subroutine get_ngfld(ng, model, ifield, ncid, piofile, nfiles, s, recordless, update, lbi, ubi, ubj, ubk, istr, iend, jrec, fout)
Definition get_ngfld.F:9
integer idodin
integer idasrf
type(t_boundary), dimension(:), allocatable boundary
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_io), dimension(:), allocatable ssf
type(t_io), dimension(:), allocatable tlf
integer, dimension(:), allocatable nclmfiles
type(t_io), dimension(:), allocatable blk
type(t_io), dimension(:,:), allocatable frc
type(t_io), dimension(:,:), allocatable bry
type(t_io), dimension(:), allocatable fwd
integer, dimension(:), allocatable nbcfiles
type(t_io), dimension(:,:), allocatable clm
integer, dimension(:), allocatable nffiles
integer iddqdt
integer idvair
integer, dimension(:,:), allocatable clmncid
type(file_desc_t), dimension(:,:), pointer frcpiofile
integer, dimension(4) idu3bc
integer idrtra
integer idubar
integer, dimension(4) idzbry
integer idvvel
integer, dimension(4) idu2bc
integer idvsms
integer idcfra
integer idwlen
integer idwdiw
integer isvvel
integer isvbar
integer idvclm
integer, dimension(:), allocatable idrtrc
integer idpair
integer idwlep
integer, dimension(:), allocatable idtbot
integer, dimension(:), allocatable idtsur
integer idempf
integer, dimension(:), allocatable idtclm
integer, dimension(:,:), allocatable frcncid
integer idfsur
integer idwbrk
integer, dimension(:), allocatable idtvar
integer, dimension(:), allocatable istvar
integer, dimension(:,:), allocatable bryncid
integer idldwn
integer isuvel
integer idsssc
integer isfsur
integer iduclm
integer idsfwf
integer, dimension(4) idv3bc
integer iduair
integer iduvel
integer idqair
integer isubar
integer idubcl
integer, dimension(4) idv2bc
integer idlrad
integer idwdip
integer idusms
integer idwamp
integer idwdir
integer idwptp
integer idwdib
type(file_desc_t), dimension(:,:), pointer clmpiofile
type(file_desc_t), dimension(:,:), pointer brypiofile
integer idrain
integer idvbcl
integer idsrad
integer idsshc
integer idwpbt
integer idworb
integer idsstc
integer, dimension(:,:), allocatable idtbry
integer idtair
integer idvbar
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer nat
Definition mod_param.F:499
integer, parameter inlm
Definition mod_param.F:662
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
type(t_lbc), dimension(:,:,:), allocatable lbc
Definition mod_param.F:375
integer, dimension(:), allocatable lm
Definition mod_param.F:455
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer, dimension(:), allocatable mm
Definition mod_param.F:456
logical, dimension(:), allocatable luvsrc
logical, dimension(:,:), allocatable ltracersrc
logical, dimension(:), allocatable lprocessobc
integer, parameter iwest
logical, dimension(:), allocatable lm3clm
logical, dimension(:), allocatable lsshclm
logical, dimension(:), allocatable frequentimpulse
logical, dimension(:), allocatable lwsrc
integer exit_flag
integer isalt
integer itemp
integer, parameter isouth
logical, dimension(:), allocatable lm2clm
integer, parameter ieast
logical, dimension(:,:), allocatable ltracerclm
integer nrun
integer, parameter inorth
integer noerror
integer outer
type(t_sources), dimension(:), allocatable sources
Definition mod_sources.F:90
integer, dimension(:), allocatable nsrc
Definition mod_sources.F:97
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3

References mod_iounits::blk, mod_boundary::boundary, mod_param::bounds, mod_iounits::bry, mod_ncparam::bryncid, mod_ncparam::brypiofile, mod_clima::clima, mod_iounits::clm, mod_ncparam::clmncid, mod_ncparam::clmpiofile, mod_scalars::exit_flag, mod_forces::forces, strings_mod::founderror(), mod_iounits::frc, mod_ncparam::frcncid, mod_ncparam::frcpiofile, mod_scalars::frequentimpulse, mod_iounits::fwd, get_2dfld(), get_3dfld(), get_ngfld(), mod_grid::grid, mod_biology::idasrf, mod_ncparam::idcfra, mod_ncparam::iddqdt, mod_ncparam::idempf, mod_ncparam::idfsur, mod_ncparam::idldwn, mod_ncparam::idlrad, mod_biology::idodin, mod_ncparam::idpair, mod_ncparam::idqair, mod_ncparam::idrain, mod_biology::idresr, mod_ncparam::idrtra, mod_ncparam::idrtrc, mod_ncparam::idsfwf, mod_ncparam::idsrad, mod_ncparam::idsshc, mod_ncparam::idsssc, mod_ncparam::idsstc, mod_ncparam::idtair, mod_ncparam::idtbot, mod_ncparam::idtbry, mod_ncparam::idtclm, mod_ncparam::idtsur, mod_ncparam::idtvar, mod_ncparam::idu2bc, mod_ncparam::idu3bc, mod_ncparam::iduair, mod_ncparam::idubar, mod_ncparam::idubcl, mod_ncparam::iduclm, mod_ncparam::idusms, mod_ncparam::iduvel, mod_ncparam::idv2bc, mod_ncparam::idv3bc, mod_ncparam::idvair, mod_ncparam::idvbar, mod_ncparam::idvbcl, mod_ncparam::idvclm, mod_ncparam::idvsms, mod_ncparam::idvvel, mod_ncparam::idwamp, mod_ncparam::idwbrk, mod_ncparam::idwdib, mod_ncparam::idwdip, mod_ncparam::idwdir, mod_ncparam::idwdiw, mod_ncparam::idwlen, mod_ncparam::idwlep, mod_ncparam::idworb, mod_ncparam::idwpbt, mod_ncparam::idwptp, mod_ncparam::idzbry, mod_scalars::ieast, mod_param::inlm, mod_scalars::inorth, mod_scalars::isalt, mod_ncparam::isfsur, mod_scalars::isouth, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::itemp, mod_scalars::iwest, mod_param::lbc, mod_param::lm, mod_scalars::lm2clm, mod_scalars::lm3clm, mod_scalars::lprocessobc, mod_scalars::lsshclm, mod_scalars::ltracerclm, mod_scalars::ltracersrc, mod_scalars::luvsrc, mod_scalars::lwsrc, mod_param::mm, mod_param::n, mod_param::nat, mod_iounits::nbcfiles, mod_iounits::nclmfiles, mod_iounits::nffiles, mod_scalars::noerror, mod_scalars::nrun, mod_sources::nsrc, mod_param::nt, mod_ocean::ocean, mod_scalars::outer, mod_sources::sources, mod_iounits::ssf, mod_iounits::tlf, wclock_off(), and wclock_on().

Referenced by initial(), main3d(), and roms_kernel_mod::nlm_initial().

Here is the call graph for this function:
Here is the caller graph for this function: