#include "cppdefs.h" SUBROUTINE get_data (ng) ! !git $Id$ !svn $Id: get_data.F 1178 2023-07-11 17:50:57Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2023 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.md ! !======================================================================= ! ! ! This routine reads in forcing, climatology and other data from ! ! NetCDF files. If there is more than one time-record, data is ! ! loaded into global two-time record arrays. The interpolation ! ! is carried elsewhere. ! ! ! ! Currently, this routine is only executed in serial mode by the ! ! main thread. ! ! ! !======================================================================= ! USE mod_param #if defined HYPOXIA_SRM || defined RED_TIDE USE mod_biology #endif USE mod_boundary USE mod_clima USE mod_forces USE mod_grid USE mod_iounits USE mod_ncparam #if defined HYPOXIA_SRM || \ defined NLM_OUTER || \ defined RBL4DVAR || \ defined RBL4DVAR_ANA_SENSITIVITY || \ defined RBL4DVAR_FCT_SENSITIVITY || \ defined RED_TIDE || \ defined SP4DVAR || \ defined TLM_CHECK || \ defined TL_RBL4DVAR USE mod_ocean #endif USE mod_scalars USE mod_sources USE mod_stepping ! USE strings_mod, ONLY : FoundError ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng ! ! Local variable declarations. ! logical :: Lprocess logical, save :: recordless = .FALSE. logical, dimension(3) :: update = (/ .FALSE., .FALSE., .FALSE. /) ! integer :: ILB, IUB, JLB, JUB integer :: LBi, UBi, LBj, UBj integer :: i, ic, my_tile ! character (len=*), parameter :: MyFile = & & __FILE__ ! ! Lower and upper bounds for nontiled (global values) boundary arrays. ! my_tile=-1 ! for global values ILB=BOUNDS(ng)%LBi(my_tile) IUB=BOUNDS(ng)%UBi(my_tile) JLB=BOUNDS(ng)%LBj(my_tile) JUB=BOUNDS(ng)%UBj(my_tile) ! ! Lower and upper bounds for tiled arrays. ! LBi=LBOUND(GRID(ng)%h,DIM=1) UBi=UBOUND(GRID(ng)%h,DIM=1) LBj=LBOUND(GRID(ng)%h,DIM=2) UBj=UBOUND(GRID(ng)%h,DIM=2) #ifdef PROFILE ! !----------------------------------------------------------------------- ! Turn on input data time wall clock. !----------------------------------------------------------------------- ! CALL wclock_on (ng, iNLM, 3, __LINE__, MyFile) #endif #ifndef ANA_PSOURCE ! !======================================================================= ! Point Sources/Sinks time dependent data. !======================================================================= ! ! Point Source/Sink vertically integrated mass transport. ! IF (LuvSrc(ng).or.LwSrc(ng)) THEN CALL get_ngfld (ng, iNLM, idRtra, SSF(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & SSF(ng)%pioFile, & # endif & 1, SSF(ng), recordless, update(1), & & 1, Nsrc(ng), 1, 2, 1, Nsrc(ng), 1, & & SOURCES(ng) % QbarG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # ifdef SOLVE3D ! ! Tracer Sources/Sinks. ! DO i=1,NT(ng) IF (LtracerSrc(i,ng)) THEN CALL get_ngfld (ng, iNLM, idRtrc(i), SSF(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & SSF(ng)%pioFile, & # endif & 1, SSF(ng), recordless, update(1), & & 1, Nsrc(ng), N(ng), 2, 1, Nsrc(ng), N(ng), & & SOURCES(ng) % TsrcG(:,:,:,i)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF END DO # endif #endif ! !======================================================================= ! Read in forcing data from FORCING NetCDF file(s). !======================================================================= ! ! Set switch to process surface atmospheric fields. ! #if defined FOUR_DVAR && \ defined BULK_FLUXES && defined PRIOR_BULK_FLUXES ! In 4D-Var data assimilation applications, the user have the option ! to fix the prior (background phase) surface fluxes in the successive ! outer loops (Nouter>1) or the final analysis phase. In such case, the ! fluxes are read from the background trajector ! IF (Nrun.eq.1) THEN Lprocess=.TRUE. ELSE Lprocess=.FALSE. END IF #else Lprocess=.TRUE. #endif #if !(defined ANA_WINDS || defined FRC_COUPLING) && \ (defined BULK_FLUXES || defined ECOSIM) ! ! Surface wind components. ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idUair, FRCncid(idUair,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idUair,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % UwindG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ! CALL get_2dfld (ng , iNLM, idVair, FRCncid(idVair,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idVair,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % VwindG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF #endif #if !(defined ANA_SMFLUX || defined FRC_COUPLING || defined BULK_FLUXES) ! ! Surface wind stress components. ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idUsms, FRCncid(idUsms,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idUsms,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % umask, & # endif & FORCES(ng) % sustrG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ! CALL get_2dfld (ng, iNLM, idVsms, FRCncid(idVsms,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idVsms,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % vmask, & # endif & FORCES(ng) % svstrG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF #endif #if !(defined ANA_PAIR || defined FRC_COUPLING) && \ (defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS) ! ! Surface air pressure. ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idPair, FRCncid(idPair,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idPair,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % PairG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF #endif #if !defined ANA_WWAVE && defined WAVE_DATA ! ! Surface wind induced wave properties. ! # ifdef WAVES_DIR CALL get_2dfld (ng, iNLM, idWdir, FRCncid(idWdir,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWdir,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % DwaveG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifdef WAVES_DIRP ! CALL get_2dfld (ng, iNLM, idWdip, FRCncid(idWdip,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWdip,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % DwavepG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifdef WAVES_HEIGHT ! CALL get_2dfld (ng, iNLM, idWamp, FRCncid(idWamp,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWamp,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % HwaveG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifdef WAVES_LENGTH ! CALL get_2dfld (ng, iNLM, idWlen, FRCncid(idWlen,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWlen,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % LwaveG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifdef WAVES_LENGTHP ! CALL get_2dfld (ng, iNLM, idWlep, FRCncid(idWlep,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWlep,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % LwavepG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifdef WAVES_TOP_PERIOD ! CALL get_2dfld (ng, iNLM, idWptp, FRCncid(idWptp,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWptp,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % Pwave_topG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifdef WAVES_BOT_PERIOD ! CALL get_2dfld (ng, iNLM, idWpbt, FRCncid(idWpbt,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWpbt,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % Pwave_botG(:,:,1)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifdef WAVES_UB ! CALL get_2dfld (ng, iNLM, idWorb, FRCncid(idWorb,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWorb,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % Uwave_rmsG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifdef WAVES_DISS ! CALL get_2dfld (ng, iNLM, idWdib, FRCncid(idWdib,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWdib,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % Dissip_breakG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ! CALL get_2dfld (ng, iNLM, idWdiw, FRCncid(idWdiw,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWdiw,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % Dissip_wcapG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifdef ROLLER_SVENDSEN ! CALL get_2dfld (ng, iNLM, idWbrk, FRCncid(idWbrk,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idWbrk,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % Wave_breakG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif #endif #ifdef SOLVE3D # if !(defined ANA_CLOUD || defined FRC_COUPLING) && defined CLOUDS ! ! Cloud fraction. ! CALL get_2dfld (ng, iNLM, idCfra, FRCncid(idCfra,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idCfra,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % cloudG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # if !(defined ANA_SRFLUX || defined FRC_COUPLING) && defined SHORTWAVE ! ! Surface solar shortwave radiation. ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idSrad, FRCncid(idSrad,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idSrad,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % srflxG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # endif # if defined RED_TIDE && defined DAILY_SHORTWAVE ! ! Daily-averaged Surface solar shortwave radiation. ! CALL get_2dfld (ng, iNLM, idAsrf, FRCncid(idAsrf,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idAsrf,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % srflxG_avg) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # if defined BULK_FLUXES && \ !(defined FRC_COUPLING || defined LONGWAVE || defined LONGWAVE_OUT) ! ! Surface net longwave radiation. ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idLrad, FRCncid(idLrad,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idLrad,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % lrflxG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # endif # if defined BULK_FLUXES && defined LONGWAVE_OUT && \ !defined FRC_COUPLING ! ! Surface downwelling longwave radiation. ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idLdwn, FRCncid(idLdwn,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idLdwn,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % lrflxG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # endif # if !(defined ANA_TAIR || defined FRC_COUPLING) && \ ( defined BULK_FLUXES || defined ECOSIM || \ (defined SHORTWAVE && defined ANA_SRFLUX && defined ALBEDO) ) ! ! Surface air temperature. ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idTair, FRCncid(idTair,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idTair,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % TairG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # endif # if !(defined ANA_HUMIDITY || defined FRC_COUPLING) && \ (defined BULK_FLUXES || defined ECOSIM) ! ! Surface air humidity. ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idQair, FRCncid(idQair,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idQair,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % HairG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # endif # if !(defined ANA_RAIN || defined FRC_COUPLING) && defined BULK_FLUXES ! ! Rain fall rate. ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idrain, FRCncid(idrain,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idrain,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % rainG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # endif # if !(defined ANA_STFLUX || defined FRC_COUPLING || \ defined BULK_FLUXES) ! ! Surface net heat flux. ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idTsur(itemp), & & FRCncid(idTsur(itemp),ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idTsur(itemp),ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % stfluxG(:,:,:,itemp)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # endif # if !defined ANA_SST && defined QCORRECTION ! ! Surface net heat flux correction fields: sea surface temperature ! (SST). ! CALL get_2dfld (ng, iNLM, idSSTc, FRCncid(idSSTc,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idSSTc,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % sstG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # if !defined ANA_DQDSST && defined QCORRECTION ! ! Surface net heat flux correction fields: heat flux sensitivity to ! SST (dQdSST). ! CALL get_2dfld (ng, iNLM, iddQdT, FRCncid(iddQdT,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(iddQdT,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % dqdtG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifndef ANA_BTFLUX ! ! Bottom net heat flux. ! CALL get_2dfld (ng, iNLM, idTbot(itemp), & & FRCncid(idTbot(itemp),ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idTbot(itemp),ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % btfluxG(:,:,:,itemp)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # if !defined ANA_SSFLUX && defined SALINITY # if !(defined BULK_FLUXES || defined EMINUSP || \ defined FRC_COUPLING) ! ! Surface net freshwater flux: E-P from NetCDF variable "swflux". ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idsfwf, FRCncid(idsfwf,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idsfwf,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % stfluxG(:,:,:,isalt)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # elif defined BULK_FLUXES && !defined EMINUSP ! ! Surface net freshwater flux (E-P) from NetCDF variable "EminusP". ! IF (Lprocess) THEN CALL get_2dfld (ng, iNLM, idEmPf, FRCncid(idEMPf,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idEMPf,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % stfluxG(:,:,:,isalt)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # endif # if !defined ANA_SSS && (defined SCORRECTION || defined SRELAXATION) ! ! Surface net freshwater flux correction field: sea surface salinity. ! CALL get_2dfld (ng, iNLM, idSSSc, FRCncid(idSSSc,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idSSSc,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % sssG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifndef ANA_BSFLUX ! ! Bottom net freshwater flux. ! CALL get_2dfld (ng, iNLM, idTbot(isalt), & & FRCncid(idTbot(isalt),ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idTbot(isalt),ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % btfluxG(:,:,:,isalt)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # endif # if defined BIOLOGY || defined SEDIMENT || defined T_PASSIVE # ifndef ANA_SPFLUX ! ! Passive tracers surface fluxes. ! DO i=NAT+1,NT(ng) CALL get_2dfld (ng, iNLM, idTsur(i), FRCncid(idTsur(i),ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idTsur(i),ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % stfluxG(:,:,:,i)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END DO # endif # ifndef ANA_BPFLUX ! ! Passive tracers bottom fluxes. ! DO i=NAT+1,NT(ng) CALL get_2dfld (ng, iNLM, idTbot(i), FRCncid(idTbot(i),ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idTbot(i),ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % btfluxG(:,:,:,i)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END DO # endif # endif #endif #if !defined ANA_RESPIRATION && defined HYPOXIA_SRM ! ! Total respiration rate for hypoxia. ! CALL get_3dfld (ng, iNLM, idResR, FRCncid(idResR,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idResR,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & OCEAN(ng) % respirationG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN #endif #ifdef RED_TIDE ! ! Red tide Observed Dissolved Inorganic Nutrient. ! CALL get_3dfld (ng, iNLM, idODIN, FRCncid(idODIN,ng), & # if defined PIO_LIB && defined DISTRIBUTE & FRCpioFile(idODIN,ng), & # endif & nFfiles(ng), FRC(1,ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & OCEAN(ng) % DIN_obsG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN #endif #if defined FOUR_DVAR && \ defined BULK_FLUXES && defined PRIOR_BULK_FLUXES ! !======================================================================= ! In 4D-Var data assimilation algorithms, the user has the option to ! impose the prior (background phase) surface fluxes in the successive ! outer loops (Nouter>1) or the final analysis phase. Such fluxes were ! computed by routine "bulk_fluxes" and stored in the NLM background ! trajectory, BLK structure. Notice that "bulk_fluxes" is called in ! "main3d" only during the prior trajectory computation. ! ! The strategy is to save the 2D surface fluxes in the quicksave ! NetCDF file to allow frequent data records to resolve the daily ! cycle while maintaining its size small. !======================================================================= ! IF (.not.Lprocess) THEN ! ! Get prior surface wind stress components. ! CALL get_2dfld (ng, iNLM, idUsms, BLK(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & BLK(ng)%pioFile, & # endif & 1, BLK(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % umask, & # endif & FORCES(ng) % sustrG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ! CALL get_2dfld (ng, iNLM, idVsms, BLK(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & BLK(ng)%pioFile, & # endif & 1, BLK(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % vmask, & # endif & FORCES(ng) % svstrG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # ifdef ATM_PRESS ! ! Get prior surface air pressure. ! CALL get_2dfld (ng, iNLM, idPair, BLK(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & BLK(ng)%pioFile, & # endif & 1, BLK(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % PairG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # ifdef SOLVE3D # ifdef SHORTWAVE ! ! Get prior shortwave radiation flux. For consistency, we need to ! process the same values using in the computation of the net heat ! flux since there is a time interpolation from snapshots. ! CALL get_2dfld (ng, iNLM, idSrad, BLK(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & BLK(ng)%pioFile, & # endif & 1, BLK(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % srflxG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif ! ! Get prior surface net heat flux. ! CALL get_2dfld (ng, iNLM, idTsur(itemp), BLK(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & BLK(ng)%pioFile, & # endif & 1, BLK(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % stfluxG(:,:,:,itemp)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # ifdef SALINITY ! ! Get prior net freshwater flux (E-P). The kinematic net freshwater ! will be computed in "set_vbc" by multiplying with surface salinity. ! CALL get_2dfld (ng, iNLM, idEmPf, BLK(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & BLK(ng)%pioFile, & # endif & 1, BLK(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & FORCES(ng) % stfluxG(:,:,:,isalt)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # endif # endif END IF #endif ! !======================================================================= ! Read in open boundary conditions from BOUNDARY NetCDF file. !======================================================================= #ifndef ANA_FSOBC ! ! Free-surface open boundary conditions. ! IF (LprocessOBC(ng)) THEN IF (LBC(iwest,isFsur,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idZbry(iwest), & & BRYncid(idZbry(iwest),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idZbry(iwest),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, 1, 2, 0, Mm(ng)+1, 1, & & BOUNDARY(ng) % zetaG_west) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(ieast,isFsur,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idZbry(ieast), & & BRYncid(idZbry(ieast),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idZbry(ieast),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, 1, 2, 0, Mm(ng)+1, 1, & & BOUNDARY(ng) % zetaG_east) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(isouth,isFsur,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idZbry(isouth), & & BRYncid(idZbry(isouth),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idZbry(isouth),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, 1, 2, 0, Lm(ng)+1, 1, & & BOUNDARY(ng) % zetaG_south) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(inorth,isFsur,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idZbry(inorth), & & BRYncid(idZbry(inorth),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idZbry(inorth),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, 1, 2, 0, Lm(ng)+1, 1, & & BOUNDARY(ng) % zetaG_north) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF END IF #endif #ifndef ANA_M2OBC ! ! 2D momentum components open boundary conditions. ! IF (LprocessOBC(ng)) THEN IF (LBC(iwest,isUbar,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idU2bc(iwest), & & BRYncid(idU2bc(iwest),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idU2bc(iwest),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, 1, 2, 0, Mm(ng)+1, 1, & & BOUNDARY(ng) % ubarG_west) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(iwest,isVbar,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idV2bc(iwest), & & BRYncid(idV2bc(iwest),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idV2bc(iwest),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, 1, 2, 1, Mm(ng)+1, 1, & & BOUNDARY(ng) % vbarG_west) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(ieast,isUbar,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idU2bc(ieast), & & BRYncid(idU2bc(ieast),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idU2bc(ieast),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, 1, 2, 0, Mm(ng)+1, 1, & & BOUNDARY(ng) % ubarG_east) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(ieast,isVbar,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idV2bc(ieast), & & BRYncid(idV2bc(ieast),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idV2bc(ieast),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, 1, 2, 1, Mm(ng)+1, 1, & & BOUNDARY(ng) % vbarG_east) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(isouth,isUbar,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idU2bc(isouth), & & BRYncid(idU2bc(isouth),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idU2bc(isouth),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, 1, 2, 1, Lm(ng)+1, 1, & & BOUNDARY(ng) % ubarG_south) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(isouth,isVbar,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idV2bc(isouth), & & BRYncid(idV2bc(isouth),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idV2bc(isouth),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, 1, 2, 0, Lm(ng)+1, 1, & & BOUNDARY(ng) % vbarG_south) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(inorth,isUbar,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idU2bc(inorth), & & BRYncid(idU2bc(inorth),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idU2bc(inorth),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, 1, 2, 1, Lm(ng)+1, 1, & & BOUNDARY(ng) % ubarG_north) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(inorth,isVbar,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idV2bc(inorth), & & BRYncid(idV2bc(inorth),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idV2bc(inorth),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, 1, 2, 0, Lm(ng)+1, 1, & & BOUNDARY(ng) % vbarG_north) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF END IF #endif #ifdef SOLVE3D # ifndef ANA_M3OBC ! ! 3D momentum components open boundary conditions. ! IF (LprocessOBC(ng)) THEN IF (LBC(iwest,isUvel,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idU3bc(iwest), & & BRYncid(idU3bc(iwest),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idU3bc(iwest),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, N(ng), 2, 0, Mm(ng)+1, N(ng), & & BOUNDARY(ng) % uG_west) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(iwest,isVvel,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idV3bc(iwest), & & BRYncid(idV3bc(iwest),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idV3bc(iwest),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, N(ng), 2, 1, Mm(ng)+1, N(ng), & & BOUNDARY(ng) % vG_west) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(ieast,isUvel,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idU3bc(ieast), & & BRYncid(idU3bc(ieast),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idU3bc(ieast),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, N(ng), 2, 0, Mm(ng)+1, N(ng), & & BOUNDARY(ng) % uG_east) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(ieast,isVvel,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idV3bc(ieast), & & BRYncid(idV3bc(ieast),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idV3bc(ieast),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, N(ng), 2, 1, Mm(ng)+1, N(ng), & & BOUNDARY(ng) % vG_east) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(isouth,isUvel,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idU3bc(isouth), & & BRYncid(idU3bc(isouth),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idU3bc(isouth),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, N(ng), 2, 1, Lm(ng)+1, N(ng), & & BOUNDARY(ng) % uG_south) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(isouth,isVvel,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idV3bc(isouth), & & BRYncid(idV3bc(isouth),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idV3bc(isouth),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, N(ng), 2, 0, Lm(ng)+1, N(ng), & & BOUNDARY(ng) % vG_south) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(inorth,isUvel,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idU3bc(inorth), & & BRYncid(idU3bc(inorth),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idU3bc(inorth),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, N(ng), 2, 1, Lm(ng)+1, N(ng), & & BOUNDARY(ng) % uG_north) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF ! IF (LBC(inorth,isVvel,ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idV3bc(inorth), & & BRYncid(idV3bc(inorth),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idV3bc(inorth),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, N(ng), 2, 0, Lm(ng)+1, N(ng), & & BOUNDARY(ng) % vG_north) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF END IF # endif # ifndef ANA_TOBC ! ! Tracer variables open boundary conditions. ! IF (LprocessOBC(ng)) THEN DO i=1,NT(ng) IF (LBC(iwest,isTvar(i),ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idTbry(iwest,i), & & BRYncid(idTbry(iwest,i),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idTbry(iwest,i),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, N(ng), 2, 0, Mm(ng)+1, N(ng), & & BOUNDARY(ng) % tG_west(:,:,:,i)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF END DO ! DO i=1,NT(ng) IF (LBC(ieast,isTvar(i),ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idTbry(ieast,i), & & BRYncid(idTbry(ieast,i),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idTbry(ieast,i),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & JLB, JUB, N(ng), 2, 0, Mm(ng)+1, N(ng), & & BOUNDARY(ng) % tG_east(:,:,:,i)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF END DO ! DO i=1,NT(ng) IF (LBC(isouth,isTvar(i),ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idTbry(isouth,i), & & BRYncid(idTbry(isouth,i),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idTbry(isouth,i),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, N(ng), 2, 0, Lm(ng)+1, N(ng), & & BOUNDARY(ng) % tG_south(:,:,:,i)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF END DO ! DO i=1,NT(ng) IF (LBC(inorth,isTvar(i),ng)%acquire) THEN CALL get_ngfld (ng, iNLM, idTbry(inorth,i), & & BRYncid(idTbry(inorth,i),ng), & # if defined PIO_LIB && defined DISTRIBUTE & BRYpioFile(idTbry(inorth,i),ng), & # endif & nBCfiles(ng), BRY(1,ng), & & recordless, update(1), & & ILB, IUB, N(ng), 2, 0, Lm(ng)+1, N(ng), & & BOUNDARY(ng) % tG_north(:,:,:,i)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF END DO END IF # endif #endif ! !======================================================================= ! Read in data from Climatology NetCDF file. !======================================================================= #ifndef ANA_SSH ! ! Free-surface climatology. ! IF (LsshCLM(ng)) THEN CALL get_2dfld (ng, iNLM, idSSHc, CLMncid(idSSHc,ng), & # if defined PIO_LIB && defined DISTRIBUTE & CLMpioFile(idSSHc,ng), & # endif & nCLMfiles(ng), CLM(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & CLIMA(ng) % sshG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF #endif #ifndef ANA_M2CLIMA ! ! 2D momentum components climatology. ! IF (Lm2CLM(ng)) THEN CALL get_2dfld (ng, iNLM, idUbcl, CLMncid(idUbcl,ng), & # if defined PIO_LIB && defined DISTRIBUTE & CLMpioFile(idUbcl,ng), & # endif & nCLMfiles(ng), CLM(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % umask, & # endif & CLIMA(ng) % ubarclmG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ! CALL get_2dfld (ng, iNLM, idVbcl, CLMncid(idVbcl,ng), & # if defined PIO_LIB && defined DISTRIBUTE & CLMpioFile(idVbcl,ng), & # endif & nCLMfiles(ng), CLM(1,ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % vmask, & # endif & CLIMA(ng) % vbarclmG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF #endif #ifdef SOLVE3D # ifndef ANA_M3CLIMA ! ! 3D momentum components climatology. ! IF (Lm3CLM(ng)) THEN CALL get_3dfld (ng, iNLM, idUclm, CLMncid(idUclm,ng), & # if defined PIO_LIB && defined DISTRIBUTE & CLMpioFile(idUclm,ng), & # endif & nCLMfiles(ng), CLM(1,ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % umask, & # endif & CLIMA(ng) % uclmG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ! CALL get_3dfld (ng, iNLM, idVclm, CLMncid(idVclm,ng), & # if defined PIO_LIB && defined DISTRIBUTE & CLMpioFile(idVclm,ng), & # endif & nCLMfiles(ng), CLM(1,ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % vmask, & # endif & CLIMA(ng) % vclmG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF # endif # ifndef ANA_TCLIMA ! ! Tracers variables climatology. ! ic=0 DO i=1,NT(ng) IF (LtracerCLM(i,ng)) THEN ic=ic+1 CALL get_3dfld (ng, iNLM, idTclm(i), & & CLMncid(idTclm(i),ng), & # if defined PIO_LIB && defined DISTRIBUTE & CLMpioFile(idTclm(i),ng), & # endif & nCLMfiles(ng), CLM(1,ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & CLIMA(ng) % tclmG(:,:,:,:,ic)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF END DO # endif #endif #ifdef TLM_CHECK ! !======================================================================= ! If tangent linear model check, read in nonlinear forward solution ! to compute dot product with perturbated nonlinear solution. Time ! interpolation between snapshot is not required (see subroutine ! "nl_dotproduct"). !======================================================================= ! IF (outer.ge.1) THEN ! ! Read in free-surface. ! CALL get_2dfld (ng, iNLM, idFsur, FWD(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & FWD(ng)%pioFile, & # endif & 1, FWD(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & OCEAN(ng) % zetaG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ! ! Read 2D momentum. ! CALL get_2dfld (ng, iNLM, idUbar, FWD(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & FWD(ng)%pioFile, & # endif & 1, FWD(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % umask, & # endif & OCEAN(ng) % ubarG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN CALL get_2dfld (ng, iNLM, idVbar, FWD(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & FWD(ng)%pioFile, & # endif & 1, FWD(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % vmask, & # endif & OCEAN(ng) % vbarG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # ifdef SOLVE3D ! ! Read in 3D momentum. ! CALL get_3dfld (ng, iNLM, idUvel, FWD(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & FWD(ng)%pioFile, & # endif & 1, FWD(ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % umask, & # endif & OCEAN(ng) % uG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN CALL get_3dfld (ng, iNLM, idVvel, FWD(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & FWD(ng)%pioFile, & # endif & 1, FWD(ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % vmask, & # endif & OCEAN(ng) % vG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ! ! Read in 3D tracers. ! DO i=1,NT(ng) CALL get_3dfld (ng, iNLM, idTvar(i), FWD(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & FWD(ng)%pioFile, & # endif & 1, FWD(ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & OCEAN(ng) % tG(:,:,:,:,i)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END DO # endif END IF #endif #if defined NLM_OUTER || \ defined RBL4DVAR || \ defined RBL4DVAR_ANA_SENSITIVITY || \ defined RBL4DVAR_FCT_SENSITIVITY || \ defined SP4DVAR || \ defined TL_RBL4DVAR ! !======================================================================= ! Read weak contraint forcing snapshots. Notice that the forward ! basic state snapshops arrays are reused here. !======================================================================= ! IF (FrequentImpulse(ng)) THEN ! ! Read in free-surface. ! CALL get_2dfld (ng, iNLM, idFsur, TLF(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & TLF(ng)%pioFile, & # endif & 1, TLF(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & OCEAN(ng) % zetaG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # ifndef SOLVE3D ! ! Read 2D momentum. ! CALL get_2dfld (ng, iNLM, idUbar, TLF(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & TLF(ng)%pioFile, & # endif & 1, TLF(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % umask, & # endif & OCEAN(ng) % ubarG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN CALL get_2dfld (ng, iNLM, idVbar, TLF(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & TLF(ng)%pioFile, & # endif & 1, TLF(ng), update(1), & & LBi, UBi, LBj, UBj, 2, 1, & # ifdef MASKING & GRID(ng) % vmask, & # endif & OCEAN(ng) % vbarG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN # else ! ! Read in 3D momentum. ! CALL get_3dfld (ng, iNLM, idUvel, TLF(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & TLF(ng)%pioFile, & # endif & 1, TLF(ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % umask, & # endif & OCEAN(ng) % uG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN CALL get_3dfld (ng, iNLM, idVvel, TLF(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & TLF(ng)%pioFile, & # endif & 1, TLF(ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % vmask, & # endif & OCEAN(ng) % vG) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ! ! Read in 3D tracers. ! DO i=1,NT(ng) CALL get_3dfld (ng, iNLM, idTvar(i), TLF(ng)%ncid, & # if defined PIO_LIB && defined DISTRIBUTE & TLF(ng)%pioFile, & # endif & 1, TLF(ng), update(1), & & LBi, UBi, LBj, UBj, 1, N(ng), 2, 1, & # ifdef MASKING & GRID(ng) % rmask, & # endif & OCEAN(ng) % tG(:,:,:,:,i)) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END DO # endif END IF #endif #ifdef PROFILE ! !----------------------------------------------------------------------- ! Turn off input data time wall clock. !----------------------------------------------------------------------- ! CALL wclock_off (ng, iNLM, 3, __LINE__, MyFile) #endif ! RETURN END SUBROUTINE get_data