#include "cppdefs.h" #ifdef TL_IOMS SUBROUTINE rp_output (ng) ! !svn $Id: rp_output.F 889 2018-02-10 03:32:52Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2019 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.txt ! !======================================================================= ! ! ! This routine manages tangent linear model output. It creates output ! ! NetCDF files and writes out data into NetCDF files. If requested, ! ! it can create several tangent history files to avoid generating too ! ! large files during a single model run. ! ! ! !======================================================================= ! USE mod_param USE mod_parallel # ifdef FOUR_DVAR USE mod_fourdvar # endif USE mod_iounits USE mod_ncparam USE mod_netcdf USE mod_scalars ! # ifdef DISTRIBUTE USE distribute_mod, ONLY : mp_bcasts # endif USE strings_mod, ONLY : FoundError ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng ! ! Local variable declarations. ! logical, save :: First = .TRUE. logical :: Ldefine, NewFile integer :: Fcount, ifile, status, tile ! SourceFile=__FILE__ # ifdef PROFILE ! !----------------------------------------------------------------------- ! Turn on output data time wall clock. !----------------------------------------------------------------------- ! CALL wclock_on (ng, iRPM, 8, __LINE__, __FILE__) # endif ! !----------------------------------------------------------------------- ! If appropriate, process tangent linear history NetCDF file. !----------------------------------------------------------------------- ! ! Turn off checking for analytical header files. ! IF (Lanafile) THEN Lanafile=.FALSE. END IF ! ! Create output tangent NetCDF file or prepare existing file to ! append new data to it. Also, notice that it is possible to ! create several files during a single model run. ! IF (LdefTLM(ng)) THEN IF (ndefTLM(ng).gt.0) THEN IF (idefTLM(ng).lt.0) THEN idefTLM(ng)=((ntstart(ng)-1)/ndefTLM(ng))*ndefTLM(ng) IF (idefTLM(ng).lt.iic(ng)-1) THEN idefTLM(ng)=idefTLM(ng)+ndefTLM(ng) END IF END IF IF ((nrrec(ng).ne.0).and.(iic(ng).eq.ntstart(ng))) THEN IF ((iic(ng)-1).eq.idefTLM(ng)) THEN Ldefine=.FALSE. ! finished file, delay ELSE ! creation of next file Ldefine=.TRUE. NewFile=.FALSE. ! unfinished file, inquire END IF ! content for appending idefTLM(ng)=idefTLM(ng)+nTLM(ng) ! restart offset ELSE IF ((iic(ng)-1).eq.idefTLM(ng)) THEN idefTLM(ng)=idefTLM(ng)+ndefTLM(ng) IF (nTLM(ng).ne.ndefTLM(ng).and.iic(ng).eq.ntstart(ng)) THEN idefTLM(ng)=idefTLM(ng)+nTLM(ng) ! multiple record offset END IF Ldefine=.TRUE. NewFile=.TRUE. ELSE Ldefine=.FALSE. END IF IF (Ldefine) THEN ! create new file or Fcount=TLM(ng)%Fcount ! inquire existing file TLM(ng)%Nrec(Fcount)=0 ifile=(iic(ng)-1)/ndefTLM(ng)+1 IF (Master) THEN WRITE (TLM(ng)%name,10) TRIM(TLM(ng)%base), ifile 10 FORMAT (a,'_',i4.4,'.nc') END IF # ifdef DISTRIBUTE CALL mp_bcasts (ng, iRPM, TLM(ng)%name) # endif IF (TLM(ng)%ncid.ne.-1) THEN CALL netcdf_close (ng, iRPM, TLM(ng)%ncid) END IF CALL tl_def_his (ng, NewFile) IF (FoundError(exit_flag, NoError, __LINE__, & & __FILE__)) RETURN END IF IF ((iic(ng).eq.ntstart(ng)).and.(nrrec(ng).ne.0)) THEN LwrtTLM(ng)=.FALSE. ! avoid writing initial ELSE ! fields during restart LwrtTLM(ng)=.TRUE. END IF ELSE IF (iic(ng).eq.ntstart(ng)) THEN CALL tl_def_his (ng, ldefout(ng)) IF (FoundError(exit_flag, NoError, __LINE__, & & __FILE__)) RETURN LwrtTLM(ng)=.TRUE. LdefTLM(ng)=.FALSE. END IF END IF END IF ! ! Write out data into tangent NetCDF file. Avoid writing initial ! conditions in perturbation mode computations. ! IF (LwrtTLM(ng)) THEN IF (LwrtPER(ng)) THEN IF ((iic(ng).gt.ntstart(ng)).and. & & (MOD(iic(ng)-1,nTLM(ng)).eq.0)) THEN CALL tl_wrt_his (ng) IF (FoundError(exit_flag, NoError, __LINE__, & & __FILE__)) RETURN END IF ELSE IF ((MOD(iic(ng)-1,nTLM(ng)).eq.0).and. & & ((nrrec(ng).eq.0).or.(iic(ng).ne.ntstart(ng)))) THEN CALL tl_wrt_his (ng) IF (FoundError(exit_flag, NoError, __LINE__, & & __FILE__)) RETURN END IF END IF END IF # ifdef RP_AVERAGES ! !----------------------------------------------------------------------- ! If appropriate, process time-averaged NetCDF file. !----------------------------------------------------------------------- ! ! Create output time-averaged NetCDF file or prepare existing file ! to append new data to it. Also, notice that it is possible to ! create several files during a single model run. ! IF (LdefAVG(ng)) THEN IF (ndefAVG(ng).gt.0) THEN IF (idefAVG(ng).lt.0) THEN idefAVG(ng)=((ntstart(ng)-1)/ndefAVG(ng))*ndefAVG(ng) IF ((ndefAVG(ng).eq.nAVG(ng)).and.(idefAVG(ng).le.0)) THEN idefAVG(ng)=ndefAVG(ng) ! one file per record ELSE IF (idefAVG(ng).lt.iic(ng)-1) THEN idefAVG(ng)=idefAVG(ng)+ndefAVG(ng) END IF END IF IF ((nrrec(ng).ne.0).and.(iic(ng).eq.ntstart(ng))) THEN IF ((iic(ng)-1).eq.idefAVG(ng)) THEN Ldefine=.FALSE. ! finished file, delay ELSE ! creation of next file NewFile=.FALSE. Ldefine=.TRUE. ! unfinished file, inquire END IF ! content for appending idefAVG(ng)=idefAVG(ng)+nAVG(ng) ! restart offset ELSE IF ((iic(ng)-1).eq.idefAVG(ng)) THEN idefAVG(ng)=idefAVG(ng)+ndefAVG(ng) IF (nAVG(ng).ne.ndefAVG(ng).and.iic(ng).eq.ntstart(ng)) THEN idefAVG(ng)=idefAVG(ng)+nAVG(ng) END IF Ldefine=.TRUE. Newfile=.TRUE. ELSE Ldefine=.FALSE. END IF IF (Ldefine) THEN Fcount=AVG(ng)%Fcount AVG(ng)%Nrec(Fcount)=0 IF (ndefAVG(ng).eq.nAVG(ng)) THEN ifile=(iic(ng)-1)/ndefAVG(ng) ELSE ifile=(iic(ng)-1)/ndefAVG(ng)+1 END IF IF (Master) THEN WRITE (AVG(ng)%name,20) TRIM(AVG(ng)%base), ifile 20 FORMAT (a,'_',i4.4,'.nc') END IF # ifdef DISTRIBUTE CALL mp_bcasts (ng, iRPM, AVG(ng)%name) # endif IF (AVG(ng)%ncid.ne.-1) THEN CALL netcdf_close (ng, iRPM, AVG(ng)%ncid) END IF CALL def_avg (ng, Newfile) IF (FoundError(exit_flag, NoError, __LINE__, & & __FILE__)) RETURN LwrtAVG(ng)=.TRUE. END IF ELSE IF (iic(ng).eq.ntstart(ng)) THEN CALL def_avg (ng, ldefout(ng)) IF (FoundError(exit_flag, NoError, __LINE__, & & __FILE__)) RETURN LwrtAVG(ng)=.TRUE. LdefAVG(ng)=.FALSE. END IF END IF END IF ! ! Write out data into time-averaged NetCDF file. ! IF (LwrtAVG(ng)) THEN IF ((iic(ng).gt.ntstart(ng)).and. & & (MOD(iic(ng)-1,nAVG(ng)).eq.0)) THEN CALL wrt_avg (ng) IF (FoundError(exit_flag, NoError, __LINE__, & & __FILE__)) RETURN END IF END IF # endif # ifdef FOUR_DVAR # ifndef WEAK_CONSTRAINT ! !----------------------------------------------------------------------- ! Create tangent linear model initial conditions file, if necessary. !----------------------------------------------------------------------- ! ! If start of descent algorithm iterations, create initial conditions ! file or prepare existing file to append new data to it. ! IF (First) THEN First=.FALSE. CALL tl_def_ini (ng) IF (FoundError(exit_flag, NoError, __LINE__, & & __FILE__)) RETURN END IF # endif # ifdef OBSERVATIONS ! !----------------------------------------------------------------------- ! If appropriate, process and write model state at observation ! locations. Compute misfit (model-observations) cost function. !----------------------------------------------------------------------- ! IF (((time(ng)-0.5_r8*dt(ng)).le.ObsTime(ng)).and. & & (ObsTime(ng).lt.(time(ng)+0.5_r8*dt(ng)))) THEN ProcessObs(ng)=.TRUE. # ifdef DISTRIBUTE tile=MyRank # else tile=-1 # endif CALL obs_read (ng, iRPM, .FALSE.) CALL obs_write (ng, tile, iRPM) # if defined TL_W4DVAR || defined W4DVAR || defined W4DVAR_SENSITIVITY CALL obs_cost (ng, iRPM) # endif ELSE ProcessObs(ng)=.FALSE. END IF # endif # endif # ifdef PROFILE ! !----------------------------------------------------------------------- ! Turn off output data time wall clock. !----------------------------------------------------------------------- ! CALL wclock_off (ng, iRPM, 8, __LINE__, __FILE__) # endif RETURN END SUBROUTINE rp_output #else SUBROUTINE rp_output RETURN END SUBROUTINE rp_output #endif