ROMS
Loading...
Searching...
No Matches
wrt_station.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#ifdef STATIONS
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This module writes STATION data into output file using the standard !
13! NetCDF library or the Parallel-IO (PIO) library. !
14! !
15!=======================================================================
16!
17 USE mod_param
18 USE mod_parallel
19# ifdef BBL_MODEL
20 USE mod_bbl
21# endif
22 USE mod_forces
23 USE mod_grid
24 USE mod_iounits
25 USE mod_mixing
26 USE mod_ncparam
27 USE mod_ocean
28 USE mod_scalars
29# if defined SEDIMENT || defined BBL_MODEL
30 USE mod_sedbed
31 USE mod_sediment
32# endif
33 USE mod_stepping
34!
35# if defined BBL_MODEL || defined WAVES_OUTPUT
37# if defined PIO_LIB && defined DISTRIBUTE
39# endif
40# endif
42# ifdef SOLVE3D
44# endif
45# ifdef ICE_MODEL
46 USE ice_output_mod, ONLY : ice_wrt_station_nf90
47# if defined PIO_LIB && defined DISTRIBUTE
48 USE ice_output_mod, ONLY : ice_wrt_station_pio
49# endif
50# endif
51# ifdef SEDIMENT
53# if defined PIO_LIB && defined DISTRIBUTE
55# endif
56# endif
57 USE uv_rotate_mod, ONLY : uv_rotate2d
58# ifdef SOLVE3D
59 USE uv_rotate_mod, ONLY : uv_rotate3d
60# endif
61 USE strings_mod, ONLY : founderror
62# if defined WEC || defined WEC_VF
63 USE wec_output_mod, ONLY : wec_wrt_station_nf90
64# if defined PIO_LIB && defined DISTRIBUTE
65 USE wec_output_mod, ONLY : wec_wrt_station_pio
66# endif
67# endif
68!
69 implicit none
70!
71 PUBLIC :: wrt_station
72 PRIVATE :: wrt_station_nf90
73# if defined PIO_LIB && defined DISTRIBUTE
74 PRIVATE :: wrt_station_pio
75# endif
76!
77 CONTAINS
78!
79!***********************************************************************
80 SUBROUTINE wrt_station (ng, tile)
81!***********************************************************************
82!
83! Imported variable declarations.
84!
85 integer, intent(in) :: ng, tile
86!
87! Local variable declarations.
88!
89 integer :: lbi, ubi, lbj, ubj
90!
91 character (len=*), parameter :: myfile = &
92 & __FILE__
93!
94!-----------------------------------------------------------------------
95! Write out history fields according to IO type.
96!-----------------------------------------------------------------------
97!
98 lbi=bounds(ng)%LBi(tile)
99 ubi=bounds(ng)%UBi(tile)
100 lbj=bounds(ng)%LBj(tile)
101 ubj=bounds(ng)%UBj(tile)
102!
103 SELECT CASE (sta(ng)%IOtype)
104 CASE (io_nf90)
105 CALL wrt_station_nf90 (ng, inlm, tile, &
106 & lbi, ubi, lbj, ubj)
107
108# if defined PIO_LIB && defined DISTRIBUTE
109 CASE (io_pio)
110 CALL wrt_station_pio (ng, inlm, tile, &
111 & lbi, ubi, lbj, ubj)
112# endif
113 CASE DEFAULT
114 IF (master) WRITE (stdout,10) sta(ng)%IOtype
115 exit_flag=3
116 END SELECT
117 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
118!
119 10 FORMAT (' WRT_STATION - Illegal output file type, io_type = ',i0, &
120 & /,15x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.')
121!
122 RETURN
123 END SUBROUTINE wrt_station
124!
125!***********************************************************************
126 SUBROUTINE wrt_station_nf90 (ng, model, tile, &
127 & LBi, UBi, LBj, UBj)
128!***********************************************************************
129!
130 USE mod_netcdf
131!
132! Imported variable declarations.
133!
134 integer, intent(in) :: ng, model, tile
135 integer, intent(in) :: lbi, ubi, lbj, ubj
136!
137! Local variable declarations.
138!
139 logical :: cgrid
140!
141 integer :: nposb, nposr, nposw
142 integer :: fcount, i, ifield, k, np, status
143!
144 real(dp) :: scale
145
146 real(r8), dimension(Nstation(ng)) :: xpos, ypos, zpos, psta
147# ifdef SOLVE3D
148# ifdef SEDIMENT
149 real(r8), dimension(Nstation(ng)*Nbed) :: xposb, yposb, zposb
150 real(r8), dimension(Nstation(ng)*Nbed) :: bsta
151# endif
152 real(r8), dimension(Nstation(ng)*(N(ng))) :: xposr, yposr, zposr
153 real(r8), dimension(Nstation(ng)*(N(ng)+1)) :: xposw, yposw, zposw
154 real(r8), dimension(Nstation(ng)*(N(ng)+1)) :: rsta
155# endif
156
157 real(r8), allocatable :: ur2d(:,:)
158 real(r8), allocatable :: vr2d(:,:)
159# ifdef SOLVE3D
160 real(r8), allocatable :: ur3d(:,:,:)
161 real(r8), allocatable :: vr3d(:,:,:)
162# endif
163!
164 character (len=*), parameter :: myfile = &
165 & __FILE__//", wrt_station_nf90"
166!
167 sourcefile=myfile
168!
169!-----------------------------------------------------------------------
170! Write out station data at RHO-points.
171!-----------------------------------------------------------------------
172!
173 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
174!
175! Set time record index.
176!
177 sta(ng)%Rindex=sta(ng)%Rindex+1
178 fcount=sta(ng)%Fcount
179 sta(ng)%Nrec(fcount)=sta(ng)%Nrec(fcount)+1
180!
181! Set switch to extract station data at native C-grid position (TRUE)
182! or at RHO-points (FALSE).
183!
184# ifdef STATIONS_CGRID
185 cgrid=.true.
186# else
187 cgrid=.false.
188# endif
189!
190! Set positions for generic extraction routine.
191!
192 nposb=nstation(ng)*nbed
193 nposr=nstation(ng)*n(ng)
194 nposw=nstation(ng)*(n(ng)+1)
195 DO i=1,nstation(ng)
196 xpos(i)=scalars(ng)%SposX(i)
197 ypos(i)=scalars(ng)%SposY(i)
198 zpos(i)=1.0_r8
199# ifdef SOLVE3D
200 DO k=1,n(ng)
201 np=k+(i-1)*n(ng)
202 xposr(np)=scalars(ng)%SposX(i)
203 yposr(np)=scalars(ng)%SposY(i)
204 zposr(np)=real(k,r8)
205 END DO
206 DO k=0,n(ng)
207 np=k+1+(i-1)*(n(ng)+1)
208 xposw(np)=scalars(ng)%SposX(i)
209 yposw(np)=scalars(ng)%SposY(i)
210 zposw(np)=real(k,r8)
211 END DO
212# ifdef SEDIMENT
213 DO k=1,nbed
214 np=k+(i-1)*nbed
215 xposb(np)=scalars(ng)%SposX(i)
216 yposb(np)=scalars(ng)%SposY(i)
217 zposb(np)=real(k,r8)
218 END DO
219# endif
220# endif
221 END DO
222!
223! Write out model time (s).
224!
225 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
226 & trim(vname(1,idtime)), time(ng:), &
227 & (/sta(ng)%Rindex/), (/1/), &
228 & ncid = sta(ng)%ncid, &
229 & varid = sta(ng)%Vid(idtime))
230 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
231!
232! Write out free-surface (m).
233!
234 IF (sout(idfsur,ng)) THEN
235 scale=1.0_dp
236 CALL extract_sta2d (ng, model, cgrid, idfsur, r2dvar, &
237 & lbi, ubi, lbj, ubj, &
238 & scale, ocean(ng)%zeta(:,:,kout), &
239 & nstation(ng), xpos, ypos, psta)
240 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
241 & trim(vname(1,idfsur)), psta, &
242 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
243 & ncid = sta(ng)%ncid, &
244 & varid = sta(ng)%Vid(idfsur))
245 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
246 END IF
247!
248! Write out 2D momentum component (m/s) in the XI-direction.
249!
250 IF (sout(idubar,ng)) THEN
251 scale=1.0_dp
252 CALL extract_sta2d (ng, model, cgrid, idubar, u2dvar, &
253 & lbi, ubi, lbj, ubj, &
254 & scale, ocean(ng)%ubar(:,:,kout), &
255 & nstation(ng), xpos, ypos, psta)
256 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
257 & trim(vname(1,idubar)), psta, &
258 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
259 & ncid = sta(ng)%ncid, &
260 & varid = sta(ng)%Vid(idubar))
261 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
262 END IF
263!
264! Write out 2D momentum component (m/s) in the ETA-direction.
265!
266 IF (sout(idvbar,ng)) THEN
267 scale=1.0_dp
268 CALL extract_sta2d (ng, model, cgrid, idvbar, v2dvar, &
269 & lbi, ubi, lbj, ubj, &
270 & scale, ocean(ng)%vbar(:,:,kout), &
271 & nstation(ng), xpos, ypos, psta)
272 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
273 & trim(vname(1,idvbar)), psta, &
274 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
275 & ncid = sta(ng)%ncid, &
276 & varid = sta(ng)%Vid(idvbar))
277 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
278 END IF
279!
280! Write out 2D Eastward and Northward momentum components (m/s) at
281! RHO-points
282!
283 IF (sout(idu2de,ng).and.sout(idv2dn,ng)) THEN
284 IF (.not.allocated(ur2d)) THEN
285 allocate (ur2d(lbi:ubi,lbj:ubj))
286 ur2d(lbi:ubi,lbj:ubj)=0.0_r8
287 END IF
288 IF (.not.allocated(vr2d)) THEN
289 allocate (vr2d(lbi:ubi,lbj:ubj))
290 vr2d(lbi:ubi,lbj:ubj)=0.0_r8
291 END IF
292 CALL uv_rotate2d (ng, tile, .false., .true., &
293 & lbi, ubi, lbj, ubj, &
294 & grid(ng) % CosAngler, &
295 & grid(ng) % SinAngler, &
296# ifdef MASKING
297 & grid(ng) % rmask_full, &
298# endif
299 & ocean(ng) % ubar(:,:,kout), &
300 & ocean(ng) % vbar(:,:,kout), &
301 & ur2d, vr2d)
302!
303 scale=1.0_dp
304 CALL extract_sta2d (ng, model, cgrid, idu2de, r2dvar, &
305 & lbi, ubi, lbj, ubj, &
306 & scale, ur2d, &
307 & nstation(ng), xpos, ypos, psta)
308 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
309 & trim(vname(1,idu2de)), psta, &
310 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
311 & ncid = sta(ng)%ncid, &
312 & varid = sta(ng)%Vid(idu2de))
313 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
314!
315 CALL extract_sta2d (ng, model, cgrid, idv2dn, r2dvar, &
316 & lbi, ubi, lbj, ubj, &
317 & scale, vr2d, &
318 & nstation(ng), xpos, ypos, psta)
319 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
320 & trim(vname(1,idv2dn)), psta, &
321 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
322 & ncid = sta(ng)%ncid, &
323 & varid = sta(ng)%Vid(idv2dn))
324 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
325
326 deallocate (ur2d)
327 deallocate (vr2d)
328 END IF
329
330# ifdef SOLVE3D
331!
332! Write out 3D momentum component (m/s) in the XI-direction.
333!
334 IF (sout(iduvel,ng)) THEN
335 scale=1.0_dp
336 CALL extract_sta3d (ng, model, cgrid, iduvel, u3dvar, &
337 & lbi, ubi, lbj, ubj, 1, n(ng), &
338 & scale, ocean(ng)%u(:,:,:,nout), &
339 & nposr, xposr, yposr, zposr, rsta)
340 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
341 & trim(vname(1,iduvel)), rsta, &
342 & (/1,1,sta(ng)%Rindex/), &
343 & (/n(ng),nstation(ng),1/), &
344 & ncid = sta(ng)%ncid, &
345 & varid = sta(ng)%Vid(iduvel))
346 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
347 END IF
348!
349! Write out 3D momentum component (m/s) in the ETA-direction.
350!
351 IF (sout(idvvel,ng)) THEN
352 scale=1.0_dp
353 CALL extract_sta3d (ng, model, cgrid, idvvel, v3dvar, &
354 & lbi, ubi, lbj, ubj, 1, n(ng), &
355 & scale, ocean(ng)%v(:,:,:,nout), &
356 & nposr, xposr, yposr, zposr, rsta)
357 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
358 & trim(vname(1,idvvel)), rsta, &
359 & (/1,1,sta(ng)%Rindex/), &
360 & (/n(ng),nstation(ng),1/), &
361 & ncid = sta(ng)%ncid, &
362 & varid = sta(ng)%Vid(idvvel))
363 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
364 END IF
365!
366! Write out 3D Eastward and Northward momentum components (m/s) at
367! RHO-points.
368!
369 IF (sout(idu3de,ng).and.sout(idv3dn,ng)) THEN
370 IF (.not.allocated(ur3d)) THEN
371 allocate (ur3d(lbi:ubi,lbj:ubj,n(ng)))
372 ur3d(lbi:ubi,lbj:ubj,1:n(ng))=0.0_r8
373 END IF
374 IF (.not.allocated(vr3d)) THEN
375 allocate (vr3d(lbi:ubi,lbj:ubj,n(ng)))
376 vr3d(lbi:ubi,lbj:ubj,1:n(ng))=0.0_r8
377 END IF
378 CALL uv_rotate3d (ng, tile, .false., .true., &
379 & lbi, ubi, lbj, ubj, 1, n(ng), &
380 & grid(ng) % CosAngler, &
381 & grid(ng) % SinAngler, &
382# ifdef MASKING
383 & grid(ng) % rmask_full, &
384# endif
385 & ocean(ng) % u(:,:,:,nout), &
386 & ocean(ng) % v(:,:,:,nout), &
387 & ur3d, vr3d)
388!
389 scale=1.0_dp
390 CALL extract_sta3d (ng, model, cgrid, idu3de, r3dvar, &
391 & lbi, ubi, lbj, ubj, 1, n(ng), &
392 & scale, ur3d, &
393 & nposr, xposr, yposr, zposr, rsta)
394 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
395 & trim(vname(1,idu3de)), rsta, &
396 & (/1,1,sta(ng)%Rindex/), &
397 & (/n(ng),nstation(ng),1/), &
398 & ncid = sta(ng)%ncid, &
399 & varid = sta(ng)%Vid(idu3de))
400 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
401!
402 CALL extract_sta3d (ng, model, cgrid, idv3dn, r3dvar, &
403 & lbi, ubi, lbj, ubj, 1, n(ng), &
404 & scale, vr3d, &
405 & nposr, xposr, yposr, zposr, rsta)
406 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
407 & trim(vname(1,idv3dn)), rsta, &
408 & (/1,1,sta(ng)%Rindex/), &
409 & (/n(ng),nstation(ng),1/), &
410 & ncid = sta(ng)%ncid, &
411 & varid = sta(ng)%Vid(idv3dn))
412 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
413
414 deallocate (ur3d)
415 deallocate (vr3d)
416 END IF
417!
418! Write out "true" vertical velocity (m/s).
419!
420 IF (sout(idwvel,ng)) THEN
421 scale=1.0_dp
422 CALL extract_sta3d (ng, model, cgrid, idwvel, w3dvar, &
423 & lbi, ubi, lbj, ubj, 0, n(ng), &
424 & scale, ocean(ng)%wvel, &
425 & nposw, xposw, yposw, zposw, rsta)
426 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
427 & trim(vname(1,idwvel)), rsta, &
428 & (/1,1,sta(ng)%Rindex/), &
429 & (/n(ng)+1,nstation(ng),1/), &
430 & ncid = sta(ng)%ncid, &
431 & varid = sta(ng)%Vid(idwvel))
432 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
433 END IF
434!
435! Write S-coordinate "omega" vertical velocity (m3/s).
436!
437 IF (sout(idovel,ng)) THEN
438 scale=1.0_dp
439 CALL extract_sta3d (ng, model, cgrid, idovel, w3dvar, &
440 & lbi, ubi, lbj, ubj, 0, n(ng), &
441 & scale, ocean(ng)%W, &
442 & nposw, xposw, yposw, zposw, rsta)
443 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
444 & trim(vname(1,idovel)), rsta, &
445 & (/1,1,sta(ng)%Rindex/), &
446 & (/n(ng)+1,nstation(ng),1/), &
447 & ncid = sta(ng)%ncid, &
448 & varid = sta(ng)%Vid(idovel))
449 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
450 END IF
451!
452! Write out tracer type variables.
453!
454 DO i=1,nt(ng)
455 ifield=idtvar(i)
456 IF (sout(ifield,ng)) THEN
457 scale=1.0_dp
458 CALL extract_sta3d (ng, model, cgrid, ifield, r3dvar, &
459 & lbi, ubi, lbj, ubj, 1, n(ng), &
460 & scale, ocean(ng)%t(:,:,:,nout,i), &
461 & nposr, xposr, yposr, zposr, rsta)
462 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
463 & trim(vname(1,idtvar(i))), rsta, &
464 & (/1,1,sta(ng)%Rindex/), &
465 & (/n(ng),nstation(ng),1/), &
466 & ncid = sta(ng)%ncid, &
467 & varid = sta(ng)%Tid(i))
468 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
469 END IF
470 END DO
471!
472! Write out density anomaly.
473!
474 IF (sout(iddano,ng)) THEN
475 scale=1.0_dp
476 CALL extract_sta3d (ng, model, cgrid, iddano, r3dvar, &
477 & lbi, ubi, lbj, ubj, 1, n(ng), &
478 & scale, ocean(ng)%rho, &
479 & nposr, xposr, yposr, zposr, rsta)
480 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
481 & trim(vname(1,iddano)), rsta, &
482 & (/1,1,sta(ng)%Rindex/), &
483 & (/n(ng),nstation(ng),1/), &
484 & ncid = sta(ng)%ncid, &
485 & varid = sta(ng)%Vid(iddano))
486 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
487 END IF
488
489# ifdef LMD_SKPP
490!
491! Write out depth of surface boundary layer.
492!
493 IF (sout(idhsbl,ng)) THEN
494 scale=1.0_dp
495 CALL extract_sta2d (ng, model, cgrid, idhsbl, r2dvar, &
496 & lbi, ubi, lbj, ubj, &
497 & scale, mixing(ng)%hsbl, &
498 & nstation(ng), xpos, ypos, psta)
499 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
500 & trim(vname(1,idhsbl)), psta, &
501 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
502 & ncid = sta(ng)%ncid, &
503 & varid = sta(ng)%Vid(idhsbl))
504 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
505 END IF
506# endif
507# ifdef LMD_BKPP
508!
509! Write out depth of bottom boundary layer.
510!
511 IF (sout(idhbbl,ng)) THEN
512 scale=1.0_dp
513 CALL extract_sta2d (ng, model, cgrid, idhbbl, r2dvar, &
514 & lbi, ubi, lbj, ubj, &
515 & scale, mixing(ng)%hbbl, &
516 & nstation(ng), xpos, ypos, psta)
517 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
518 & trim(vname(1,idhbbl)), psta, &
519 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
520 & ncid = sta(ng)%ncid, &
521 & varid = sta(ng)%Vid(idhbbl))
522 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
523 END IF
524# endif
525!
526! Write out vertical viscosity coefficient.
527!
528 IF (sout(idvvis,ng)) THEN
529 scale=1.0_dp
530 CALL extract_sta3d (ng, model, cgrid, idvvis, w3dvar, &
531 & lbi, ubi, lbj, ubj, 0, n(ng), &
532 & scale, mixing(ng)%Akv, &
533 & nposw, xposw, yposw, zposw, rsta)
534 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
535 & trim(vname(1,idvvis)), rsta, &
536 & (/1,1,sta(ng)%Rindex/), &
537 & (/n(ng)+1,nstation(ng),1/), &
538 & ncid = sta(ng)%ncid, &
539 & varid = sta(ng)%Vid(idvvis))
540 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
541 END IF
542!
543! Write out vertical diffusion coefficient for potential temperature.
544!
545 IF (sout(idtdif,ng)) THEN
546 scale=1.0_dp
547 CALL extract_sta3d (ng, model, cgrid, idtdif, w3dvar, &
548 & lbi, ubi, lbj, ubj, 0, n(ng), &
549 & scale, mixing(ng)%Akt(:,:,:,itemp), &
550 & nposw, xposw, yposw, zposw, rsta)
551 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
552 & trim(vname(1,idtdif)), rsta, &
553 & (/1,1,sta(ng)%Rindex/), &
554 & (/n(ng)+1,nstation(ng),1/), &
555 & ncid = sta(ng)%ncid, &
556 & varid = sta(ng)%Vid(idtdif))
557 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
558 END IF
559
560# ifdef SALINITY
561!
562! Write out vertical diffusion coefficient for salinity.
563!
564 IF (sout(idsdif,ng)) THEN
565 scale=1.0_dp
566 CALL extract_sta3d (ng, model, cgrid, idsdif, w3dvar, &
567 & lbi, ubi, lbj, ubj, 0, n(ng), &
568 & scale, mixing(ng)%Akt(:,:,:,isalt), &
569 & nposw, xposw, yposw, zposw, rsta)
570 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
571 & trim(vname(1,idsdif)), rsta, &
572 & (/1,1,sta(ng)%Rindex/), &
573 & (/n(ng)+1,nstation(ng),1/), &
574 & ncid = sta(ng)%ncid, &
575 & varid = sta(ng)%Vid(idsdif))
576 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
577 END IF
578# endif
579# if defined GLS_MIXING || defined MY25_MIXING
580!
581! Write out turbulent kinetic energy.
582!
583 IF (sout(idmtke,ng)) THEN
584 scale=1.0_dp
585 CALL extract_sta3d (ng, model, cgrid, idmtke, w3dvar, &
586 & lbi, ubi, lbj, ubj, 0, n(ng), &
587 & scale, mixing(ng)%tke(:,:,:,nout), &
588 & nposw, xposw, yposw, zposw, rsta)
589 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
590 & trim(vname(1,idmtke)), rsta, &
591 & (/1,1,sta(ng)%Rindex/), &
592 & (/n(ng)+1,nstation(ng),1/), &
593 & ncid = sta(ng)%ncid, &
594 & varid = sta(ng)%Vid(idmtke))
595 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
596 END IF
597!
598! Write out turbulent kinetic energy times length scale.
599!
600 IF (sout(idmtls,ng)) THEN
601 scale=1.0_dp
602 CALL extract_sta3d (ng, model, cgrid, idmtls, w3dvar, &
603 & lbi, ubi, lbj, ubj, 0, n(ng), &
604 & scale, mixing(ng)%gls(:,:,:,nout), &
605 & nposw, xposw, yposw, zposw, rsta)
606 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
607 & trim(vname(1,idmtls)), rsta, &
608 & (/1,1,sta(ng)%Rindex/), &
609 & (/n(ng)+1,nstation(ng),1/), &
610 & ncid = sta(ng)%ncid, &
611 & varid = sta(ng)%Vid(idmtls))
612 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
613 END IF
614# endif
615# if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS
616!
617! Write out surface air pressure.
618!
619 IF (sout(idpair,ng)) THEN
620 scale=1.0_dp
621 CALL extract_sta2d (ng, model, cgrid, idpair, r2dvar, &
622 & lbi, ubi, lbj, ubj, &
623 & scale, forces(ng)%Pair, &
624 & nstation(ng), xpos, ypos, psta)
625 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
626 & trim(vname(1,idpair)), psta, &
627 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
628 & ncid = sta(ng)%ncid, &
629 & varid = sta(ng)%Vid(idpair))
630 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
631 END IF
632# endif
633# if defined BULK_FLUXES || defined ECOSIM
634!
635! Write out surface winds.
636!
637 IF (sout(iduair,ng)) THEN
638 scale=1.0_dp
639 CALL extract_sta2d (ng, model, cgrid, iduair, r2dvar, &
640 & lbi, ubi, lbj, ubj, &
641 & scale, forces(ng)%Uwind, &
642 & nstation(ng), xpos, ypos, psta)
643 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
644 & trim(vname(1,iduair)), psta, &
645 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
646 & ncid = sta(ng)%ncid, &
647 & varid = sta(ng)%Vid(iduair))
648 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
649 END IF
650!
651 IF (sout(idvair,ng)) THEN
652 scale=1.0_dp
653 CALL extract_sta2d (ng, model, cgrid, idvair, r2dvar, &
654 & lbi, ubi, lbj, ubj, &
655 & scale, forces(ng)%Vwind, &
656 & nstation(ng), xpos, ypos, psta)
657 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
658 & trim(vname(1,idvair)), psta, &
659 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
660 & ncid = sta(ng)%ncid, &
661 & varid = sta(ng)%Vid(idvair))
662 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
663 END IF
664!
665! Write out 2D Eastward and Northward surface winds (m/s) at
666! RHO-points
667!
668 IF (sout(iduaie,ng).and.sout(idvain,ng)) THEN
669 IF (.not.allocated(ur2d)) THEN
670 allocate (ur2d(lbi:ubi,lbj:ubj))
671 ur2d(lbi:ubi,lbj:ubj)=0.0_r8
672 END IF
673 IF (.not.allocated(vr2d)) THEN
674 allocate (vr2d(lbi:ubi,lbj:ubj))
675 vr2d(lbi:ubi,lbj:ubj)=0.0_r8
676 END IF
677 CALL uv_rotate2d (ng, tile, .false., .true., &
678 & lbi, ubi, lbj, ubj, &
679 & grid(ng) % CosAngler, &
680 & grid(ng) % SinAngler, &
681# ifdef MASKING
682 & grid(ng) % rmask_full, &
683# endif
684 & forces(ng) % Uwind, &
685 & forces(ng) % Vwind, &
686 & ur2d, vr2d)
687!
688 scale=1.0_dp
689 CALL extract_sta2d (ng, model, cgrid, iduaie, r2dvar, &
690 & lbi, ubi, lbj, ubj, &
691 & scale, ur2d, &
692 & nstation(ng), xpos, ypos, psta)
693 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
694 & trim(vname(1,iduaie)), psta, &
695 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
696 & ncid = sta(ng)%ncid, &
697 & varid = sta(ng)%Vid(iduaie))
698 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
699!
700 CALL extract_sta2d (ng, model, cgrid, idvain, r2dvar, &
701 & lbi, ubi, lbj, ubj, &
702 & scale, vr2d, &
703 & nstation(ng), xpos, ypos, psta)
704 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
705 & trim(vname(1,idvain)), psta, &
706 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
707 & ncid = sta(ng)%ncid, &
708 & varid = sta(ng)%Vid(idvain))
709 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
710
711 deallocate (ur2d)
712 deallocate (vr2d)
713 END IF
714# endif
715!
716! Write out surface net heat flux.
717!
718 IF (sout(idtsur(itemp),ng)) THEN
719 ifield=idtsur(itemp)
720 scale=rho0*cp
721 CALL extract_sta2d (ng, model, cgrid, ifield, r2dvar, &
722 & lbi, ubi, lbj, ubj, &
723 & scale, forces(ng)%stflx(:,:,itemp), &
724 & nstation(ng), xpos, ypos, psta)
725 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
726 & trim(vname(1,ifield)), psta, &
727 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
728 & ncid = sta(ng)%ncid, &
729 & varid = sta(ng)%Vid(ifield))
730 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
731 END IF
732
733# ifdef SALINITY
734!
735! Write out surface salt flux.
736!
737 IF (sout(idtsur(isalt),ng)) THEN
738 ifield=idtsur(isalt)
739 scale=1.0_dp
740 CALL extract_sta2d (ng, model, cgrid, ifield, r2dvar, &
741 & lbi, ubi, lbj, ubj, &
742 & scale, forces(ng)%stflx(:,:,isalt), &
743 & nstation(ng), xpos, ypos, psta)
744 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
745 & trim(vname(1,ifield)), psta, &
746 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
747 & ncid = sta(ng)%ncid, &
748 & varid = sta(ng)%Vid(ifield))
749 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
750 END IF
751# endif
752
753# ifdef BULK_FLUXES
754!
755! Write out latent heat flux.
756!
757 IF (sout(idlhea,ng)) THEN
758 scale=rho0*cp
759 CALL extract_sta2d (ng, model, cgrid, idlhea, r2dvar, &
760 & lbi, ubi, lbj, ubj, &
761 & scale, forces(ng)%lhflx, &
762 & nstation(ng), xpos, ypos, psta)
763 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
764 & trim(vname(1,idlhea)), psta, &
765 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
766 & ncid = sta(ng)%ncid, &
767 & varid = sta(ng)%Vid(idlhea))
768 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
769 END IF
770!
771! Write out sensible heat flux.
772!
773 IF (sout(idshea,ng)) THEN
774 scale=rho0*cp
775 CALL extract_sta2d (ng, model, cgrid, idshea, r2dvar, &
776 & lbi, ubi, lbj, ubj, &
777 & scale, forces(ng)%shflx, &
778 & nstation(ng), xpos, ypos, psta)
779 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
780 & trim(vname(1,idshea)), psta, &
781 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
782 & ncid = sta(ng)%ncid, &
783 & varid = sta(ng)%Vid(idshea))
784 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
785 END IF
786!
787! Write out longwave radiation flux.
788!
789 IF (sout(idlrad,ng)) THEN
790 scale=rho0*cp
791 CALL extract_sta2d (ng, model, cgrid, idlrad, r2dvar, &
792 & lbi, ubi, lbj, ubj, &
793 & scale, forces(ng)%lrflx, &
794 & nstation(ng), xpos, ypos, psta)
795 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
796 & trim(vname(1,idlrad)), psta, &
797 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
798 & ncid = sta(ng)%ncid, &
799 & varid = sta(ng)%Vid(idlrad))
800 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
801 END IF
802
803# endif
804# ifdef SHORTWAVE
805!
806! Write out shortwave radiation flux.
807!
808 IF (sout(idsrad,ng)) THEN
809 scale=rho0*cp
810 CALL extract_sta2d (ng, model, cgrid, idsrad, r2dvar, &
811 & lbi, ubi, lbj, ubj, &
812 & scale, forces(ng)%srflx, &
813 & nstation(ng), xpos, ypos, psta)
814 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
815 & trim(vname(1,idsrad)), psta, &
816 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
817 & ncid = sta(ng)%ncid, &
818 & varid = sta(ng)%Vid(idsrad))
819 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
820 END IF
821# endif
822# if defined EMINUSP && defined BULK_FLUXES
823!
824! Write out E-P (m/s).
825!
826 IF (sout(idempf,ng)) THEN
827 scale=1.0_dp
828 CALL extract_sta2d (ng, model, cgrid, idempf, r2dvar, &
829 & lbi, ubi, lbj, ubj, &
830 & scale, forces(ng)%stflux(:,:,isalt), &
831 & nstation(ng), xpos, ypos, psta)
832 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
833 & trim(vname(1,idempf)), psta, &
834 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
835 & ncid = sta(ng)%ncid, &
836 & varid = sta(ng)%Vid(idempf))
837 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
838 END IF
839!
840! Write out evaporation rate (kg/m2/s).
841!
842 IF (sout(idevap,ng)) THEN
843 scale=1.0_dp
844 CALL extract_sta2d (ng, model, cgrid, idevap, r2dvar, &
845 & lbi, ubi, lbj, ubj, &
846 & scale, forces(ng)%evap, &
847 & nstation(ng), xpos, ypos, psta)
848 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
849 & trim(vname(1,idevap)), psta, &
850 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
851 & ncid = sta(ng)%ncid, &
852 & varid = sta(ng)%Vid(idevap))
853 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
854 END IF
855!
856! Write out precipitation rate (kg/m2/s).
857!
858 IF (sout(idrain,ng)) THEN
859 scale=1.0_dp
860 CALL extract_sta2d (ng, model, cgrid, idrain, r2dvar, &
861 & lbi, ubi, lbj, ubj, &
862 & scale, forces(ng)%rain, &
863 & nstation(ng), xpos, ypos, psta)
864 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
865 & trim(vname(1,idrain)), psta, &
866 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
867 & ncid = sta(ng)%ncid, &
868 & varid = sta(ng)%Vid(idrain))
869 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
870 END IF
871# endif
872# endif
873!
874! Write out surface U-momentum stress.
875!
876 IF (sout(idusms,ng)) THEN
877 scale=rho0
878 CALL extract_sta2d (ng, model, cgrid, idusms, u2dvar, &
879 & lbi, ubi, lbj, ubj, &
880 & scale, forces(ng)%sustr, &
881 & nstation(ng), xpos, ypos, psta)
882 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
883 & trim(vname(1,idusms)), psta, &
884 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
885 & ncid = sta(ng)%ncid, &
886 & varid = sta(ng)%Vid(idusms))
887 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
888 END IF
889!
890! Write out surface V-momentum stress.
891!
892 IF (sout(idvsms,ng)) THEN
893 scale=rho0
894 CALL extract_sta2d (ng, model, cgrid, idvsms, v2dvar, &
895 & lbi, ubi, lbj, ubj, &
896 & scale, forces(ng)%svstr, &
897 & nstation(ng), xpos, ypos, psta)
898 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
899 & trim(vname(1,idvsms)), psta, &
900 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
901 & ncid = sta(ng)%ncid, &
902 & varid = sta(ng)%Vid(idvsms))
903 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
904 END IF
905!
906! Write out bottom U-momentum stress.
907!
908 IF (sout(idubms,ng)) THEN
909 scale=-rho0
910 CALL extract_sta2d (ng, model, cgrid, idubms, u2dvar, &
911 & lbi, ubi, lbj, ubj, &
912 & scale, forces(ng)%bustr, &
913 & nstation(ng), xpos, ypos, psta)
914 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
915 & trim(vname(1,idubms)), psta, &
916 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
917 & ncid = sta(ng)%ncid, &
918 & varid = sta(ng)%Vid(idubms))
919 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
920 END IF
921!
922! Write out bottom V-momentum stress.
923!
924 IF (sout(idvbms,ng)) THEN
925 scale=-rho0
926 CALL extract_sta2d (ng, model, cgrid, idvbms, v2dvar, &
927 & lbi, ubi, lbj, ubj, &
928 & scale, forces(ng)%bvstr, &
929 & nstation(ng), xpos, ypos, psta)
930 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
931 & trim(vname(1,idvbms)), psta, &
932 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
933 & ncid = sta(ng)%ncid, &
934 & varid = sta(ng)%Vid(idvbms))
935 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
936 END IF
937
938# ifdef WET_DRY
939!
940! Write out wet/dry mask at RHO-points.
941!
942 IF (sout(idrwet,ng)) THEN
943 scale=1.0_dp
944 CALL extract_sta2d (ng, model, cgrid, idrwet, r2dvar, &
945 & lbi, ubi, lbj, ubj, &
946 & scale, grid(ng)%rmask_wet, &
947 & nstation(ng), xpos, ypos, psta)
948 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
949 & trim(vname(1,idrwet)), psta, &
950 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
951 & ncid = sta(ng)%ncid, &
952 & varid = sta(ng)%Vid(idrwet))
953 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
954 END IF
955!
956! Write out wet/dry mask at U-points.
957!
958 IF (sout(iduwet,ng)) THEN
959 scale=1.0_dp
960 CALL extract_sta2d (ng, model, cgrid, iduwet, u2dvar, &
961 & lbi, ubi, lbj, ubj, &
962 & scale, grid(ng)%umask_wet, &
963 & nstation(ng), xpos, ypos, psta)
964 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
965 & trim(vname(1,iduwet)), psta, &
966 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
967 & ncid = sta(ng)%ncid, &
968 & varid = sta(ng)%Vid(iduwet))
969 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
970 END IF
971!
972! Write out wet/dry mask at V-points.
973!
974 IF (sout(idvwet,ng)) THEN
975 scale=1.0_dp
976 CALL extract_sta2d (ng, model, cgrid, idvwet, v2dvar, &
977 & lbi, ubi, lbj, ubj, &
978 & scale, grid(ng)%vmask_wet, &
979 & nstation(ng), xpos, ypos, psta)
980 CALL netcdf_put_fvar (ng, model, sta(ng)%name, &
981 & trim(vname(1,idvwet)), psta, &
982 & (/1,sta(ng)%Rindex/), (/nstation(ng),1/), &
983 & ncid = sta(ng)%ncid, &
984 & varid = sta(ng)%Vid(idvwet))
985 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
986 END IF
987# endif
988
989# if defined BBL_MODEL || defined WAVES_OUTPUT
990!
991!-----------------------------------------------------------------------
992! Write out the bottom boundary layer model or waves variables.
993!-----------------------------------------------------------------------
994!
995 CALL bbl_wrt_station_nf90 (ng, model, tile, &
996 & lbi, ubi, lbj, ubj, &
997 & sout, sta)
998 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
999# endif
1000
1001# ifdef ICE_MODEL
1002!
1003!-----------------------------------------------------------------------
1004! Write out the ice model variables.
1005!-----------------------------------------------------------------------
1006!
1007 CALL ice_wrt_station_nf90 (ng, model, tile, &
1008 & lbi, ubi, lbj, ubj, &
1009 & sout, sta)
1010 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1011# endif
1012
1013# ifdef SEDIMENT
1014!
1015!-----------------------------------------------------------------------
1016! Write out the sediment model variables.
1017!-----------------------------------------------------------------------
1018!
1019 CALL sediment_wrt_station_nf90 (ng, model, tile, &
1020 & lbi, ubi, lbj, ubj, &
1021 & sout, sta)
1022 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1023# endif
1024
1025# if defined WEC || defined WEC_VF
1026!
1027!-----------------------------------------------------------------------
1028! Write out the Waves Effect on Currents station variables.
1029!-----------------------------------------------------------------------
1030!
1031 CALL wec_wrt_station_nf90 (ng, model, tile, &
1032 & lbi, ubi, lbj, ubj, &
1033 & sout, sta)
1034 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1035# endif
1036!
1037!-----------------------------------------------------------------------
1038! Synchronize stations NetCDF file to disk.
1039!-----------------------------------------------------------------------
1040!
1041 CALL netcdf_sync (ng, model, sta(ng)%name, sta(ng)%ncid)
1042
1043 RETURN
1044 END SUBROUTINE wrt_station_nf90
1045
1046# if defined PIO_LIB && defined DISTRIBUTE
1047!
1048!***********************************************************************
1049 SUBROUTINE wrt_station_pio (ng, model, tile, &
1050 & LBi, UBi, LBj, UBj)
1051!***********************************************************************
1052!
1053 USE mod_pio_netcdf
1054!
1055! Imported variable declarations.
1056!
1057 integer, intent(in) :: ng, tile
1058 integer, intent(in) :: lbi, ubi, lbj, ubj
1059!
1060! Local variable declarations.
1061!
1062 logical :: cgrid
1063!
1064 integer :: nposb, nposr, nposw
1065 integer :: fcount, i, ifield, k, np, status
1066!
1067 real(dp) :: scale
1068
1069 real(r8), dimension(Nstation(ng)) :: xpos, ypos, zpos, psta
1070# ifdef SOLVE3D
1071# ifdef SEDIMENT
1072 real(r8), dimension(Nstation(ng)*Nbed) :: xposb, yposb, zposb
1073 real(r8), dimension(Nstation(ng)*Nbed) :: bsta
1074# endif
1075 real(r8), dimension(Nstation(ng)*(N(ng))) :: xposr, yposr, zposr
1076 real(r8), dimension(Nstation(ng)*(N(ng)+1)) :: xposw, yposw, zposw
1077 real(r8), dimension(Nstation(ng)*(N(ng)+1)) :: rsta
1078# endif
1079
1080 real(r8), allocatable :: ur2d(:,:)
1081 real(r8), allocatable :: vr2d(:,:)
1082# ifdef SOLVE3D
1083 real(r8), allocatable :: ur3d(:,:,:)
1084 real(r8), allocatable :: vr3d(:,:,:)
1085# endif
1086!
1087 character (len=*), parameter :: myfile = &
1088 & __FILE__//", wrt_station_pio"
1089!
1090 sourcefile=myfile
1091!
1092!-----------------------------------------------------------------------
1093! Write out station data at RHO-points.
1094!-----------------------------------------------------------------------
1095!
1096 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1097!
1098! Set time record index.
1099!
1100 sta(ng)%Rindex=sta(ng)%Rindex+1
1101 fcount=sta(ng)%Fcount
1102 sta(ng)%Nrec(fcount)=sta(ng)%Nrec(fcount)+1
1103!
1104! Set switch to extract station data at native C-grid position (TRUE)
1105! or at RHO-points (FALSE).
1106!
1107# ifdef STATIONS_CGRID
1108 cgrid=.true.
1109# else
1110 cgrid=.false.
1111# endif
1112!
1113! Set positions for generic extraction routine.
1114!
1115 nposb=nstation(ng)*nbed
1116 nposr=nstation(ng)*n(ng)
1117 nposw=nstation(ng)*(n(ng)+1)
1118 DO i=1,nstation(ng)
1119 xpos(i)=scalars(ng)%SposX(i)
1120 ypos(i)=scalars(ng)%SposY(i)
1121 zpos(i)=1.0_r8
1122# ifdef SOLVE3D
1123 DO k=1,n(ng)
1124 np=k+(i-1)*n(ng)
1125 xposr(np)=scalars(ng)%SposX(i)
1126 yposr(np)=scalars(ng)%SposY(i)
1127 zposr(np)=real(k,r8)
1128 END DO
1129 DO k=0,n(ng)
1130 np=k+1+(i-1)*(n(ng)+1)
1131 xposw(np)=scalars(ng)%SposX(i)
1132 yposw(np)=scalars(ng)%SposY(i)
1133 zposw(np)=real(k,r8)
1134 END DO
1135# ifdef SEDIMENT
1136 DO k=1,nbed
1137 np=k+(i-1)*nbed
1138 xposb(np)=scalars(ng)%SposX(i)
1139 yposb(np)=scalars(ng)%SposY(i)
1140 zposb(np)=real(k,r8)
1141 END DO
1142# endif
1143# endif
1144 END DO
1145!
1146! Write out model time (s).
1147!
1148 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1149 & trim(vname(1,idtime)), time(ng:), &
1150 & (/sta(ng)%Rindex/), (/1/), &
1151 & piofile = sta(ng)%pioFile, &
1152 & piovar = sta(ng)%pioVar(idtime)%vd)
1153 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1154!
1155! Write out free-surface (m).
1156!
1157 IF (sout(idfsur,ng)) THEN
1158 scale=1.0_dp
1159 CALL extract_sta2d (ng, model, cgrid, idfsur, r2dvar, &
1160 & lbi, ubi, lbj, ubj, &
1161 & scale, ocean(ng)%zeta(:,:,kout), &
1162 & nstation(ng), xpos, ypos, psta)
1163 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1164 & trim(vname(1,idfsur)), psta, &
1165 & (/1,sta(ng)%Rindex/), &
1166 & (/nstation(ng),1/), &
1167 & piofile = sta(ng)%pioFile, &
1168 & piovar = sta(ng)%pioVar(idfsur)%vd)
1169 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1170 END IF
1171
1172!
1173! Write out 2D momentum component (m/s) in the XI-direction.
1174!
1175 IF (sout(idubar,ng)) THEN
1176 scale=1.0_dp
1177 CALL extract_sta2d (ng, model, cgrid, idubar, u2dvar, &
1178 & lbi, ubi, lbj, ubj, &
1179 & scale, ocean(ng)%ubar(:,:,kout), &
1180 & nstation(ng), xpos, ypos, psta)
1181 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1182 & trim(vname(1,idubar)), psta, &
1183 & (/1,sta(ng)%Rindex/), &
1184 & (/nstation(ng),1/), &
1185 & piofile = sta(ng)%pioFile, &
1186 & piovar = sta(ng)%pioVar(idubar)%vd)
1187 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1188 END IF
1189!
1190! Write out 2D momentum component (m/s) in the ETA-direction.
1191!
1192 IF (sout(idvbar,ng)) THEN
1193 scale=1.0_dp
1194 CALL extract_sta2d (ng, model, cgrid, idvbar, v2dvar, &
1195 & lbi, ubi, lbj, ubj, &
1196 & scale, ocean(ng)%vbar(:,:,kout), &
1197 & nstation(ng), xpos, ypos, psta)
1198 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1199 & trim(vname(1,idvbar)), psta, &
1200 & (/1,sta(ng)%Rindex/), &
1201 & (/nstation(ng),1/), &
1202 & piofile = sta(ng)%pioFile, &
1203 & piovar = sta(ng)%pioVar(idvbar)%vd)
1204 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1205 END IF
1206!
1207! Write out 2D Eastward and Northward momentum components (m/s) at
1208! RHO-points
1209!
1210 IF (sout(idu2de,ng).and.sout(idv2dn,ng)) THEN
1211 IF (.not.allocated(ur2d)) THEN
1212 allocate (ur2d(lbi:ubi,lbj:ubj))
1213 ur2d(lbi:ubi,lbj:ubj)=0.0_r8
1214 END IF
1215 IF (.not.allocated(vr2d)) THEN
1216 allocate (vr2d(lbi:ubi,lbj:ubj))
1217 vr2d(lbi:ubi,lbj:ubj)=0.0_r8
1218 END IF
1219 CALL uv_rotate2d (ng, tile, .false., .true., &
1220 & lbi, ubi, lbj, ubj, &
1221 & grid(ng) % CosAngler, &
1222 & grid(ng) % SinAngler, &
1223# ifdef MASKING
1224 & grid(ng) % rmask_full, &
1225# endif
1226 & ocean(ng) % ubar(:,:,kout), &
1227 & ocean(ng) % vbar(:,:,kout), &
1228 & ur2d, vr2d)
1229!
1230 scale=1.0_dp
1231 CALL extract_sta2d (ng, model, cgrid, idu2de, r2dvar, &
1232 & lbi, ubi, lbj, ubj, &
1233 & scale, ur2d, &
1234 & nstation(ng), xpos, ypos, psta)
1235 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1236 & trim(vname(1,idu2de)), psta, &
1237 & (/1,sta(ng)%Rindex/), &
1238 & (/nstation(ng),1/), &
1239 & piofile = sta(ng)%pioFile, &
1240 & piovar = sta(ng)%pioVar(idu2de)%vd)
1241 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1242!
1243 CALL extract_sta2d (ng, model, cgrid, idv2dn, r2dvar, &
1244 & lbi, ubi, lbj, ubj, &
1245 & scale, vr2d, &
1246 & nstation(ng), xpos, ypos, psta)
1247 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1248 & trim(vname(1,idv2dn)), psta, &
1249 & (/1,sta(ng)%Rindex/), &
1250 & (/nstation(ng),1/), &
1251 & piofile = sta(ng)%pioFile, &
1252 & piovar = sta(ng)%pioVar(idv2dn)%vd)
1253 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1254
1255 deallocate (ur2d)
1256 deallocate (vr2d)
1257 END IF
1258
1259# ifdef SOLVE3D
1260!
1261! Write out 3D momentum component (m/s) in the XI-direction.
1262!
1263 IF (sout(iduvel,ng)) THEN
1264 scale=1.0_dp
1265 CALL extract_sta3d (ng, model, cgrid, iduvel, u3dvar, &
1266 & lbi, ubi, lbj, ubj, 1, n(ng), &
1267 & scale, ocean(ng)%u(:,:,:,nout), &
1268 & nposr, xposr, yposr, zposr, rsta)
1269 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1270 & trim(vname(1,iduvel)), rsta, &
1271 & (/1,1,sta(ng)%Rindex/), &
1272 & (/n(ng),nstation(ng),1/), &
1273 & piofile = sta(ng)%pioFile, &
1274 & piovar = sta(ng)%pioVar(iduvel)%vd)
1275 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1276 END IF
1277!
1278! Write out 3D momentum component (m/s) in the ETA-direction.
1279!
1280 IF (sout(idvvel,ng)) THEN
1281 scale=1.0_dp
1282 CALL extract_sta3d (ng, model, cgrid, idvvel, v3dvar, &
1283 & lbi, ubi, lbj, ubj, 1, n(ng), &
1284 & scale, ocean(ng)%v(:,:,:,nout), &
1285 & nposr, xposr, yposr, zposr, rsta)
1286 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1287 & trim(vname(1,idvvel)), rsta, &
1288 & (/1,1,sta(ng)%Rindex/), &
1289 & (/n(ng),nstation(ng),1/), &
1290 & piofile = sta(ng)%pioFile, &
1291 & piovar = sta(ng)%pioVar(idvvel)%vd)
1292 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1293 END IF
1294!
1295! Write out 3D Eastward and Northward momentum components (m/s) at
1296! RHO-points.
1297!
1298 IF (sout(idu3de,ng).and.sout(idv3dn,ng)) THEN
1299 IF (.not.allocated(ur3d)) THEN
1300 allocate (ur3d(lbi:ubi,lbj:ubj,n(ng)))
1301 ur3d(lbi:ubi,lbj:ubj,1:n(ng))=0.0_r8
1302 END IF
1303 IF (.not.allocated(vr3d)) THEN
1304 allocate (vr3d(lbi:ubi,lbj:ubj,n(ng)))
1305 vr3d(lbi:ubi,lbj:ubj,1:n(ng))=0.0_r8
1306 END IF
1307 CALL uv_rotate3d (ng, tile, .false., .true., &
1308 & lbi, ubi, lbj, ubj, 1, n(ng), &
1309 & grid(ng) % CosAngler, &
1310 & grid(ng) % SinAngler, &
1311# ifdef MASKING
1312 & grid(ng) % rmask_full, &
1313# endif
1314 & ocean(ng) % u(:,:,:,nout), &
1315 & ocean(ng) % v(:,:,:,nout), &
1316 & ur3d, vr3d)
1317!
1318 scale=1.0_dp
1319 CALL extract_sta3d (ng, model, cgrid, idu3de, r3dvar, &
1320 & lbi, ubi, lbj, ubj, 1, n(ng), &
1321 & scale, ur3d, &
1322 & nposr, xposr, yposr, zposr, rsta)
1323 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1324 & trim(vname(1,idu3de)), rsta, &
1325 & (/1,1,sta(ng)%Rindex/), &
1326 & (/n(ng),nstation(ng),1/), &
1327 & piofile = sta(ng)%pioFile, &
1328 & piovar = sta(ng)%pioVar(idu3de)%vd)
1329 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1330!
1331 CALL extract_sta3d (ng, model, cgrid, idv3dn, r3dvar, &
1332 & lbi, ubi, lbj, ubj, 1, n(ng), &
1333 & scale, vr3d, &
1334 & nposr, xposr, yposr, zposr, rsta)
1335 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1336 & trim(vname(1,idv3dn)), rsta, &
1337 & (/1,1,sta(ng)%Rindex/), &
1338 & (/n(ng),nstation(ng),1/), &
1339 & piofile = sta(ng)%pioFile, &
1340 & piovar = sta(ng)%pioVar(idv3dn)%vd)
1341 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1342
1343 deallocate (ur3d)
1344 deallocate (vr3d)
1345 END IF
1346!
1347! Write out "true" vertical velocity (m/s).
1348!
1349 IF (sout(idwvel,ng)) THEN
1350 scale=1.0_dp
1351 CALL extract_sta3d (ng, model, cgrid, idwvel, w3dvar, &
1352 & lbi, ubi, lbj, ubj, 0, n(ng), &
1353 & scale, ocean(ng)%wvel, &
1354 & nposw, xposw, yposw, zposw, rsta)
1355 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1356 & trim(vname(1,idwvel)), rsta, &
1357 & (/1,1,sta(ng)%Rindex/), &
1358 & (/n(ng)+1,nstation(ng),1/), &
1359 & piofile = sta(ng)%pioFile, &
1360 & piovar = sta(ng)%pioVar(idwvel)%vd)
1361 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1362 END IF
1363!
1364! Write S-coordinate "omega" vertical velocity (m3/s).
1365!
1366 IF (sout(idovel,ng)) THEN
1367 scale=1.0_dp
1368 CALL extract_sta3d (ng, model, cgrid, idovel, w3dvar, &
1369 & lbi, ubi, lbj, ubj, 0, n(ng), &
1370 & scale, ocean(ng)%W, &
1371 & nposw, xposw, yposw, zposw, rsta)
1372 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1373 & trim(vname(1,idovel)), rsta, &
1374 & (/1,1,sta(ng)%Rindex/), &
1375 & (/n(ng)+1,nstation(ng),1/), &
1376 & piofile = sta(ng)%pioFile, &
1377 & piovar = sta(ng)%pioVar(idovel)%vd)
1378 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1379 END IF
1380!
1381! Write out tracer type variables.
1382!
1383 DO i=1,nt(ng)
1384 ifield=idtvar(i)
1385 IF (sout(ifield,ng)) THEN
1386 scale=1.0_dp
1387 CALL extract_sta3d (ng, model, cgrid, ifield, r3dvar, &
1388 & lbi, ubi, lbj, ubj, 1, n(ng), &
1389 & scale, ocean(ng)%t(:,:,:,nout,i), &
1390 & nposr, xposr, yposr, zposr, rsta)
1391 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1392 & trim(vname(1,idtvar(i))), rsta, &
1393 & (/1,1,sta(ng)%Rindex/), &
1394 & (/n(ng),nstation(ng),1/), &
1395 & piofile = sta(ng)%pioFile, &
1396 & piovar = sta(ng)%pioTrc(i)%vd)
1397 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1398 END IF
1399 END DO
1400!
1401! Write out density anomaly.
1402!
1403 IF (sout(iddano,ng)) THEN
1404 scale=1.0_dp
1405 CALL extract_sta3d (ng, model, cgrid, iddano, r3dvar, &
1406 & lbi, ubi, lbj, ubj, 1, n(ng), &
1407 & scale, ocean(ng)%rho, &
1408 & nposr, xposr, yposr, zposr, rsta)
1409 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1410 & trim(vname(1,iddano)), rsta, &
1411 & (/1,1,sta(ng)%Rindex/), &
1412 & (/n(ng),nstation(ng),1/), &
1413 & piofile = sta(ng)%pioFile, &
1414 & piovar = sta(ng)%pioVar(iddano)%vd)
1415 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1416 END IF
1417
1418# ifdef LMD_SKPP
1419!
1420! Write out depth of surface boundary layer.
1421!
1422 IF (sout(idhsbl,ng)) THEN
1423 scale=1.0_dp
1424 CALL extract_sta2d (ng, model, cgrid, idhsbl, r2dvar, &
1425 & lbi, ubi, lbj, ubj, &
1426 & scale, mixing(ng)%hsbl, &
1427 & nstation(ng), xpos, ypos, psta)
1428 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1429 & trim(vname(1,idhsbl)), psta, &
1430 & (/1,sta(ng)%Rindex/), &
1431 & (/nstation(ng),1/), &
1432 & piofile = sta(ng)%pioFile, &
1433 & piovar = sta(ng)%pioVar(idhsbl)%vd)
1434 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1435 END IF
1436# endif
1437# ifdef LMD_BKPP
1438!
1439! Write out depth of bottom boundary layer.
1440!
1441 IF (sout(idhbbl,ng)) THEN
1442 scale=1.0_dp
1443 CALL extract_sta2d (ng, model, cgrid, idhbbl, r2dvar, &
1444 & lbi, ubi, lbj, ubj, &
1445 & scale, mixing(ng)%hbbl, &
1446 & nstation(ng), xpos, ypos, psta)
1447 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1448 & trim(vname(1,idhbbl)), psta, &
1449 & (/1,sta(ng)%Rindex/), &
1450 & (/nstation(ng),1/), &
1451 & piofile = sta(ng)%pioFile, &
1452 & piovar = sta(ng)%pioVar(idhbbl)%vd)
1453 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1454 END IF
1455# endif
1456!
1457! Write out vertical viscosity coefficient.
1458!
1459 IF (sout(idvvis,ng)) THEN
1460 scale=1.0_dp
1461 CALL extract_sta3d (ng, model, cgrid, idvvis, w3dvar, &
1462 & lbi, ubi, lbj, ubj, 0, n(ng), &
1463 & scale, mixing(ng)%Akv, &
1464 & nposw, xposw, yposw, zposw, rsta)
1465 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1466 & trim(vname(1,idvvis)), rsta, &
1467 & (/1,1,sta(ng)%Rindex/), &
1468 & (/n(ng)+1,nstation(ng),1/), &
1469 & piofile = sta(ng)%pioFile, &
1470 & piovar = sta(ng)%pioVar(idvvis)%vd)
1471 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1472 END IF
1473!
1474! Write out vertical diffusion coefficient for potential temperature.
1475!
1476 IF (sout(idtdif,ng)) THEN
1477 scale=1.0_dp
1478 CALL extract_sta3d (ng, model, cgrid, idtdif, w3dvar, &
1479 & lbi, ubi, lbj, ubj, 0, n(ng), &
1480 & scale, mixing(ng)%Akt(:,:,:,itemp), &
1481 & nposw, xposw, yposw, zposw, rsta)
1482 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1483 & trim(vname(1,idtdif)), rsta, &
1484 & (/1,1,sta(ng)%Rindex/), &
1485 & (/n(ng)+1,nstation(ng),1/), &
1486 & piofile = sta(ng)%pioFile, &
1487 & piovar = sta(ng)%pioVar(idtdif)%vd)
1488 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1489 END IF
1490
1491# ifdef SALINITY
1492!
1493! Write out vertical diffusion coefficient for salinity.
1494!
1495 IF (sout(idsdif,ng)) THEN
1496 scale=1.0_dp
1497 CALL extract_sta3d (ng, model, cgrid, idsdif, w3dvar, &
1498 & lbi, ubi, lbj, ubj, 0, n(ng), &
1499 & scale, mixing(ng)%Akt(:,:,:,isalt), &
1500 & nposw, xposw, yposw, zposw, rsta)
1501 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1502 & trim(vname(1,idsdif)), rsta, &
1503 & (/1,1,sta(ng)%Rindex/), &
1504 & (/n(ng)+1,nstation(ng),1/), &
1505 & piofile = sta(ng)%pioFile, &
1506 & piovar = sta(ng)%pioVar(idsdif)%vd)
1507 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1508 END IF
1509# endif
1510# if defined GLS_MIXING || defined MY25_MIXING
1511!
1512! Write out turbulent kinetic energy.
1513!
1514 IF (sout(idmtke,ng)) THEN
1515 scale=1.0_dp
1516 CALL extract_sta3d (ng, model, cgrid, idmtke, w3dvar, &
1517 & lbi, ubi, lbj, ubj, 0, n(ng), &
1518 & scale, mixing(ng)%tke(:,:,:,nout), &
1519 & nposw, xposw, yposw, zposw, rsta)
1520 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1521 & trim(vname(1,idmtke)), rsta, &
1522 & (/1,1,sta(ng)%Rindex/), &
1523 & (/n(ng)+1,nstation(ng),1/), &
1524 & piofile = sta(ng)%pioFile, &
1525 & piovar = sta(ng)%pioVar(idmtke)%vd)
1526 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1527 END IF
1528!
1529! Write out turbulent kinetic energy times length scale.
1530!
1531 IF (sout(idmtls,ng)) THEN
1532 scale=1.0_dp
1533 CALL extract_sta3d (ng, model, cgrid, idmtls, w3dvar, &
1534 & lbi, ubi, lbj, ubj, 0, n(ng), &
1535 & scale, mixing(ng)%gls(:,:,:,nout), &
1536 & nposw, xposw, yposw, zposw, rsta)
1537 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1538 & trim(vname(1,idmtls)), rsta, &
1539 & (/1,1,sta(ng)%Rindex/), &
1540 & (/n(ng)+1,nstation(ng),1/), &
1541 & piofile = sta(ng)%pioFile, &
1542 & piovar = sta(ng)%pioVar(idmtls)%vd)
1543 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1544 END IF
1545# endif
1546# if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS
1547!
1548! Write out surface air pressure.
1549!
1550 IF (sout(idpair,ng)) THEN
1551 scale=1.0_dp
1552 CALL extract_sta2d (ng, model, cgrid, idpair, r2dvar, &
1553 & lbi, ubi, lbj, ubj, &
1554 & scale, forces(ng)%Pair, &
1555 & nstation(ng), xpos, ypos, psta)
1556 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1557 & trim(vname(1,idpair)), psta, &
1558 & (/1,sta(ng)%Rindex/), &
1559 & (/nstation(ng),1/), &
1560 & piofile = sta(ng)%pioFile, &
1561 & piovar = sta(ng)%pioVar(idpair)%vd)
1562 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1563 END IF
1564# endif
1565# if defined BULK_FLUXES || defined ECOSIM
1566!
1567! Write out surface winds.
1568!
1569 IF (sout(iduair,ng)) THEN
1570 scale=1.0_dp
1571 CALL extract_sta2d (ng, model, cgrid, iduair, r2dvar, &
1572 & lbi, ubi, lbj, ubj, &
1573 & scale, forces(ng)%Uwind, &
1574 & nstation(ng), xpos, ypos, psta)
1575 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1576 & trim(vname(1,iduair)), psta, &
1577 & (/1,sta(ng)%Rindex/), &
1578 & (/nstation(ng),1/), &
1579 & piofile = sta(ng)%pioFile, &
1580 & piovar = sta(ng)%pioVar(iduair)%vd)
1581 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1582 END IF
1583!
1584 IF (sout(idvair,ng)) THEN
1585 scale=1.0_dp
1586 CALL extract_sta2d (ng, model, cgrid, idvair, r2dvar, &
1587 & lbi, ubi, lbj, ubj, &
1588 & scale, forces(ng)%Vwind, &
1589 & nstation(ng), xpos, ypos, psta)
1590 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1591 & trim(vname(1,idvair)), psta, &
1592 & (/1,sta(ng)%Rindex/), &
1593 & (/nstation(ng),1/), &
1594 & piofile = sta(ng)%pioFile, &
1595 & piovar = sta(ng)%pioVar(idvair)%vd)
1596 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1597 END IF
1598!
1599! Write out 2D Eastward and Northward surface winds (m/s) at
1600! RHO-points
1601!
1602 IF (sout(iduaie,ng).and.sout(idvain,ng)) THEN
1603 IF (.not.allocated(ur2d)) THEN
1604 allocate (ur2d(lbi:ubi,lbj:ubj))
1605 ur2d(lbi:ubi,lbj:ubj)=0.0_r8
1606 END IF
1607 IF (.not.allocated(vr2d)) THEN
1608 allocate (vr2d(lbi:ubi,lbj:ubj))
1609 vr2d(lbi:ubi,lbj:ubj)=0.0_r8
1610 END IF
1611 CALL uv_rotate2d (ng, tile, .false., .true., &
1612 & lbi, ubi, lbj, ubj, &
1613 & grid(ng) % CosAngler, &
1614 & grid(ng) % SinAngler, &
1615# ifdef MASKING
1616 & grid(ng) % rmask_full, &
1617# endif
1618 & forces(ng) % Uwind, &
1619 & forces(ng) % Vwind, &
1620 & ur2d, vr2d)
1621!
1622 scale=1.0_dp
1623 CALL extract_sta2d (ng, model, cgrid, iduaie, r2dvar, &
1624 & lbi, ubi, lbj, ubj, &
1625 & scale, ur2d, &
1626 & nstation(ng), xpos, ypos, psta)
1627 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1628 & trim(vname(1,iduaie)), psta, &
1629 & (/1,sta(ng)%Rindex/), &
1630 & (/nstation(ng),1/), &
1631 & piofile = sta(ng)%pioFile, &
1632 & piovar = sta(ng)%pioVar(iduaie)%vd)
1633 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1634!
1635 CALL extract_sta2d (ng, model, cgrid, idvain, r2dvar, &
1636 & lbi, ubi, lbj, ubj, &
1637 & scale, vr2d, &
1638 & nstation(ng), xpos, ypos, psta)
1639 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1640 & trim(vname(1,idvain)), psta, &
1641 & (/1,sta(ng)%Rindex/), &
1642 & (/nstation(ng),1/), &
1643 & piofile = sta(ng)%pioFile, &
1644 & piovar = sta(ng)%pioVar(idvain)%vd)
1645 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1646
1647 deallocate (ur2d)
1648 deallocate (vr2d)
1649 END IF
1650# endif
1651!
1652! Write out surface net heat flux.
1653!
1654 IF (sout(idtsur(itemp),ng)) THEN
1655 ifield=idtsur(itemp)
1656 scale=rho0*cp
1657 CALL extract_sta2d (ng, model, cgrid, ifield, r2dvar, &
1658 & lbi, ubi, lbj, ubj, &
1659 & scale, forces(ng)%stflx(:,:,itemp), &
1660 & nstation(ng), xpos, ypos, psta)
1661 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1662 & trim(vname(1,ifield)), psta, &
1663 & (/1,sta(ng)%Rindex/), &
1664 & (/nstation(ng),1/), &
1665 & piofile = sta(ng)%pioFile, &
1666 & piovar = sta(ng)%pioVar(ifield)%vd)
1667 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1668 END IF
1669
1670# ifdef SALINITY
1671!
1672! Write out surface salt flux.
1673!
1674 IF (sout(idtsur(isalt),ng)) THEN
1675 ifield=idtsur(isalt)
1676 scale=1.0_dp
1677 CALL extract_sta2d (ng, model, cgrid, ifield, r2dvar, &
1678 & lbi, ubi, lbj, ubj, &
1679 & scale, forces(ng)%stflx(:,:,isalt), &
1680 & nstation(ng), xpos, ypos, psta)
1681 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1682 & trim(vname(1,ifield)), psta, &
1683 & (/1,sta(ng)%Rindex/), &
1684 & (/nstation(ng),1/), &
1685 & piofile = sta(ng)%pioFile, &
1686 & piovar = sta(ng)%pioVar(ifield)%vd)
1687 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1688 END IF
1689# endif
1690
1691# ifdef BULK_FLUXES
1692!
1693! Write out latent heat flux.
1694!
1695 IF (sout(idlhea,ng)) THEN
1696 scale=rho0*cp
1697 CALL extract_sta2d (ng, model, cgrid, idlhea, r2dvar, &
1698 & lbi, ubi, lbj, ubj, &
1699 & scale, forces(ng)%lhflx, &
1700 & nstation(ng), xpos, ypos, psta)
1701 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1702 & trim(vname(1,idlhea)), psta, &
1703 & (/1,sta(ng)%Rindex/), &
1704 & (/nstation(ng),1/), &
1705 & piofile = sta(ng)%pioFile, &
1706 & piovar = sta(ng)%pioVar(idlhea)%vd)
1707 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1708 END IF
1709!
1710! Write out sensible heat flux.
1711!
1712 IF (sout(idshea,ng)) THEN
1713 scale=rho0*cp
1714 CALL extract_sta2d (ng, model, cgrid, idshea, r2dvar, &
1715 & lbi, ubi, lbj, ubj, &
1716 & scale, forces(ng)%shflx, &
1717 & nstation(ng), xpos, ypos, psta)
1718 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1719 & trim(vname(1,idshea)), psta, &
1720 & (/1,sta(ng)%Rindex/), &
1721 & (/nstation(ng),1/), &
1722 & piofile = sta(ng)%pioFile, &
1723 & piovar = sta(ng)%pioVar(idshea)%vd)
1724 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1725 END IF
1726!
1727! Write out longwave radiation flux.
1728!
1729 IF (sout(idlrad,ng)) THEN
1730 scale=rho0*cp
1731 CALL extract_sta2d (ng, model, cgrid, idlrad, r2dvar, &
1732 & lbi, ubi, lbj, ubj, &
1733 & scale, forces(ng)%lrflx, &
1734 & nstation(ng), xpos, ypos, psta)
1735 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1736 & trim(vname(1,idlrad)), psta, &
1737 & (/1,sta(ng)%Rindex/), &
1738 & (/nstation(ng),1/), &
1739 & piofile = sta(ng)%pioFile, &
1740 & piovar = sta(ng)%pioVar(idlrad)%vd)
1741 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1742 END IF
1743
1744# endif
1745# ifdef SHORTWAVE
1746!
1747! Write out shortwave radiation flux.
1748!
1749 IF (sout(idsrad,ng)) THEN
1750 scale=rho0*cp
1751 CALL extract_sta2d (ng, model, cgrid, idsrad, r2dvar, &
1752 & lbi, ubi, lbj, ubj, &
1753 & scale, forces(ng)%srflx, &
1754 & nstation(ng), xpos, ypos, psta)
1755 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1756 & trim(vname(1,idsrad)), psta, &
1757 & (/1,sta(ng)%Rindex/), &
1758 & (/nstation(ng),1/), &
1759 & piofile = sta(ng)%pioFile, &
1760 & piovar = sta(ng)%pioVar(idsrad)%vd)
1761 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1762 END IF
1763# endif
1764# if defined EMINUSP && defined BULK_FLUXES
1765!
1766! Write out E-P (m/s).
1767!
1768 IF (sout(idempf,ng)) THEN
1769 scale=1.0_dp
1770 CALL extract_sta2d (ng, model, cgrid, idempf, r2dvar, &
1771 & lbi, ubi, lbj, ubj, &
1772 & scale, forces(ng)%stflux(:,:,isalt), &
1773 & nstation(ng), xpos, ypos, psta)
1774 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1775 & trim(vname(1,idempf)), psta, &
1776 & (/1,sta(ng)%Rindex/), &
1777 & (/nstation(ng),1/), &
1778 & piofile = sta(ng)%pioFile, &
1779 & piovar = sta(ng)%pioVar(idempf)%vd)
1780 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1781 END IF
1782!
1783! Write out evaportaion rate (kg/m2/s).
1784!
1785 IF (sout(idevap,ng)) THEN
1786 scale=1.0_dp
1787 CALL extract_sta2d (ng, model, cgrid, idevap, r2dvar, &
1788 & lbi, ubi, lbj, ubj, &
1789 & scale, forces(ng)%evap, &
1790 & nstation(ng), xpos, ypos, psta)
1791 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1792 & trim(vname(1,idevap)), psta, &
1793 & (/1,sta(ng)%Rindex/), &
1794 & (/nstation(ng),1/), &
1795 & piofile = sta(ng)%pioFile, &
1796 & piovar = sta(ng)%pioVar(idevap)%vd)
1797 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1798 END IF
1799!
1800! Write out precipitation rate (kg/m2/s).
1801!
1802 IF (sout(idrain,ng)) THEN
1803 scale=1.0_dp
1804 CALL extract_sta2d (ng, model, cgrid, idrain, r2dvar, &
1805 & lbi, ubi, lbj, ubj, &
1806 & scale, forces(ng)%rain, &
1807 & nstation(ng), xpos, ypos, psta)
1808 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1809 & trim(vname(1,idrain)), psta, &
1810 & (/1,sta(ng)%Rindex/), &
1811 & (/nstation(ng),1/), &
1812 & piofile = sta(ng)%pioFile, &
1813 & piovar = sta(ng)%pioVar(idrain)%vd)
1814 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1815 END IF
1816# endif
1817# endif
1818!
1819! Write out surface U-momentum stress.
1820!
1821 IF (sout(idusms,ng)) THEN
1822 scale=rho0
1823 CALL extract_sta2d (ng, model, cgrid, idusms, u2dvar, &
1824 & lbi, ubi, lbj, ubj, &
1825 & scale, forces(ng)%sustr, &
1826 & nstation(ng), xpos, ypos, psta)
1827 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1828 & trim(vname(1,idusms)), psta, &
1829 & (/1,sta(ng)%Rindex/), &
1830 & (/nstation(ng),1/), &
1831 & piofile = sta(ng)%pioFile, &
1832 & piovar = sta(ng)%pioVar(idusms)%vd)
1833 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1834 END IF
1835!
1836! Write out surface V-momentum stress.
1837!
1838 IF (sout(idvsms,ng)) THEN
1839 scale=rho0
1840 CALL extract_sta2d (ng, model, cgrid, idvsms, v2dvar, &
1841 & lbi, ubi, lbj, ubj, &
1842 & scale, forces(ng)%svstr, &
1843 & nstation(ng), xpos, ypos, psta)
1844 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1845 & trim(vname(1,idvsms)), psta, &
1846 & (/1,sta(ng)%Rindex/), &
1847 & (/nstation(ng),1/), &
1848 & piofile = sta(ng)%pioFile, &
1849 & piovar = sta(ng)%pioVar(idvsms)%vd)
1850 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1851 END IF
1852!
1853! Write out bottom U-momentum stress.
1854!
1855 IF (sout(idubms,ng)) THEN
1856 scale=-rho0
1857 CALL extract_sta2d (ng, model, cgrid, idubms, u2dvar, &
1858 & lbi, ubi, lbj, ubj, &
1859 & scale, forces(ng)%bustr, &
1860 & nstation(ng), xpos, ypos, psta)
1861 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1862 & trim(vname(1,idubms)), psta, &
1863 & (/1,sta(ng)%Rindex/), &
1864 & (/nstation(ng),1/), &
1865 & piofile = sta(ng)%pioFile, &
1866 & piovar = sta(ng)%pioVar(idubms)%vd)
1867 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1868 END IF
1869!
1870! Write out bottom V-momentum stress.
1871!
1872 IF (sout(idvbms,ng)) THEN
1873 scale=-rho0
1874 CALL extract_sta2d (ng, model, cgrid, idvbms, v2dvar, &
1875 & lbi, ubi, lbj, ubj, &
1876 & scale, forces(ng)%bvstr, &
1877 & nstation(ng), xpos, ypos, psta)
1878 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1879 & trim(vname(1,idvbms)), psta, &
1880 & (/1,sta(ng)%Rindex/), &
1881 & (/nstation(ng),1/), &
1882 & piofile = sta(ng)%pioFile, &
1883 & piovar = sta(ng)%pioVar(idvbms)%vd)
1884 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1885 END IF
1886
1887# ifdef WET_DRY
1888!
1889! Write out wet/dry mask at RHO-points.
1890!
1891 IF (sout(idrwet,ng)) THEN
1892 scale=1.0_dp
1893 CALL extract_sta2d (ng, model, cgrid, idrwet, r2dvar, &
1894 & lbi, ubi, lbj, ubj, &
1895 & scale, grid(ng)%rmask_wet, &
1896 & nstation(ng), xpos, ypos, psta)
1897 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1898 & trim(vname(1,idrwet)), psta, &
1899 & (/1,sta(ng)%Rindex/), &
1900 & (/nstation(ng),1/), &
1901 & piofile = sta(ng)%pioFile, &
1902 & piovar = sta(ng)%pioVar(idrwet)%vd)
1903 END IF
1904!
1905! Write out wet/dry mask at U-points.
1906!
1907 IF (sout(iduwet,ng)) THEN
1908 scale=1.0_dp
1909 CALL extract_sta2d (ng, model, cgrid, iduwet, u2dvar, &
1910 & lbi, ubi, lbj, ubj, &
1911 & scale, grid(ng)%umask_wet, &
1912 & nstation(ng), xpos, ypos, psta)
1913 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1914 & trim(vname(1,iduwet)), psta, &
1915 & (/1,sta(ng)%Rindex/), &
1916 & (/nstation(ng),1/), &
1917 & piofile = sta(ng)%pioFile, &
1918 & piovar = sta(ng)%pioVar(iduwet)%vd)
1919 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1920 END IF
1921!
1922! Write out wet/dry mask at V-points.
1923!
1924 IF (sout(idvwet,ng)) THEN
1925 scale=1.0_dp
1926 CALL extract_sta2d (ng, model, cgrid, idvwet, v2dvar, &
1927 & lbi, ubi, lbj, ubj, &
1928 & scale, grid(ng)%vmask_wet, &
1929 & nstation(ng), xpos, ypos, psta)
1930 CALL pio_netcdf_put_fvar (ng, model, sta(ng)%name, &
1931 & trim(vname(1,idvwet)), psta, &
1932 & (/1,sta(ng)%Rindex/), &
1933 & (/nstation(ng),1/), &
1934 & piofile = sta(ng)%pioFile, &
1935 & piovar = sta(ng)%pioVar(idvwet)%vd)
1936 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1937 END IF
1938# endif
1939
1940# if defined BBL_MODEL || defined WAVES_OUTPUT
1941!
1942!-----------------------------------------------------------------------
1943! Write out the bottom boundary layer model or waves variables.
1944!-----------------------------------------------------------------------
1945!
1946 CALL bbl_wrt_station_pio (ng, model, tile, &
1947 & lbi, ubi, lbj, ubj, &
1948 & sout, sta)
1949 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1950# endif
1951
1952# ifdef ICE_MODEL
1953!
1954!-----------------------------------------------------------------------
1955! Write out the ice model variables.
1956!-----------------------------------------------------------------------
1957!
1958 CALL ice_wrt_station_pio (ng, model, tile, &
1959 & lbi, ubi, lbj, ubj, &
1960 & sout, sta)
1961 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1962# endif
1963
1964# ifdef SEDIMENT
1965!
1966!-----------------------------------------------------------------------
1967! Write out the sediment model variables.
1968!-----------------------------------------------------------------------
1969!
1970 CALL sediment_wrt_station_pio (ng, model, tile, &
1971 & lbi, ubi, lbj, ubj, &
1972 & sout, sta)
1973 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1974# endif
1975
1976# if defined WEC || defined WEC_VF
1977!
1978!-----------------------------------------------------------------------
1979! Write out the Waves Effect on Currents station variables.
1980!-----------------------------------------------------------------------
1981!
1982 CALL wec_wrt_station_pio (ng, model, tile, &
1983 & lbi, ubi, lbj, ubj, &
1984 & sout, sta)
1985 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1986# endif
1987!
1988!-----------------------------------------------------------------------
1989! Synchronize stations NetCDF file to disk.
1990!-----------------------------------------------------------------------
1991!
1992 CALL pio_netcdf_sync (ng, model, sta(ng)%name, sta(ng)%pioFile)
1993
1994 RETURN
1995 END SUBROUTINE wrt_station_pio
1996# endif
1997#endif
1998 END MODULE wrt_station_mod
1999
subroutine, public bbl_wrt_station_pio(ng, model, tile, lbi, ubi, lbj, ubj, varout, s)
subroutine, public bbl_wrt_station_nf90(ng, model, tile, lbi, ubi, lbj, ubj, varout, s)
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_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_io), dimension(:), allocatable sta
integer stdout
character(len=256) sourcefile
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
integer iddano
integer idvair
integer idevap
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer idubar
integer idwvel
integer idvvel
integer idhsbl
integer idvsms
integer, parameter io_pio
Definition mod_ncparam.F:96
integer idpair
integer idrwet
integer idv2dn
integer idsdif
integer, dimension(:), allocatable idtsur
integer idempf
integer idvain
integer idtdif
integer idfsur
integer, dimension(:), allocatable idtvar
integer idhbbl
integer idvbms
integer iduair
integer idmtke
integer iduvel
integer idv3dn
integer idovel
integer iduwet
character(len=maxlen), dimension(6, 0:nv) vname
integer idtime
integer idshea
integer idlrad
integer idusms
logical, dimension(:,:), allocatable sout
integer idvvis
integer idu3de
integer idu2de
integer idlhea
integer idrain
integer idubms
integer idvwet
integer idsrad
integer idmtls
integer iduaie
integer idvbar
subroutine, public netcdf_sync(ng, model, ncname, ncid)
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
logical master
integer, parameter inlm
Definition mod_param.F:662
integer nbed
Definition mod_param.F:517
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer, parameter r3dvar
Definition mod_param.F:721
integer, dimension(:), allocatable nstation
Definition mod_param.F:557
integer, parameter u3dvar
Definition mod_param.F:722
integer, parameter u2dvar
Definition mod_param.F:718
integer, parameter w3dvar
Definition mod_param.F:724
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer, parameter r2dvar
Definition mod_param.F:717
integer, parameter v2dvar
Definition mod_param.F:719
integer, parameter v3dvar
Definition mod_param.F:723
subroutine, public pio_netcdf_sync(ng, model, ncname, piofile)
real(dp) cp
integer exit_flag
integer isalt
integer itemp
type(t_scalars), dimension(:), allocatable scalars
Definition mod_scalars.F:65
real(dp), dimension(:), allocatable time
real(dp) rho0
integer noerror
subroutine, public sediment_wrt_station_nf90(ng, model, tile, lbi, ubi, lbj, ubj, varout, s)
subroutine, public sediment_wrt_station_pio(ng, model, tile, lbi, ubi, lbj, ubj, varout, s)
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine, public uv_rotate3d(ng, tile, add, lboundary, lbi, ubi, lbj, ubj, lbk, ubk, cosangler, sinangler, rmask_full, uinp, vinp, uout, vout)
Definition uv_rotate.F:155
subroutine, public uv_rotate2d(ng, tile, add, lboundary, lbi, ubi, lbj, ubj, cosangler, sinangler, rmask_full, uinp, vinp, uout, vout)
Definition uv_rotate.F:35
subroutine, private wrt_station_nf90(ng, model, tile, lbi, ubi, lbj, ubj)
subroutine, private wrt_station_pio(ng, model, tile, lbi, ubi, lbj, ubj)
subroutine, public wrt_station(ng, tile)
Definition wrt_station.F:81