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

Go to the source code of this file.

Functions/Subroutines

subroutine read_stapar (model, inp, out, lwrite)
 

Function/Subroutine Documentation

◆ read_stapar()

subroutine read_stapar ( integer, intent(in) model,
integer, intent(in) inp,
integer, intent(in) out,
logical, intent(in) lwrite )

Definition at line 3 of file read_stapar.F.

4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine reads and reports stations input parameters. !
13! !
14!=======================================================================
15!
16 USE mod_param
17 USE mod_parallel
18# ifdef ICE_MODEL
19 USE mod_ice
20# endif
21 USE mod_iounits
22 USE mod_ncparam
23# if defined SEDIMENT || defined BBL_MODEL
24 USE mod_sediment
25# endif
26 USE mod_scalars
27!
29!
30 implicit none
31!
32! Imported variable declarations
33!
34 logical, intent(in) :: Lwrite
35 integer, intent(in) :: model, inp, out
36!
37! Local variable declarations.
38!
39 integer :: Mstation, Npts, Nval
40 integer :: flag, i, igrid, ista, itrc, ng, status
41
42 real(r8) :: Xpos, Ypos
43
44 logical, dimension(Ngrids) :: Lswitch
45 logical, dimension(MT,Ngrids) :: Ltracer
46
47 integer, dimension(Ngrids) :: is
48
49 real(dp), dimension(nRval) :: Rval
50
51 character (len=40 ) :: KeyWord
52 character (len=256) :: line
53 character (len=256), dimension(nCval) :: Cval
54!
55!-----------------------------------------------------------------------
56! Read in stations parameters.
57!-----------------------------------------------------------------------
58!
59 DO WHILE (.true.)
60 READ (inp,'(a)',err=20,END=30) line
61 status=decode_line(line, keyword, nval, cval, rval)
62 IF (status.gt.0) THEN
63 SELECT CASE (trim(keyword))
64 CASE ('Lstations')
65 npts=load_l(nval, cval, ngrids, lstations)
66 CASE ('Sout(idUvel)')
67 npts=load_l(nval, cval, ngrids, lswitch)
68 sout(iduvel,1:ngrids)=lswitch(1:ngrids)
69 CASE ('Sout(idVvel)')
70 npts=load_l(nval, cval, ngrids, lswitch)
71 sout(idvvel,1:ngrids)=lswitch(1:ngrids)
72 CASE ('Sout(idu3dE)')
73 npts=load_l(nval, cval, ngrids, lswitch)
74 sout(idu3de,1:ngrids)=lswitch(1:ngrids)
75 CASE ('Sout(idv3dN)')
76 npts=load_l(nval, cval, ngrids, lswitch)
77 sout(idv3dn,1:ngrids)=lswitch(1:ngrids)
78 CASE ('Sout(idWvel)')
79 npts=load_l(nval, cval, ngrids, lswitch)
80 sout(idwvel,1:ngrids)=lswitch(1:ngrids)
81 CASE ('Sout(idOvel)')
82 npts=load_l(nval, cval, ngrids, lswitch)
83 sout(idovel,1:ngrids)=lswitch(1:ngrids)
84 CASE ('Sout(idUbar)')
85 npts=load_l(nval, cval, ngrids, lswitch)
86 sout(idubar,1:ngrids)=lswitch(1:ngrids)
87 CASE ('Sout(idVbar)')
88 npts=load_l(nval, cval, ngrids, lswitch)
89 sout(idvbar,1:ngrids)=lswitch(1:ngrids)
90 CASE ('Sout(idu2dE)')
91 npts=load_l(nval, cval, ngrids, lswitch)
92 sout(idu2de,1:ngrids)=lswitch(1:ngrids)
93 CASE ('Sout(idv2dN)')
94 npts=load_l(nval, cval, ngrids, lswitch)
95 sout(idv2dn,1:ngrids)=lswitch(1:ngrids)
96 CASE ('Sout(idFsur)')
97 npts=load_l(nval, cval, ngrids, lswitch)
98 sout(idfsur,1:ngrids)=lswitch(1:ngrids)
99# if defined SEDIMENT && defined SED_MORPH
100 CASE ('Sout(idBath)')
101 npts=load_l(nval, cval, ngrids, lswitch)
102 sout(idbath,1:ngrids)=lswitch(1:ngrids)
103# endif
104 CASE ('Sout(idTvar)')
105 npts=load_l(nval, cval, mt, ngrids, ltracer)
106 DO ng=1,ngrids
107 DO itrc=1,nt(ng)
108 sout(idtvar(itrc),ng)=ltracer(itrc,ng)
109 END DO
110 END DO
111 CASE ('Sout(idUsms)')
112 npts=load_l(nval, cval, ngrids, lswitch)
113 sout(idusms,1:ngrids)=lswitch(1:ngrids)
114 CASE ('Sout(idVsms)')
115 npts=load_l(nval, cval, ngrids, lswitch)
116 sout(idvsms,1:ngrids)=lswitch(1:ngrids)
117 CASE ('Sout(idUbms)')
118 npts=load_l(nval, cval, ngrids, lswitch)
119 sout(idubms,1:ngrids)=lswitch(1:ngrids)
120 CASE ('Sout(idVbms)')
121 npts=load_l(nval, cval, ngrids, lswitch)
122 sout(idvbms,1:ngrids)=lswitch(1:ngrids)
123# ifdef BBL_MODEL
124 CASE ('Sout(idUbrs)')
125 npts=load_l(nval, cval, ngrids, lswitch)
126 sout(idubrs,1:ngrids)=lswitch(1:ngrids)
127 CASE ('Sout(idVbrs)')
128 npts=load_l(nval, cval, ngrids, lswitch)
129 sout(idvbrs,1:ngrids)=lswitch(1:ngrids)
130 CASE ('Sout(idUbws)')
131 npts=load_l(nval, cval, ngrids, lswitch)
132 sout(idubws,1:ngrids)=lswitch(1:ngrids)
133 CASE ('Sout(idVbws)')
134 npts=load_l(nval, cval, ngrids, lswitch)
135 sout(idvbws,1:ngrids)=lswitch(1:ngrids)
136 CASE ('Sout(idUbcs)')
137 npts=load_l(nval, cval, ngrids, lswitch)
138 sout(idubcs,1:ngrids)=lswitch(1:ngrids)
139 CASE ('Sout(idVbcs)')
140 npts=load_l(nval, cval, ngrids, lswitch)
141 sout(idvbcs,1:ngrids)=lswitch(1:ngrids)
142 CASE ('Sout(idUbot)')
143 npts=load_l(nval, cval, ngrids, lswitch)
144 sout(idubot,1:ngrids)=lswitch(1:ngrids)
145 CASE ('Sout(idVbot)')
146 npts=load_l(nval, cval, ngrids, lswitch)
147 sout(idvbot,1:ngrids)=lswitch(1:ngrids)
148 CASE ('Sout(idUbur)')
149 npts=load_l(nval, cval, ngrids, lswitch)
150 sout(idubur,1:ngrids)=lswitch(1:ngrids)
151 CASE ('Sout(idVbvr)')
152 npts=load_l(nval, cval, ngrids, lswitch)
153 sout(idvbvr,1:ngrids)=lswitch(1:ngrids)
154# endif
155# ifdef ICE_MODEL
156 CASE ('Sout(idUice)')
157 npts=load_l(nval, cval, ngrids, lswitch)
158 sout(iduice,1:ngrids)=lswitch(1:ngrids)
159 CASE ('Sout(idVice)')
160 npts=load_l(nval, cval, ngrids, lswitch)
161 sout(idvice,1:ngrids)=lswitch(1:ngrids)
162 CASE ('Sout(idUiER)')
163 npts=load_l(nval, cval, ngrids, lswitch)
164 sout(iduier,1:ngrids)=lswitch(1:ngrids)
165 CASE ('Sout(idViNR)')
166 npts=load_l(nval, cval, ngrids, lswitch)
167 sout(idvinr,1:ngrids)=lswitch(1:ngrids)
168 CASE ('Sout(idAice)')
169 npts=load_l(nval, cval, ngrids, lswitch)
170 sout(idaice,1:ngrids)=lswitch(1:ngrids)
171 CASE ('Sout(idIage)')
172 npts=load_l(nval, cval, ngrids, lswitch)
173 sout(idiage,1:ngrids)=lswitch(1:ngrids)
174 CASE ('Sout(idHice)')
175 npts=load_l(nval, cval, ngrids, lswitch)
176 sout(idhice,1:ngrids)=lswitch(1:ngrids)
177 CASE ('Sout(idHmel)')
178 npts=load_l(nval, cval, ngrids, lswitch)
179 sout(idhmel,1:ngrids)=lswitch(1:ngrids)
180 CASE ('Sout(idHsno)')
181 npts=load_l(nval, cval, ngrids, lswitch)
182 sout(idhsno,1:ngrids)=lswitch(1:ngrids)
183 CASE ('Sout(idTice)')
184 npts=load_l(nval, cval, ngrids, lswitch)
185 sout(idtice,1:ngrids)=lswitch(1:ngrids)
186 CASE ('Sout(idISxx)')
187 npts=load_l(nval, cval, ngrids, lswitch)
188 sout(idisxx,1:ngrids)=lswitch(1:ngrids)
189 CASE ('Sout(idISxy)')
190 npts=load_l(nval, cval, ngrids, lswitch)
191 sout(idisxy,1:ngrids)=lswitch(1:ngrids)
192 CASE ('Sout(idISyy)')
193 npts=load_l(nval, cval, ngrids, lswitch)
194 sout(idisyy,1:ngrids)=lswitch(1:ngrids)
195 CASE ('Sout(idIsst)')
196 npts=load_l(nval, cval, ngrids, lswitch)
197 sout(idisst,1:ngrids)=lswitch(1:ngrids)
198 CASE ('Sout(idIOmf)')
199 npts=load_l(nval, cval, ngrids, lswitch)
200 sout(idiomf,1:ngrids)=lswitch(1:ngrids)
201 CASE ('Sout(idIOfv)')
202 npts=load_l(nval, cval, ngrids, lswitch)
203 sout(idiofv,1:ngrids)=lswitch(1:ngrids)
204 CASE ('Sout(idIOmt)')
205 npts=load_l(nval, cval, ngrids, lswitch)
206 sout(idiomt,1:ngrids)=lswitch(1:ngrids)
207 CASE ('Sout(idS0mk)')
208 npts=load_l(nval, cval, ngrids, lswitch)
209 sout(ids0mk,1:ngrids)=lswitch(1:ngrids)
210 CASE ('Sout(idT0mk)')
211 npts=load_l(nval, cval, ngrids, lswitch)
212 sout(idt0mk,1:ngrids)=lswitch(1:ngrids)
213 CASE ('Sout(idWdiv)')
214 npts=load_l(nval, cval, ngrids, lswitch)
215 sout(idwdiv,1:ngrids)=lswitch(1:ngrids)
216 CASE ('Sout(idW_fr)')
217 npts=load_l(nval, cval, ngrids, lswitch)
218 sout(idw_fr,1:ngrids)=lswitch(1:ngrids)
219 CASE ('Sout(idW_ai)')
220 npts=load_l(nval, cval, ngrids, lswitch)
221 sout(idw_ai,1:ngrids)=lswitch(1:ngrids)
222 CASE ('Sout(idW_ao)')
223 npts=load_l(nval, cval, ngrids, lswitch)
224 sout(idw_ao,1:ngrids)=lswitch(1:ngrids)
225 CASE ('Sout(idW_io)')
226 npts=load_l(nval, cval, ngrids, lswitch)
227 sout(idw_io,1:ngrids)=lswitch(1:ngrids)
228 CASE ('Sout(idW_ro)')
229 npts=load_l(nval, cval, ngrids, lswitch)
230 sout(idw_ro,1:ngrids)=lswitch(1:ngrids)
231# endif
232# ifdef WEC_VF
233 CASE ('Sout(idWztw)')
234 npts=load_l(nval, cval, ngrids, lswitch)
235 sout(idwztw,1:ngrids)=lswitch(1:ngrids)
236 CASE ('Sout(idWqsp)')
237 npts=load_l(nval, cval, ngrids, lswitch)
238 sout(idwqsp,1:ngrids)=lswitch(1:ngrids)
239 CASE ('Sout(idWbeh)')
240 npts=load_l(nval, cval, ngrids, lswitch)
241 sout(idwbeh,1:ngrids)=lswitch(1:ngrids)
242# endif
243# ifdef WEC
244 CASE ('Sout(idU2rs)')
245 npts=load_l(nval, cval, ngrids, lswitch)
246 sout(idu2rs,1:ngrids)=lswitch(1:ngrids)
247 CASE ('Sout(idV2rs)')
248 npts=load_l(nval, cval, ngrids, lswitch)
249 sout(idv2rs,1:ngrids)=lswitch(1:ngrids)
250 CASE ('Sout(idU3rs)')
251 npts=load_l(nval, cval, ngrids, lswitch)
252 sout(idu3rs,1:ngrids)=lswitch(1:ngrids)
253 CASE ('Sout(idV3rs)')
254 npts=load_l(nval, cval, ngrids, lswitch)
255 sout(idv3rs,1:ngrids)=lswitch(1:ngrids)
256 CASE ('Sout(idU2Sd)')
257 npts=load_l(nval, cval, ngrids, lswitch)
258 sout(idu2sd,1:ngrids)=lswitch(1:ngrids)
259 CASE ('Sout(idV2Sd)')
260 npts=load_l(nval, cval, ngrids, lswitch)
261 sout(idv2sd,1:ngrids)=lswitch(1:ngrids)
262 CASE ('Sout(idU3Sd)')
263 npts=load_l(nval, cval, ngrids, lswitch)
264 sout(idu3sd,1:ngrids)=lswitch(1:ngrids)
265 CASE ('Sout(idV3Sd)')
266 npts=load_l(nval, cval, ngrids, lswitch)
267 sout(idv3sd,1:ngrids)=lswitch(1:ngrids)
268 CASE ('Sout(idW3Sd)')
269 npts=load_l(nval, cval, ngrids, lswitch)
270 sout(idw3sd,1:ngrids)=lswitch(1:ngrids)
271 CASE ('Sout(idW3St)')
272 npts=load_l(nval, cval, ngrids, lswitch)
273 sout(idw3st,1:ngrids)=lswitch(1:ngrids)
274# endif
275# ifdef WAVES_HEIGHT
276 CASE ('Sout(idWamp)')
277 npts=load_l(nval, cval, ngrids, lswitch)
278 sout(idwamp,1:ngrids)=lswitch(1:ngrids)
279# endif
280# ifdef WAVES_LENGTH
281 CASE ('Sout(idWlen)')
282 npts=load_l(nval, cval, ngrids, lswitch)
283 sout(idwlen,1:ngrids)=lswitch(1:ngrids)
284# endif
285# ifdef WAVES_LENGTHP
286 CASE ('Sout(idWlep)')
287 npts=load_l(nval, cval, ngrids, lswitch)
288 sout(idwlep,1:ngrids)=lswitch(1:ngrids)
289# endif
290# ifdef WAVES_DIR
291 CASE ('Sout(idWdir)')
292 npts=load_l(nval, cval, ngrids, lswitch)
293 sout(idwdir,1:ngrids)=lswitch(1:ngrids)
294# endif
295# ifdef WAVES_TOP_PERIOD
296 CASE ('Sout(idWptp)')
297 npts=load_l(nval, cval, ngrids, lswitch)
298 sout(idwptp,1:ngrids)=lswitch(1:ngrids)
299# endif
300# ifdef WAVES_DIRP
301 CASE ('Sout(idWdip)')
302 npts=load_l(nval, cval, ngrids, lswitch)
303 sout(idwdip,1:ngrids)=lswitch(1:ngrids)
304# endif
305# ifdef WAVES_BOT_PERIOD
306 CASE ('Sout(idWpbt)')
307 npts=load_l(nval, cval, ngrids, lswitch)
308 sout(idwpbt,1:ngrids)=lswitch(1:ngrids)
309# endif
310# if defined BBL_MODEL || defined BEDLOAD_SOULSBY || \
311 defined bedload_vandera || defined wav_coupling
312 CASE ('Sout(idWorb)')
313 npts=load_l(nval, cval, ngrids, lswitch)
314 sout(idworb,1:ngrids)=lswitch(1:ngrids)
315# endif
316# if defined WAV_COUPLING || (defined WEC_VF && defined BOTTOM_STREAMING)
317 CASE ('Sout(idWdif)')
318 npts=load_l(nval, cval, ngrids, lswitch)
319 sout(idwdif,1:ngrids)=lswitch(1:ngrids)
320
321# endif
322# if defined TKE_WAVEDISS || defined WAV_COUPLING || \
323 defined wdiss_thorguza || defined wdiss_churthor || \
324 defined waves_diss || defined wdiss_inwave
325 CASE ('Sout(idWdib)')
326 npts=load_l(nval, cval, ngrids, lswitch)
327 sout(idwdib,1:ngrids)=lswitch(1:ngrids)
328 CASE ('Sout(idWdiw)')
329 npts=load_l(nval, cval, ngrids, lswitch)
330 sout(idwdiw,1:ngrids)=lswitch(1:ngrids)
331# endif
332# if defined ROLLER_SVENDSEN
333 CASE ('Sout(idWbrk)')
334 npts=load_l(nval, cval, ngrids, lswitch)
335 sout(idwbrk,1:ngrids)=lswitch(1:ngrids)
336# endif
337# if defined WEC_ROLLER
338 CASE ('Sout(idWdis)')
339 npts=load_l(nval, cval, ngrids, lswitch)
340 sout(idwdis,1:ngrids)=lswitch(1:ngrids)
341# endif
342# if defined ROLLER_RENIERS
343 CASE ('Sout(idWrol)')
344 npts=load_l(nval, cval, ngrids, lswitch)
345 sout(idwrol,1:ngrids)=lswitch(1:ngrids)
346# endif
347# if defined WAVES_DSPR
348 CASE ('Sout(idWvds)')
349 npts=load_l(nval, cval, ngrids, lswitch)
350 sout(idwvds,1:ngrids)=lswitch(1:ngrids)
351 CASE ('Sout(idWvqp)')
352 npts=load_l(nval, cval, ngrids, lswitch)
353 sout(idwvqp,1:ngrids)=lswitch(1:ngrids)
354# endif
355# ifdef SOLVE3D
356# if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS
357 CASE ('Sout(idPair)')
358 npts=load_l(nval, cval, ngrids, lswitch)
359 sout(idpair,1:ngrids)=lswitch(1:ngrids)
360# endif
361# if defined BULK_FLUXES || defined ECOSIM
362 CASE ('Sout(idUair)')
363 npts=load_l(nval, cval, ngrids, lswitch)
364 sout(iduair,1:ngrids)=lswitch(1:ngrids)
365 CASE ('Sout(idVair)')
366 npts=load_l(nval, cval, ngrids, lswitch)
367 sout(idvair,1:ngrids)=lswitch(1:ngrids)
368# endif
369 CASE ('Sout(idTsur)')
370 npts=load_l(nval, cval, mt, ngrids, ltracer)
371 DO ng=1,ngrids
372 DO itrc=1,nat
373 sout(idtsur(itrc),ng)=ltracer(itrc,ng)
374 END DO
375 END DO
376 CASE ('Sout(idLhea)')
377 npts=load_l(nval, cval, ngrids, lswitch)
378 sout(idlhea,1:ngrids)=lswitch(1:ngrids)
379 CASE ('Sout(idShea)')
380 npts=load_l(nval, cval, ngrids, lswitch)
381 sout(idshea,1:ngrids)=lswitch(1:ngrids)
382 CASE ('Sout(idLrad)')
383 npts=load_l(nval, cval, ngrids, lswitch)
384 sout(idlrad,1:ngrids)=lswitch(1:ngrids)
385 CASE ('Sout(idSrad)')
386 npts=load_l(nval, cval, ngrids, lswitch)
387 sout(idsrad,1:ngrids)=lswitch(1:ngrids)
388 CASE ('Sout(idEmPf)')
389 npts=load_l(nval, cval, ngrids, lswitch)
390 sout(idempf,1:ngrids)=lswitch(1:ngrids)
391 CASE ('Sout(idevap)')
392 npts=load_l(nval, cval, ngrids, lswitch)
393 sout(idevap,1:ngrids)=lswitch(1:ngrids)
394 CASE ('Sout(idrain)')
395 npts=load_l(nval, cval, ngrids, lswitch)
396 sout(idrain,1:ngrids)=lswitch(1:ngrids)
397 CASE ('Sout(idDano)')
398 npts=load_l(nval, cval, ngrids, lswitch)
399 sout(iddano,1:ngrids)=lswitch(1:ngrids)
400 CASE ('Sout(idVvis)')
401 npts=load_l(nval, cval, ngrids, lswitch)
402 sout(idvvis,1:ngrids)=lswitch(1:ngrids)
403 CASE ('Sout(idTdif)')
404 npts=load_l(nval, cval, ngrids, lswitch)
405 sout(idtdif,1:ngrids)=lswitch(1:ngrids)
406 CASE ('Sout(idSdif)')
407 npts=load_l(nval, cval, ngrids, lswitch)
408 sout(idsdif,1:ngrids)=lswitch(1:ngrids)
409 CASE ('Sout(idHsbl)')
410 npts=load_l(nval, cval, ngrids, lswitch)
411 sout(idhsbl,1:ngrids)=lswitch(1:ngrids)
412 CASE ('Sout(idHbbl)')
413 npts=load_l(nval, cval, ngrids, lswitch)
414 sout(idhbbl,1:ngrids)=lswitch(1:ngrids)
415 CASE ('Sout(idMtke)')
416 npts=load_l(nval, cval, ngrids, lswitch)
417 sout(idmtke,1:ngrids)=lswitch(1:ngrids)
418 CASE ('Sout(idMtls)')
419 npts=load_l(nval, cval, ngrids, lswitch)
420 sout(idmtls,1:ngrids)=lswitch(1:ngrids)
421# if defined BBL_MODEL || defined SEDIMENT
422 CASE ('Sout(isd50)')
423 npts=load_l(nval, cval, ngrids, lswitch)
424 i=idbott(isd50)
425 DO ng=1,ngrids
426 sout(i,ng)=lswitch(ng)
427 END DO
428 CASE ('Sout(idens)')
429 npts=load_l(nval, cval, ngrids, lswitch)
430 i=idbott(idens)
431 DO ng=1,ngrids
432 sout(i,ng)=lswitch(ng)
433 END DO
434 CASE ('Sout(iwsed)')
435 npts=load_l(nval, cval, ngrids, lswitch)
436 i=idbott(iwsed)
437 DO ng=1,ngrids
438 sout(i,ng)=lswitch(ng)
439 END DO
440 CASE ('Sout(itauc)')
441 npts=load_l(nval, cval, ngrids, lswitch)
442 i=idbott(itauc)
443 DO ng=1,ngrids
444 sout(i,ng)=lswitch(ng)
445 END DO
446 CASE ('Sout(irlen)')
447 npts=load_l(nval, cval, ngrids, lswitch)
448 i=idbott(irlen)
449 DO ng=1,ngrids
450 sout(i,ng)=lswitch(ng)
451 END DO
452 CASE ('Sout(irhgt)')
453 npts=load_l(nval, cval, ngrids, lswitch)
454 i=idbott(irhgt)
455 DO ng=1,ngrids
456 sout(i,ng)=lswitch(ng)
457 END DO
458 CASE ('Sout(ibwav)')
459 npts=load_l(nval, cval, ngrids, lswitch)
460 i=idbott(ibwav)
461 DO ng=1,ngrids
462 sout(i,ng)=lswitch(ng)
463 END DO
464 CASE ('Sout(izdef)')
465 npts=load_l(nval, cval, ngrids, lswitch)
466 i=idbott(izdef)
467 DO ng=1,ngrids
468 sout(i,ng)=lswitch(ng)
469 END DO
470 CASE ('Sout(izapp)')
471 npts=load_l(nval, cval, ngrids, lswitch)
472 i=idbott(izapp)
473 DO ng=1,ngrids
474 sout(i,ng)=lswitch(ng)
475 END DO
476 CASE ('Sout(izNik)')
477 npts=load_l(nval, cval, ngrids, lswitch)
478 i=idbott(iznik)
479 DO ng=1,ngrids
480 sout(i,ng)=lswitch(ng)
481 END DO
482 CASE ('Sout(izbio)')
483 npts=load_l(nval, cval, ngrids, lswitch)
484 i=idbott(izbio)
485 DO ng=1,ngrids
486 sout(i,ng)=lswitch(ng)
487 END DO
488 CASE ('Sout(izbfm)')
489 npts=load_l(nval, cval, ngrids, lswitch)
490 i=idbott(izbfm)
491 DO ng=1,ngrids
492 sout(i,ng)=lswitch(ng)
493 END DO
494 CASE ('Sout(izbld)')
495 npts=load_l(nval, cval, ngrids, lswitch)
496 i=idbott(izbld)
497 DO ng=1,ngrids
498 sout(i,ng)=lswitch(ng)
499 END DO
500 CASE ('Sout(izwbl)')
501 npts=load_l(nval, cval, ngrids, lswitch)
502 i=idbott(izwbl)
503 DO ng=1,ngrids
504 sout(i,ng)=lswitch(ng)
505 END DO
506 CASE ('Sout(iactv)')
507 npts=load_l(nval, cval, ngrids, lswitch)
508 i=idbott(iactv)
509 DO ng=1,ngrids
510 sout(i,ng)=lswitch(ng)
511 END DO
512 CASE ('Sout(ishgt)')
513 npts=load_l(nval, cval, ngrids, lswitch)
514 i=idbott(ishgt)
515 DO ng=1,ngrids
516 sout(i,ng)=lswitch(ng)
517 END DO
518# endif
519# endif
520 CASE ('NSTATION')
521 npts=load_i(nval, rval, ngrids, nstation)
522 DO ng=1,ngrids
523 IF (.not.lstations(ng)) THEN
524 nstation(ng)=0
525 ELSE
526 IF (nstation(ng).le.0) THEN
527 IF (master) WRITE (out,40) ng, 'Nstation', &
528 & nstation(ng), &
529 & 'Must be positive and greater than zero.'
530 exit_flag=4
531 RETURN
532 END IF
533 END IF
534 END DO
535 CASE ('POS')
536 DO ng=1,ngrids
537 IF (lstations(ng)) THEN
538 allocate ( scalars(ng) % Sflag(nstation(ng)) )
539 allocate ( scalars(ng) % SposX(nstation(ng)) )
540 allocate ( scalars(ng) % SposY(nstation(ng)) )
541 dmem(ng)=dmem(ng)+3.0_r8*real(nstation(ng),r8)
542 END IF
543 END DO
544 is(1:ngrids)=0
545 DO WHILE (.true.)
546 READ (inp,*,err=10,END=10) igrid, flag, Xpos, ypos
547 ng=max(1,min(abs(igrid),ngrids))
548 IF (lstations(ng)) THEN
549 is(ng)=is(ng)+1
550 scalars(ng)%Sflag(is(ng))=flag
551 scalars(ng)%SposX(is(ng))=xpos
552 scalars(ng)%SposY(is(ng))=ypos
553 END IF
554 END DO
555 10 DO ng=1,ngrids
556 IF (lstations(ng).and.(nstation(ng).ne.is(ng))) THEN
557 IF (master) WRITE (out,50) nstation(ng), is(ng)
558 exit_flag=4
559 RETURN
560 END IF
561 END DO
562 END SELECT
563 END IF
564 END DO
565 20 IF (master) WRITE (out,60) line
566 exit_flag=4
567 RETURN
568 30 CONTINUE
569!
570!-----------------------------------------------------------------------
571! Process input parameters.
572!-----------------------------------------------------------------------
573!
574! Turn off the processing of stations if not running long enough to
575! create a stations file (LdefSTA=.FALSE. because nSTA < ntimes or
576! nSTA = 0 when nrrec = 0).
577!
578 DO ng=1,ngrids
579 IF (.not.ldefsta(ng).and.lstations(ng)) THEN
580 lstations(ng)=.false.
581 END IF
582 END DO
583!
584! Make sure that both component switches are activated when processing
585! (Eastward,Northward) momentum components at RHO-points.
586!
587 DO ng=1,ngrids
588 IF (.not.sout(idu2de,ng).and.sout(idv2dn,ng)) THEN
589 sout(idu2de,ng)=.true.
590 END IF
591 IF (.not.sout(idv2dn,ng).and.sout(idu2de,ng)) THEN
592 sout(idv2dn,ng)=.true.
593 END IF
594# ifdef SOLVE3D
595 IF (.not.sout(idu3de,ng).and.sout(idv3dn,ng)) THEN
596 sout(idu3de,ng)=.true.
597 END IF
598 IF (.not.sout(idv3dn,ng).and.sout(idu3de,ng)) THEN
599 sout(idv3dn,ng)=.true.
600 END IF
601# endif
602 END DO
603!
604!-----------------------------------------------------------------------
605! Report input parameters.
606!-----------------------------------------------------------------------
607!
608 IF (master.and.lwrite) THEN
609 DO ng=1,ngrids
610 IF (lstations(ng)) THEN
611 WRITE (out,70) ng
612 WRITE (out,80) nstation(ng), 'Nstation', &
613 & 'Number of stations to write out into stations file.'
614#if defined SEDIMENT && defined SED_MORPH
615 IF (sout(idbath,ng)) WRITE (out,90) sout(idbath,ng), &
616 & 'Sout(idbath)', &
617 & 'Write out free-surface.'
618#endif
619 IF (sout(idfsur,ng)) WRITE (out,90) sout(idfsur,ng), &
620 & 'Sout(idFsur)', &
621 & 'Write out free-surface.'
622 IF (sout(idubar,ng)) WRITE (out,90) sout(idubar,ng), &
623 & 'Sout(idUbar)', &
624 & 'Write out 2D U-momentum component.'
625 IF (sout(idvbar,ng)) WRITE (out,90) sout(idvbar,ng), &
626 & 'Sout(idVbar)', &
627 & 'Write out 2D V-momentum component.'
628 IF (sout(idu2de,ng)) WRITE (out,90) sout(idu2de,ng), &
629 & 'Sout(idu2dE)', &
630 & 'Write out 2D U-eastward at RHO-points.'
631 IF (sout(idv2dn,ng)) WRITE (out,90) sout(idv2dn,ng), &
632 & 'Sout(idv2dN)', &
633 & 'Write out 2D V-northward at RHO-points.'
634# ifdef SOLVE3D
635 IF (sout(iduvel,ng)) WRITE (out,90) sout(iduvel,ng), &
636 & 'Sout(idUvel)', &
637 & 'Write out 3D U-momentum component.'
638 IF (sout(idvvel,ng)) WRITE (out,90) sout(idvvel,ng), &
639 & 'Sout(idVvel)', &
640 & 'Write out 3D V-momentum component.'
641 IF (sout(idu3de,ng)) WRITE (out,90) sout(idu3de,ng), &
642 & 'Sout(idu3dE)', &
643 & 'Write out 3D U-eastward at RHO-points.'
644 IF (sout(idv3dn,ng)) WRITE (out,90) sout(idv3dn,ng), &
645 & 'Sout(idv3dN)', &
646 & 'Write out 3D V-northward at RHO-points.'
647 IF (sout(idwvel,ng)) WRITE (out,90) sout(idwvel,ng), &
648 & 'Sout(idWvel)', &
649 & 'Write out W-momentum component.'
650 IF (sout(idovel,ng)) WRITE (out,90) sout(idovel,ng), &
651 & 'Sout(idOvel)', &
652 & 'Write out omega vertical velocity.'
653 DO itrc=1,nt(ng)
654 IF (sout(idtvar(itrc),ng)) WRITE (out,100) &
655 & sout(idtvar(itrc),ng), 'Sout(idTvar)', &
656 & 'Write out tracer ', itrc, trim(vname(1,idtvar(itrc)))
657 END DO
658# endif
659 IF (sout(idusms,ng)) WRITE (out,90) sout(idusms,ng), &
660 & 'Sout(idUsms)', &
661 & 'Write out surface U-momentum stress.'
662 IF (sout(idvsms,ng)) WRITE (out,90) sout(idvsms,ng), &
663 & 'Sout(idVsms)', &
664 & 'Write out surface V-momentum stress.'
665 IF (sout(idubms,ng)) WRITE (out,90) sout(idubms,ng), &
666 & 'Sout(idUbms)', &
667 & 'Write out bottom U-momentum stress.'
668 IF (sout(idvbms,ng)) WRITE (out,90) sout(idvbms,ng), &
669 & 'Sout(idVbms)', &
670 & 'Write out bottom V-momentum stress.'
671# ifdef BBL_MODEL
672 IF (sout(idubrs,ng)) WRITE (out,90) sout(idubrs,ng), &
673 & 'Sout(idUbrs)', &
674 & 'Write out bottom U-current stress.'
675 IF (sout(idvbrs,ng)) WRITE (out,90) sout(idvbrs,ng), &
676 & 'Sout(idVbrs)', &
677 & 'Write out bottom V-current stress.'
678 IF (sout(idubws,ng)) WRITE (out,90) sout(idubws,ng), &
679 & 'Sout(idUbws)', &
680 & 'Write out wind-induced, bottom U-wave stress.'
681 IF (sout(idvbws,ng)) WRITE (out,90) sout(idvbws,ng), &
682 & 'Sout(idVbws)', &
683 & 'Write out wind-induced, bottom V-wave stress.'
684 IF (sout(idubcs,ng)) WRITE (out,90) sout(idubcs,ng), &
685 & 'Sout(idUbcs)', &
686 & 'Write out max wind + current, bottom U-wave stress.'
687 IF (sout(idvbcs,ng)) WRITE (out,90) sout(idvbcs,ng), &
688 & 'Sout(idVbcs)', &
689 & 'Write out max wind + current, bottom V-wave stress.'
690 IF (sout(idubot,ng)) WRITE (out,90) sout(idubot,ng), &
691 & 'Sout(idUbot)', &
692 & 'Write out bed wave orbital U-velocity.'
693 IF (sout(idvbot,ng)) WRITE (out,90) sout(idvbot,ng), &
694 & 'Sout(idVbot)', &
695 & 'Write out bed wave orbital V-velocity.'
696 IF (sout(idubur,ng)) WRITE (out,90) sout(idubur,ng), &
697 & 'Sout(idUbur)', &
698 & 'Write out bottom U-velocity above bed.'
699 IF (sout(idvbvr,ng)) WRITE (out,90) sout(idvbvr,ng), &
700 & 'Sout(idVbvr)', &
701 & 'Write out bottom V-velocity above bed.'
702# endif
703# ifdef ICE_MODEL
704 IF (sout(iduice,ng)) WRITE (out,90) sout(iduice,ng), &
705 & 'Sout(idUice)', &
706 & 'Write out ice U-velocity component.'
707 IF (sout(idvice,ng)) WRITE (out,90) sout(idvice,ng), &
708 & 'Sout(idVice)', &
709 & 'Write out ice V-velocity component.'
710 IF (sout(iduier,ng)) WRITE (out,90) sout(iduier,ng), &
711 & 'Sout(idUiER)', &
712 & 'Write out ice Eastward velocity component at RHO-points.'
713 IF (sout(idvinr,ng)) WRITE (out,90) sout(idvinr,ng), &
714 & 'Sout(idViNR)', &
715 & 'Write out ice Northward velocity component at RHO-points.'
716 IF (sout(idaice,ng)) WRITE (out,90) sout(idaice,ng), &
717 & 'Sout(idAice)', &
718 & 'Write out ice concentration (fractional area coverage).'
719 IF (sout(idiage,ng)) WRITE (out,90) sout(idiage,ng), &
720 & 'Sout(idIage)', &
721 & 'Write out age of ice.'
722 IF (sout(idhice,ng)) WRITE (out,90) sout(idhice,ng), &
723 & 'Sout(idHice)', &
724 & 'Write out average ice thickness (ice mass per area).'
725 IF (sout(idhmel,ng)) WRITE (out,90) sout(idhmel,ng), &
726 & 'Sout(idHmel)', &
727 & 'Write out surface melt water thickness on ice.'
728 IF (sout(idhsno,ng)) WRITE (out,90) sout(idhsno,ng), &
729 & 'Sout(idHsno)', &
730 & 'Write out average snow coverage thickness.'
731 IF (sout(idtice,ng)) WRITE (out,90) sout(idtice,ng), &
732 & 'Sout(idTice)', &
733 & 'Write out interior ice temperature.'
734 IF (sout(idisxx,ng)) WRITE (out,90) sout(idisxx,ng), &
735 & 'Sout(idISxx)', &
736 & 'Write out internal ice stress tensor, xx-component.'
737 IF (sout(idisxy,ng)) WRITE (out,90) sout(idisxy,ng), &
738 & 'Sout(idISxy)', &
739 & 'Write out internal ice stress tensor, xy-component.'
740 IF (sout(idisyy,ng)) WRITE (out,90) sout(idisyy,ng), &
741 & 'Sout(idISyy)', &
742 & 'Write out internal ice stress tensor, yy-component.'
743 IF (sout(idisst,ng)) WRITE (out,90) sout(idisst,ng), &
744 & 'Sout(idIsst)', &
745 & 'Write out ice/snow surface temperature.'
746 IF (sout(idiomf,ng)) WRITE (out,90) sout(idiomf,ng), &
747 & 'Sout(idIOmf)', &
748 & 'Write out ice-ocean mass flux.'
749 IF (sout(idiofv,ng)) WRITE (out,90) sout(idiofv,ng), &
750 & 'Sout(idIOfv)', &
751 & 'Write out ice-ocean friction velocity.'
752 IF (sout(idiomt,ng)) WRITE (out,90) sout(idiomt,ng), &
753 & 'Sout(idIOmt)', &
754 & 'Write out ice-ocean momentum transfer coefficient.'
755 IF (sout(ids0mk,ng)) WRITE (out,90) sout(ids0mk,ng), &
756 & 'Sout(idS0mk)', &
757 & 'Write out salinity of molecular sublayer under ice.'
758 IF (sout(idt0mk,ng)) WRITE (out,90) sout(idt0mk,ng), &
759 & 'Sout(idT0mk)', &
760 & 'Write out temperature of molecular sublayer under ice.'
761 IF (sout(idwdiv,ng)) WRITE (out,90) sout(idwdiv,ng), &
762 & 'Sout(idWdiv)', &
763 & 'Write out rate of ice divergence.'
764 IF (sout(idw_fr,ng)) WRITE (out,90) sout(idw_fr,ng), &
765 & 'Sout(idW_fr)', &
766 & 'Write out rate of ice accretion by frazil ice growth.'
767 IF (sout(idw_ai,ng)) WRITE (out,90) sout(idw_ai,ng), &
768 & 'Sout(idW_ai)', &
769 & 'Write out rate of melt/freeze at air/ice interface.'
770 IF (sout(idw_ao,ng)) WRITE (out,90) sout(idw_ao,ng), &
771 & 'Sout(idW_ao)', &
772 & 'Write out rate of melt/freeze at air/ocean interface.'
773 IF (sout(idw_io,ng)) WRITE (out,90) sout(idw_io,ng), &
774 & 'Sout(idW_io)', &
775 & 'Write out rate of melt/freeze at ice/ocean interface.'
776 IF (sout(idw_ro,ng)) WRITE (out,90) sout(idw_ro,ng), &
777 & 'Sout(idW_ro)', &
778 & 'Write out rate of melt/freeze runoff into ocean.'
779# endif
780# ifdef WEC
781 IF (sout(idu2rs,ng)) WRITE (out,90) sout(idu2rs,ng), &
782 & 'Sout(idU2rs)', &
783 & 'Write out 2D barotropic wec u-stress.'
784 IF (sout(idv2rs,ng)) WRITE (out,90) sout(idv2rs,ng), &
785 & 'Sout(idV2rs)', &
786 & 'Write out 2D barotropic wec v-stress.'
787 IF (sout(idu2sd,ng)) WRITE (out,90) sout(idu2sd,ng), &
788 & 'Sout(idU2Sd)', &
789 & 'Write out 2D barotropic Stokes u-velocity.'
790 IF (sout(idv2sd,ng)) WRITE (out,90) sout(idv2sd,ng), &
791 & 'Sout(idV2Sd)', &
792 & 'Write out 2D barotropic Stokes v-velocity.'
793# endif
794# ifdef SOLVE3D
795# ifdef WEC
796 IF (sout(idu3rs,ng)) WRITE (out,90) sout(idu3rs,ng), &
797 & 'Sout(idU3rs)', &
798 & 'Write out 3D total wec u-stress.'
799 IF (sout(idv3rs,ng)) WRITE (out,90) sout(idv3rs,ng), &
800 & 'Sout(idV3rs)', &
801 & 'Write out 3D total wec v-stress.'
802 IF (sout(idu3sd,ng)) WRITE (out,90) sout(idu3sd,ng), &
803 & 'Sout(idU3Sd)', &
804 & 'Write out 3D total wec Stokes u-velocity.'
805 IF (sout(idv3sd,ng)) WRITE (out,90) sout(idv3sd,ng), &
806 & 'Sout(idV3Sd)', &
807 & 'Write out 3D total wec Stokes v-velocity.'
808 IF (sout(idw3sd,ng)) WRITE (out,90) sout(idw3sd,ng), &
809 & 'Sout(idW3Sd)', &
810 & 'Write out 3D wec omega Stokes vertical velocity.'
811 IF (sout(idw3st,ng)) WRITE (out,90) sout(idw3st,ng), &
812 & 'Sout(idW3St)', &
813 & 'Write out 3D wec Stokes vertical velocity.'
814# endif
815# ifdef WEC_VF
816 IF (sout(idwztw,ng)) WRITE (out,90) sout(idwztw,ng), &
817 & 'Sout(idWztw)', &
818 & 'Write out wec quasi-static sea level adjustment.'
819 IF (sout(idwqsp,ng)) WRITE (out,90) sout(idwqsp,ng), &
820 & 'Sout(idWqsp)', &
821 & 'Write out wec quasi-static sea pressure adjustment.'
822 IF (sout(idwbeh,ng)) WRITE (out,90) sout(idwbeh,ng), &
823 & 'Sout(idWbeh)', &
824 & 'Write out wec Bernoulli head sea level adjustment.'
825# endif
826# endif
827# ifdef WAVES_HEIGHT
828 IF (sout(idwamp,ng)) WRITE (out,90) sout(idwamp,ng), &
829 & 'Sout(idWamp)', &
830 & 'Write out wave height.'
831# endif
832# ifdef WAVES_LENGTH
833 IF (sout(idwlen,ng)) WRITE (out,90) sout(idwlen,ng), &
834 & 'Sout(idWlen)', &
835 & 'Write out waves mean wavelength.'
836# endif
837# ifdef WAVES_LENGTHP
838 IF (sout(idwlep,ng)) WRITE (out,90) sout(idwlep,ng), &
839 & 'Sout(idWlep)', &
840 & 'Write out waves peak wavelength.'
841# endif
842# ifdef WAVES_DIR
843 IF (sout(idwdir,ng)) WRITE (out,90) sout(idwdir,ng), &
844 & 'Sout(idWdir)', &
845 & 'Write out waves mean direction.'
846# endif
847# ifdef WAVES_DIRP
848 IF (hout(idwdip,ng)) WRITE (out,90) sout(idwdip,ng), &
849 & 'Sout(idWdip)', &
850 & 'Write out peak waves direction.'
851# endif
852# ifdef WAVES_TOP_PERIOD
853 IF (sout(idwptp,ng)) WRITE (out,90) sout(idwptp,ng), &
854 & 'Sout(idWptp)', &
855 & 'Write out wave surface period.'
856# endif
857# ifdef WAVES_BOT_PERIOD
858 IF (sout(idwpbt,ng)) WRITE (out,90) sout(idwpbt,ng), &
859 & 'Sout(idWpbt)', &
860 & 'Write out wave bottom period.'
861# endif
862# if defined BBL_MODEL || defined BEDLOAD_SOULSBY || \
863 defined bedload_vandera || defined wav_coupling
864 IF (sout(idworb,ng)) WRITE (out,90) sout(idworb,ng), &
865 & 'Sout(idWorb)', &
866 & 'Write out wave bottom orbital velocity.'
867# endif
868# if defined WAV_COUPLING || (defined WEC_VF && defined BOTTOM_STREAMING)
869 IF (sout(idwdif,ng)) WRITE (out,90) sout(idwdif,ng), &
870 & 'Sout(idWdif)', &
871 & 'Write out wave dissipation due to bottom friction.'
872# endif
873# if defined WAV_COUPLING || defined TKE_WAVEDISS || \
874 defined wdiss_thorguza || defined wdiss_churthor
875 IF (sout(idwdib,ng)) WRITE (out,90) sout(idwdib,ng), &
876 & 'Sout(idWdib)', &
877 & 'Write out wave dissipation due to breaking.'
878 IF (sout(idwdiw,ng)) WRITE (out,90) sout(idwdiw,ng), &
879 & 'Sout(idWdiw)', &
880 & 'Write out wave dissipation due to whitecapping.'
881# endif
882# ifdef ROLLER_SVENDSEN
883 IF (sout(idwbrk,ng)) WRITE (out,90) sout(idwbrk,ng), &
884 & 'Sout(idWbrk)', &
885 & 'Write out percent wave breaking.'
886# endif
887# ifdef WEC_ROLLER
888 IF (sout(idwdis,ng)) WRITE (out,90) sout(idwdis,ng), &
889 & 'Sout(idWdis)', &
890 & 'Write out wave roller dissipation.'
891# endif
892# ifdef ROLLER_RENIERS
893 IF (sout(idwrol,ng)) WRITE (out,170) sout(idwrol,ng), &
894 & 'Sout(idWrol)', &
895 & 'Write out wave roller action density.'
896# endif
897# ifdef WAVES_DSPR
898 IF (sout(idwvds,ng)) WRITE (out,170) sout(idwvds,ng), &
899 & 'Sout(idWvds)', &
900 & 'Write out wave directional spread.'
901 IF (sout(idwvqp,ng)) WRITE (out,170) sout(idwvqp,ng), &
902 & 'Sout(idWvqp)', &
903 & 'Write out wave spectrum peakedness.'
904# endif
905# if defined SOLVE3D && (defined BBL_MODEL || defined SEDIMENT)
906 DO itrc=1,mbotp
907 IF (sout(idbott(itrc),ng)) WRITE (out,100) &
908 & sout(idbott(itrc),ng), 'Sout(idBott)', &
909 & 'Write out bottom property ', itrc, &
910 & trim(vname(1,idbott(itrc)))
911 END DO
912# endif
913# ifdef SOLVE3D
914# if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS
915 IF (sout(idpair,ng)) WRITE (out,90) sout(idpair,ng), &
916 & 'Sout(idPair)', &
917 & 'Write out surface air pressure.'
918# endif
919# if defined BULK_FLUXES || defined ECOSIM
920 IF (sout(iduair,ng)) WRITE (out,90) sout(iduair,ng), &
921 & 'Sout(idUair)', &
922 & 'Write out surface U-wind component.'
923 IF (sout(idvair,ng)) WRITE (out,90) sout(idvair,ng), &
924 & 'Sout(idVair)', &
925 & 'Write out surface V-wind component.'
926# endif
927 IF (sout(idtsur(itemp),ng)) WRITE (out,90) &
928 & sout(idtsur(itemp),ng), 'Sout(idTsur)', &
929 & 'Write out surface net heat flux.'
930# ifdef SALINITY
931 IF (sout(idtsur(isalt),ng)) WRITE (out,90) &
932 & sout(idtsur(isalt),ng), 'Sout(idTsur)', &
933 & 'Write out surface net salt flux.'
934# endif
935# ifdef SHORTWAVE
936 IF (sout(idsrad,ng)) WRITE (out,90) sout(idsrad,ng), &
937 & 'Sout(idSrad)', &
938 & 'Write out shortwave radiation flux.'
939# endif
940# ifdef BULK_FLUXES
941 IF (sout(idlrad,ng)) WRITE (out,90) sout(idlrad,ng), &
942 & 'Sout(idLrad)', &
943 & 'Write out longwave radiation flux.'
944 IF (sout(idlhea,ng)) WRITE (out,90) sout(idlhea,ng), &
945 & 'Sout(idLhea)', &
946 & 'Write out latent heat flux.'
947 IF (sout(idshea,ng)) WRITE (out,90) sout(idshea,ng), &
948 & 'Sout(idShea)', &
949 & 'Write out sensible heat flux.'
950# ifdef EMINUSP
951 IF (sout(idempf,ng)) WRITE (out,90) sout(idempf,ng), &
952 & 'Sout(idEmPf)', &
953 & 'Write out E-P flux.'
954 IF (sout(idevap,ng)) WRITE (out,90) sout(idevap,ng), &
955 & 'Sout(idevap)', &
956 & 'Write out evaporation rate.'
957 IF (sout(idrain,ng)) WRITE (out,90) sout(idrain,ng), &
958 & 'Sout(idrain)', &
959 & 'Write out rain rate.'
960# endif
961# endif
962 IF (sout(iddano,ng)) WRITE (out,90) sout(iddano,ng), &
963 & 'Sout(idDano)', &
964 & 'Write out density anomaly.'
965 IF (sout(idvvis,ng)) WRITE (out,90) sout(idvvis,ng), &
966 & 'Sout(idVvis)', &
967 & 'Write out vertical viscosity coefficient.'
968 IF (sout(idtdif,ng)) WRITE (out,90) sout(idtdif,ng), &
969 & 'Sout(idTdif)', &
970 & 'Write out vertical T-diffusion coefficient.'
971 IF (sout(idsdif,ng)) WRITE (out,90) sout(idsdif,ng), &
972 & 'Sout(idSdif)', &
973 & 'Write out vertical S-diffusion coefficient.'
974# ifdef LMD_SKPP
975 IF (sout(idhsbl,ng)) WRITE (out,90) sout(idhsbl,ng), &
976 & 'Sout(idHsbl)', &
977 & 'Write out depth of surface boundary layer.'
978# endif
979# ifdef LMD_BKPP
980 IF (sout(idhbbl,ng)) WRITE (out,90) sout(idhbbl,ng), &
981 & 'Sout(idHbbl)', &
982 & 'Write out depth of bottom boundary layer.'
983# endif
984# if defined GLS_MIXING || defined MY25_MIXING
985 IF (sout(idmtke,ng)) WRITE (out,90) sout(idmtke,ng), &
986 & 'Sout(idMtke)', &
987 & 'Write out turbulent kinetic energy.'
988 IF (sout(idmtls,ng)) WRITE (out,90) sout(idmtls,ng), &
989 & 'Sout(idMtls)', &
990 & 'Write out turbulent generic length-scale.'
991# endif
992# endif
993 WRITE (out,*)
994 DO i=1,nstation(ng)
995 WRITE (out,110) i, scalars(ng)%Sflag(i), &
996 & scalars(ng)%SposX(i), &
997 & scalars(ng)%SposY(i)
998 END DO
999 END IF
1000 END DO
1001 END IF
1002
1003 40 FORMAT (/,' READ_StaPar - Grid = ',i2.2,',',3x, &
1004 & 'Illegal value for ',a,' = ', i8,/,15x,a)
1005 50 FORMAT (/,' READ_StaPar - Inconsistent number of stations, ', &
1006 & 'Nstation = ',2i8,/,15x, &
1007 & 'change stations input script values.')
1008 60 FORMAT (/,' READ_StaPar - Error while processing line: ',/,a)
1009
1010 70 FORMAT (/,/,' Stations Parameters, Grid: ',i2.2, &
1011 & /, ' =============================',/)
1012 80 FORMAT (1x,i10,2x,a,t30,a)
1013 90 FORMAT (10x,l1,2x,a,t30,a)
1014 100 FORMAT (10x,l1,2x,a,t30,a,i2.2,':',1x,a)
1015 110 FORMAT (13x,'Flag and positions for station ',i4.4,':', &
1016 & i3,1x,2f10.4)
1017
1018 RETURN
integer function decode_line(line_text, keyword, nval, cval, rval)
Definition inp_decode.F:97
integer idisxy
Definition ice_mod.h:109
integer idvinr
Definition ice_mod.h:119
integer idhmel
Definition ice_mod.h:101
integer idw_ai
Definition ice_mod.h:121
integer idw_ro
Definition ice_mod.h:125
integer idiomf
Definition ice_mod.h:105
integer idvice
Definition ice_mod.h:117
integer ids0mk
Definition ice_mod.h:111
integer idiofv
Definition ice_mod.h:104
integer idwdiv
Definition ice_mod.h:120
integer idt0mk
Definition ice_mod.h:112
integer idw_ao
Definition ice_mod.h:122
integer idw_fr
Definition ice_mod.h:123
integer idisst
Definition ice_mod.h:107
integer idiage
Definition ice_mod.h:103
integer idiomt
Definition ice_mod.h:106
integer idhice
Definition ice_mod.h:99
integer idisxx
Definition ice_mod.h:108
integer idisyy
Definition ice_mod.h:110
integer idhsno
Definition ice_mod.h:102
integer idaice
Definition ice_mod.h:96
integer idw_io
Definition ice_mod.h:124
integer iduice
Definition ice_mod.h:114
integer idtice
Definition ice_mod.h:113
integer iduier
Definition ice_mod.h:116
type(t_io), dimension(:), allocatable err
integer iddano
integer idvair
integer idwqsp
logical, dimension(:,:), allocatable hout
integer idevap
integer idvbrs
integer idv2rs
integer idwdis
integer idwdif
integer idv3sd
integer idubar
integer idwvel
integer idvvel
integer idubur
integer idhsbl
integer idvsms
integer idwlen
integer idwdiw
integer idvbvr
integer idpair
integer idwlep
integer idw3st
integer idu2rs
integer idv2dn
integer idvbot
integer idsdif
integer idwrol
integer, dimension(:), allocatable idtsur
integer idempf
integer idtdif
integer idfsur
integer idwbrk
integer, dimension(:), allocatable idtvar
integer idhbbl
integer idvbws
integer idbath
integer idvbms
integer iduair
integer idmtke
integer iduvel
integer idv3dn
integer idovel
integer idv3rs
character(len=maxlen), dimension(6, 0:nv) vname
integer idu3rs
integer idshea
integer idlrad
integer idwdip
integer idusms
integer idvbcs
integer idwbeh
logical, dimension(:,:), allocatable sout
integer idwamp
integer idwdir
integer idvvis
integer idwptp
integer idu3de
integer idwdib
integer idv2sd
integer idubcs
integer idu2de
integer idlhea
integer idubot
integer idrain
integer idubms
integer idubrs
integer idsrad
integer idmtls
integer idubws
integer idwpbt
integer idworb
integer idu3sd
integer idu2sd
integer idwvds
integer idw3sd
integer idwvqp
integer idwztw
integer idvbar
logical master
integer nat
Definition mod_param.F:499
real(r8), dimension(:), allocatable dmem
Definition mod_param.F:137
integer, dimension(:), allocatable nstation
Definition mod_param.F:557
integer ngrids
Definition mod_param.F:113
integer mt
Definition mod_param.F:490
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer exit_flag
integer isalt
integer itemp
type(t_scalars), dimension(:), allocatable scalars
Definition mod_scalars.F:65
logical, dimension(:), allocatable ldefsta
logical, dimension(:), allocatable lstations
integer, parameter iwsed
integer, parameter irlen
integer, parameter irhgt
integer, parameter mbotp
integer, parameter iactv
integer, parameter iznik
integer, parameter ibwav
integer, parameter ishgt
integer, parameter idens
integer, parameter isd50
integer, dimension(mbotp) idbott
integer, parameter izbld
integer, parameter izapp
integer, parameter izbio
integer, parameter izdef
integer, parameter itauc
integer, parameter izwbl
integer, parameter izbfm

References mod_iounits::err.

Referenced by inp_par_mod::inp_par().

Here is the caller graph for this function: