ROMS
Loading...
Searching...
No Matches
sediment_output.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#if (defined SEDIMENT || defined BBL_MODEL) && defined SOLVE3D
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 defines/writes Waves Effect on Currents variables into !
14! output NetCDF files. !
15! !
16!=======================================================================
17!
18 USE mod_param
19 USE mod_parallel
20# ifdef AVERAGES
21 USE mod_average
22# endif
23# ifdef BBL_MODEL
24 USE mod_bbl
25# endif
26 USE mod_forces
27 USE mod_grid
28 USE mod_iounits
29 USE mod_mixing
30 USE mod_ncparam
31 USE mod_ocean
32 USE mod_scalars
33 USE mod_sedbed
34# ifdef SEDIMENT
35 USE mod_sediment
36# endif
37 USE mod_stepping
38!
39 USE def_var_mod, ONLY : def_var
40# ifdef STATIONS
42# ifdef SOLVE3D
44# endif
45# endif
47# ifdef SOLVE3D
49 USE omega_mod, ONLY : scale_omega
50# endif
51 USE strings_mod, ONLY : founderror
52!
53 implicit none
54!
55 PUBLIC :: sediment_def_nf90
56# if defined PIO_LIB && defined DISTRIBUTE
57 PUBLIC :: sediment_def_pio
58# endif
59# ifdef STATIONS
61# if defined PIO_LIB && defined DISTRIBUTE
63# endif
64# endif
65 PUBLIC :: sediment_wrt_nf90
66# if defined PIO_LIB && defined DISTRIBUTE
67 PUBLIC :: sediment_wrt_pio
68# endif
69# ifdef STATIONS
71# if defined PIO_LIB && defined DISTRIBUTE
73# endif
74# endif
75!
76 CONTAINS
77!
78!***********************************************************************
79 SUBROUTINE sediment_def_nf90 (ng, model, ldef, VarOut, S, &
80 & t2dgrd, u2dgrd, v2dgrd, &
81 & b3dgrd)
82!***********************************************************************
83!
84 USE mod_netcdf
85!
86! Imported variable declarations.
87!
88 logical, intent(in) :: ldef, varout(nv,ngrids)
89!
90 integer, intent(in) :: ng, model
91 integer, intent(in), optional :: t2dgrd(:), u2dgrd(:), v2dgrd(:)
92 integer, intent(in), optional :: b3dgrd(:)
93!
94 TYPE(t_io), intent(inout) :: s(ngrids)
95!
96! Local variable declarations.
97!
98 logical :: got_var(nv)
99!
100 integer, parameter :: natt = 25
101
102 integer :: i, itrc, j, nvd3, nvd4, status
103!
104 real(r8) :: aval(6)
105!
106# ifdef ADJOINT
107 character (len=21) :: prefix
108# else
109 character (len=13) :: prefix
110# endif
111 character (len=120) :: vinfo(natt)
112 character (len=256) :: ncname
113!
114 character (len=*), parameter :: myfile = &
115 & __FILE__//", sediment_def_nf90"
116!
117 sourcefile=myfile
118!
119!-----------------------------------------------------------------------
120! Define sediment output variables.
121!-----------------------------------------------------------------------
122!
123 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
124 ncname=s(ng)%name
125!
126 define : IF (ldef) THEN
127!
128! Set number of dimensions for output variables.
129!
130# if defined WRITE_WATER && defined MASKING
131 nvd3=2
132 nvd4=2
133# else
134 nvd3=3
135 nvd4=4
136# endif
137!
138! Set long name prefix string.
139!
140# ifdef ADJOINT
141!! Prefix='time-averaged adjoint'
142 prefix='adjoint'
143# else
144!! Prefix='time-averaged'
145 prefix=char(32) ! blank space
146# endif
147!
148! Initialize local information variable arrays.
149!
150 DO i=1,natt
151 DO j=1,len(vinfo(1))
152 vinfo(i)(j:j)=' '
153 END DO
154 END DO
155 DO i=1,6
156 aval(i)=0.0_r8
157 END DO
158
159# if defined SEDIMENT && defined SED_MORPH
160!
161! Define time-varying bathymetry.
162!
163 IF (varout(idbath,ng)) THEN
164 vinfo( 1)=vname(1,idbath)
165 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
166 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idbath))
167 ELSE
168 vinfo( 3)=vname(3,idbath)
169 END IF
170 vinfo(14)=vname(4,idbath)
171 vinfo(16)=vname(1,idtime)
172 vinfo(21)=vname(6,idbath)
173 vinfo(22)='coordinates'
174 aval(5)=real(iinfo(1,idbath,ng),r8)
175 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idbath), &
176 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname, &
177 & setfillval = .false.)
178 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
179 END IF
180# endif
181
182# if defined SEDIMENT && defined BEDLOAD
183!
184! Define Bedload transport U-direction.
185!
186 DO i=1,nst
187 IF (varout(idubld(i),ng)) THEN
188 vinfo( 1)=vname(1,idubld(i))
189 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
190 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
191 & trim(vname(2,idubld(i)))
192 ELSE
193 vinfo( 2)=vname(2,idubld(i))
194 END IF
195 vinfo( 3)=vname(3,idubld(i))
196 vinfo(14)=vname(4,idubld(i))
197 vinfo(16)=vname(1,idtime)
198# if defined WRITE_WATER && defined MASKING
199 vinfo(20)='mask_u'
200# endif
201 vinfo(21)=vname(6,idubld(i))
202 vinfo(22)='coordinates'
203 aval(5)=real(iinfo(1,idubld(i),ng),r8)
204 status=def_var(ng, model, s(ng)%ncid, &
205 & s(ng)%Vid(idubld(i)), nf_fout, &
206 & nvd3, u2dgrd, aval, vinfo, ncname)
207 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
208 END IF
209!
210! Define Bedload transport V-direction.
211!
212 IF (varout(idvbld(i),ng)) THEN
213 vinfo( 1)=vname(1,idvbld(i))
214 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
215 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
216 & trim(vname(2,idvbld(i)))
217 ELSE
218 vinfo( 2)=vname(2,idvbld(i))
219 END IF
220 vinfo( 3)=vname(3,idvbld(i))
221 vinfo(14)=vname(4,idvbld(i))
222 vinfo(16)=vname(1,idtime)
223# if defined WRITE_WATER && defined MASKING
224 vinfo(20)='mask_v'
225# endif
226 vinfo(21)=vname(6,idvbld(i))
227 vinfo(22)='coordinates'
228 aval(5)=real(iinfo(1,idvbld(i),ng),r8)
229 status=def_var(ng, model, s(ng)%ncid, &
230 & s(ng)%Vid(idvbld(i)), nf_fout, &
231 & nvd3, v2dgrd, aval, vinfo, ncname)
232 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
233 END IF
234 END DO
235# endif
236
237# ifdef SEDIMENT
238!
239! Define sediment fraction of each size class in each bed layer.
240!
241 DO i=1,nst
242 IF (varout(idfrac(i),ng)) THEN
243 vinfo( 1)=vname(1,idfrac(i))
244 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
245 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
246 & trim(vname(2,idfrac(i)))
247 ELSE
248 vinfo( 2)=vname(2,idfrac(i))
249 END IF
250 vinfo( 3)=vname(3,idfrac(i))
251 vinfo(14)=vname(4,idfrac(i))
252 vinfo(16)=vname(1,idtime)
253 WRITE (vinfo(19),10) 1000.0_r8*sd50(i,ng)
254# if defined WRITE_WATER && defined MASKING
255 vinfo(20)='mask_rho'
256# endif
257 vinfo(21)=vname(6,idfrac(i))
258 vinfo(22)='coordinates'
259 aval(5)=real(iinfo(1,idfrac(i),ng))
260 status=def_var(ng, model, s(ng)%ncid, &
261 & s(ng)%Vid(idfrac(i)), nf_fout, &
262 & nvd4, b3dgrd, aval, vinfo, ncname)
263 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
264 END IF
265 END DO
266!
267! Define sediment mass of each size class in each bed layer.
268!
269 DO i=1,nst
270 IF (varout(idbmas(i),ng)) THEN
271 vinfo( 1)=vname(1,idbmas(i))
272 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
273 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
274 & trim(vname(2,idbmas(i)))
275 ELSE
276 vinfo( 2)=vname(2,idbmas(i))
277 END IF
278 vinfo( 3)=vname(3,idbmas(i))
279 vinfo(14)=vname(4,idbmas(i))
280 vinfo(16)=vname(1,idtime)
281 WRITE (vinfo(19),10) 1000.0_r8*sd50(i,ng)
282# if defined WRITE_WATER && defined MASKING
283 vinfo(20)='mask_rho'
284# endif
285 vinfo(21)=vname(6,idbmas(i))
286 vinfo(22)='coordinates'
287 aval(5)=real(iinfo(1,idbmas(i),ng),r8)
288 status=def_var(ng, model, s(ng)%ncid, &
289 & s(ng)%Vid(idbmas(i)), nf_fout, &
290 & nvd4, b3dgrd, aval, vinfo, ncname)
291 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
292 END IF
293 END DO
294!
295! Define sediment properties in each bed layer.
296!
297 DO i=1,mbedp
298 IF (varout(idsbed(i),ng)) THEN
299 vinfo( 1)=vname(1,idsbed(i))
300 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
301 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
302 & trim(vname(2,idsbed(i)))
303 ELSE
304 vinfo( 2)=vname(2,idsbed(i))
305 END IF
306 vinfo( 3)=vname(3,idsbed(i))
307 vinfo(14)=vname(4,idsbed(i))
308 vinfo(16)=vname(1,idtime)
309# if defined WRITE_WATER && defined MASKING
310 vinfo(20)='mask_rho'
311# endif
312 vinfo(21)=vname(6,idsbed(i))
313 vinfo(22)='coordinates'
314 aval(5)=real(iinfo(1,idsbed(i),ng),r8)
315 status=def_var(ng, model, s(ng)%ncid, &
316 & s(ng)%Vid(idsbed(i)), nf_fout, &
317 & nvd4, b3dgrd, aval, vinfo, ncname)
318 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
319 END IF
320 END DO
321# endif
322
323# if defined SEDIMENT || defined BBL_MODEL
324!
325! Define exposed sediment layer properties.
326!
327 DO i=1,mbotp
328 IF (varout(idbott(i),ng)) THEN
329 vinfo( 1)=vname(1,idbott(i))
330 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
331 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
332 & trim(vname(2,idbott(i)))
333 ELSE
334 vinfo( 2)=vname(2,idbott(i))
335 END IF
336 vinfo( 3)=vname(3,idbott(i))
337 vinfo(14)=vname(4,idbott(i))
338 vinfo(16)=vname(1,idtime)
339# if defined WRITE_WATER && defined MASKING
340 vinfo(20)='mask_rho'
341# endif
342 vinfo(21)=vname(6,idbott(i))
343 vinfo(22)='coordinates'
344 aval(5)=real(iinfo(1,idbott(i),ng),r8)
345 status=def_var(ng, model, s(ng)%ncid, &
346 & s(ng)%Vid(idbott(i)), nf_fout, &
347 & nvd3, t2dgrd, aval, vinfo, ncname)
348 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
349 END IF
350 END DO
351# endif
352
353# if defined SEDIMENT && defined BEDLOAD && defined BEDLOAD_VANDERA
354!
355! Define Ursell number of the asymmetric wave form.
356!
357 IF (varout(idsurs,ng)) THEN
358 vinfo( 1)=vname(1,idsurs)
359 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
360 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsurs))
361 ELSE
362 vinfo( 2)=vname(2,idsurs)
363 END IF
364 vinfo( 3)=vname(3,idsurs)
365 vinfo(14)=vname(4,idsurs)
366 vinfo(16)=vname(1,idsurs)
367# if defined WRITE_WATER && defined MASKING
368 vinfo(20)='mask_rho'
369# endif
370 vinfo(21)=vname(6,idsurs)
371 vinfo(22)='coordinates'
372 aval(5)=real(iinfo(1,idsurs,ng),r8)
373 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idsurs), &
374 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
375 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
376 END IF
377!
378! Define velocity skewness parameter of the asymmetric wave form.
379!
380 IF (varout(idsrrw,ng)) THEN
381 vinfo( 1)=vname(1,idsrrw)
382 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
383 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsrrw))
384 ELSE
385 vinfo( 2)=vname(2,idsrrw)
386 END IF
387 vinfo( 3)=vname(3,idsrrw)
388 vinfo(14)=vname(4,idsrrw)
389 vinfo(16)=vname(1,idsrrw)
390# if defined WRITE_WATER && defined MASKING
391 vinfo(20)='mask_rho'
392# endif
393 vinfo(21)=vname(6,idsrrw)
394 vinfo(22)='coordinates'
395 aval(5)=real(iinfo(1,idsrrw,ng),r8)
396 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idsrrw), &
397 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
398 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
399 END IF
400!
401! Define acceleration asymmetry parameter of the asymmetric wave form.
402!
403 IF (varout(idsbtw,ng)) THEN
404 vinfo( 1)=vname(1,idsbtw)
405 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
406 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsbtw))
407 ELSE
408 vinfo( 2)=vname(2,idsbtw)
409 END IF
410 vinfo( 3)=vname(3,idsbtw)
411 vinfo(14)=vname(4,idsbtw)
412 vinfo(16)=vname(1,idsbtw)
413# if defined WRITE_WATER && defined MASKING
414 vinfo(20)='mask_rho'
415# endif
416 vinfo(21)=vname(6,idsbtw)
417 vinfo(22)='coordinates'
418 aval(5)=real(iinfo(1,idsbtw,ng),r8)
419 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idsbtw), &
420 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
421 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
422 END IF
423!
424! Define crest velocity of the asymmetric wave form.
425!
426 IF (varout(idsucr,ng)) THEN
427 vinfo( 1)=vname(1,idsucr)
428 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
429 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsucr))
430 ELSE
431 vinfo( 2)=vname(2,idsucr)
432 END IF
433 vinfo( 3)=vname(3,idsucr)
434 vinfo(14)=vname(4,idsucr)
435 vinfo(16)=vname(1,idsucr)
436# if defined WRITE_WATER && defined MASKING
437 vinfo(20)='mask_rho'
438# endif
439 vinfo(21)=vname(6,idsucr)
440 vinfo(22)='coordinates'
441 aval(5)=real(iinfo(1,idsucr,ng),r8)
442 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idsucr), &
443 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
444 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
445 END IF
446!
447! Define trough velocity of the asymmetric wave form.
448!
449 IF (varout(idsutr,ng)) THEN
450 vinfo( 1)=vname(1,idsutr)
451 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
452 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsutr))
453 ELSE
454 vinfo( 2)=vname(2,idsutr)
455 END IF
456 vinfo( 3)=vname(3,idsutr)
457 vinfo(14)=vname(4,idsutr)
458 vinfo(16)=vname(1,idsutr)
459# if defined WRITE_WATER && defined MASKING
460 vinfo(20)='mask_rho'
461# endif
462 vinfo(21)=vname(6,idsutr)
463 vinfo(22)='coordinates'
464 aval(5)=real(iinfo(1,idsutr,ng),r8)
465 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idsutr), &
466 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
467 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
468 END IF
469!
470! Define crest time period of the asymmetric wave form.
471!
472 IF (varout(idstcr,ng)) THEN
473 vinfo( 1)=vname(1,idstcr)
474 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
475 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idstcr))
476 ELSE
477 vinfo( 2)=vname(2,idstcr)
478 END IF
479 vinfo( 3)=vname(3,idstcr)
480 vinfo(14)=vname(4,idstcr)
481 vinfo(16)=vname(1,idstcr)
482# if defined WRITE_WATER && defined MASKING
483 vinfo(20)='mask_rho'
484# endif
485 vinfo(21)=vname(6,idstcr)
486 vinfo(22)='coordinates'
487 aval(5)=real(iinfo(1,idstcr,ng),r8)
488 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idstcr), &
489 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
490 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
491 END IF
492!
493! Trough time period of the asymmetric wave form.
494!
495 IF (varout(idsttr,ng)) THEN
496 vinfo( 1)=vname(1,idsttr)
497 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
498 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsttr))
499 ELSE
500 vinfo( 2)=vname(2,idsttr)
501 END IF
502 vinfo( 3)=vname(3,idsttr)
503 vinfo(14)=vname(4,idsttr)
504 vinfo(16)=vname(1,idsttr)
505# if defined WRITE_WATER && defined MASKING
506 vinfo(20)='mask_rho'
507# endif
508 vinfo(21)=vname(6,idsttr)
509 vinfo(22)='coordinates'
510 aval(5)=real(iinfo(1,idsttr,ng),r8)
511 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idsttr), &
512 & nf_fout, nvd3, t2dgrd, aval, vinfo, ncname)
513 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
514 END IF
515# endif
516
517 END IF define
518!
519!-----------------------------------------------------------------------
520! Otherwise, check existing output file and prepare for appending
521! data.
522!-----------------------------------------------------------------------
523!
524 query : IF (.not.ldef) THEN
525!
526! Initialize local logical switches.
527!
528 DO i=1,nv
529 got_var(i)=.false.
530 END DO
531!
532! Scan variable list from input NetCDF and activate switches for
533! Waves Effect on Currents variables. Get variable IDs.
534!
535 DO i=1,n_var
536 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
537 got_var(idtime)=.true.
538 s(ng)%Vid(idtime)=var_id(i)
539# if defined SEDIMENT && defined SED_MORPH
540 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idbath))) THEN
541 got_var(idbath)=.true.
542 s(ng)%Vid(idbath)=var_id(i)
543# endif
544# if defined SEDIMENT && defined BEDLOAD && defined BEDLOAD_VANDERA
545 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsurs))) THEN
546 got_var(idsurs)=.true.
547 s(ng)%Vid(idsurs)=var_id(i)
548 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsrrw))) THEN
549 got_var(idsrrw)=.true.
550 s(ng)%Vid(idsrrw)=var_id(i)
551 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsbtw))) THEN
552 got_var(idsbtw)=.true.
553 s(ng)%Vid(idsbtw)=var_id(i)
554 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsucr))) THEN
555 got_var(idsucr)=.true.
556 s(ng)%Vid(idsucr)=var_id(i)
557 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsutr))) THEN
558 got_var(idsutr)=.true.
559 s(ng)%Vid(idsutr)=var_id(i)
560 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idstcr))) THEN
561 got_var(idstcr)=.true.
562 s(ng)%Vid(idstcr)=var_id(i)
563 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsttr))) THEN
564 got_var(idsttr)=.true.
565 s(ng)%Vid(idsttr)=var_id(i)
566# endif
567 END IF
568# ifdef SEDIMENT
569 DO itrc=1,nst
570 IF (trim(var_name(i)).eq. &
571 & trim(vname(1,idfrac(itrc)))) THEN
572 got_var(idfrac(itrc))=.true.
573 s(ng)%Vid(idfrac(itrc))=var_id(i)
574 ELSE IF (trim(var_name(i)).eq. &
575 & trim(vname(1,idbmas(itrc)))) THEN
576 got_var(idbmas(itrc))=.true.
577 s(ng)%Vid(idbmas(itrc))=var_id(i)
578# ifdef BEDLOAD
579 ELSE IF (trim(var_name(i)).eq. &
580 & trim(vname(1,idubld(itrc)))) THEN
581 got_var(idubld(itrc))=.true.
582 s(ng)%Vid(idubld(itrc))=var_id(i)
583 ELSE IF (trim(var_name(i)).eq. &
584 & trim(vname(1,idvbld(itrc)))) THEN
585 got_var(idvbld(itrc))=.true.
586 s(ng)%Vid(idvbld(itrc))=var_id(i)
587# endif
588 END IF
589 END DO
590 DO itrc=1,mbedp
591 IF (trim(var_name(i)).eq.trim(vname(1,idsbed(itrc)))) THEN
592 got_var(idsbed(itrc))=.true.
593 s(ng)%Vid(idsbed(itrc))=var_id(i)
594 END IF
595 END DO
596# endif
597# if defined SEDIMENT || defined BBL_MODEL
598 DO itrc=1,mbotp
599 IF (trim(var_name(i)).eq.trim(vname(1,idbott(itrc)))) THEN
600 got_var(idbott(itrc))=.true.
601 s(ng)%Vid(idbott(itrc))=var_id(i)
602 END IF
603 END DO
604# endif
605 END DO
606!
607! Check if output variables are available in input NetCDF file.
608!
609 IF (.not.got_var(idtime)) THEN
610 IF (master) WRITE (stdout,20) trim(vname(1,idtime)), &
611 & trim(ncname)
612 exit_flag=3
613 RETURN
614 END IF
615# if defined SEDIMENT && defined SED_MORPH
616 IF (.not.got_var(idbath).and.varout(idbath,ng)) THEN
617 IF (master) WRITE (stdout,20) trim(vname(1,idbath)), &
618 & trim(ncname)
619 exit_flag=3
620 RETURN
621 END IF
622# endif
623# if defined SEDIMENT && defined BEDLOAD && defined BEDLOAD_VANDERA
624 IF (.not.got_var(idsurs).and.varout(idsurs,ng)) THEN
625 IF (master) WRITE (stdout,20) trim(vname(1,idsurs)), &
626 & trim(ncname)
627 exit_flag=3
628 RETURN
629 END IF
630 IF (.not.got_var(idsrrw).and.varout(idsrrw,ng)) THEN
631 IF (master) WRITE (stdout,20) trim(vname(1,idsrrw)), &
632 & trim(ncname)
633 exit_flag=3
634 RETURN
635 END IF
636 IF (.not.got_var(idsbtw).and.varout(idsbtw,ng)) THEN
637 IF (master) WRITE (stdout,20) trim(vname(1,idsbtw)), &
638 & trim(ncname)
639 exit_flag=3
640 RETURN
641 END IF
642 IF (.not.got_var(idsucr).and.varout(idsucr,ng)) THEN
643 IF (master) WRITE (stdout,20) trim(vname(1,idsucr)), &
644 & trim(ncname)
645 exit_flag=3
646 RETURN
647 END IF
648 IF (.not.got_var(idsutr).and.varout(idsutr,ng)) THEN
649 IF (master) WRITE (stdout,20) trim(vname(1,idsutr)), &
650 & trim(ncname)
651 exit_flag=3
652 RETURN
653 END IF
654 IF (.not.got_var(idstcr).and.varout(idstcr,ng)) THEN
655 IF (master) WRITE (stdout,20) trim(vname(1,idstcr)), &
656 & trim(ncname)
657 exit_flag=3
658 RETURN
659 END IF
660 IF (.not.got_var(idsttr).and.varout(idsttr,ng)) THEN
661 IF (master) WRITE (stdout,20) trim(vname(1,idsttr)), &
662 & trim(ncname)
663 exit_flag=3
664 RETURN
665 END IF
666# endif
667# ifdef SEDIMENT
668 DO i=1,nst
669 IF (.not.got_var(idfrac(i)).and.varout(idfrac(i),ng)) THEN
670 IF (master) WRITE (stdout,20) trim(vname(1,idfrac(i))), &
671 & trim(ncname)
672 exit_flag=3
673 RETURN
674 END IF
675
676 IF(.not.got_var(idbmas(i)).and.varout(idbmas(i),ng)) THEN
677 IF (master) WRITE (stdout,20) trim(vname(1,idbmas(i))), &
678 & trim(ncname)
679 exit_flag=3
680 RETURN
681 END IF
682# ifdef BEDLOAD
683 IF (.not.got_var(idubld(i)).and.varout(idubld(i),ng)) THEN
684 IF (master) WRITE (stdout,20) trim(vname(1,idubld(i))), &
685 & trim(ncname)
686 exit_flag=3
687 RETURN
688 END IF
689 IF (.not.got_var(idvbld(i)).and.varout(idvbld(i),ng)) THEN
690 IF (master) WRITE (stdout,20) trim(vname(1,idvbld(i))), &
691 & trim(ncname)
692 exit_flag=3
693 RETURN
694 END IF
695# endif
696 END DO
697 DO i=1,mbedp
698 IF (.not.got_var(idsbed(i)).and.varout(idsbed(i),ng)) THEN
699 IF (master) WRITE (stdout,20) trim(vname(1,idsbed(i))), &
700 & trim(ncname)
701 exit_flag=3
702 RETURN
703 END IF
704 END DO
705# endif
706# if defined SEDIMENT || defined BBL_MODEL
707 DO i=1,mbotp
708 IF (.not.got_var(idbott(i)).and.varout(idbott(i),ng)) THEN
709 IF (master) WRITE (stdout,20) trim(vname(1,idbott(i))), &
710 & trim(ncname)
711 exit_flag=3
712 RETURN
713 END IF
714 END DO
715# endif
716 END IF query
717!
718 10 FORMAT (1pe11.4,1x,'millimeter')
719 20 FORMAT (/,' SEDIMENT_DEF_NF90 - unable to find variable: ', &
720 & a,2x,' in output NetCDF file: ',a)
721!
722 RETURN
723 END SUBROUTINE sediment_def_nf90
724
725# ifdef STATIONS
726!
727!***********************************************************************
728 SUBROUTINE sediment_def_station_nf90 (ng, model, ldef, VarOut, S, &
729 & bgrd, pgrd, rgrd)
730!***********************************************************************
731!
732 USE mod_netcdf
733!
734! Imported variable declarations.
735!
736 logical, intent(in) :: ldef, varout(nv,ngrids)
737!
738 integer, intent(in) :: ng, model
739 integer, intent(in), optional :: bgrd(:), pgrd(:), rgrd(:)
740!
741 TYPE(t_io), intent(inout) :: s(ngrids)
742!
743! Local variable declarations.
744!
745 logical :: got_var(nv)
746!
747 integer, parameter :: natt = 25
748
749 integer :: i, itrc, j, status
750!
751 real(r8) :: aval(6)
752!
753 character (len=120) :: vinfo(natt)
754 character (len=256) :: ncname
755!
756 character (len=*), parameter :: myfile = &
757 & __FILE__//", sediment_def_station_nf90"
758!
759 sourcefile=myfile
760!
761!-----------------------------------------------------------------------
762! Define sediment output stations variables.
763!-----------------------------------------------------------------------
764!
765 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
766 ncname=s(ng)%name
767!
768 define : IF (ldef) THEN
769!
770! Initialize local information variable arrays.
771!
772 DO i=1,natt
773 DO j=1,len(vinfo(1))
774 vinfo(i)(j:j)=' '
775 END DO
776 END DO
777 DO i=1,6
778 aval(i)=0.0_r8
779 END DO
780
781# if defined SEDIMENT && defined SED_MORPH
782!
783! Define time-varying bathymetry.
784!
785 IF (varout(idbath,ng)) THEN
786 vinfo( 1)=vname(1,idbath)
787 vinfo( 2)=vname(2,idbath)
788 vinfo( 3)=vname(3,idbath)
789 vinfo(14)=vname(4,idbath)
790 vinfo(16)=vname(1,idtime)
791 status=def_var(ng, model, s(ng)%ncid, s(ng)%Vid(idbath), &
792 & nf_fout, 2, pgrd, aval, vinfo, ncname, &
793 & setfillval = .false., &
794 & setparaccess = .true.)
795 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
796 END IF
797# endif
798
799# ifdef SEDIMENT
800!
801! Define sediment fraction of each size class in each bed layer.
802!
803 DO i=1,nst
804 IF (varout(idfrac(i),ng)) THEN
805 vinfo( 1)=vname(1,idfrac(i))
806 vinfo( 2)=vname(2,idfrac(i))
807 vinfo( 3)=vname(3,idfrac(i))
808 vinfo(14)=vname(4,idfrac(i))
809 vinfo(16)=vname(1,idtime)
810 WRITE (vinfo(19),10) 1000.0_r8*sd50(i,ng)
811 status=def_var(ng, model, s(ng)%ncid, &
812 & s(ng)%Vid(idfrac(i)), nf_fout, &
813 & 3, bgrd, aval, vinfo, ncname, &
814 & setfillval = .true., &
815 & setparaccess = .true.)
816 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
817 END IF
818!
819! Define sediment mass of each size class in each bed layer.
820!
821 IF (varout(idbmas(i),ng)) THEN
822 vinfo( 1)=vname(1,idbmas(i))
823 vinfo( 2)=vname(2,idbmas(i))
824 vinfo( 3)=vname(3,idbmas(i))
825 vinfo(14)=vname(4,idbmas(i))
826 vinfo(16)=vname(1,idtime)
827 WRITE (vinfo(19),10) 1000.0_r8*sd50(i,ng)
828 status=def_var(ng, model, s(ng)%ncid, &
829 & s(ng)%Vid(idbmas(i)), nf_fout, &
830 & 3, bgrd, aval, vinfo, ncname, &
831 & setfillval = .true., &
832 & setparaccess = .true.)
833 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
834 END IF
835 END DO
836!
837! Define sediment properties in each bed layer.
838!
839 DO i=1,mbedp
840 IF (varout(idsbed(i),ng)) THEN
841 vinfo( 1)=vname(1,idsbed(i))
842 vinfo( 2)=vname(2,idsbed(i))
843 vinfo( 3)=vname(3,idsbed(i))
844 vinfo(14)=vname(4,idsbed(i))
845 vinfo(16)=vname(1,idtime)
846 status=def_var(ng, model, s(ng)%ncid, &
847 & s(ng)%Vid(idsbed(i)), nf_fout, &
848 & 3, bgrd, aval, vinfo, ncname, &
849 & setfillval = .true., &
850 & setparaccess = .true.)
851 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
852 END IF
853 END DO
854# endif
855
856# if defined SEDIMENT || defined BBL_MODEL
857!
858! Define exposed sediment layer properties.
859!
860 DO i=1,mbotp
861 IF (varout(idbott(i),ng)) THEN
862 vinfo( 1)=vname(1,idbott(i))
863 vinfo( 2)=vname(2,idbott(i))
864 vinfo( 3)=vname(3,idbott(i))
865 vinfo(14)=vname(4,idbott(i))
866 vinfo(16)=vname(1,idtime)
867 status=def_var(ng, model, s(ng)%ncid, &
868 & s(ng)%Vid(idbott(i)), nf_fout, &
869 & 2, pgrd, aval, vinfo, ncname, &
870 & setfillval = .true., &
871 & setparaccess = .true.)
872 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
873 END IF
874 END DO
875# endif
876
877 END IF define
878!
879!-----------------------------------------------------------------------
880! Otherwise, check existing output file and prepare for appending
881! data.
882!-----------------------------------------------------------------------
883!
884 query : IF (.not.ldef) THEN
885!
886! Initialize locallogical switches.
887!
888 DO i=1,nv
889 got_var(i)=.false.
890 END DO
891!
892! Scan variable list from input NetCDF and activate switches for
893! Waves Effect on Currents variables. Get variable IDs.
894!
895 DO i=1,n_var
896 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
897 got_var(idtime)=.true.
898 s(ng)%Vid(idtime)=var_id(i)
899# if defined SEDIMENT && defined SED_MORPH
900 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idbath))) THEN
901 got_var(idbath)=.true.
902 s(ng)%Vid(idbath)=var_id(i)
903# endif
904 END IF
905# ifdef SEDIMENT
906 DO itrc=1,nst
907 IF (trim(var_name(i)).eq.trim(vname(1,idfrac(itrc)))) THEN
908 got_var(idfrac(itrc))=.true.
909 s(ng)%Vid(idfrac(itrc))=var_id(i)
910 ELSE IF (trim(var_name(i)).eq. &
911 & trim(vname(1,idbmas(itrc)))) THEN
912 got_var(idbmas(itrc))=.true.
913 s(ng)%Vid(idbmas(itrc))=var_id(i)
914 END IF
915 END DO
916 DO itrc=1,mbedp
917 IF (trim(var_name(i)).eq.trim(vname(1,idsbed(itrc)))) THEN
918 got_var(idsbed(itrc))=.true.
919 s(ng)%Vid(idsbed(itrc))=var_id(i)
920 END IF
921 END DO
922# endif
923# if defined SEDIMENT || defined BBL_MODEL
924 DO itrc=1,mbotp
925 IF (trim(var_name(i)).eq.trim(vname(1,idbott(itrc)))) THEN
926 got_var(idbott(itrc))=.true.
927 s(ng)%Vid(idbott(itrc))=var_id(i)
928 END IF
929 END DO
930# endif
931 END DO
932!
933! Check if output variables are available in input NetCDF file.
934!
935 IF (.not.got_var(idtime)) THEN
936 IF (master) WRITE (stdout,20) trim(vname(1,idtime)), &
937 & trim(ncname)
938 exit_flag=3
939 RETURN
940 END IF
941# if defined SEDIMENT && defined SED_MORPH
942 IF (.not.got_var(idbath).and.varout(idbath,ng)) THEN
943 IF (master) WRITE (stdout,20) trim(vname(1,idbath)), &
944 & trim(ncname)
945 exit_flag=3
946 RETURN
947 END IF
948# endif
949# ifdef SEDIMENT
950 DO i=1,nst
951 IF (.not.got_var(idfrac(i)).and.varout(idfrac(i),ng)) THEN
952 IF (master) WRITE (stdout,20) trim(vname(1,idfrac(i))), &
953 & trim(ncname)
954 exit_flag=3
955 RETURN
956 END IF
957 IF (.not.got_var(idbmas(i)).and.varout(idbmas(i),ng)) THEN
958 IF (master) WRITE (stdout,20) trim(vname(1,idbmas(i))), &
959 & trim(ncname)
960 exit_flag=3
961 RETURN
962 END IF
963 END DO
964 DO i=1,mbedp
965 IF (.not.got_var(idsbed(i)).and.varout(idsbed(i),ng)) THEN
966 IF (master) WRITE (stdout,20) trim(vname(1,idsbed(i))), &
967 & trim(ncname)
968 exit_flag=3
969 RETURN
970 END IF
971 END DO
972# endif
973# if defined SEDIMENT || defined BBL_MODEL
974 DO i=1,mbotp
975 IF (.not.got_var(idbott(i)).and.varout(idbott(i),ng)) THEN
976 IF (master) WRITE (stdout,20) trim(vname(1,idbott(i))), &
977 & trim(ncname)
978 exit_flag=3
979 RETURN
980 END IF
981 END DO
982# endif
983
984 END IF query
985!
986 10 FORMAT (1pe11.4,1x,'millimeter')
987 20 FORMAT (/,' SEDIMENT_DEF_STATION_NF90 - unable to find variable:',&
988 & 1x,a,2x,' in output NetCDF file: ',a)
989!
990 RETURN
991 END SUBROUTINE sediment_def_station_nf90
992# endif
993!
994!***********************************************************************
995 SUBROUTINE sediment_wrt_nf90 (ng, model, tile, &
996 & LBi, UBi, LBj, UBj, &
997 & VarOut, S)
998!***********************************************************************
999!
1000 USE mod_netcdf
1001!
1002! Imported variable declarations.
1003!
1004 logical, intent(in) :: varout(nv,ngrids)
1005!
1006 integer, intent(in) :: ng, model, tile
1007 integer, intent(in) :: lbi, ubi, lbj, ubj
1008!
1009 TYPE(t_io), intent(inout) :: s(ngrids)
1010!
1011! Local variable declarations.
1012!
1013 logical :: linstataneous
1014!
1015 integer :: gfactor, gtype, i, status
1016!
1017 real(dp) :: scale
1018!
1019 character (len=*), parameter :: myfile = &
1020 & __FILE__//", sediment_wrt_nf90"
1021!
1022 sourcefile=myfile
1023!
1024!-----------------------------------------------------------------------
1025! Write out Waves Effect on Currents output variables into specified
1026! output NetCDF file.
1027!-----------------------------------------------------------------------
1028!
1029 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1030!
1031! Set grid type factor to write full (gfactor=1) fields or water
1032! points (gfactor=-1) fields only.
1033!
1034# if defined WRITE_WATER && defined MASKING
1035 gfactor=-1
1036# else
1037 gfactor=1
1038# endif
1039!
1040! Set instantaneous fields.
1041!
1042 IF ((s(ng)%ncid.eq.s(ng)%ncid).or. &
1043 & (s(ng)%ncid.eq.qck(ng)%ncid)) THEN
1044 linstataneous=.true.
1045 ELSE
1046 linstataneous=.false. ! time-averged fiels
1047 END IF
1048
1049# if defined SEDIMENT && defined SED_MORPH
1050!
1051! Write out time-dependent bathymetry (m)
1052!
1053 IF (varout(idbath,ng)) THEN
1054 scale=1.0_dp
1055 gtype=gfactor*r2dvar
1056 IF (linstataneous) THEN
1057 status=nf_fwrite2d(ng, model, s(ng)%ncid, idbath, &
1058 & s(ng)%Vid(idbath), &
1059 & s(ng)%Rindex, gtype, &
1060 & lbi, ubi, lbj, ubj, scale, &
1061# ifdef MASKING
1062 & grid(ng) % rmask, &
1063# endif
1064 & grid(ng) % h, &
1065 & setfillval = .false.)
1066 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1067 IF (master) THEN
1068 WRITE (stdout,10) trim(vname(1,idbath)), s(ng)%Rindex
1069 END IF
1070 exit_flag=3
1071 ioerror=status
1072 RETURN
1073 END IF
1074 END IF
1075 END IF
1076# endif
1077
1078# if defined SEDIMENT && defined BEDLOAD
1079!
1080! Write out bed load transport in U-direction.
1081!
1082 DO i=1,nst
1083 IF (varout(idubld(i),ng)) THEN
1084 scale=1.0_dp
1085 gtype=gfactor*u2dvar
1086 IF (linstataneous) THEN
1087 status=nf_fwrite2d(ng, model, s(ng)%ncid, idubld(i), &
1088 & s(ng)%Vid(idubld(i)), &
1089 & s(ng)%Rindex, gtype, &
1090 & lbi, ubi, lbj, ubj, scale, &
1091# ifdef MASKING
1092 & grid(ng) % umask, &
1093# endif
1094 & sedbed(ng) % bedldu(:,:,i))
1095# ifdef AVERAGES
1096 ELSE
1097 status=nf_fwrite2d(ng, model, s(ng)%ncid, idubld(i), &
1098 & s(ng)%Vid(idubld(i)), &
1099 & s(ng)%Rindex, gtype, &
1100 & lbi, ubi, lbj, ubj, scale, &
1101# ifdef MASKING
1102 & grid(ng) % umask, &
1103# endif
1104 & sedbed(ng) % avgbedldu(:,:,i))
1105# endif
1106 END IF
1107 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1108 IF (master) THEN
1109 WRITE (stdout,10) trim(vname(1,idubld(i))), s(ng)%Rindex
1110 END IF
1111 exit_flag=3
1112 ioerror=status
1113 RETURN
1114 END IF
1115 END IF
1116!
1117! Write out bed load transport in V-direction.
1118!
1119 IF (varout(idvbld(i),ng)) THEN
1120 scale=1.0_dp
1121 gtype=gfactor*v2dvar
1122 IF (linstataneous) THEN
1123 status=nf_fwrite2d(ng, model, s(ng)%ncid, idvbld(i), &
1124 & s(ng)%Vid(idvbld(i)), &
1125 & s(ng)%Rindex, gtype, &
1126 & lbi, ubi, lbj, ubj, scale, &
1127# ifdef MASKING
1128 & grid(ng) % vmask, &
1129# endif
1130 & sedbed(ng) % bedldv(:,:,i))
1131# ifdef AVERAGES
1132 ELSE
1133 status=nf_fwrite2d(ng, model, s(ng)%ncid, idvbld(i), &
1134 & s(ng)%Vid(idvbld(i)), &
1135 & s(ng)%Rindex, gtype, &
1136 & lbi, ubi, lbj, ubj, scale, &
1137# ifdef MASKING
1138 & grid(ng) % vmask, &
1139# endif
1140 & sedbed(ng) % avgbedldv(:,:,i))
1141# endif
1142 END IF
1143 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1144 IF (master) THEN
1145 WRITE (stdout,10) trim(vname(1,idvbld(i))), s(ng)%Rindex
1146 END IF
1147 exit_flag=3
1148 ioerror=status
1149 RETURN
1150 END IF
1151 END IF
1152 END DO
1153# endif
1154
1155# ifdef SEDIMENT
1156!
1157! Write out sediment fraction of each size class in each bed layer.
1158!
1159 DO i=1,nst
1160 IF (varout(idfrac(i),ng)) THEN
1161 scale=1.0_dp
1162 gtype=gfactor*b3dvar
1163 status=nf_fwrite3d(ng, model, s(ng)%ncid, idfrac(i), &
1164 & s(ng)%Vid(idfrac(i)), &
1165 & s(ng)%Rindex, gtype, &
1166 & lbi, ubi, lbj, ubj, 1, nbed, scale, &
1167# ifdef MASKING
1168 & grid(ng) % rmask, &
1169# endif
1170 & sedbed(ng) % bed_frac(:,:,:,i))
1171 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1172 IF (master) THEN
1173 WRITE (stdout,10) trim(vname(1,idfrac(i))), s(ng)%Rindex
1174 END IF
1175 exit_flag=3
1176 ioerror=status
1177 RETURN
1178 END IF
1179 END IF
1180 END DO
1181!
1182! Write out sediment mass of each size class in each bed layer.
1183!
1184 DO i=1,nst
1185 IF (varout(idbmas(i),ng)) THEN
1186 scale=1.0_dp
1187 gtype=gfactor*b3dvar
1188 status=nf_fwrite3d(ng, model, s(ng)%ncid, idbmas(i), &
1189 & s(ng)%Vid(idbmas(i)), &
1190 & s(ng)%Rindex, gtype, &
1191 & lbi, ubi, lbj, ubj, 1, nbed, scale, &
1192# ifdef MASKING
1193 & grid(ng) % rmask, &
1194# endif
1195 & sedbed(ng) % bed_mass(:,:,:,nout,i))
1196 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1197 IF (master) THEN
1198 WRITE (stdout,10) trim(vname(1,idbmas(i))), s(ng)%Rindex
1199 END IF
1200 exit_flag=3
1201 ioerror=status
1202 RETURN
1203 END IF
1204 END IF
1205 END DO
1206!
1207! Write out sediment properties in each bed layer.
1208!
1209 DO i=1,mbedp
1210 IF (varout(idsbed(i),ng)) THEN
1211 scale=1.0_dp
1212 gtype=gfactor*b3dvar
1213 status=nf_fwrite3d(ng, model, s(ng)%ncid, idsbed(i), &
1214 & s(ng)%Vid(idsbed(i)), &
1215 & s(ng)%Rindex, gtype, &
1216 & lbi, ubi, lbj, ubj, 1, nbed, scale, &
1217# ifdef MASKING
1218 & grid(ng) % rmask, &
1219# endif
1220 & sedbed(ng) % bed(:,:,:,i))
1221 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1222 IF (master) THEN
1223 WRITE (stdout,10) trim(vname(1,idsbed(i))), s(ng)%Rindex
1224 END IF
1225 exit_flag=3
1226 ioerror=status
1227 RETURN
1228 END IF
1229 END IF
1230 END DO
1231# endif
1232
1233# if defined SEDIMENT || defined BBL_MODEL
1234!
1235! Write out exposed sediment layer properties.
1236!
1237 DO i=1,mbotp
1238 IF (varout(idbott(i),ng)) THEN
1239 IF (i.eq.itauc) THEN
1240 scale=rho0
1241 ELSE
1242 scale=1.0_dp
1243 END IF
1244 gtype=gfactor*r2dvar
1245 status=nf_fwrite2d(ng, model, s(ng)%ncid, idbott(i), &
1246 & s(ng)%Vid(idbott(i)), &
1247 & s(ng)%Rindex, gtype, &
1248 & lbi, ubi, lbj, ubj, scale, &
1249# ifdef MASKING
1250 & grid(ng) % rmask, &
1251# endif
1252 & sedbed(ng) % bottom(:,:,i))
1253 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1254 IF (master) THEN
1255 WRITE (stdout,10) trim(vname(1,idbott(i))), s(ng)%Rindex
1256 END IF
1257 exit_flag=3
1258 ioerror=status
1259 RETURN
1260 END IF
1261 END IF
1262 END DO
1263# endif
1264
1265# if defined SEDIMENT && defined BEDLOAD && defined SED_BEDLOAD_VANDERA
1266!
1267! Write out Ursell number of the asymmetric wave form.
1268!
1269 IF (varout(idsurs,ng)) THEN
1270 scale=1.0_dp
1271 gtype=gfactor*r2dvar
1272 status=nf_fwrite2d(ng, model, s(ng)%ncid, idsurs, &
1273 & s(ng)%Vid(idsurs), &
1274 & s(ng)%Rindex, gtype, &
1275 & lbi, ubi, lbj, ubj, scale, &
1276# ifdef MASKING
1277 & grid(ng) % rmask, &
1278# endif
1279 & sedbed(ng) % ursell_no)
1280 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1281 IF (master) THEN
1282 WRITE (stdout,10) trim(vname(1,idsurs)), s(ng)%Rindex
1283 END IF
1284 exit_flag=3
1285 ioerror=status
1286 RETURN
1287 END IF
1288 END IF
1289!
1290! Write out velocity skewness parameter of the asymmetric wave form.
1291!
1292 IF (varout(idsrrw,ng)) THEN
1293 scale=1.0_dp
1294 gtype=gfactor*r2dvar
1295 status=nf_fwrite2d(ng, model, s(ng)%ncid, idsrrw, &
1296 & s(ng)%Vid(idsrrw), &
1297 & s(ng)%Rindex, gtype, &
1298 & lbi, ubi, lbj, ubj, scale, &
1299# ifdef MASKING
1300 & grid(ng) % rmask, &
1301# endif
1302 & sedbed(ng) % RR_asymwave)
1303 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1304 IF (master) THEN
1305 WRITE (stdout,10) trim(vname(1,idsrrw)), s(ng)%Rindex
1306 END IF
1307 exit_flag=3
1308 ioerror=status
1309 RETURN
1310 END IF
1311 END IF
1312!
1313! Write out acceleration asymmetry parameter of the asymmetric wave form.
1314!
1315 IF (varout(idsbtw,ng)) THEN
1316 scale=1.0_dp
1317 gtype=gfactor*r2dvar
1318 status=nf_fwrite2d(ng, model, s(ng)%ncid, idsbtw, &
1319 & s(ng)%Vid(idsbtw), &
1320 & s(ng)%Rindex, gtype, &
1321 & lbi, ubi, lbj, ubj, scale, &
1322# ifdef MASKING
1323 & grid(ng) % rmask, &
1324# endif
1325 & sedbed(ng) % beta_asymwave)
1326 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1327 IF (master) THEN
1328 WRITE (stdout,10) trim(vname(1,idsbtw)), s(ng)%Rindex
1329 END IF
1330 exit_flag=3
1331 ioerror=status
1332 RETURN
1333 END IF
1334 END IF
1335!
1336! Write out crest velocity of the asymmetric wave form.
1337!
1338 IF (varout(idsucr,ng)) THEN
1339 scale=1.0_dp
1340 gtype=gfactor*r2dvar
1341 status=nf_fwrite2d(ng, model, s(ng)%ncid, idsucr, &
1342 & s(ng)%Vid(idsucr), &
1343 & s(ng)%Rindex, gtype, &
1344 & lbi, ubi, lbj, ubj, scale, &
1345# ifdef MASKING
1346 & grid(ng) % rmask, &
1347# endif
1348 & sedbed(ng) % ucrest_r)
1349 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1350 IF (master) THEN
1351 WRITE (stdout,10) trim(vname(1,idsucr)), s(ng)%Rindex
1352 END IF
1353 exit_flag=3
1354 ioerror=status
1355 RETURN
1356 END IF
1357 END IF
1358!
1359! Write out trough velocity of the asymmetric wave form.
1360!
1361 IF (varout(idsutr,ng)) THEN
1362 scale=1.0_dp
1363 gtype=gfactor*r2dvar
1364 status=nf_fwrite2d(ng, model, s(ng)%ncid, idsutr, &
1365 & s(ng)%Vid(idsutr), &
1366 & s(ng)%Rindex, gtype, &
1367 & lbi, ubi, lbj, ubj, scale, &
1368# ifdef MASKING
1369 & grid(ng) % rmask, &
1370# endif
1371 & sedbed(ng) % utrough_r)
1372 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1373 IF (master) THEN
1374 WRITE (stdout,10) trim(vname(1,idsutr)), s(ng)%Rindex
1375 END IF
1376 exit_flag=3
1377 ioerror=status
1378 RETURN
1379 END IF
1380 END IF
1381!
1382! Write out crest time period of the asymmetric wave form.
1383!
1384 IF (varout(idstcr,ng)) THEN
1385 scale=1.0_dp
1386 gtype=gfactor*r2dvar
1387 status=nf_fwrite2d(ng, model, s(ng)%ncid, idstcr, &
1388 & s(ng)%Vid(idstcr), &
1389 & s(ng)%Rindex, gtype, &
1390 & lbi, ubi, lbj, ubj, scale, &
1391# ifdef MASKING
1392 & grid(ng) % rmask, &
1393# endif
1394 & sedbed(ng) % T_crest)
1395 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1396 IF (master) THEN
1397 WRITE (stdout,10) trim(vname(1,idstcr)), s(ng)%Rindex
1398 END IF
1399 exit_flag=3
1400 ioerror=status
1401 RETURN
1402 END IF
1403 END IF
1404!
1405! Write out trough time period of the asymmetric wave form.
1406!
1407 IF (varout(idsttr,ng)) THEN
1408 scale=1.0_dp
1409 gtype=gfactor*r2dvar
1410 status=nf_fwrite2d(ng, model, s(ng)%ncid, idsttr, &
1411 & s(ng)%Vid(idsttr), &
1412 & s(ng)%Rindex, gtype, &
1413 & lbi, ubi, lbj, ubj, scale, &
1414# ifdef MASKING
1415 & grid(ng) % rmask, &
1416# endif
1417 & sedbed(ng) % T_trough)
1418 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
1419 IF (master) THEN
1420 WRITE (stdout,10) trim(vname(1,idsttr)), s(ng)%Rindex
1421 END IF
1422 exit_flag=3
1423 ioerror=status
1424 RETURN
1425 END IF
1426 END IF
1427# endif
1428!
1429 10 FORMAT (/," SEDIMENT_WRT_NF90 - error while writing variable '", &
1430 & a,"', time record = ",i0,/,11x,'into file: ',a)
1431!
1432 RETURN
1433 END SUBROUTINE sediment_wrt_nf90
1434
1435# ifdef STATIONS
1436!
1437!***********************************************************************
1438 SUBROUTINE sediment_wrt_station_nf90 (ng, model, tile, &
1439 & LBi, UBi, LBj, UBj, &
1440 & VarOut, S)
1441!***********************************************************************
1442!
1443 USE mod_netcdf
1444!
1445! Imported variable declarations.
1446!
1447 logical, intent(in) :: varout(nv,ngrids)
1448!
1449 integer, intent(in) :: ng, model, tile
1450 integer, intent(in) :: lbi, ubi, lbj, ubj
1451!
1452 TYPE(t_io), intent(inout) :: s(ngrids)
1453!
1454! Local variable declarations.
1455!
1456 logical :: cgrid
1457!
1458 integer :: nposb
1459 integer :: i, k, np, status
1460!
1461 real(dp) :: scale
1462!
1463 real(r8), dimension(Nstation(ng)) :: xpos, ypos, zpos, psta
1464# ifdef SEDIMENT
1465 real(r8), dimension(Nstation(ng)*Nbed) :: xposb, yposb, zposb
1466 real(r8), dimension(Nstation(ng)*Nbed) :: bsta
1467# endif
1468!
1469 character (len=*), parameter :: myfile = &
1470 & __FILE__//", sediment_wrt_station_nf90"
1471!
1472 sourcefile=myfile
1473!
1474!-----------------------------------------------------------------------
1475! Write out sediment output variables into specified stations NetCDF
1476! file.
1477!-----------------------------------------------------------------------
1478!
1479 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1480!
1481! Set switch to extract station data at native C-grid position (TRUE)
1482! or at RHO-points (FALSE).
1483!
1484# ifdef STATIONS_CGRID
1485 cgrid=.true.
1486# else
1487 cgrid=.false.
1488# endif
1489!
1490! Set positions for generic extraction routine.
1491!
1492 nposb=nstation(ng)*nbed
1493 DO i=1,nstation(ng)
1494 xpos(i)=scalars(ng)%SposX(i)
1495 ypos(i)=scalars(ng)%SposY(i)
1496 zpos(i)=1.0_r8
1497# ifdef SEDIMENT
1498 DO k=1,nbed
1499 np=k+(i-1)*nbed
1500 xposb(np)=scalars(ng)%SposX(i)
1501 yposb(np)=scalars(ng)%SposY(i)
1502 zposb(np)=real(k,r8)
1503 END DO
1504# endif
1505 END DO
1506
1507# if defined SEDIMENT && defined SED_MORPH
1508!
1509! Write out time-varying bathymetry.
1510!
1511 IF (varout(idbath,ng)) THEN
1512 scale=1.0_dp
1513 CALL extract_sta2d (ng, model, cgrid, idbath, r2dvar, &
1514 & lbi, ubi, lbj, ubj, &
1515 & scale, grid(ng)%h, &
1516 & nstation(ng), xpos, ypos, psta)
1517 CALL netcdf_put_fvar (ng, model, s(ng)%name, &
1518 & trim(vname(1,idbath)), psta, &
1519 & (/1,s(ng)%Rindex/), (/nstation(ng),1/), &
1520 & ncid = s(ng)%ncid, &
1521 & varid = s(ng)%Vid(idbath))
1522 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1523 END IF
1524# endif
1525
1526# ifdef SEDIMENT
1527!
1528! Write out sediment fraction of each size class in each bed layer.
1529!
1530 DO i=1,nst
1531 IF (varout(idfrac(i),ng)) THEN
1532 scale=1.0_dp
1533 CALL extract_sta3d (ng, model, cgrid, idfrac(i), b3dvar, &
1534 & lbi, ubi, lbj, ubj, 1, nbed, &
1535 & scale, sedbed(ng)%bed_frac(:,:,:,i), &
1536 & nposb, xposb, yposb, zposb, bsta)
1537 CALL netcdf_put_fvar (ng, model, s(ng)%name, &
1538 & trim(vname(1,idfrac(i))), bsta, &
1539 & (/1,1,s(ng)%Rindex/), &
1540 & (/nbed,nstation(ng),1/), &
1541 & ncid = s(ng)%ncid, &
1542 & varid = s(ng)%Vid(idfrac(i)))
1543 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1544 END IF
1545!
1546! Write out sediment mass of each size class in each bed layer.
1547!
1548 IF (varout(idbmas(i),ng)) THEN
1549 scale=1.0_dp
1550 CALL extract_sta3d (ng, model, cgrid, idbmas(i), b3dvar, &
1551 & lbi, ubi, lbj, ubj, 1, nbed, &
1552 & scale, &
1553 & sedbed(ng)%bed_mass(:,:,:,nout,i), &
1554 & nposb, xposb, yposb, zposb, bsta)
1555 CALL netcdf_put_fvar (ng, model, s(ng)%name, &
1556 & trim(vname(1,idbmas(i))), bsta, &
1557 & (/1,1,s(ng)%Rindex/), &
1558 & (/nbed,nstation(ng),1/), &
1559 & ncid = s(ng)%ncid, &
1560 & varid = s(ng)%Vid(idbmas(i)))
1561 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1562 END IF
1563 END DO
1564!
1565! Write out sediment properties in each bed layer.
1566!
1567 DO i=1,mbedp
1568 IF (varout(idsbed(i),ng)) THEN
1569 scale=1.0_dp
1570 CALL extract_sta3d (ng, model, cgrid, idsbed(i), b3dvar, &
1571 & lbi, ubi, lbj, ubj, 1, nbed, &
1572 & scale, sedbed(ng)%bed(:,:,:,i), &
1573 & nposb, xposb, yposb, zposb, bsta)
1574 CALL netcdf_put_fvar (ng, model, s(ng)%name, &
1575 & trim(vname(1,idsbed(i))), bsta, &
1576 & (/1,1,s(ng)%Rindex/), &
1577 & (/nbed,nstation(ng),1/), &
1578 & ncid = s(ng)%ncid, &
1579 & varid = s(ng)%Vid(idsbed(i)))
1580 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1581 END IF
1582 END DO
1583# endif
1584
1585# if defined SEDIMENT || defined BBL_MODEL
1586!
1587! Write out exposed sediment layer properties.
1588!
1589 DO i=1,mbedp
1590 IF (varout(idbott(i),ng)) THEN
1591 scale=1.0_dp
1592 CALL extract_sta2d (ng, model, cgrid, idbott(i), r2dvar, &
1593 & lbi, ubi, lbj, ubj, &
1594 & scale, sedbed(ng)%bottom(:,:,i), &
1595 & nstation(ng), xpos, ypos, psta)
1596 CALL netcdf_put_fvar (ng, model, s(ng)%name, &
1597 & trim(vname(1,idbott(i))), bsta, &
1598 & (/1,s(ng)%Rindex/), &
1599 & (/nstation(ng),1/), &
1600 & ncid = s(ng)%ncid, &
1601 & varid = s(ng)%Vid(idbott(i)))
1602 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1603 END IF
1604 END DO
1605# endif
1606!
1607 RETURN
1608 END SUBROUTINE sediment_wrt_station_nf90
1609# endif
1610
1611# if defined PIO_LIB && defined DISTRIBUTE
1612!
1613!***********************************************************************
1614 SUBROUTINE sediment_def_pio (ng, model, ldef, VarOut, S, &
1615 & t2dgrd, u2dgrd, v2dgrd, &
1616 & t3dgrd, u3dgrd, v3dgrd, w3dgrd)
1617!***********************************************************************
1618!
1619 USE mod_pio_netcdf
1620!
1621! Imported variable declarations.
1622!
1623 logical, intent(in) :: ldef, varout(nv,ngrids)
1624!
1625 integer, intent(in) :: ng, model
1626 integer, intent(in), optional :: t2dgrd(:), u2dgrd(:), v2dgrd(:)
1627 integer, intent(in), optional :: t3dgrd(:), u3dgrd(:), v3dgrd(:)
1628 integer, intent(in), optional :: w3dgrd(:)
1629!
1630 TYPE(t_io), intent(inout) :: s(ngrids)
1631!
1632! Local variable declarations.
1633!
1634 logical :: got_var(nv)
1635!
1636 integer, parameter :: natt = 25
1637
1638 integer :: i, itrc, j, nvd3, nvd4, status
1639!
1640 real(r8) :: aval(6)
1641!
1642# ifdef ADJOINT
1643 character (len=21) :: prefix
1644# else
1645 character (len=13) :: prefix
1646# endif
1647 character (len=120) :: vinfo(natt)
1648 character (len=256) :: ncname
1649!
1650 character (len=*), parameter :: myfile = &
1651 & __FILE__//", sediment_def_pio"
1652!
1653 sourcefile=myfile
1654!
1655!-----------------------------------------------------------------------
1656! Define Waves Effect on Currents output variables.
1657!-----------------------------------------------------------------------
1658!
1659 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1660 ncname=s(ng)%name
1661!
1662 define : IF (ldef) THEN
1663!
1664! Set number of dimensions for output variables.
1665!
1666# if defined WRITE_WATER && defined MASKING
1667 nvd3=2
1668 nvd4=2
1669# else
1670 nvd3=3
1671 nvd4=4
1672# endif
1673!
1674! Set long name prefix string.
1675!
1676# ifdef ADJOINT
1677!! Prefix='time-averaged adjoint'
1678 prefix='adjoint'
1679# else
1680!! Prefix='time-averaged'
1681 prefix=char(32) ! blank space
1682# endif
1683!
1684! Initialize local information variable arrays.
1685!
1686 DO i=1,natt
1687 DO j=1,len(vinfo(1))
1688 vinfo(i)(j:j)=' '
1689 END DO
1690 END DO
1691 DO i=1,6
1692 aval(i)=0.0_r8
1693 END DO
1694
1695# if defined SEDIMENT && defined SED_MORPH
1696!
1697! Define time-varying bathymetry.
1698!
1699 IF (varout(idbath,ng)) THEN
1700 vinfo( 1)=vname(1,idbath)
1701 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1702 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idbath))
1703 ELSE
1704 vinfo( 2)=vname(2,idbath)
1705 END IF
1706 vinfo( 3)=vname(3,idbath)
1707 vinfo(14)=vname(4,idbath)
1708 vinfo(16)=vname(1,idtime)
1709 vinfo(21)=vname(6,idbath)
1710 vinfo(22)='coordinates'
1711 aval(5)=real(iinfo(1,idbath,ng),r8)
1712 s(ng)%pioVar(idbath)%dkind=pio_fout
1713 s(ng)%pioVar(idbath)%gtype=r2dvar
1714!
1715 status=def_var(ng, model, s(ng)%pioFile, &
1716 & s(ng)%pioVar(idbath)%vd, &
1717 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname, &
1718 & setfillval = .false.)
1719 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1720 END IF
1721# endif
1722
1723# if defined SEDIMENT && defined BEDLOAD
1724!
1725! Define Bedload transbort U-direction.
1726!
1727 DO i=1,nst
1728 IF (varout(idubld(i),ng)) THEN
1729 vinfo( 1)=vname(1,idubld(i))
1730 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1731 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
1732 & trim(vname(2,idubld(i)))
1733 ELSE
1734 vinfo( 2)=vname(2,idubld(i))
1735 END IF
1736 vinfo( 3)=vname(3,idubld(i))
1737 vinfo(14)=vname(4,idubld(i))
1738 vinfo(16)=vname(1,idtime)
1739# if defined WRITE_WATER && defined MASKING
1740 vinfo(20)='mask_u'
1741# endif
1742 vinfo(21)=vname(6,idubld(i))
1743 vinfo(22)='coordinates'
1744 aval(5)=real(iinfo(1,idubld(i),ng),r8)
1745 s(ng)%pioVar(idubld(i))%dkind=pio_fout
1746 s(ng)%pioVar(idubld(i))%gtype=u2dvar
1747!
1748 status=def_var(ng, model, s(ng)%pioFile, &
1749 & s(ng)%pioVar(idubld(i))%vd, &
1750 & pio_fout, nvd3, u2dgrd, aval, vinfo, ncname)
1751 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1752 END IF
1753!
1754! Define Bedload transport V-direction.
1755!
1756 IF (varout(idvbld(i),ng)) THEN
1757 vinfo( 1)=vname(1,idvbld(i))
1758 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1759 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
1760 & trim(vname(2,idvbld(i)))
1761 ELSE
1762 vinfo( 2)=vname(2,idvbld(i))
1763 END IF
1764 vinfo( 3)=vname(3,idvbld(i))
1765 vinfo(14)=vname(4,idvbld(i))
1766 vinfo(16)=vname(1,idtime)
1767# if defined WRITE_WATER && defined MASKING
1768 vinfo(20)='mask_v'
1769# endif
1770 vinfo(21)=vname(6,idvbld(i))
1771 vinfo(22)='coordinates'
1772 aval(5)=real(iinfo(1,idvbld(i),ng),r8)
1773 s(ng)%pioVar(idvbld(i))%dkind=pio_fout
1774 s(ng)%pioVar(idvbld(i))%gtype=v2dvar
1775!
1776 status=def_var(ng, model, s(ng)%pioFile, &
1777 & s(ng)%pioVar(idvbld(i))%vd, &
1778 & pio_fout, nvd3, v2dgrd, aval, vinfo, ncname)
1779 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1780 END IF
1781 END DO
1782# endif
1783
1784# ifdef SEDIMENT
1785!
1786! Define sediment fraction of each size class in each bed layer.
1787!
1788 DO i=1,nst
1789 IF (varout(idfrac(i),ng)) THEN
1790 vinfo( 1)=vname(1,idfrac(i))
1791 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1792 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
1793 & trim(vname(2,idfrac(i)))
1794 ELSE
1795 vinfo( 2)=vname(2,idfrac(i))
1796 END IF
1797 vinfo( 3)=vname(3,idfrac(i))
1798 vinfo(14)=vname(4,idfrac(i))
1799 vinfo(16)=vname(1,idtime)
1800 WRITE (vinfo(19),10) 1000.0_r8*sd50(i,ng)
1801# if defined WRITE_WATER && defined MASKING
1802 vinfo(20)='mask_rho'
1803# endif
1804 vinfo(21)=vname(6,idfrac(i))
1805 vinfo(22)='coordinates'
1806 aval(5)=real(iinfo(1,idfrac(i),ng))
1807 s(ng)%pioVar(idfrac(i))%dkind=pio_fout
1808 s(ng)%pioVar(idfrac(i))%gtype=b3dvar
1809!
1810 status=def_var(ng, model, s(ng)%pioFile, &
1811 & s(ng)%pioVar(idfrac(i))%vd, &
1812 & pio_fout, nvd4, b3dgrd, aval, vinfo, ncname)
1813 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1814 END IF
1815 END DO
1816!
1817! Define sediment mass of each size class in each bed layer.
1818!
1819 DO i=1,nst
1820 IF (varout(idbmas(i),ng)) THEN
1821 vinfo( 1)=vname(1,idbmas(i))
1822 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1823 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
1824 & trim(vname(2,idbmas(i)))
1825 ELSE
1826 vinfo( 2)=vname(2,idbmas(i))
1827 END IF
1828 vinfo( 3)=vname(3,idbmas(i))
1829 vinfo(14)=vname(4,idbmas(i))
1830 vinfo(16)=vname(1,idtime)
1831 WRITE (vinfo(19),10) 1000.0_r8*sd50(i,ng)
1832# if defined WRITE_WATER && defined MASKING
1833 vinfo(20)='mask_rho'
1834# endif
1835 vinfo(21)=vname(6,idbmas(i))
1836 vinfo(22)='coordinates'
1837 aval(5)=real(iinfo(1,idbmas(i),ng),r8)
1838 s(ng)%pioVar(idbmas(i))%dkind=pio_fout
1839 s(ng)%pioVar(idbmas(i))%gtype=b3dvar
1840!
1841 status=def_var(ng, model, s(ng)%pioFile, &
1842 & s(ng)%pioVar(idbmas(i))%vd, &
1843 & pio_fout, nvd4, b3dgrd, aval, vinfo, ncname)
1844 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1845 END IF
1846 END DO
1847!
1848! Define sediment properties in each bed layer.
1849!
1850 DO i=1,mbedp
1851 IF (varout(idsbed(i),ng)) THEN
1852 vinfo( 1)=vname(1,idsbed(i))
1853 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1854 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
1855 & trim(vname(2,idsbed(i)))
1856 ELSE
1857 vinfo( 2)=vname(2,idsbed(i))
1858 END IF
1859 vinfo( 3)=vname(3,idsbed(i))
1860 vinfo(14)=vname(4,idsbed(i))
1861 vinfo(16)=vname(1,idtime)
1862# if defined WRITE_WATER && defined MASKING
1863 vinfo(20)='mask_rho'
1864# endif
1865 vinfo(21)=vname(6,idsbed(i))
1866 vinfo(22)='coordinates'
1867 aval(5)=real(iinfo(1,idsbed(i),ng),r8)
1868 s(ng)%pioVar(idsbed(i))%dkind=pio_fout
1869 s(ng)%pioVar(idsbed(i))%gtype=b3dvar
1870!
1871 status=def_var(ng, model, s(ng)%pioFile, &
1872 & s(ng)%pioVar(idsbed(i))%vd, &
1873 pio_fout, nvd4, b3dgrd, aval, vinfo, ncname)
1874 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1875 END IF
1876 END DO
1877# endif
1878
1879# if defined SEDIMENT || defined BBL_MODEL
1880!
1881! Define exposed sediment layer properties.
1882!
1883 DO i=1,mbotp
1884 IF (varout(idbott(i),ng)) THEN
1885 vinfo( 1)=vname(1,idbott(i))
1886 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1887 WRITE (vinfo( 2),'(a,1x,a)') prefix, &
1888 & trim(vname(2,idbott(i)))
1889 ELSE
1890 vinfo( 2)=vname(2,idbott(i))
1891 END IF
1892 vinfo( 3)=vname(3,idbott(i))
1893 vinfo(14)=vname(4,idbott(i))
1894 vinfo(16)=vname(1,idtime)
1895# if defined WRITE_WATER && defined MASKING
1896 vinfo(20)='mask_rho'
1897# endif
1898 vinfo(21)=vname(6,idbott(i))
1899 vinfo(22)='coordinates'
1900 aval(5)=real(iinfo(1,idbott(i),ng),r8)
1901 s(ng)%pioVar(idbott(i))%dkind=pio_fout
1902 s(ng)%pioVar(idbott(i))%gtype=r2dvar
1903!
1904 status=def_var(ng, model, s(ng)%pioFile, &
1905 & s(ng)%pioVar(idbott(i))%vd, &
1906 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
1907 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1908 END IF
1909 END DO
1910# endif
1911
1912# if defined SEDIMENT && defined BEDLOAD && defined BEDLOAD_VANDERA
1913!
1914! Define Ursell number of the asymmetric wave form.
1915!
1916 IF (varout(idsurs,ng)) THEN
1917 vinfo( 1)=vname(1,idsurs)
1918 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1919 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsurs))
1920 ELSE
1921 vinfo( 2)=vname(2,idsurs)
1922 END IF
1923 vinfo( 3)=vname(3,idsurs)
1924 vinfo(14)=vname(4,idsurs)
1925 vinfo(16)=vname(1,idsurs)
1926# if defined WRITE_WATER && defined MASKING
1927 vinfo(20)='mask_rho'
1928# endif
1929 vinfo(21)=vname(6,idsurs)
1930 vinfo(22)='coordinates'
1931 aval(5)=real(iinfo(1,idsurs,ng),r8)
1932 s(ng)%pioVar(idsurs)%dkind=pio_fout
1933 s(ng)%pioVar(idsurs)%gtype=r2dvar
1934!
1935 status=def_var(ng, model, s(ng)%pioFile, &
1936 & s(ng)%pioVar(idsurs)%vd, &
1937 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
1938 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1939 END IF
1940!
1941! Define velocity skewness parameter of the asymmetric wave form.
1942!
1943 IF (varout(idsrrw,ng)) THEN
1944 vinfo( 1)=vname(1,idsrrw)
1945 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1946 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsrrw))
1947 ELSE
1948 vinfo( 2)=vname(2,idsrrw)
1949 END IF
1950 vinfo( 3)=vname(3,idsrrw)
1951 vinfo(14)=vname(4,idsrrw)
1952 vinfo(16)=vname(1,idsrrw)
1953# if defined WRITE_WATER && defined MASKING
1954 vinfo(20)='mask_rho'
1955# endif
1956 vinfo(21)=vname(6,idsrrw)
1957 vinfo(22)='coordinates'
1958 aval(5)=real(iinfo(1,idsrrw,ng),r8)
1959 s(ng)%pioVar(idsrrw)%dkind=pio_fout
1960 s(ng)%pioVar(idsrrw)%gtype=r2dvar
1961!
1962 status=def_var(ng, model, s(ng)%pioFile, &
1963 & s(ng)%pioVar(idsrrw)%vd, &
1964 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
1965 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1966 END IF
1967!
1968! Define acceleration asymmetry parameter of the asymmetric wave form.
1969!
1970 IF (varout(idsbtw,ng)) THEN
1971 vinfo( 1)=vname(1,idsbtw)
1972 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1973 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsbtw))
1974 ELSE
1975 vinfo( 2)=vname(2,idsbtw)
1976 END IF
1977 vinfo( 3)=vname(3,idsbtw)
1978 vinfo(14)=vname(4,idsbtw)
1979 vinfo(16)=vname(1,idsbtw)
1980# if defined WRITE_WATER && defined MASKING
1981 vinfo(20)='mask_rho'
1982# endif
1983 vinfo(21)=vname(6,idsbtw)
1984 vinfo(22)='coordinates'
1985 aval(5)=real(iinfo(1,idsbtw,ng),r8)
1986 s(ng)%pioVar(idsbtw)%dkind=pio_fout
1987 s(ng)%pioVar(idsbtw)%gtype=r2dvar
1988!
1989 status=def_var(ng, model, s(ng)%pioFile, &
1990 & s(ng)%pioVar(idsbtw)%vd, &
1991 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
1992 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1993 END IF
1994!
1995! Define crest velocity of the asymmetric wave form.
1996!
1997 IF (varout(idsucr,ng)) THEN
1998 vinfo( 1)=vname(1,idsucr)
1999 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
2000 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsucr))
2001 ELSE
2002 vinfo( 2)=vname(2,idsucr)
2003 END IF
2004 vinfo( 3)=vname(3,idsucr)
2005 vinfo(14)=vname(4,idsucr)
2006 vinfo(16)=vname(1,idsucr)
2007# if defined WRITE_WATER && defined MASKING
2008 vinfo(20)='mask_rho'
2009# endif
2010 vinfo(21)=vname(6,idsucr)
2011 vinfo(22)='coordinates'
2012 aval(5)=real(iinfo(1,idsucr,ng),r8)
2013 s(ng)%pioVar(idsucr)%dkind=pio_fout
2014 s(ng)%pioVar(idsucr)%gtype=r2dvar
2015!
2016 status=def_var(ng, model, s(ng)%pioFile, &
2017 & s(ng)%pioVar(idsucr)%vd, &
2018 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
2019 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2020 END IF
2021!
2022! Define trough velocity of the asymmetric wave form.
2023!
2024 IF (varout(idsutr,ng)) THEN
2025 vinfo( 1)=vname(1,idsutr)
2026 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
2027 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsutr))
2028 ELSE
2029 vinfo( 2)=vname(2,idsutr)
2030 END IF
2031 vinfo( 3)=vname(3,idsutr)
2032 vinfo(14)=vname(4,idsutr)
2033 vinfo(16)=vname(1,idsutr)
2034# if defined WRITE_WATER && defined MASKING
2035 vinfo(20)='mask_rho'
2036# endif
2037 vinfo(21)=vname(6,idsutr)
2038 vinfo(22)='coordinates'
2039 aval(5)=real(iinfo(1,idsutr,ng),r8)
2040 s(ng)%pioVar(idsutr)%dkind=pio_fout
2041 s(ng)%pioVar(idsutr)%gtype=r2dvar
2042!
2043 status=def_var(ng, model, s(ng)%pioFile, &
2044 & s(ng)%pioVar(idsutr)%vd, &
2045 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
2046 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2047 END IF
2048!
2049! Define crest time period of the asymmetric wave form.
2050!
2051 IF (varout(idstcr,ng)) THEN
2052 vinfo( 1)=vname(1,idstcr)
2053 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
2054 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idstcr))
2055 ELSE
2056 vinfo( 2)=vname(2,idstcr)
2057 END IF
2058 vinfo( 3)=vname(3,idstcr)
2059 vinfo(14)=vname(4,idstcr)
2060 vinfo(16)=vname(1,idstcr)
2061# if defined WRITE_WATER && defined MASKING
2062 vinfo(20)='mask_rho'
2063# endif
2064 vinfo(21)=vname(6,idstcr)
2065 vinfo(22)='coordinates'
2066 aval(5)=real(iinfo(1,idstcr,ng),r8)
2067 s(ng)%pioVar(idstcr)%dkind=pio_fout
2068 s(ng)%pioVar(idstcr)%gtype=r2dvar
2069!
2070 status=def_var(ng, model, s(ng)%pioFile, &
2071 & s(ng)%pioVar(idstcr)%vd, &
2072 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
2073 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2074 END IF
2075!
2076! Define trough time period of the asymmetric wave form.
2077!
2078 IF (varout(idsttr,ng)) THEN
2079 vinfo( 1)=vname(1,idsttr)
2080 IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
2081 WRITE (vinfo( 2),'(a,1x,a)') prefix, trim(vname(2,idsttr))
2082 ELSE
2083 vinfo( 2)=vname(2,idsttr)
2084 END IF
2085 vinfo( 3)=vname(3,idsttr)
2086 vinfo(14)=vname(4,idsttr)
2087 vinfo(16)=vname(1,idsttr)
2088# if defined WRITE_WATER && defined MASKING
2089 vinfo(20)='mask_rho'
2090# endif
2091 vinfo(21)=vname(6,idsttr)
2092 vinfo(22)='coordinates'
2093 aval(5)=real(iinfo(1,idsttr,ng),r8)
2094 s(ng)%pioVar(idsttr)%dkind=pio_fout
2095 s(ng)%pioVar(idsttr)%gtype=r2dvar
2096!
2097 status=def_var(ng, model, s(ng)%pioFile, &
2098 & s(ng)%pioVar(idsttr)%vd, &
2099 & pio_fout, nvd3, t2dgrd, aval, vinfo, ncname)
2100 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2101 END IF
2102# endif
2103
2104 END IF define
2105!
2106!-----------------------------------------------------------------------
2107! Otherwise, check existing output file and prepare for appending
2108! data.
2109!-----------------------------------------------------------------------
2110!
2111 query : IF (.not.ldef) THEN
2112!
2113! Initialize locallogical switches.
2114!
2115 DO i=1,nv
2116 got_var(i)=.false.
2117 END DO
2118!
2119! Scan variable list from input NetCDF and activate switches for
2120! Waves Effect on Currents variables. Get variable IDs.
2121!
2122 DO i=1,n_var
2123 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
2124 got_var(idtime)=.true.
2125 s(ng)%pioVar(idtime)%vd=var_desc(i)
2126 s(ng)%pioVar(idtime)%dkind=pio_tout
2127 s(ng)%pioVar(idtime)%gtype=0
2128# if defined SEDIMENT && defined SED_MORPH
2129 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idbath))) THEN
2130 got_var(idbath)=.true.
2131 s(ng)%pioVar(idbath)%vd=var_desc(i)
2132 s(ng)%pioVar(idbath)%dkind=pio_fout
2133 s(ng)%pioVar(idbath)%gtype=r2dvar
2134# endif
2135# if defined SEDIMENT && defined BEDLOAD && defined BEDLOAD_VANDERA
2136 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsurs))) THEN
2137 got_var(idsurs)=.true.
2138 s(ng)%pioVar(idsurs)%vd=var_desc(i)
2139 s(ng)%pioVar(idsurs)%dkind=pio_fout
2140 s(ng)%pioVar(idsurs)%gtype=r2dvar
2141 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsrrw))) THEN
2142 got_var(idsrrw)=.true.
2143 s(ng)%pioVar(idsrrw)%vd=var_desc(i)
2144 s(ng)%pioVar(idsrrw)%dkind=pio_fout
2145 s(ng)%pioVar(idsrrw)%gtype=r2dvar
2146 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsbtw))) THEN
2147 got_var(idsbtw)=.true.
2148 s(ng)%pioVar(idsbtw)%vd=var_desc(i)
2149 s(ng)%pioVar(idsbtw)%dkind=pio_fout
2150 s(ng)%pioVar(idsbtw)%gtype=r2dvar
2151 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsucr))) THEN
2152 got_var(idsucr)=.true.
2153 s(ng)%pioVar(idsucr)%vd=var_desc(i)
2154 s(ng)%pioVar(idsucr)%dkind=pio_fout
2155 s(ng)%pioVar(idsucr)%gtype=r2dvar
2156 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsutr))) THEN
2157 got_var(idsutr)=.true.
2158 s(ng)%pioVar(idsutr)%vd=var_desc(i)
2159 s(ng)%pioVar(idsutr)%dkind=pio_fout
2160 s(ng)%pioVar(idsutr)%gtype=r2dvar
2161 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idstcr))) THEN
2162 got_var(idstcr)=.true.
2163 s(ng)%pioVar(idstcr)%vd=var_desc(i)
2164 s(ng)%pioVar(idstcr)%dkind=pio_fout
2165 s(ng)%pioVar(idstcr)%gtype=r2dvar
2166 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsttr))) THEN
2167 got_var(idsttr)=.true.
2168 s(ng)%pioVar(idsttr)%vd=var_desc(i)
2169 s(ng)%pioVar(idsttr)%dkind=pio_fout
2170 s(ng)%pioVar(idsttr)%gtype=r2dvar
2171# endif
2172# ifdef SEDIMENT
2173 DO itrc=1,nst
2174 IF (trim(var_name(i)).eq. &
2175 & trim(vname(1,idfrac(itrc)))) THEN
2176 got_var(idfrac(itrc))=.true.
2177 s(ng)%pioVar(idfrac(itrc))%vd=var_desc(i)
2178 s(ng)%pioVar(idfrac(itrc))%dkind=pio_fout
2179 s(ng)%pioVar(idfrac(itrc))%gtype=b3dvar
2180 ELSE IF (trim(var_name(i)).eq. &
2181 & trim(vname(1,idbmas(itrc)))) THEN
2182 got_var(idbmas(itrc))=.true.
2183 s(ng)%pioVar(idbmas(itrc))%vd=var_desc(i)
2184 s(ng)%pioVar(idbmas(itrc))%dkind=pio_fout
2185 s(ng)%pioVar(idbmas(itrc))%gtype=b3dvar
2186# ifdef BEDLOAD
2187 ELSE IF (trim(var_name(i)).eq. &
2188 & trim(vname(1,idubld(itrc)))) THEN
2189 got_var(idubld(itrc))=.true.
2190 s(ng)%pioVar(idubld(itrc))%vd=var_desc(i)
2191 s(ng)%pioVar(idubld(itrc))%dkind=pio_fout
2192 s(ng)%pioVar(idubld(itrc))%gtype=u2dvar
2193 ELSE IF (trim(var_name(i)).eq. &
2194 & trim(vname(1,idvbld(itrc)))) THEN
2195 got_var(idvbld(itrc))=.true.
2196 s(ng)%pioVar(idvbld(itrc))%vd=var_desc(i)
2197 s(ng)%pioVar(idvbld(itrc))%dkind=pio_fout
2198 s(ng)%pioVar(idvbld(itrc))%gtype=v2dvar
2199# endif
2200 END IF
2201 END DO
2202 DO itrc=1,mbedp
2203 IF (trim(var_name(i)).eq.trim(vname(1,idsbed(itrc)))) THEN
2204 got_var(idsbed(itrc))=.true.
2205 s(ng)%pioVar(idsbed(itrc))%vd=var_desc(i)
2206 s(ng)%pioVar(idsbed(itrc))%dkind=pio_fout
2207 s(ng)%pioVar(idsbed(itrc))%gtype=b3dvar
2208 END IF
2209 END DO
2210# endif
2211# if defined SEDIMENT || defined BBL_MODEL
2212 DO itrc=1,mbotp
2213 IF (trim(var_name(i)).eq.trim(vname(1,idbott(itrc)))) THEN
2214 got_var(idbott(itrc))=.true.
2215 s(ng)%pioVar(idbott(itrc))%vd=var_desc(i)
2216 s(ng)%pioVar(idbott(itrc))%dkind=pio_fout
2217 s(ng)%pioVar(idbott(itrc))%gtype=r2dvar
2218 END IF
2219 END DO
2220# endif
2221 END DO
2222!
2223! Check if output variables are available in input NetCDF file.
2224!
2225 IF (.not.got_var(idtime)) THEN
2226 IF (master) WRITE (stdout,20) trim(vname(1,idtime)), &
2227 & trim(ncname)
2228 exit_flag=3
2229 RETURN
2230 END IF
2231# if defined SEDIMENT && defined SED_MORPH
2232 IF (.not.got_var(idbath).and.varout(idbath,ng)) THEN
2233 IF (master) WRITE (stdout,20) trim(vname(1,idbath)), &
2234 & trim(ncname)
2235 exit_flag=3
2236 RETURN
2237 END IF
2238# endif
2239# if defined SEDIMENT && defined BEDLOAD && defined BEDLOAD_VANDERA
2240 IF (.not.got_var(idsurs).and.varout(idsurs,ng)) THEN
2241 IF (master) WRITE (stdout,20) trim(vname(1,idsurs)), &
2242 & trim(ncname)
2243 exit_flag=3
2244 RETURN
2245 END IF
2246 IF (.not.got_var(idsrrw).and.varout(idsrrw,ng)) THEN
2247 IF (master) WRITE (stdout,20) trim(vname(1,idsrrw)), &
2248 & trim(ncname)
2249 exit_flag=3
2250 RETURN
2251 END IF
2252 IF (.not.got_var(idsbtw).and.varout(idsbtw,ng)) THEN
2253 IF (master) WRITE (stdout,20) trim(vname(1,idsbtw)), &
2254 & trim(ncname)
2255 exit_flag=3
2256 RETURN
2257 END IF
2258 IF (.not.got_var(idsucr).and.varout(idsucr,ng)) THEN
2259 IF (master) WRITE (stdout,20) trim(vname(1,idsucr)), &
2260 & trim(ncname)
2261 exit_flag=3
2262 RETURN
2263 END IF
2264 IF (.not.got_var(idsutr).and.varout(idsutr,ng)) THEN
2265 IF (master) WRITE (stdout,20) trim(vname(1,idsutr)), &
2266 & trim(ncname)
2267 exit_flag=3
2268 RETURN
2269 END IF
2270 IF (.not.got_var(idstcr).and.varout(idstcr,ng)) THEN
2271 IF (master) WRITE (stdout,20) trim(vname(1,idstcr)), &
2272 & trim(ncname)
2273 exit_flag=3
2274 RETURN
2275 END IF
2276 IF (.not.got_var(idsttr).and.varout(idsttr,ng)) THEN
2277 IF (master) WRITE (stdout,20) trim(vname(1,idsttr)), &
2278 & trim(ncname)
2279 exit_flag=3
2280 RETURN
2281 END IF
2282# endif
2283# ifdef SEDIMENT
2284 DO i=1,nst
2285 IF (.not.got_var(idfrac(i)).and.varout(idfrac(i),ng)) THEN
2286 IF (master) WRITE (stdout,20) trim(vname(1,idfrac(i))), &
2287 & trim(ncname)
2288 exit_flag=3
2289 RETURN
2290 END IF
2291
2292 IF(.not.got_var(idbmas(i)).and.varout(idbmas(i),ng)) THEN
2293 IF (master) WRITE (stdout,20) trim(vname(1,idbmas(i))), &
2294 & trim(ncname)
2295 exit_flag=3
2296 RETURN
2297 END IF
2298# ifdef BEDLOAD
2299 IF (.not.got_var(idubld(i)).and.varout(idubld(i),ng)) THEN
2300 IF (master) WRITE (stdout,20) trim(vname(1,idubld(i))), &
2301 & trim(ncname)
2302 exit_flag=3
2303 RETURN
2304 END IF
2305 IF (.not.got_var(idvbld(i)).and.varout(idvbld(i),ng)) THEN
2306 IF (master) WRITE (stdout,20) trim(vname(1,idvbld(i))), &
2307 & trim(ncname)
2308 exit_flag=3
2309 RETURN
2310 END IF
2311# endif
2312 END DO
2313 DO i=1,mbedp
2314 IF (.not.got_var(idsbed(i)).and.varout(idsbed(i),ng)) THEN
2315 IF (master) WRITE (stdout,20) trim(vname(1,idsbed(i))), &
2316 & trim(ncname)
2317 exit_flag=3
2318 RETURN
2319 END IF
2320 END DO
2321# endif
2322# if defined SEDIMENT || defined BBL_MODEL
2323 DO i=1,mbotp
2324 IF (.not.got_var(idbott(i)).and.varout(idbott(i),ng)) THEN
2325 IF (master) WRITE (stdout,20) trim(vname(1,idbott(i))), &
2326 & trim(ncname)
2327 exit_flag=3
2328 RETURN
2329 END IF
2330 END DO
2331# endif
2332 END IF query
2333!
2334 10 FORMAT (1pe11.4,1x,'millimeter')
2335 20 FORMAT (/,' SEDIMENT_DEF_PIO - unable to find variable: ', &
2336 & a,2x,' in output NetCDF file: ',a)
2337!
2338 RETURN
2339 END SUBROUTINE sediment_def_pio
2340!
2341!***********************************************************************
2342 SUBROUTINE sediment_wrt_pio (ng, model, tile, &
2343 & LBi, UBi, LBj, UBj, &
2344 & VarOut, S)
2345!***********************************************************************
2346!
2347 USE mod_pio_netcdf
2348!
2349! Imported variable declarations.
2350!
2351 logical, intent(in) :: varout(nv,ngrids)
2352!
2353 integer, intent(in) :: ng, model, tile
2354 integer, intent(in) :: lbi, ubi, lbj, ubj
2355!
2356 TYPE(t_io), intent(inout) :: s(ngrids)
2357!
2358! Local variable declarations.
2359!
2360 logical :: linstataneous
2361!
2362 integer :: status
2363!
2364 real(dp) :: scale
2365!
2366 character (len=*), parameter :: myfile = &
2367 & __FILE__//", sediment_wrt_pio"
2368!
2369 TYPE (io_desc_t), pointer :: iodesc
2370!
2371 sourcefile=myfile
2372!
2373!-----------------------------------------------------------------------
2374! Write out sediment output variables into specified NetCDF file.
2375!-----------------------------------------------------------------------
2376!
2377 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2378!
2379! Set instantaneous fields.
2380!
2381 IF ((s(ng)%ncid.eq.s(ng)%ncid).or. &
2382 & (s(ng)%ncid.eq.qck(ng)%ncid)) THEN
2383 linstataneous=.true.
2384 ELSE
2385 linstataneous=.false. ! time-averged fiels
2386 END IF
2387
2388# if defined SEDIMENT && defined SED_MORPH
2389!
2390! Write out time-dependent bathymetry (m)
2391!
2392 IF (varout(idbath,ng)) THEN
2393 scale=1.0_dp
2394 IF (s(ng)%pioVar(idbath)%dkind.eq.pio_double) THEN
2395 iodesc => iodesc_dp_r2dvar(ng)
2396 ELSE
2397 iodesc => iodesc_sp_r2dvar(ng)
2398 END IF
2399 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idbath, &
2400 & s(ng)%pioVar(idbath), &
2401 & s(ng)%Rindex, &
2402 & iodesc, &
2403 & lbi, ubi, lbj, ubj, scale, &
2404# ifdef MASKING
2405 & grid(ng) % rmask, &
2406# endif
2407 & grid(ng) % h, &
2408 & setfillval = .false.)
2409 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2410 IF (master) THEN
2411 WRITE (stdout,10) trim(vname(1,idbath)), s(ng)%Rindex
2412 END IF
2413 exit_flag=3
2414 ioerror=status
2415 RETURN
2416 END IF
2417 END IF
2418# endif
2419
2420# if defined SEDIMENT && bedload BEDLOAD
2421!
2422! Write out bed load transport in U-direction.
2423!
2424 DO i=1,nst
2425 IF (varout(idubld(i),ng)) THEN
2426 scale=1.0_dp
2427 IF (s(ng)%pioVar(idubld(i))%dkind.eq.pio_double) THEN
2428 iodesc => iodesc_dp_u2dvar(ng)
2429 ELSE
2430 iodesc => iodesc_sp_u2dvar(ng)
2431 END IF
2432 IF (linstataneous) THEN
2433 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idubld(i), &
2434 & s(ng)%pioVar(idubld(i)), &
2435 & s(ng)%Rindex, &
2436 & iodesc, &
2437 & lbi, ubi, lbj, ubj, scale, &
2438# ifdef MASKING
2439 & grid(ng) % umask, &
2440# endif
2441 & sedbed(ng) % bedldu(:,:,i))
2442# ifdef AVERAGES
2443 ELSE
2444 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idubld(i), &
2445 & s(ng)%pioVar(idubld(i)), &
2446 & s(ng)%Rindex, &
2447 & iodesc, &
2448 & lbi, ubi, lbj, ubj, scale, &
2449# ifdef MASKING
2450 & grid(ng) % umask, &
2451# endif
2452 & sedbed(ng) % avgbedldu(:,:,i))
2453# endif
2454 END IF
2455 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2456 IF (master) THEN
2457 WRITE (stdout,10) trim(vname(1,idubld(i))), s(ng)%Rindex
2458 END IF
2459 exit_flag=3
2460 ioerror=status
2461 RETURN
2462 END IF
2463 END IF
2464!
2465! Write out bed load transport in V-direction.
2466!
2467 IF (varout(idvbld(i),ng)) THEN
2468 scale=1.0_dp
2469 IF (s(ng)%pioVar(idvbld(i))%dkind.eq.pio_double) THEN
2470 iodesc => iodesc_dp_v2dvar(ng)
2471 ELSE
2472 iodesc => iodesc_sp_v2dvar(ng)
2473 END IF
2474 IF (linstataneous) THEN
2475 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idvbld(i), &
2476 & s(ng)%pioVar(idvbld(i)), &
2477 & s(ng)%Rindex, &
2478 & iodesc, &
2479 & lbi, ubi, lbj, ubj, scale, &
2480# ifdef MASKING
2481 & grid(ng) % vmask, &
2482# endif
2483 & sedbed(ng) % bedldv(:,:,i))
2484# ifdef AVERAGES
2485 ELSE
2486 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idvbld(i), &
2487 & s(ng)%pioVar(idvbld(i)), &
2488 & s(ng)%Rindex, &
2489 & iodesc, &
2490 & lbi, ubi, lbj, ubj, scale, &
2491# ifdef MASKING
2492 & grid(ng) % vmask, &
2493# endif
2494 & sedbed(ng) % avgbedldv(:,:,i))
2495# endif
2496 END IF
2497 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2498 IF (master) THEN
2499 WRITE (stdout,10) trim(vname(1,idvbld(i))), s(ng)%Rindex
2500 END IF
2501 exit_flag=3
2502 ioerror=status
2503 RETURN
2504 END IF
2505 END IF
2506 END DO
2507# endif
2508
2509# ifdef SEDIMENT
2510!
2511! Write out sediment fraction of each size class in each bed layer.
2512!
2513 DO i=1,nst
2514 IF (varout(idfrac(i),ng)) THEN
2515 scale=1.0_dp
2516 IF (s(ng)%pioVar(idfrac(i))%dkind.eq.pio_double) THEN
2517 iodesc => iodesc_dp_b3dvar(ng)
2518 ELSE
2519 iodesc => iodesc_sp_b3dvar(ng)
2520 END IF
2521 status=nf_fwrite3d(ng, model, s(ng)%pioFile, idfrac(i), &
2522 & s(ng)%pioVar(idfrac(i)), &
2523 & s(ng)%Rindex, &
2524 & iodesc, &
2525 & lbi, ubi, lbj, ubj, 1, nbed, scale, &
2526# ifdef MASKING
2527 & grid(ng) % rmask, &
2528# endif
2529 & sedbed(ng) % bed_frac(:,:,:,i))
2530 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2531 IF (master) THEN
2532 WRITE (stdout,10) trim(vname(1,idfrac(i))), s(ng)%Rindex
2533 END IF
2534 exit_flag=3
2535 ioerror=status
2536 RETURN
2537 END IF
2538 END IF
2539 END DO
2540!
2541! Write out sediment mass of each size class in each bed layer.
2542!
2543 DO i=1,nst
2544 IF (varout(idbmas(i),ng)) THEN
2545 scale=1.0_dp
2546 IF (s(ng)%pioVar(idbmas(i))%dkind.eq.pio_double) THEN
2547 iodesc => iodesc_dp_b3dvar(ng)
2548 ELSE
2549 iodesc => iodesc_sp_b3dvar(ng)
2550 END IF
2551 status=nf_fwrite3d(ng, model, s(ng)%pioFile, idbmas(i), &
2552 & s(ng)%pioVar(idbmas(i)), &
2553 & s(ng)%Rindex, &
2554 & iodesc, &
2555 & lbi, ubi, lbj, ubj, 1, nbed, scale, &
2556# ifdef MASKING
2557 & grid(ng) % rmask, &
2558# endif
2559 & sedbed(ng) % bed_mass(:,:,:,nout,i))
2560 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2561 IF (master) THEN
2562 WRITE (stdout,10) trim(vname(1,idbmas(i))), s(ng)%Rindex
2563 END IF
2564 exit_flag=3
2565 ioerror=status
2566 RETURN
2567 END IF
2568 END IF
2569 END DO
2570!
2571! Write out sediment properties in each bed layer.
2572!
2573 DO i=1,mbedp
2574 IF (varout(idsbed(i),ng)) THEN
2575 scale=1.0_dp
2576 IF (s(ng)%pioVar(idsbed(i))%dkind.eq.pio_double) THEN
2577 iodesc => iodesc_dp_b3dvar(ng)
2578 ELSE
2579 iodesc => iodesc_sp_b3dvar(ng)
2580 END IF
2581 status=nf_fwrite3d(ng, model, s(ng)%pioFile, idsbed(i), &
2582 & s(ng)%pioVar(idsbed(i)), &
2583 & s(ng)%Rindex, &
2584 & iodesc, &
2585 & lbi, ubi, lbj, ubj, 1, nbed, scale, &
2586# ifdef MASKING
2587 & grid(ng) % rmask, &
2588# endif
2589 & sedbed(ng) % bed(:,:,:,i))
2590 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2591 IF (master) THEN
2592 WRITE (stdout,10) trim(vname(1,idsbed(i))), s(ng)%Rindex
2593 END IF
2594 exit_flag=3
2595 ioerror=status
2596 RETURN
2597 END IF
2598 END IF
2599 END DO
2600# endif
2601
2602# if defined SEDIMENT || defined BBL_MODEL
2603!
2604! Write out exposed sediment layer properties.
2605!
2606 DO i=1,mbotp
2607 IF (varout(idbott(i),ng)) THEN
2608 IF (i.eq.itauc) THEN
2609 scale=rho0
2610 ELSE
2611 scale=1.0_dp
2612 END IF
2613 IF (s(ng)%pioVar(idbott(i))%dkind.eq.pio_double) THEN
2614 iodesc => iodesc_dp_r2dvar(ng)
2615 ELSE
2616 iodesc => iodesc_sp_r2dvar(ng)
2617 END IF
2618 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idbott(i), &
2619 & s(ng)%pioVar(idbott(i)), &
2620 & s(ng)%Rindex, &
2621 & iodesc, &
2622 & lbi, ubi, lbj, ubj, scale, &
2623# ifdef MASKING
2624 & grid(ng) % rmask, &
2625# endif
2626 & sedbed(ng) % bottom(:,:,i))
2627 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2628 IF (master) THEN
2629 WRITE (stdout,10) trim(vname(1,idbott(i))), s(ng)%Rindex
2630 END IF
2631 exit_flag=3
2632 ioerror=status
2633 RETURN
2634 END IF
2635 END IF
2636 END DO
2637# endif
2638
2639# if defined SEDIMENT && defined BEDLOAD && defined SED_BEDLOAD_VANDERA
2640!
2641! Write out Ursell number of the asymmetric wave form.
2642!
2643 IF (varout(idsurs,ng)) THEN
2644 scale=1.0_dp
2645 IF (s(ng)%pioVar(idsurs)%dkind.eq.pio_double) THEN
2646 iodesc => iodesc_dp_r2dvar(ng)
2647 ELSE
2648 iodesc => iodesc_sp_r2dvar(ng)
2649 END IF
2650 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idsurs, &
2651 & s(ng)%pioVar(idsurs), &
2652 & s(ng)%Rindex, &
2653 & iodesc, &
2654 & lbi, ubi, lbj, ubj, scale, &
2655# ifdef MASKING
2656 & grid(ng) % rmask, &
2657# endif
2658 & sedbed(ng) % ursell_no)
2659 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2660 IF (master) THEN
2661 WRITE (stdout,10) trim(vname(1,idsurs)), s(ng)%Rindex
2662 END IF
2663 exit_flag=3
2664 ioerror=status
2665 RETURN
2666 END IF
2667 END IF
2668!
2669! Write out velocity skewness parameter of the asymmetric wave form.
2670!
2671 IF (varout(idsrrw,ng)) THEN
2672 scale=1.0_dp
2673 IF (s(ng)%pioVar(idsrrw)%dkind.eq.pio_double) THEN
2674 iodesc => iodesc_dp_r2dvar(ng)
2675 ELSE
2676 iodesc => iodesc_sp_r2dvar(ng)
2677 END IF
2678 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idsrrw, &
2679 & s(ng)%pioVar(idsrrw), &
2680 & s(ng)%Rindex, &
2681 & iodesc, &
2682 & lbi, ubi, lbj, ubj, scale, &
2683# ifdef MASKING
2684 & grid(ng) % rmask, &
2685# endif
2686 & sedbed(ng) % RR_asymwave)
2687 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2688 IF (master) THEN
2689 WRITE (stdout,10) trim(vname(1,idsrrw)), s(ng)%Rindex
2690 END IF
2691 exit_flag=3
2692 ioerror=status
2693 RETURN
2694 END IF
2695 END IF
2696!
2697! Write out acceleration asymmetry parameter of the asymmetric wave form.
2698!
2699 IF (varout(idsbtw,ng)) THEN
2700 scale=1.0_dp
2701 IF (s(ng)%pioVar(idsbtw)%dkind.eq.pio_double) THEN
2702 iodesc => iodesc_dp_r2dvar(ng)
2703 ELSE
2704 iodesc => iodesc_sp_r2dvar(ng)
2705 END IF
2706 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idsbtw, &
2707 & s(ng)%pioVar(idsbtw), &
2708 & s(ng)%Rindex, &
2709 & iodesc, &
2710 & lbi, ubi, lbj, ubj, scale, &
2711# ifdef MASKING
2712 & grid(ng) % rmask, &
2713# endif
2714 & sedbed(ng) % beta_asymwave)
2715 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2716 IF (master) THEN
2717 WRITE (stdout,10) trim(vname(1,idsbtw)), s(ng)%Rindex
2718 END IF
2719 exit_flag=3
2720 ioerror=status
2721 RETURN
2722 END IF
2723 END IF
2724!
2725! Write out crest velocity of the asymmetric wave form.
2726!
2727 IF (varout(idsucr,ng)) THEN
2728 scale=1.0_dp
2729 IF (s(ng)%pioVar(idsucr)%dkind.eq.pio_double) THEN
2730 iodesc => iodesc_dp_r2dvar(ng)
2731 ELSE
2732 iodesc => iodesc_sp_r2dvar(ng)
2733 END IF
2734 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idsucr, &
2735 & s(ng)%pioVar(idsucr), &
2736 & s(ng)%Rindex, &
2737 & iodesc, &
2738 & lbi, ubi, lbj, ubj, scale, &
2739# ifdef MASKING
2740 & grid(ng) % rmask, &
2741# endif
2742 & sedbed(ng) % ucrest_r)
2743 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2744 IF (master) THEN
2745 WRITE (stdout,10) trim(vname(1,idsucr)), s(ng)%Rindex
2746 END IF
2747 exit_flag=3
2748 ioerror=status
2749 RETURN
2750 END IF
2751 END IF
2752!
2753! Write out trough velocity of the asymmetric wave form.
2754!
2755 IF (varout(idsutr,ng)) THEN
2756 scale=1.0_dp
2757 IF (s(ng)%pioVar(idsutr)%dkind.eq.pio_double) THEN
2758 iodesc => iodesc_dp_r2dvar(ng)
2759 ELSE
2760 iodesc => iodesc_sp_r2dvar(ng)
2761 END IF
2762 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idsutr, &
2763 & s(ng)%pioVar(idsutr), &
2764 & s(ng)%Rindex, &
2765 & iodesc, &
2766 & lbi, ubi, lbj, ubj, scale, &
2767# ifdef MASKING
2768 & grid(ng) % rmask, &
2769# endif
2770 & sedbed(ng) % utrough_r)
2771 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2772 IF (master) THEN
2773 WRITE (stdout,10) trim(vname(1,idsutr)), s(ng)%Rindex
2774 END IF
2775 exit_flag=3
2776 ioerror=status
2777 RETURN
2778 END IF
2779 END IF
2780!
2781! Write out crest time period of the asymmetric wave form.
2782!
2783 IF (varout(idstcr,ng)) THEN
2784 scale=1.0_dp
2785 IF (s(ng)%pioVar(idstcr)%dkind.eq.pio_double) THEN
2786 iodesc => iodesc_dp_r2dvar(ng)
2787 ELSE
2788 iodesc => iodesc_sp_r2dvar(ng)
2789 END IF
2790 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idstcr, &
2791 & s(ng)%pioVar(idstcr), &
2792 & s(ng)%Rindex, &
2793 & iodesc, &
2794 & lbi, ubi, lbj, ubj, scale, &
2795# ifdef MASKING
2796 & grid(ng) % rmask, &
2797# endif
2798 & sedbed(ng) % T_crest)
2799 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2800 IF (master) THEN
2801 WRITE (stdout,10) trim(vname(1,idstcr)), s(ng)%Rindex
2802 END IF
2803 exit_flag=3
2804 ioerror=status
2805 RETURN
2806 END IF
2807 END IF
2808!
2809! Write out trough time period of the asymmetric wave form.
2810!
2811 IF (varout(idsttr,ng)) THEN
2812 scale=1.0_dp
2813 IF (s(ng)%pioVar(idsttr)%dkind.eq.pio_double) THEN
2814 iodesc => iodesc_dp_r2dvar(ng)
2815 ELSE
2816 iodesc => iodesc_sp_r2dvar(ng)
2817 END IF
2818 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idsttr, &
2819 & s(ng)%pioVar(idsttr), &
2820 & s(ng)%Rindex, &
2821 & iodesc, &
2822 & lbi, ubi, lbj, ubj, scale, &
2823# ifdef MASKING
2824 & grid(ng) % rmask, &
2825# endif
2826 & sedbed(ng) % T_trough)
2827 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
2828 IF (master) THEN
2829 WRITE (stdout,10) trim(vname(1,idsttr)), s(ng)%Rindex
2830 END IF
2831 exit_flag=3
2832 ioerror=status
2833 RETURN
2834 END IF
2835 END IF
2836# endif
2837!
2838 10 FORMAT (/," SEDIMENT_WRT_PIO - error while writing variable '", &
2839 & a,"', time record = ",i0,/,11x,'into file: ',a)
2840!
2841 RETURN
2842 END SUBROUTINE sediement_wrt_pio
2843
2844# ifdef STATIONS
2845!
2846!***********************************************************************
2847 SUBROUTINE sediment_def_station_pio (ng, model, ldef, VarOut, S, &
2848 & bgrd, pgrd, rgrd)
2849!***********************************************************************
2850!
2851 USE mod_netcdf
2852!
2853! Imported variable declarations.
2854!
2855 logical, intent(in) :: ldef, varout(nv,ngrids)
2856!
2857 integer, intent(in) :: ng, model
2858 integer, intent(in), optional :: bgrd(:), pgrd(:), rgrd(:)
2859!
2860 TYPE(t_io), intent(inout) :: s(ngrids)
2861!
2862! Local variable declarations.
2863!
2864 logical :: got_var(nv)
2865!
2866 integer, parameter :: natt = 25
2867
2868 integer :: i, itrc, j, status
2869!
2870 real(r8) :: aval(6)
2871!
2872 character (len=120) :: vinfo(natt)
2873 character (len=256) :: ncname
2874!
2875 character (len=*), parameter :: myfile = &
2876 & __FILE__//", sediment_def_station_nf90"
2877!
2878 sourcefile=myfile
2879!
2880!-----------------------------------------------------------------------
2881! Define sediment output stations variables.
2882!-----------------------------------------------------------------------
2883!
2884 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2885 ncname=s(ng)%name
2886!
2887 define : IF (ldef) THEN
2888!
2889! Initialize local information variable arrays.
2890!
2891 DO i=1,natt
2892 DO j=1,len(vinfo(1))
2893 vinfo(i)(j:j)=' '
2894 END DO
2895 END DO
2896 DO i=1,6
2897 aval(i)=0.0_r8
2898 END DO
2899
2900# if defined SEDIMENT && defined SED_MORPH
2901!
2902! Define time-varying bathymetry.
2903!
2904 IF (varout(idbath,ng)) THEN
2905 vinfo( 1)=vname(1,idbath)
2906 vinfo( 2)=vname(2,idbath)
2907 vinfo( 3)=vname(3,idbath)
2908 vinfo(14)=vname(4,idbath)
2909 vinfo(16)=vname(1,idtime)
2910 s(ng)%pioVar(idbath)%dkind=pio_fout
2911 s(ng)%pioVar(idbath)%gtype=0
2912!
2913 status=def_var(ng, model, s(ng)%pioFile, &
2914 & s(ng)%pioVar(idbath)%vd, &
2915 & pio_fout, 2, pgrd, aval, vinfo, ncname, &
2916 & setfillval = .false., &
2917 & setparaccess = .true.)
2918 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2919 END IF
2920# endif
2921
2922# ifdef SEDIMENT
2923!
2924! Define sediment fraction of each size class in each bed layer.
2925!
2926 DO i=1,nst
2927 IF (varout(idfrac(i),ng)) THEN
2928 vinfo( 1)=vname(1,idfrac(i))
2929 vinfo( 2)=vname(2,idfrac(i))
2930 vinfo( 3)=vname(3,idfrac(i))
2931 vinfo(14)=vname(4,idfrac(i))
2932 vinfo(16)=vname(1,idtime)
2933 WRITE (vinfo(19),10) 1000.0_r8*sd50(i,ng)
2934 s(ng)%pioVar(idfrac(i))%dkind=pio_fout
2935 s(ng)%pioVar(idfrac(i))%gtype=0
2936!
2937 status=def_var(ng, model, s(ng)%pioFile, &
2938 & s(ng)%pioVar(idfrac(i))%vd, &
2939 & pio_fout, 3, bgrd, aval, vinfo, ncname, &
2940 & setfillval = .true., &
2941 & setparaccess = .true.)
2942 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2943 END IF
2944!
2945! Define sediment mass of each size class in each bed layer.
2946!
2947 IF (varout(idbmas(i),ng)) THEN
2948 vinfo( 1)=vname(1,idbmas(i))
2949 vinfo( 2)=vname(2,idbmas(i))
2950 vinfo( 3)=vname(3,idbmas(i))
2951 vinfo(14)=vname(4,idbmas(i))
2952 vinfo(16)=vname(1,idtime)
2953 WRITE (vinfo(19),10) 1000.0_r8*sd50(i,ng)
2954 s(ng)%pioVar(idbmas(i))%dkind=pio_fout
2955 s(ng)%pioVar(idbmas(i))%gtype=0
2956!
2957 status=def_var(ng, model, s(ng)%pioFile, &
2958 & s(ng)%pioVar(idbmas(i))%vd, &
2959 & pio_fout, 3, bgrd, aval, vinfo, ncname, &
2960 & setfillval = .true., &
2961 & setparaccess = .true.)
2962 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2963 END IF
2964 END DO
2965!
2966! Define sediment properties in each bed layer.
2967!
2968 DO i=1,mbedp
2969 IF (varout(idsbed(i),ng)) THEN
2970 vinfo( 1)=vname(1,idsbed(i))
2971 vinfo( 2)=vname(2,idsbed(i))
2972 vinfo( 3)=vname(3,idsbed(i))
2973 vinfo(14)=vname(4,idsbed(i))
2974 vinfo(16)=vname(1,idtime)
2975 s(ng)%pioVar(idsbed(i))%dkind=pio_fout
2976 s(ng)%pioVar(idsbed(i))%gtype=0
2977!
2978 status=def_var(ng, model, s(ng)%pioFile, &
2979 & s(ng)%pioVar(idsbed(i))%vd, &
2980 & pio_fout, 3, bgrd, aval, vinfo, ncname, &
2981 & setfillval = .true., &
2982 & setparaccess = .true.)
2983 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
2984 END IF
2985 END DO
2986# endif
2987
2988# if defined SEDIMENT || defined BBL_MODEL
2989!
2990! Define exposed sediment layer properties.
2991!
2992 DO i=1,mbotp
2993 IF (varout(idbott(i),ng)) THEN
2994 vinfo( 1)=vname(1,idbott(i))
2995 vinfo( 2)=vname(2,idbott(i))
2996 vinfo( 3)=vname(3,idbott(i))
2997 vinfo(14)=vname(4,idbott(i))
2998 vinfo(16)=vname(1,idtime)
2999 s(ng)%pioVar(idbott(i))%dkind=pio_fout
3000 s(ng)%pioVar(idbott(i))%gtype=0
3001!
3002 status=def_var(ng, model, s(ng)%pioFile, &
3003 & s(ng)%pioVar(idbott(i))%vd, &
3004 & pio_fout, 2, pgrd, aval, vinfo, ncname, &
3005 & setfillval = .true., &
3006 & setparaccess = .true.)
3007 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
3008 END IF
3009 END DO
3010# endif
3011
3012 END IF define
3013!
3014!-----------------------------------------------------------------------
3015! Open an existing stations file, check its contents, and prepare for
3016! appending data.
3017!-----------------------------------------------------------------------
3018!
3019 query : IF (.not.ldef) THEN
3020!
3021! Initialize logical switches.
3022!
3023 DO i=1,nv
3024 got_var(i)=.false.
3025 END DO
3026!
3027! Scan variable list from input NetCDF and activate switches for
3028! stations variables. Get variable IDs.
3029!
3030 DO i=1,n_var
3031 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
3032 got_var(idtime)=.true.
3033 s(ng)%pioVar(idtime)%vd=var_desc(i)
3034 s(ng)%pioVar(idtime)%dkind=pio_tout
3035 s(ng)%pioVar(idtime)%gtype=0
3036# if defined SEDIMENT && defined SED_MORPH
3037 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idbath))) THEN
3038 got_var(idbath)=.true.
3039 s(ng)%pioVar(idbath)%vd=var_desc(i)
3040 s(ng)%pioVar(idbath)%dkind=pio_fout
3041 s(ng)%pioVar(idbath)%gtype=0
3042# endif
3043 END IF
3044# ifdef SEDIMENT
3045 DO itrc=1,nst
3046 IF (trim(var_name(i)).eq.trim(vname(1,idfrac(itrc)))) THEN
3047 got_var(idfrac(itrc))=.true.
3048 s(ng)%pioVar(idfrac(itrc))%vd=var_desc(i)
3049 s(ng)%pioVar(idfrac(itrc))%dkind=pio_fout
3050 s(ng)%pioVar(idfrac(itrc))%gtype=0
3051 ELSE IF (trim(var_name(i)).eq. &
3052 & trim(vname(1,idbmas(itrc)))) THEN
3053 got_var(idbmas(itrc))=.true.
3054 s(ng)%pioVar(idbmas(itrc))%vd=var_desc(i)
3055 s(ng)%pioVar(idbmas(itrc))%dkind=pio_fout
3056 s(ng)%pioVar(idbmas(itrc))%gtype=0
3057 END IF
3058 END DO
3059 DO itrc=1,mbedp
3060 IF (trim(var_name(i)).eq.trim(vname(1,idsbed(itrc)))) THEN
3061 got_var(idsbed(itrc))=.true.
3062 s(ng)%pioVar(idsbed(itrc))%vd=var_desc(i)
3063 s(ng)%pioVar(idsbed(itrc))%dkind=pio_fout
3064 s(ng)%pioVar(idsbed(itrc))%gtype=0
3065 END IF
3066 END DO
3067# endif
3068# if defined SEDIMENT || defined BBL_MODEL
3069 DO itrc=1,mbotp
3070 IF (trim(var_name(i)).eq.trim(vname(1,idbott(itrc)))) THEN
3071 got_var(idbott(itrc))=.true.
3072 s(ng)%pioVar(idbott(itrc))%vd=var_desc(i)
3073 s(ng)%pioVar(idbott(itrc))%dkind=pio_fout
3074 s(ng)%pioVar(idbott(itrc))%gtype=0
3075 END IF
3076 END DO
3077# endif
3078 END DO
3079!
3080! Check if station variables are available in input NetCDF file.
3081!
3082 IF (.not.got_var(idtime)) THEN
3083 IF (master) WRITE (stdout,20) trim(vname(1,idtime)), &
3084 & trim(ncname)
3085 exit_flag=3
3086 RETURN
3087 END IF
3088# if defined SEDIMENT && defined SED_MORPH
3089 IF (.not.got_var(idbath).and.varout(idbath,ng)) THEN
3090 IF (master) WRITE (stdout,20) trim(vname(1,idbath)), &
3091 & trim(ncname)
3092 exit_flag=3
3093 RETURN
3094 END IF
3095# endif
3096# ifdef SEDIMENT
3097 DO i=1,nst
3098 IF (.not.got_var(idfrac(i)).and.varout(idfrac(i),ng)) THEN
3099 IF (master) WRITE (stdout,20) trim(vname(1,idfrac(i))), &
3100 & trim(ncname)
3101 exit_flag=3
3102 RETURN
3103 END IF
3104 IF (.not.got_var(idbmas(i)).and.varout(idbmas(i),ng)) THEN
3105 IF (master) WRITE (stdout,20) trim(vname(1,idbmas(i))), &
3106 & trim(ncname)
3107 exit_flag=3
3108 RETURN
3109 END IF
3110 END DO
3111 DO i=1,mbedp
3112 IF (.not.got_var(idsbed(i)).and.varout(idsbed(i),ng)) THEN
3113 IF (master) WRITE (stdout,20) trim(vname(1,idsbed(i))), &
3114 & trim(ncname)
3115 exit_flag=3
3116 RETURN
3117 END IF
3118 END DO
3119# endif
3120# if defined SEDIMENT || defined BBL_MODEL
3121 DO i=1,mbotp
3122 IF (.not.got_var(idbott(i)).and.varout(idbott(i),ng)) THEN
3123 IF (master) WRITE (stdout,20) trim(vname(1,idbott(i))), &
3124 & trim(ncname)
3125 exit_flag=3
3126 RETURN
3127 END IF
3128 END DO
3129# endif
3130!
3131 10 FORMAT (1pe11.4,1x,'millimeter')
3132 20 FORMAT (/,' SEDIMENT_DEF_STATION_PIO - unable to find variable:', &
3133 & 1x,a,2x,' in stations NetCDF file: ',a)
3134!
3135 RETURN
3136 END SUBROUTINE sediment_def_station_pio
3137!
3138!***********************************************************************
3139 SUBROUTINE sediment_wrt_station_pio (ng, model, tile, &
3140 & LBi, UBi, LBj, UBj, &
3141 & VarOut, S)
3142!***********************************************************************
3143!
3144 USE mod_pio_netcdf
3145!
3146! Imported variable declarations.
3147!
3148 logical, intent(in) :: varout(nv,ngrids)
3149!
3150 integer, intent(in) :: ng, model, tile
3151 integer, intent(in) :: lbi, ubi, lbj, ubj
3152!
3153 TYPE(t_io), intent(inout) :: s(ngrids)
3154!
3155! Local variable declarations.
3156!
3157 logical :: cgrid
3158!
3159 integer :: nposb
3160 integer :: i, k, np, status
3161!
3162 real(dp) :: scale
3163!
3164 real(r8), dimension(Nstation(ng)) :: xpos, ypos, zpos, psta
3165
3166# ifdef SEDIMENT
3167 real(r8), dimension(Nstation(ng)*Nbed) :: xposb, yposb, zposb
3168 real(r8), dimension(Nstation(ng)*Nbed) :: bsta
3169# endif
3170!
3171 character (len=*), parameter :: myfile = &
3172 & __FILE__//", sediment_wrt_station_pio"
3173!
3174 sourcefile=myfile
3175!
3176!-----------------------------------------------------------------------
3177! Write out sediment output variables into specified stations NetCDF
3178! file.
3179!-----------------------------------------------------------------------
3180!
3181 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
3182!
3183! Set switch to extract station data at native C-grid position (TRUE)
3184! or at RHO-points (FALSE).
3185!
3186# ifdef STATIONS_CGRID
3187 cgrid=.true.
3188# else
3189 cgrid=.false.
3190# endif
3191!
3192! Set positions for generic extraction routine.
3193!
3194 nposb=nstation(ng)*nbed
3195 DO i=1,nstation(ng)
3196 xpos(i)=scalars(ng)%SposX(i)
3197 ypos(i)=scalars(ng)%SposY(i)
3198 zpos(i)=1.0_r8
3199# ifdef SEDIMENT
3200 DO k=1,nbed
3201 np=k+(i-1)*nbed
3202 xposb(np)=scalars(ng)%SposX(i)
3203 yposb(np)=scalars(ng)%SposY(i)
3204 zposb(np)=real(k,r8)
3205 END DO
3206# endif
3207 END DO
3208
3209# if defined SEDIMENT && defined SED_MORPH
3210!
3211! Write out time-varying bathymetry.
3212!
3213 IF (varout(idbath,ng)) THEN
3214 scale=1.0_dp
3215 CALL extract_sta2d (ng, model, cgrid, idbath, r2dvar, &
3216 & lbi, ubi, lbj, ubj, &
3217 & scale, grid(ng)%h, &
3218 & nstation(ng), xpos, ypos, psta)
3219 CALL pio_netcdf_put_fvar (ng, model, s(ng)%name, &
3220 & trim(vname(1,idbath)), psta, &
3221 & (/1,s(ng)%Rindex/), &
3222 & (/nstation(ng),1/), &
3223 & piofile = s(ng)%pioFile, &
3224 & piovar = s(ng)%pioVar(idbath)%vd)
3225 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
3226 END IF
3227# endif
3228
3229# ifdef SEDIMENT
3230!
3231! Write out sediment fraction of each size class in each bed layer.
3232!
3233 DO i=1,nst
3234 IF (varout(idfrac(i),ng)) THEN
3235 scale=1.0_dp
3236 CALL extract_sta3d (ng, model, cgrid, idfrac(i), b3dvar, &
3237 & lbi, ubi, lbj, ubj, 1, nbed, &
3238 & scale, sedbed(ng)%bed_frac(:,:,:,i), &
3239 & nposb, xposb, yposb, zposb, bsta)
3240 CALL pio_netcdf_put_fvar (ng, model, s(ng)%name, &
3241 & trim(vname(1,idfrac(i))), bsta, &
3242 & (/1,1,s(ng)%Rindex/), &
3243 & (/nbed,nstation(ng),1/), &
3244 & piofile = s(ng)%pioFile, &
3245 & piovar = s(ng)%pioVar(idfrac(i))%vd)
3246 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
3247 END IF
3248!
3249! Write out sediment mass of each size class in each bed layer.
3250!
3251 IF (varout(idbmas(i),ng)) THEN
3252 scale=1.0_dp
3253 CALL extract_sta3d (ng, model, cgrid, idbmas(i), b3dvar, &
3254 & lbi, ubi, lbj, ubj, 1, nbed, &
3255 & scale, &
3256 & sedbed(ng)%bed_mass(:,:,:,nout,i), &
3257 & nposb, xposb, yposb, zposb, bsta)
3258 CALL pio_netcdf_put_fvar (ng, model, s(ng)%name, &
3259 & trim(vname(1,idbmas(i))), bsta, &
3260 & (/1,1,s(ng)%Rindex/), &
3261 & (/nbed,nstation(ng),1/), &
3262 & piofile = s(ng)%pioFile, &
3263 & piovar = s(ng)%pioVar(idbmas(i))%vd)
3264 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
3265 END IF
3266 END DO
3267!
3268! Write out sediment properties in each bed layer.
3269!
3270 DO i=1,mbedp
3271 IF (varout(idsbed(i),ng)) THEN
3272 scale=1.0_dp
3273 CALL extract_sta3d (ng, model, cgrid, idsbed(i), b3dvar, &
3274 & lbi, ubi, lbj, ubj, 1, nbed, &
3275 & scale, sedbed(ng)%bed(:,:,:,i), &
3276 & nposb, xposb, yposb, zposb, bsta)
3277 CALL pio_netcdf_put_fvar (ng, model, s(ng)%name, &
3278 & trim(vname(1,idsbed(i))), bsta, &
3279 & (/1,1,s(ng)%Rindex/), &
3280 & (/nbed,nstation(ng),1/), &
3281 & piofile = s(ng)%pioFile, &
3282 & piovar = s(ng)%pioVar(idsbed(i))%vd)
3283 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
3284 END IF
3285 END DO
3286# endif
3287
3288# if defined SEDIMENT || defined BBL_MODEL
3289!
3290! Write out exposed sediment layer properties.
3291!
3292 DO i=1,mbedp
3293 IF (varout(idbott(i),ng)) THEN
3294 scale=1.0_dp
3295 CALL extract_sta2d (ng, model, cgrid, idbott(i), r2dvar, &
3296 & lbi, ubi, lbj, ubj, &
3297 & scale, sedbed(ng)%bottom(:,:,i), &
3298 & nstation(ng), xpos, ypos, psta)
3299 CALL pio_netcdf_put_fvar (ng, model, s(ng)%name, &
3300 & trim(vname(1,idbott(i))), bsta, &
3301 & (/1,s(ng)%Rindex/), &
3302 & (/nstation(ng),1/), &
3303 & piofile = s(ng)%pioFile, &
3304 & piovar = s(ng)%pioVar(idbott(i))%vd)
3305 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
3306 END IF
3307 END DO
3308# endif
3309!
3310 RETURN
3311 END SUBROUTINE sediment_wrt_station_pio
3312# endif
3313# endif
3314#endif
3315!
3316 END MODULE sediment_output_mod
subroutine extract_sta2d(ng, model, cgrid, ifield, gtype, lbi, ubi, lbj, ubj, ascl, a, npos, xpos, ypos, apos)
Definition extract_sta.F:65
subroutine extract_sta3d(ng, model, cgrid, ifield, gtype, lbi, ubi, lbj, ubj, lbk, ubk, ascl, a, npos, xpos, ypos, zpos, apos)
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer ioerror
type(t_io), dimension(:), allocatable qck
type(t_io), dimension(:), allocatable avg
integer stdout
character(len=256) sourcefile
integer, parameter nv
integer idbath
character(len=maxlen), dimension(6, 0:nv) vname
integer idtime
integer, dimension(:,:,:), allocatable iinfo
integer, parameter nf_fout
Definition mod_netcdf.F:188
character(len=100), dimension(mvars) var_name
Definition mod_netcdf.F:169
integer, dimension(mvars) var_id
Definition mod_netcdf.F:160
integer n_var
Definition mod_netcdf.F:152
logical master
integer nbed
Definition mod_param.F:517
integer, parameter b3dvar
Definition mod_param.F:725
integer, dimension(:), allocatable nstation
Definition mod_param.F:557
integer, parameter u2dvar
Definition mod_param.F:718
integer ngrids
Definition mod_param.F:113
integer, parameter r2dvar
Definition mod_param.F:717
integer, parameter v2dvar
Definition mod_param.F:719
integer nst
Definition mod_param.F:521
type(io_desc_t), dimension(:), pointer iodesc_dp_b3dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dvar
integer, parameter pio_fout
type(var_desc_t), dimension(:), pointer var_desc
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_b3dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dvar
integer, parameter pio_tout
type(io_desc_t), dimension(:), pointer iodesc_dp_v2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dvar
integer exit_flag
type(t_scalars), dimension(:), allocatable scalars
Definition mod_scalars.F:65
real(dp) rho0
integer noerror
type(t_sedbed), dimension(:), allocatable sedbed
Definition sedbed_mod.h:157
integer, dimension(:), allocatable idubld
integer, parameter mbotp
integer, dimension(mbedp) idsbed
integer, dimension(:), allocatable idfrac
integer, dimension(mbotp) idbott
real(r8), dimension(:,:), allocatable sd50
integer, dimension(:), allocatable idvbld
integer, dimension(:), allocatable idbmas
integer, parameter itauc
integer, parameter mbedp
subroutine, public scale_omega(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, pm, pn, w, wscl)
Definition omega.F:382
subroutine, public sediment_def_station_pio(ng, model, ldef, varout, s, bgrd, pgrd, rgrd)
subroutine, public sediment_wrt_nf90(ng, model, tile, lbi, ubi, lbj, ubj, varout, s)
subroutine, public sediment_wrt_station_nf90(ng, model, tile, lbi, ubi, lbj, ubj, varout, s)
subroutine, public sediment_def_station_nf90(ng, model, ldef, varout, s, bgrd, pgrd, rgrd)
subroutine, public sediment_wrt_station_pio(ng, model, tile, lbi, ubi, lbj, ubj, varout, s)
subroutine, public sediment_def_pio(ng, model, ldef, varout, s, t2dgrd, u2dgrd, v2dgrd, t3dgrd, u3dgrd, v3dgrd, w3dgrd)
subroutine, public sediment_wrt_pio(ng, model, tile, lbi, ubi, lbj, ubj, varout, s)
subroutine, public sediment_def_nf90(ng, model, ldef, varout, s, t2dgrd, u2dgrd, v2dgrd, b3dgrd)
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52