ROMS
Loading...
Searching...
No Matches
wrt_tides.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined AVERAGES && defined AVERAGES_DETIDE && \
4 (defined ssh_tides || defined uv_tides)
5!
6!git $Id$
7!================================================== Hernan G. Arango ===
8! Copyright (c) 2002-2025 The ROMS Group !
9! Licensed under a MIT/X style license !
10! See License_ROMS.md !
11!=======================================================================
12! !
13! This module writes time-accumulated tide harmonic fields used for !
14! detiding into harmonics file using the standard NetCDF library or !
15! the Parallel-IO (PIO) library. !
16! !
17!=======================================================================
18!
19 USE mod_param
20 USE mod_parallel
21 USE mod_grid
22 USE mod_iounits
23 USE mod_ncparam
24 USE mod_scalars
25 USE mod_stepping
26 USE mod_tides
27!
29# ifdef SOLVE3D
31# endif
32 USE strings_mod, ONLY : founderror
33!
34 implicit none
35!
36 PUBLIC :: wrt_tides
37 PRIVATE :: wrt_tides_nf90
38# if defined PIO_LIB && defined DISTRIBUTE
39 PRIVATE :: wrt_tides_pio
40# endif
41!
42 CONTAINS
43!
44!***********************************************************************
45 SUBROUTINE wrt_tides (ng, tile)
46!***********************************************************************
47!
48! Imported variable declarations.
49!
50 integer, intent(in) :: ng, tile
51!
52! Local variable declarations.
53!
54 integer :: lbi, ubi, lbj, ubj
55!
56 character (len=*), parameter :: myfile = &
57 & __FILE__
58!
59!-----------------------------------------------------------------------
60! Write out time-averaged fields according to IO type.
61!-----------------------------------------------------------------------
62!
63 lbi=bounds(ng)%LBi(tile)
64 ubi=bounds(ng)%UBi(tile)
65 lbj=bounds(ng)%LBj(tile)
66 ubj=bounds(ng)%UBj(tile)
67!
68 SELECT CASE (har(ng)%IOtype)
69 CASE (io_nf90)
70 CALL wrt_tides_nf90 (ng, tile, &
71 & lbi, ubi, lbj, ubj)
72
73# if defined PIO_LIB && defined DISTRIBUTE
74 CASE (io_pio)
75 CALL wrt_tides_pio (ng, tile, &
76 & lbi, ubi, lbj, ubj)
77# endif
78 CASE DEFAULT
79 IF (master) THEN
80 WRITE (stdout,10) har(ng)%IOtype
81 END IF
82 exit_flag=3
83 END SELECT
84 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
85!
86 10 FORMAT (' WRT_TIDES - Illegal output file type, io_type = ',i0, &
87 & /,13x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.' )
88!
89 RETURN
90 END SUBROUTINE wrt_tides
91!
92!***********************************************************************
93 SUBROUTINE wrt_tides_nf90 (ng, tile, &
94 & LBi, UBi, LBj, UBj)
95!***********************************************************************
96!
97 USE mod_netcdf
98!
99! Imported variable declarations.
100!
101 integer, intent(in) :: ng, tile
102!
103! Local variable declarations.
104!
105 integer :: lbi, ubi, lbj, ubj
106 integer :: gtype, itide, itrc, status, varid
107!
108 real(dp) :: scale
109 real(r8) :: work(ntc(ng))
110!
111 character (len=*), parameter :: myfile = &
112 & __FILE__//", wrt_tides_nf90"
113!
114 sourcefile=myfile
115!
116!-----------------------------------------------------------------------
117! Write out time-accumulated harmonic fields.
118!-----------------------------------------------------------------------
119!
120 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
121!
122! Report.
123!
124 IF (master) WRITE (stdout,20) ng
125!
126! Write out number of time-accumulated harmonics.
127!
128 CALL netcdf_put_ivar (ng, inlm, har(ng)%name, 'Hcount', &
129 & hcount(ng), (/0/), (/0/), &
130 & ncid = har(ng)%ncid)
131 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
132!
133! Write out model time for current time-accumulated harmonics.
134!
135 CALL netcdf_put_fvar (ng, inlm, har(ng)%name, &
136 & trim(vname(1,idtime)), time(ng), &
137 & (/0/), (/0/), &
138 & ncid = har(ng)%ncid, &
139 & varid = har(ng)%Vid(idtime))
140 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
141!
142! Write out tidal period (hours).
143!
144 DO itide=1,ntc(ng)
145 work(itide)=tides(ng)%Tperiod(itide)/3600.0_r8
146 END DO
147
148 CALL netcdf_put_fvar (ng, inlm, har(ng)%name, &
149 & trim(vname(1,idtper)), &
150 & work, &
151 & (/1/), (/ntc(ng)/), &
152 & ncid = har(ng)%ncid, &
153 & varid = har(ng)%Vid(idtper))
154 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
155!
156! Write out time-accumulated COS(omega(k)*t) harmonics.
157!
158 CALL netcdf_put_fvar (ng, inlm, har(ng)%name, &
159 & trim(vname(1,idcosw)), &
160 & tides(ng) % CosW_sum, &
161 & (/1/), (/ntc(ng)/), &
162 & ncid = har(ng)%ncid, &
163 & varid = har(ng)%Vid(idcosw))
164 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
165!
166! Write out time-accumulated SIN(omega(k)*t) harmonics.
167!
168 CALL netcdf_put_fvar (ng, inlm, har(ng)%name, &
169 & trim(vname(1,idsinw)), &
170 & tides(ng) % SinW_sum, &
171 & (/1/), (/ntc(ng)/), &
172 & ncid = har(ng)%ncid, &
173 & varid = har(ng)%Vid(idsinw))
174 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
175!
176! Write out time-accumulated COS(omega(k)*t)*COS(omega(l)*t) harmonics.
177!
178 CALL netcdf_put_fvar (ng, inlm, har(ng)%name, &
179 & trim(vname(1,idcos2)), &
180 & tides(ng) % CosWCosW, &
181 & (/1,1/), (/ntc(ng),ntc(ng)/), &
182 & ncid = har(ng)%ncid, &
183 & varid = har(ng)%Vid(idcos2))
184 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
185!
186! Write out time-accumulated SIN(omega(k)*t)*SIN(omega(l)*t) harmonics.
187!
188 CALL netcdf_put_fvar (ng, inlm, har(ng)%name, &
189 & trim(vname(1,idsin2)), &
190 & tides(ng) % SinWSinW, &
191 & (/1,1/), (/ntc(ng),ntc(ng)/), &
192 & ncid = har(ng)%ncid, &
193 & varid = har(ng)%Vid(idsin2))
194 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
195!
196! Write out time-accumulated SIN(omega(k)*t)*COS(omega(l)*t) harmonics.
197!
198 CALL netcdf_put_fvar (ng, inlm, har(ng)%name, &
199 & trim(vname(1,idswcw)), &
200 & tides(ng) % SinWCosW, &
201 & (/1,1/), (/ntc(ng),ntc(ng)/), &
202 & ncid = har(ng)%ncid, &
203 & varid = har(ng)%Vid(idswcw))
204 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
205!
206! Write out free-surface time-accumulated tide harmonics (m).
207!
208 IF (aout(idfsud,ng)) THEN
209 scale=1.0_dp
210 gtype=r3dvar
211 status=nf_fwrite3d(ng, inlm, har(ng)%ncid, idfsuh, &
212 & har(ng)%Vid(idfsuh), 0, gtype, &
213 & lbi, ubi, lbj, ubj, 0, 2*ntc(ng), scale, &
214# ifdef MASKING
215 & grid(ng) % rmask, &
216# endif
217 & tides(ng) % zeta_tide, &
218 & setfillval = .false.)
219 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
220 IF (master) THEN
221 WRITE (stdout,20) trim(vname(1,idfsuh)), trim(har(ng)%name)
222 END IF
223 exit_flag=3
224 ioerror=status
225 RETURN
226 END IF
227 END IF
228!
229! Write out 2D u-momentum time-accumulated tide harmonics (m/s).
230!
231 IF (aout(idu2dd,ng)) THEN
232 scale=1.0_dp
233 gtype=u3dvar
234 status=nf_fwrite3d(ng, inlm, har(ng)%ncid, idu2dh, &
235 & har(ng)%Vid(idu2dh), 0, gtype, &
236 & lbi, ubi, lbj, ubj, 0, 2*ntc(ng), scale, &
237# ifdef MASKING
238 & grid(ng) % umask, &
239# endif
240 & tides(ng) % ubar_tide, &
241 & setfillval = .false.)
242 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
243 IF (master) THEN
244 WRITE (stdout,20) trim(vname(1,idu2dh)), trim(har(ng)%name)
245 END IF
246 exit_flag=3
247 ioerror=status
248 RETURN
249 END IF
250 END IF
251!
252! Write out 2D v-momentum time-accumulated tide harmonics (m/s).
253!
254 IF (aout(idv2dd,ng)) THEN
255 scale=1.0_dp
256 gtype=v3dvar
257 status=nf_fwrite3d(ng, inlm, har(ng)%ncid, idv2dh, &
258 & har(ng)%Vid(idv2dh), 0, gtype, &
259 & lbi, ubi, lbj, ubj, 0, 2*ntc(ng), scale, &
260# ifdef MASKING
261 & grid(ng) % vmask, &
262# endif
263 & tides(ng) % vbar_tide, &
264 & setfillval = .false.)
265 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
266 IF (master) THEN
267 WRITE (stdout,20) trim(vname(1,idv2dh)), trim(har(ng)%name)
268 END IF
269 exit_flag=3
270 ioerror=status
271 RETURN
272 END IF
273 END IF
274
275# ifdef SOLVE3D
276!
277! Write out 3D u-momentum time-accumulated tide harmonics (m/s).
278!
279 IF (aout(idu3dd,ng)) THEN
280 scale=1.0_dp
281 gtype=u3dvar
282 status=nf_fwrite4d(ng, inlm, har(ng)%ncid, idu3dh, &
283 & har(ng)%Vid(idu3dh), 0, gtype, &
284 & lbi, ubi, lbj, ubj, 1, n(ng), 0, 2*ntc(ng), &
285 & scale, &
286# ifdef MASKING
287 & grid(ng) % umask, &
288# endif
289 & tides(ng) % u_tide, &
290 & setfillval = .false.)
291 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
292 IF (master) THEN
293 WRITE (stdout,20) trim(vname(1,idu3dh)), trim(har(ng)%name)
294 END IF
295 exit_flag=3
296 ioerror=status
297 RETURN
298 END IF
299 END IF
300!
301! Write out 3D v-momentum time-accumulated tide harmonics (m/s).
302!
303 IF (aout(idv3dd,ng)) THEN
304 scale=1.0_dp
305 gtype=v3dvar
306 status=nf_fwrite4d(ng, inlm, har(ng)%ncid, idv3dh, &
307 & har(ng)%Vid(idv3dh), 0, gtype, &
308 & lbi, ubi, lbj, ubj, 1, n(ng), 0, 2*ntc(ng), &
309 & scale, &
310# ifdef MASKING
311 & grid(ng) % vmask, &
312# endif
313 & tides(ng) % v_tide, &
314 & setfillval = .false.)
315 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
316 IF (master) THEN
317 WRITE (stdout,20) trim(vname(1,idv3dh)), trim(har(ng)%name)
318 END IF
319 exit_flag=3
320 ioerror=status
321 RETURN
322 END IF
323 END IF
324!
325! Write out temperature and salinity time-accumulated tide harmonics.
326!
327 DO itrc=1,nat
328 IF (aout(idtrcd(itrc),ng)) THEN
329 scale=1.0_dp
330 gtype=r3dvar
331 status=nf_fwrite4d(ng, inlm, har(ng)%ncid, idtrch(itrc), &
332 & har(ng)%Vid(idtrch(itrc)), 0, gtype, &
333 & lbi, ubi, lbj, ubj, 1, n(ng), 0, 2*ntc(ng),&
334 & scale, &
335# ifdef MASKING
336 & grid(ng) % vmask, &
337# endif
338 & tides(ng) % t_tide(:,:,:,:,itrc), &
339 & setfillval = .false.)
340 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
341 IF (master) THEN
342 WRITE (stdout,20) trim(vname(1,idtrch(itrc))), &
343 & trim(har(ng)%name)
344 END IF
345 exit_flag=3
346 ioerror=status
347 RETURN
348 END IF
349 END IF
350 END DO
351# endif
352!
353!-----------------------------------------------------------------------
354! Synchronize tide forcing NetCDF file to disk to allow other processes
355! to access data immediately after it is written.
356!-----------------------------------------------------------------------
357!
358 CALL netcdf_sync (ng, inlm, har(ng)%name, har(ng)%ncid)
359 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
360!
361 10 FORMAT (2x,'WRT_TIDES_NF90 - writing time-accumulated tide ', &
362 & 'harmonics, Grid ',i2.2)
363 20 FORMAT (/,' WRT_TIDES_NF90 - error while writing variable: ',a, &
364 & /,13x,'into detide harmonics NetCDF file: ',/,13x,a)
365!
366 RETURN
367 END SUBROUTINE wrt_tides_nf90
368
369# if defined PIO_LIB && defined DISTRIBUTE
370!
371!***********************************************************************
372 SUBROUTINE wrt_tides_pio (ng, tile, &
373 & LBi, UBi, LBj, UBj)
374!***********************************************************************
375!
377!
378! Imported variable declarations.
379!
380 integer, intent(in) :: ng, tile
381!
382! Local variable declarations.
383!
384 integer :: lbi, ubi, lbj, ubj
385 integer :: itide, itrc, status
386!
387 real(dp) :: scale
388 real(r8) :: work(ntc(ng))
389!
390 character (len=*), parameter :: myfile = &
391 & __FILE__//", wrt_tides_pio"
392!
393 sourcefile=myfile
394!
395!-----------------------------------------------------------------------
396! Write out time-accumulated harmonic fields.
397!-----------------------------------------------------------------------
398!
399 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
400!
401! Report.
402!
403 IF (master) WRITE (stdout,10) ng
404!
405! Write out number of time-accumulated harmonics.
406!
407 CALL pio_netcdf_put_ivar (ng, inlm, har(ng)%name, 'Hcount', &
408 & hcount(ng), (/0/), (/0/), &
409 & piofile = har(ng)%pioFile)
410 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
411!
412! Write out model time for current time-accumulated harmonics.
413!
414 CALL pio_netcdf_put_fvar (ng, inlm, har(ng)%name, &
415 & trim(vname(1,idtime)), time(ng), &
416 & (/0/), (/0/), &
417 & piofile = har(ng)%pioFile, &
418 & piovar = har(ng)%pioVar(idtime)%vd)
419 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
420!
421! Write out tidal period (hours).
422!
423 DO itide=1,ntc(ng)
424 work(itide)=tides(ng)%Tperiod(itide)/3600.0_r8
425 END DO
426
427 CALL pio_netcdf_put_fvar (ng, inlm, har(ng)%name, &
428 & trim(vname(1,idtper)), &
429 & work, &
430 & (/1/), (/ntc(ng)/), &
431 & piofile = har(ng)%pioFile, &
432 & piovar = har(ng)%pioVar(idtper)%vd)
433 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
434!
435! Write out time-accumulated COS(omega(k)*t) harmonics.
436!
437 CALL pio_netcdf_put_fvar (ng, inlm, har(ng)%name, &
438 & trim(vname(1,idcosw)), &
439 & tides(ng) % CosW_sum, &
440 & (/1/), (/ntc(ng)/), &
441 & piofile = har(ng)%pioFile, &
442 & piovar = har(ng)%pioVar(idcosw)%vd)
443 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
444!
445! Write out time-accumulated SIN(omega(k)*t) harmonics.
446!
447 CALL pio_netcdf_put_fvar (ng, inlm, har(ng)%name, &
448 & trim(vname(1,idsinw)), &
449 & tides(ng) % SinW_sum, &
450 & (/1/), (/ntc(ng)/), &
451 & piofile = har(ng)%pioFile, &
452 & piovar = har(ng)%pioVar(idsinw)%vd)
453 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
454!
455! Write out time-accumulated COS(omega(k)*t)*COS(omega(l)*t) harmonics.
456!
457 CALL pio_netcdf_put_fvar (ng, inlm, har(ng)%name, &
458 & trim(vname(1,idcos2)), &
459 & tides(ng) % CosWCosW, &
460 & (/1,1/), (/ntc(ng),ntc(ng)/), &
461 & piofile = har(ng)%pioFile, &
462 & piovar = har(ng)%pioVar(idcos2)%vd)
463 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
464!
465! Write out time-accumulated SIN(omega(k)*t)*SIN(omega(l)*t) harmonics.
466!
467 CALL pio_netcdf_put_fvar (ng, inlm, har(ng)%name, &
468 & trim(vname(1,idsin2)), &
469 & tides(ng) % SinWSinW, &
470 & (/1,1/), (/ntc(ng),ntc(ng)/), &
471 & piofile = har(ng)%pioFile, &
472 & piovar = har(ng)%pioVar(idsin2)%vd)
473 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
474!
475! Write out time-accumulated SIN(omega(k)*t)*COS(omega(l)*t) harmonics.
476!
477 CALL pio_netcdf_put_fvar (ng, inlm, har(ng)%name, &
478 & trim(vname(1,idswcw)), &
479 & tides(ng) % SinWCosW, &
480 & (/1,1/), (/ntc(ng),ntc(ng)/), &
481 & piofile = har(ng)%pioFile, &
482 & piovar = har(ng)%pioVar(idswcw)%vd)
483 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
484!
485! Write out free-surface time-accumulated tide harmonics (m).
486!
487 IF (aout(idfsud,ng)) THEN
488 scale=1.0_dp
489 IF (har(ng)%pioVar(idfsuh)%dkind.eq.pio_double) THEN
490 iodesc => iodesc_dp_r2dhar(ng)
491 ELSE
492 iodesc => iodesc_sp_r2dhar(ng)
493 END IF
494!
495 status=nf_fwrite3d(ng, inlm, har(ng)%pioFile, idfsuh, &
496 & har(ng)%pioVar(idfsuh), 0, &
497 & iodesc, &
498 & lbi, ubi, lbj, ubj, 0, 2*ntc(ng), scale, &
499# ifdef MASKING
500 & grid(ng) % rmask, &
501# endif
502 & tides(ng) % zeta_tide, &
503 & setfillval = .false.)
504 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
505 IF (master) THEN
506 WRITE (stdout,20) trim(vname(1,idfsuh)), trim(har(ng)%name)
507 END IF
508 exit_flag=3
509 ioerror=status
510 RETURN
511 END IF
512 END IF
513!
514! Write out 2D u-momentum time-accumulated tide harmonics (m/s).
515!
516 IF (aout(idu2dd,ng)) THEN
517 scale=1.0_dp
518 IF (har(ng)%pioVar(idu2dh)%dkind.eq.pio_double) THEN
519 iodesc => iodesc_dp_u2dhar(ng)
520 ELSE
521 iodesc => iodesc_sp_u2dhar(ng)
522 END IF
523!
524 status=nf_fwrite3d(ng, inlm, har(ng)%pioFile, idu2dh, &
525 & har(ng)%pioVar(idu2dh), 0, &
526 & iodesc, &
527 & lbi, ubi, lbj, ubj, 0, 2*ntc(ng), scale, &
528# ifdef MASKING
529 & grid(ng) % umask, &
530# endif
531 & tides(ng) % ubar_tide, &
532 & setfillval = .false.)
533 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
534 IF (master) THEN
535 WRITE (stdout,20) trim(vname(1,idu2dh)), trim(har(ng)%name)
536 END IF
537 exit_flag=3
538 ioerror=status
539 RETURN
540 END IF
541 END IF
542!
543! Write out 2D v-momentum time-accumulated tide harmonics (m/s).
544!
545 IF (aout(idv2dd,ng)) THEN
546 scale=1.0_dp
547 IF (har(ng)%pioVar(idv2dh)%dkind.eq.pio_double) THEN
548 iodesc => iodesc_dp_v2dhar(ng)
549 ELSE
550 iodesc => iodesc_sp_v2dhar(ng)
551 END IF
552!
553 status=nf_fwrite3d(ng, inlm, har(ng)%pioFile, idv2dh, &
554 & har(ng)%pioVar(idv2dh), 0, &
555 & iodesc, &
556 & lbi, ubi, lbj, ubj, 0, 2*ntc(ng), scale, &
557# ifdef MASKING
558 & grid(ng) % vmask, &
559# endif
560 & tides(ng) % vbar_tide, &
561 & setfillval = .false.)
562 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
563 IF (master) THEN
564 WRITE (stdout,20) trim(vname(1,idv2dh)), trim(har(ng)%name)
565 END IF
566 exit_flag=3
567 ioerror=status
568 RETURN
569 END IF
570 END IF
571
572# ifdef SOLVE3D
573!
574! Write out 3D u-momentum time-accumulated tide harmonics (m/s).
575!
576 IF (aout(idu3dd,ng)) THEN
577 scale=1.0_dp
578 IF (har(ng)%pioVar(idu3dh)%dkind.eq.pio_double) THEN
579 iodesc => iodesc_dp_u3dhar(ng)
580 ELSE
581 iodesc => iodesc_sp_u3dhar(ng)
582 END IF
583!
584 status=nf_fwrite4d(ng, inlm, har(ng)%pioFile, idu3dh, &
585 & har(ng)%pioVar(idu3dh), 0, &
586 & iodesc, &
587 & lbi, ubi, lbj, ubj, 1, n(ng), 0, 2*ntc(ng), &
588 & scale, &
589# ifdef MASKING
590 & grid(ng) % umask, &
591# endif
592 & tides(ng) % u_tide, &
593 & setfillval = .false.)
594 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
595 IF (master) THEN
596 WRITE (stdout,20) trim(vname(1,idu3dh)), trim(har(ng)%name)
597 END IF
598 exit_flag=3
599 ioerror=status
600 RETURN
601 END IF
602 END IF
603!
604! Write out 3D v-momentum time-accumulated tide harmonics (m/s).
605!
606 IF (aout(idv3dd,ng)) THEN
607 scale=1.0_dp
608 IF (har(ng)%pioVar(idfsuh)%dkind.eq.pio_double) THEN
609 iodesc => iodesc_dp_v3dhar(ng)
610 ELSE
611 iodesc => iodesc_sp_v3dhar(ng)
612 END IF
613!
614 status=nf_fwrite4d(ng, inlm, har(ng)%pioFile, idv3dh, &
615 & har(ng)%pioVar(idv3dh), 0, &
616 & iodesc, &
617 & lbi, ubi, lbj, ubj, 1, n(ng), 0, 2*ntc(ng), &
618 & scale, &
619# ifdef MASKING
620 & grid(ng) % vmask, &
621# endif
622 & tides(ng) % v_tide, &
623 & setfillval = .false.)
624 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
625 IF (master) THEN
626 WRITE (stdout,20) trim(vname(1,idv3dh)), trim(har(ng)%name)
627 END IF
628 exit_flag=3
629 ioerror=status
630 RETURN
631 END IF
632 END IF
633!
634! Write out temperature and salinity time-accumulated tide harmonics.
635!
636 DO itrc=1,nat
637 IF (aout(idtrcd(itrc),ng)) THEN
638 scale=1.0_dp
639 IF (har(ng)%pioVar(idtrch)%dkind.eq.pio_double) THEN
640 iodesc => iodesc_dp_r3dhar(ng)
641 ELSE
642 iodesc => iodesc_sp_r3dhar(ng)
643 END IF
644!
645 status=nf_fwrite4d(ng, inlm, har(ng)%pioFile, idtrch(itrc), &
646 & har(ng)%pioVar(idtrch(itrc)), 0, &
647 & iodesc, &
648 & lbi, ubi, lbj, ubj, 1, n(ng), 0, 2*ntc(ng),&
649 & scale, &
650# ifdef MASKING
651 & grid(ng) % vmask, &
652# endif
653 & tides(ng) % t_tide(:,:,:,:,itrc), &
654 & setfillval = .false.)
655 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
656 IF (master) THEN
657 WRITE (stdout,20) trim(vname(1,idtrch(itrc))), &
658 & trim(har(ng)%name)
659 END IF
660 exit_flag=3
661 ioerror=status
662 RETURN
663 END IF
664 END IF
665 END DO
666# endif
667!
668!-----------------------------------------------------------------------
669! Synchronize tide forcing NetCDF file to disk to allow other processes
670! to access data immediately after it is written.
671!-----------------------------------------------------------------------
672!
673 CALL pio_netcdf_sync (ng, inlm, har(ng)%name, har(ng)%pioFile)
674 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
675!
676 10 FORMAT (2x,'WRT_TIDES_PIO - writing time-accumulated tide ', &
677 & 'harmonics, Grid ',i2.2)
678 20 FORMAT (/,' WRT_TIDES_PIO - error while writing variable: ',a, &
679 & /,13x,'into detide harmonics NetCDF file: ',/,13x,a)
680!
681 RETURN
682 END SUBROUTINE wrt_tides_pio
683# endif
684#endif
685 END MODULE wrt_tides_mod
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer ioerror
type(t_io), dimension(:), allocatable har
integer stdout
character(len=256) sourcefile
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer idu2dh
integer idu3dd
integer idv3dd
integer idsin2
integer idtper
integer, parameter io_pio
Definition mod_ncparam.F:96
integer idu3dh
integer idv3dh
integer idcos2
integer idfsud
integer idswcw
integer, dimension(:), allocatable idtrcd
integer, dimension(:), allocatable idtrch
integer idcosw
integer idv2dd
character(len=maxlen), dimension(6, 0:nv) vname
integer idtime
integer idv2dh
integer idu2dd
integer idfsuh
integer idsinw
logical, dimension(:,:), allocatable aout
subroutine, public netcdf_sync(ng, model, ncname, ncid)
logical master
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
integer, parameter r3dvar
Definition mod_param.F:721
integer, parameter u3dvar
Definition mod_param.F:722
integer, parameter v3dvar
Definition mod_param.F:723
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dhar
subroutine, public pio_netcdf_sync(ng, model, ncname, piofile)
type(io_desc_t), dimension(:), pointer iodesc_dp_u3dhar
type(io_desc_t), dimension(:), pointer iodesc_dp_v3dhar
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dhar
type(io_desc_t), dimension(:), pointer iodesc_dp_v2dhar
type(io_desc_t), dimension(:), pointer iodesc_sp_r3dhar
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dhar
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dhar
type(io_desc_t), dimension(:), pointer iodesc_sp_v3dhar
type(io_desc_t), dimension(:), pointer iodesc_sp_u3dhar
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dhar
type(io_desc_t), dimension(:), pointer iodesc_dp_r3dhar
integer, dimension(:), allocatable hcount
integer exit_flag
real(dp), dimension(:), allocatable time
integer noerror
integer, dimension(:), allocatable ntc
type(t_tides), dimension(:), allocatable tides
Definition mod_tides.F:133
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine, private wrt_tides_pio(ng, tile, lbi, ubi, lbj, ubj)
Definition wrt_tides.F:374
subroutine, private wrt_tides_nf90(ng, tile, lbi, ubi, lbj, ubj)
Definition wrt_tides.F:95
subroutine, public wrt_tides(ng, tile)
Definition wrt_tides.F:46