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

Go to the source code of this file.

Functions/Subroutines

subroutine tl_output (ng)
 

Function/Subroutine Documentation

◆ tl_output()

subroutine tl_output ( integer, intent(in) ng)

Definition at line 3 of file tl_output.F.

4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine manages tangent linear model output. It creates output !
13! NetCDF files and writes out data into NetCDF files. If requested, !
14! it can create several tangent history files to avoid generating too !
15! large files during a single model run. !
16! !
17!=======================================================================
18!
19 USE mod_param
20 USE mod_parallel
21# ifdef FOUR_DVAR
22 USE mod_fourdvar
23# endif
24 USE mod_iounits
25 USE mod_ncparam
26 USE mod_scalars
27!
28 USE close_io_mod, ONLY : close_file
29# ifdef TL_AVERAGES
30 USE def_avg_mod, ONLY : def_avg
31# endif
32# ifdef DISTRIBUTE
33 USE distribute_mod, ONLY : mp_bcasts
34# endif
35# ifdef OBSERVATIONS
36 USE obs_read_mod, ONLY : obs_read
37 USE obs_write_mod, ONLY : obs_write
38# endif
39 USE strings_mod, ONLY : founderror
40 USE tl_def_his_mod, ONLY : tl_def_his
41 USE tl_wrt_his_mod, ONLY : tl_wrt_his
42# ifdef TL_AVERAGES
43 USE wrt_avg_mod, ONLY : def_avg
44# endif
45!
46 implicit none
47!
48! Imported variable declarations.
49!
50 integer, intent(in) :: ng
51!
52! Local variable declarations.
53!
54 logical :: Ldefine, Lupdate, NewFile
55!
56 integer :: Fcount, ifile, status, tile
57!
58 character (len=*), parameter :: MyFile = &
59 & __FILE__
60!
61 sourcefile=myfile
62
63# ifdef PROFILE
64!
65!-----------------------------------------------------------------------
66! Turn on output data time wall clock.
67!-----------------------------------------------------------------------
68!
69 CALL wclock_on (ng, itlm, 8, __line__, myfile)
70# endif
71!
72!-----------------------------------------------------------------------
73! If appropriate, process tangent linear history NetCDF file.
74!-----------------------------------------------------------------------
75!
76! Turn off checking for analytical header files.
77!
78 IF (lanafile) THEN
79 lanafile=.false.
80 END IF
81!
82! If appropriate, set switch for updating biology header file global
83! attribute in output NetCDF files.
84!
85#ifdef BIOLOGY
86 lupdate=.true.
87#else
88 lupdate=.false.
89#endif
90!
91! Create output tangent NetCDF file or prepare existing file to
92! append new data to it. Also, notice that it is possible to
93! create several files during a single model run.
94!
95 IF (ldeftlm(ng)) THEN
96 IF (ndeftlm(ng).gt.0) THEN
97 IF (ideftlm(ng).lt.0) THEN
98 ideftlm(ng)=((ntstart(ng)-1)/ndeftlm(ng))*ndeftlm(ng)
99 IF (ideftlm(ng).lt.iic(ng)-1) THEN
100 ideftlm(ng)=ideftlm(ng)+ndeftlm(ng)
101 END IF
102 END IF
103 IF ((nrrec(ng).ne.0).and.(iic(ng).eq.ntstart(ng))) THEN
104 IF ((iic(ng)-1).eq.ideftlm(ng)) THEN
105 tlm(ng)%load=0 ! restart, reset counter
106 ldefine=.false. ! finished file, delay
107 ELSE ! creation of next file
108 ldefine=.true.
109 newfile=.false. ! unfinished file, inquire
110 END IF ! content for appending
111 ideftlm(ng)=ideftlm(ng)+ntlm(ng) ! restart offset
112 ELSE IF ((iic(ng)-1).eq.ideftlm(ng)) THEN
113 ideftlm(ng)=ideftlm(ng)+ndeftlm(ng)
114 IF (ntlm(ng).ne.ndeftlm(ng).and.iic(ng).eq.ntstart(ng)) THEN
115 ideftlm(ng)=ideftlm(ng)+ntlm(ng) ! multiple record offset
116 END IF
117 ldefine=.true.
118 newfile=.true.
119 ELSE
120 ldefine=.false.
121 END IF
122 IF (ldefine) THEN ! create new file or
123 IF (iic(ng).eq.ntstart(ng)) THEN ! inquire existing file
124 tlm(ng)%load=0 ! reset filename counter
125 END IF
126 ifile=(iic(ng)-1)/ndeftlm(ng)+1 ! next filename suffix
127 tlm(ng)%load=tlm(ng)%load+1
128 IF (tlm(ng)%load.gt.tlm(ng)%Nfiles) THEN
129 IF (master) THEN
130 WRITE (stdout,10) 'TLM(ng)%load = ', tlm(ng)%load, &
131 & tlm(ng)%Nfiles, trim(tlm(ng)%base), &
132 & ifile
133 END IF
134 exit_flag=4
136 & __line__, myfile)) RETURN
137 END IF
138 fcount=tlm(ng)%load
139 tlm(ng)%Nrec(fcount)=0
140 IF (master) THEN
141 WRITE (tlm(ng)%name,20) trim(tlm(ng)%base), ifile
142 END IF
143# ifdef DISTRIBUTE
144 CALL mp_bcasts (ng, itlm, tlm(ng)%name)
145# endif
146 tlm(ng)%files(fcount)=trim(tlm(ng)%name)
147 CALL close_file (ng, itlm, tlm(ng), tlm(ng)%name, lupdate)
148 CALL tl_def_his (ng, newfile)
149 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
150 END IF
151 IF ((iic(ng).eq.ntstart(ng)).and.(nrrec(ng).ne.0)) THEN
152 lwrttlm(ng)=.false. ! avoid writing initial
153 ELSE ! fields during restart
154 lwrttlm(ng)=.true.
155 END IF
156 ELSE
157 IF (iic(ng).eq.ntstart(ng)) THEN
158 CALL tl_def_his (ng, ldefout(ng))
159 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
160 lwrttlm(ng)=.true.
161 ldeftlm(ng)=.false.
162 END IF
163 END IF
164 END IF
165!
166! Write out data into tangent NetCDF file. Avoid writing initial
167! conditions in perturbation mode computations.
168!
169 IF (lwrttlm(ng)) THEN
170# if defined STOCHASTIC_OPT || defined FORCING_SV || \
171 defined hessian_so || defined hessian_fsv
172!
173! Write only on first time step.
174!
175 IF (iic(ng).eq.1) THEN
176# ifdef DISTRIBUTE
177 CALL tl_wrt_his (ng, myrank)
178# else
179 CALL tl_wrt_his (ng, -1)
180# endif
181 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
182 END IF
183# else
184 IF (lwrtper(ng)) THEN
185 IF ((iic(ng).gt.ntstart(ng)).and. &
186 & (mod(iic(ng)-1,ntlm(ng)).eq.0)) THEN
187# ifdef DISTRIBUTE
188 CALL tl_wrt_his (ng, myrank)
189# else
190 CALL tl_wrt_his (ng, -1)
191# endif
192 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
193 END IF
194 ELSE
195 IF ((mod(iic(ng)-1,ntlm(ng)).eq.0).and. &
196 & ((nrrec(ng).eq.0).or.(iic(ng).ne.ntstart(ng)))) THEN
197# ifdef DISTRIBUTE
198 CALL tl_wrt_his (ng, myrank)
199# else
200 CALL tl_wrt_his (ng, -1)
201# endif
202 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
203 END IF
204 END IF
205# endif
206 END IF
207
208# ifdef TL_AVERAGES
209!
210!-----------------------------------------------------------------------
211! If appropriate, process time-averaged NetCDF file.
212!-----------------------------------------------------------------------
213!
214! Create output time-averaged NetCDF file or prepare existing file
215! to append new data to it. Also, notice that it is possible to
216! create several files during a single model run.
217!
218 IF (ldefavg(ng)) THEN
219 IF (ndefavg(ng).gt.0) THEN
220 IF (idefavg(ng).lt.0) THEN
221 idefavg(ng)=((ntstart(ng)-1)/ndefavg(ng))*ndefavg(ng)
222 IF ((ndefavg(ng).eq.navg(ng)).and.(idefavg(ng).le.0)) THEN
223 idefavg(ng)=ndefavg(ng) ! one file per record
224 ELSE IF (idefavg(ng).lt.iic(ng)-1) THEN
225 idefavg(ng)=idefavg(ng)+ndefavg(ng)
226 END IF
227 END IF
228 IF ((nrrec(ng).ne.0).and.(iic(ng).eq.ntstart(ng))) THEN
229 IF ((iic(ng)-1).eq.idefavg(ng)) THEN
230 avg(ng)%load=0 ! restart, reset counter
231 ldefine=.false. ! finished file, delay
232 ELSE ! creation of next file
233 newfile=.false.
234 ldefine=.true. ! unfinished file, inquire
235 END IF ! content for appending
236 idefavg(ng)=idefavg(ng)+navg(ng) ! restart offset
237 ELSE IF ((iic(ng)-1).eq.idefavg(ng)) THEN
238 idefavg(ng)=idefavg(ng)+ndefavg(ng)
239 IF (navg(ng).ne.ndefavg(ng).and.iic(ng).eq.ntstart(ng)) THEN
240 idefavg(ng)=idefavg(ng)+navg(ng)
241 END IF
242 ldefine=.true.
243 newfile=.true.
244 ELSE
245 ldefine=.false.
246 END IF
247 IF (ldefine) THEN
248 IF (iic(ng).eq.ntstart(ng)) THEN
249 avg(ng)%load=0 ! reset filename counter
250 END IF
251 avg(ng)%load=avg(ng)%load+1
252 IF (ndefavg(ng).eq.navg(ng)) THEN ! next filename suffix
253 ifile=(iic(ng)-1)/ndefavg(ng)
254 ELSE
255 ifile=(iic(ng)-1)/ndefavg(ng)+1
256 END IF
257 IF (avg(ng)%load.gt.avg(ng)%Nfiles) THEN
258 IF (master) THEN
259 WRITE (stdout,10) 'AVG(ng)%load = ', avg(ng)%load, &
260 & avg(ng)%Nfiles, trim(avg(ng)%base), &
261 & ifile
262 END IF
263 exit_flag=4
265 & __line__, myfile)) RETURN
266 END IF
267 fcount=avg(ng)%load
268 avg(ng)%Nrec(fcount)=0
269 IF (master) THEN
270 WRITE (avg(ng)%name,20) trim(avg(ng)%base), ifile
271 END IF
272# ifdef DISTRIBUTE
273 CALL mp_bcasts (ng, itlm, avg(ng)%name)
274# endif
275 avg(ng)%files(fcount)=trim(avg(ng)%name)
276 CALL close_file (ng, itlm, avg(ng), avg(ng)%name, lupdate)
277 CALL def_avg (ng, newfile)
278 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
279 lwrtavg(ng)=.true.
280 END IF
281 ELSE
282 IF (iic(ng).eq.ntstart(ng)) THEN
283 CALL def_avg (ng, ldefout(ng))
284 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
285 lwrtavg(ng)=.true.
286 ldefavg(ng)=.false.
287 END IF
288 END IF
289 END IF
290!
291! Write out data into time-averaged NetCDF file.
292!
293 IF (lwrtavg(ng)) THEN
294 IF ((iic(ng).gt.ntstart(ng)).and. &
295 & (mod(iic(ng)-1,navg(ng)).eq.0)) THEN
296# ifdef DISTRIBUTE
297 CALL wrt_avg (ng, myrank)
298# else
299 CALL wrt_avg (ng, -1)
300# endif
301 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
302 END IF
303 END IF
304# endif
305
306# if defined FOUR_DVAR
307# ifdef OBSERVATIONS
308!
309!-----------------------------------------------------------------------
310! If appropriate, process and write model state at observation
311! locations. Compute misfit (model-observations) cost function.
312!-----------------------------------------------------------------------
313!
314 IF (((time(ng)-0.5_r8*dt(ng)).le.obstime(ng)).and. &
315 & (obstime(ng).lt.(time(ng)+0.5_r8*dt(ng)))) THEN
316 processobs(ng)=.true.
317# ifdef DISTRIBUTE
318 tile=myrank
319# else
320 tile=-1
321# endif
322 CALL obs_read (ng, itlm, .false.)
323# ifdef SP4DVAR
324!
325! Do not assimilate obs on the last time step unless inter=Nsaddle.
326!
327 IF (iic(ng).ne.ntend(ng)+1.or.lsadd(ng)) THEN
328 CALL obs_write (ng, tile, itlm)
329 END IF
330# else
331 CALL obs_write (ng, tile, itlm)
332# endif
333# ifndef WEAK_CONSTRAINT
334 CALL obs_cost (ng, itlm)
335# endif
336 ELSE
337 processobs(ng)=.false.
338 END IF
339# endif
340# endif
341# ifdef PROFILE
342!
343!-----------------------------------------------------------------------
344! Turn off output data time wall clock.
345!-----------------------------------------------------------------------
346!
347 CALL wclock_off (ng, itlm, 8, __line__, myfile)
348# endif
349!
350 10 FORMAT (/,' TL_OUTPUT - multi-file counter ',a,i0, &
351 & ', is greater than Nfiles = ',i0,1x,'dimension', &
352 & /,13x,'in structure when creating next file: ', &
353 & a,'_',i4.4,'.nc', &
354 & /,13x,'Incorrect OutFiles logic in ''read_phypar''.')
355 20 FORMAT (a,'_',i4.4,'.nc')
356!
357 RETURN
subroutine, public close_file(ng, model, s, ncname, lupdate)
Definition close_io.F:43
subroutine, public def_avg(ng, ldef)
Definition def_avg.F:79
logical, dimension(:), allocatable lsadd
logical, dimension(:), allocatable processobs
type(t_io), dimension(:), allocatable tlm
type(t_io), dimension(:), allocatable avg
integer stdout
character(len=256) sourcefile
integer, dimension(:), allocatable idefavg
integer, dimension(:), allocatable ideftlm
logical lanafile
logical master
integer, parameter itlm
Definition mod_param.F:663
integer, dimension(:), allocatable nrrec
real(dp), dimension(:), allocatable obstime
integer, dimension(:), allocatable iic
integer, dimension(:), allocatable ntlm
real(dp), dimension(:), allocatable dt
integer, dimension(:), allocatable ndeftlm
logical, dimension(:), allocatable ldefavg
integer, dimension(:), allocatable navg
logical, dimension(:), allocatable lwrtavg
integer, dimension(:), allocatable ntend
integer exit_flag
integer, dimension(:), allocatable ndefavg
logical, dimension(:), allocatable lwrtper
logical, dimension(:), allocatable lwrttlm
logical, dimension(:), allocatable ldefout
real(dp), dimension(:), allocatable time
logical, dimension(:), allocatable ldeftlm
integer, dimension(:), allocatable ntstart
integer noerror
subroutine, public obs_read(ng, model, backward)
Definition obs_read.F:42
subroutine, public obs_write(ng, tile, model)
Definition obs_write.F:56
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine, public tl_def_his(ng, ldef)
Definition tl_def_his.F:51
subroutine, public tl_wrt_his(ng, tile)
Definition tl_wrt_his.F:68
subroutine, public wrt_avg(ng, tile)
Definition wrt_avg.F:83
subroutine obs_cost(ng, model)
Definition obs_cost.F:4
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::avg, close_io_mod::close_file(), def_avg_mod::def_avg(), mod_scalars::dt, mod_scalars::exit_flag, strings_mod::founderror(), mod_ncparam::idefavg, mod_ncparam::ideftlm, mod_scalars::iic, mod_param::itlm, mod_ncparam::lanafile, mod_scalars::ldefavg, mod_scalars::ldefout, mod_scalars::ldeftlm, mod_fourdvar::lsadd, mod_scalars::lwrtavg, mod_scalars::lwrtper, mod_scalars::lwrttlm, mod_parallel::master, mod_parallel::myrank, mod_scalars::navg, mod_scalars::ndefavg, mod_scalars::ndeftlm, mod_scalars::noerror, mod_scalars::nrrec, mod_scalars::ntend, mod_scalars::ntlm, mod_scalars::ntstart, obs_cost(), obs_read_mod::obs_read(), obs_write_mod::obs_write(), mod_scalars::obstime, mod_fourdvar::processobs, mod_iounits::sourcefile, mod_iounits::stdout, mod_scalars::time, tl_def_his_mod::tl_def_his(), tl_wrt_his_mod::tl_wrt_his(), mod_iounits::tlm, wclock_off(), wclock_on(), and wrt_avg_mod::wrt_avg().

Referenced by tl_main3d().

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