ROMS
Loading...
Searching...
No Matches
rp_def_ini_mod Module Reference

Functions/Subroutines

subroutine, public rp_def_ini (ng)
 
subroutine, private rp_def_ini_nf90 (ng)
 
subroutine, private rp_def_ini_pio (ng)
 

Function/Subroutine Documentation

◆ rp_def_ini()

subroutine, public rp_def_ini_mod::rp_def_ini ( integer, intent(in) ng)

Definition at line 40 of file rp_def_ini.F.

41!***********************************************************************
42!
43! Imported variable declarations.
44!
45 integer, intent(in) :: ng
46!
47! Local variable declarations.
48!
49 character (len=*), parameter :: MyFile = &
50 & __FILE__
51!
52!-----------------------------------------------------------------------
53! Create a new history file according to IO type.
54!-----------------------------------------------------------------------
55!
56 SELECT CASE (irp(ng)%IOtype)
57 CASE (io_nf90)
58 CALL rp_def_ini_nf90 (ng)
59
60# if defined PIO_LIB && defined DISTRIBUTE
61 CASE (io_pio)
62 CALL rp_def_ini_pio (ng)
63# endif
64 CASE DEFAULT
65 IF (master) WRITE (stdout,10) irp(ng)%IOtype
66 exit_flag=3
67 END SELECT
68 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
69!
70 10 FORMAT (' RP_DEF_INI - Illegal output type, io_type = ',i0, &
71 & /,14x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.')
72!
73 RETURN

References mod_scalars::exit_flag, strings_mod::founderror(), mod_ncparam::io_nf90, mod_ncparam::io_pio, mod_iounits::irp, mod_parallel::master, mod_scalars::noerror, rp_def_ini_nf90(), rp_def_ini_pio(), and mod_iounits::stdout.

Referenced by r4dvar_mod::increment(), and rp_initial().

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

◆ rp_def_ini_nf90()

subroutine, private rp_def_ini_mod::rp_def_ini_nf90 ( integer, intent(in) ng)
private

Definition at line 77 of file rp_def_ini.F.

78!***********************************************************************
79!
80 USE mod_netcdf
81!
82! Imported variable declarations.
83!
84 integer, intent(in) :: ng
85!
86! Local variable declarations.
87!
88 logical :: Ldefine = .false.
89# ifdef ADJUST_BOUNDARY
90 logical :: got_IorJ, got_boundary, got_obc_adjust
91# endif
92# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
93 logical :: got_frc_adjust
94# endif
95 logical :: got_var(NV)
96!
97 integer, parameter :: Natt = 25
98
99 integer :: i, j, ifield, itrc, status
100 integer :: Fcount
101# ifdef ADJUST_BOUNDARY
102 integer :: IorJdim, brecdim
103# endif
104# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
105 integer :: frecdim
106# endif
107 integer :: DimIDs(nDimID)
108# ifdef ADJUST_BOUNDARY
109 integer :: t2dobc(4)
110# ifdef SOLVE3D
111 integer :: t3dobc(5)
112# endif
113# endif
114# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
115 integer :: t4dfrc(4), u4dfrc(4), v4dfrc(4)
116# endif
117!
118 real(r8) :: Aval(6)
119!
120 character (len=256) :: ncname
121 character (len=MaxLen) :: Vinfo(Natt)
122
123 character (len=*), parameter :: MyFile = &
124 & __FILE__//", rp_def_ini_nf90"
125!
126 sourcefile=myfile
127
128# if defined ADJUST_BOUNDARY || \
129 defined adjust_stflux || defined adjust_wstress
130!
131!=======================================================================
132! Open existing representer model initial conditions and define new
133! variables.
134!=======================================================================
135!
136 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
137 ncname=irp(ng)%name
138!
139! Open representer model initialization file for read/write.
140!
141 CALL netcdf_open (ng, irpm, ncname, 1, irp(ng)%ncid)
142 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
143 IF (master) WRITE (stdout,10) trim(ncname)
144 RETURN
145 END IF
146!
147! Inquire about the dimensions and check for consistency.
148!
149 CALL netcdf_check_dim (ng, irpm, ncname, &
150 & ncid = irp(ng)%ncid)
151 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
152!
153! Inquire about the variables.
154!
155 CALL netcdf_inq_var (ng, irpm, ncname, &
156 & ncid = irp(ng)%ncid)
157 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
158!
159! Check if surface forcing variables have been already defined.
160!
161 DO i=1,nv
162 got_var(i)=.false.
163 END DO
164!
165 DO i=1,n_var
166# ifdef ADJUST_BOUNDARY
167 IF (trim(var_name(i)).eq. &
168 & trim(vname(1,idsbry(isfsur)))) THEN
169 got_var(idsbry(isfsur))=.true.
170 irp(ng)%Vid(idsbry(isfsur))=var_id(i)
171 ELSE IF (trim(var_name(i)).eq. &
172 & trim(vname(1,idsbry(isubar)))) THEN
173 got_var(idsbry(isubar))=.true.
174 irp(ng)%Vid(idsbry(isubar))=var_id(i)
175 ELSE IF (trim(var_name(i)).eq. &
176 & trim(vname(1,idsbry(isvbar)))) THEN
177 got_var(idsbry(isvbar))=.true.
178 irp(ng)%Vid(idsbry(isvbar))=var_id(i)
179# ifdef SOLVE3D
180 ELSE IF (trim(var_name(i)).eq. &
181 & trim(vname(1,idsbry(isuvel)))) THEN
182 got_var(idsbry(isuvel))=.true.
183 irp(ng)%Vid(idsbry(isuvel))=var_id(i)
184 ELSE IF (trim(var_name(i)).eq. &
185 & trim(vname(1,idsbry(isvvel)))) THEN
186 got_var(idsbry(isvvel))=.true.
187 irp(ng)%Vid(idsbry(isvvel))=var_id(i)
188# endif
189 END IF
190# ifdef SOLVE3D
191 DO itrc=1,nt(ng)
192 IF (trim(var_name(i)).eq. &
193 & trim(vname(1,idsbry(istvar(itrc))))) THEN
194 got_var(idsbry(istvar(itrc)))=.true.
195 irp(ng)%Vid(idsbry(istvar(itrc)))=var_id(i)
196 END IF
197 END DO
198# endif
199# endif
200# ifdef ADJUST_WSTRESS
201 IF (trim(var_name(i)).eq.trim(vname(1,idusms))) THEN
202 got_var(idusms)=.true.
203 irp(ng)%Vid(idusms)=var_id(i)
204 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvsms))) THEN
205 got_var(idvsms)=.true.
206 irp(ng)%Vid(idvsms)=var_id(i)
207 END IF
208# endif
209# if defined ADJUST_STFLUX && defined SOLVE3D
210 DO itrc=1,nt(ng)
211 IF (lstflux(itrc,ng)) THEN
212 IF (trim(var_name(i)).eq.trim(vname(1,idtsur(itrc)))) THEN
213 got_var(idtsur(itrc))=.true.
214 irp(ng)%Vid(idtsur(itrc))=var_id(i)
215 END IF
216 END IF
217 END DO
218# endif
219 END DO
220
221# ifdef ADJUST_BOUNDARY
222 IF (.not.got_var(idsbry(isfsur))) ldefine=.true.
223 IF (.not.got_var(idsbry(isubar))) ldefine=.true.
224 IF (.not.got_var(idsbry(isvbar))) ldefine=.true.
225# ifdef SOLVE3D
226 IF (.not.got_var(idsbry(isuvel))) ldefine=.true.
227 IF (.not.got_var(idsbry(isvvel))) ldefine=.true.
228 DO itrc=1,nt(ng)
229 IF (.not.got_var(idsbry(istvar(itrc)))) ldefine=.true.
230 END DO
231# endif
232# endif
233# ifdef ADJUST_WSTRESS
234 IF (.not.got_var(idusms)) ldefine=.true.
235 IF (.not.got_var(idvsms)) ldefine=.true.
236# endif
237# if defined ADJUST_STFLUX && defined SOLVE3D
238 DO itrc=1,nt(ng)
239 IF (lstflux(itrc,ng)) THEN
240 IF (.not.got_var(idtsur(itrc))) ldefine=.true.
241 END IF
242 END DO
243# endif
244!
245! Put existing file into define mode so new variables can be added.
246!
247 IF (ldefine) THEN
248 CALL netcdf_redef (ng, irpm, ncname, irp(ng)%ncid)
249 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
250 IF (master) WRITE (stdout,20) trim(ncname)
251 RETURN
252 END IF
253 END IF
254!
255!-----------------------------------------------------------------------
256! Define the dimensions of staggered fields.
257!-----------------------------------------------------------------------
258!
259 define: IF (ldefine) THEN
260
261# ifdef ADJUST_BOUNDARY
262 got_iorj=.false.
263 got_boundary=.false.
264 got_obc_adjust=.false.
265# endif
266# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
267 got_frc_adjust=.false.
268# endif
269 DO i=1,n_dim
270 SELECT CASE (trim(adjustl(dim_name(i))))
271 CASE ('xi_rho')
272 dimids( 1)=dim_id(i)
273 CASE ('xi_u')
274 dimids( 2)=dim_id(i)
275 CASE ('xi_v')
276 dimids( 3)=dim_id(i)
277 CASE ('eta_rho')
278 dimids( 5)=dim_id(i)
279 CASE ('eta_u')
280 dimids( 6)=dim_id(i)
281 CASE ('eta_v')
282 dimids( 7)=dim_id(i)
283# ifdef SOLVE3D
284 CASE ('s_rho')
285 dimids( 9)=dim_id(i)
286 CASE ('s_w')
287 dimids(10)=dim_id(i)
288# endif
289# ifdef ADJUST_BOUNDARY
290 CASE ('boundary')
291 dimids(14)=dim_id(i)
292 got_boundary=.true.
293 CASE ('IorJ')
294 iorjdim=dim_id(i)
295 got_iorj=.true.
296# endif
297# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
298 CASE ('frc_adjust')
299 frecdim=dim_id(i)
300 got_frc_adjust=.true.
301# endif
302# ifdef ADJUST_BOUNDARY
303 CASE ('obc_adjust')
304 brecdim=dim_id(i)
305 got_obc_adjust=.true.
306# endif
307 END SELECT
308 END DO
309
310 dimids(12)=rec_id
311# ifdef ADJUST_BOUNDARY
312 IF (.not.got_boundary) THEN
313 status=def_dim(ng, irpm, irp(ng)%ncid, ncname, 'boundary', &
314 & 4, dimids(14))
315 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
316 END IF
317 IF (.not.got_iorj) THEN
318 status=def_dim(ng, irpm, irp(ng)%ncid, ncname, 'IorJ', &
319 & iobounds(ng)%IorJ, iorjdim)
320 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
321 END IF
322# endif
323# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
324 IF (.not.got_frc_adjust) THEN
325 status=def_dim(ng, irpm, irp(ng)%ncid, ncname, 'frc_adjust', &
326 & nfrec(ng), frecdim)
327 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
328 END IF
329# endif
330# ifdef ADJUST_BOUNDARY
331 IF (.not.got_obc_adjust) THEN
332 status=def_dim(ng, irpm, irp(ng)%ncid, ncname, 'obc_adjust', &
333 & nbrec(ng), brecdim)
334 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
335 END IF
336# endif
337!
338! Define dimension vectors for staggered tracer type variables.
339!
340# ifdef ADJUST_BOUNDARY
341 t2dobc(1)=iorjdim
342 t2dobc(2)=dimids(14)
343 t2dobc(3)=brecdim
344 t2dobc(4)=dimids(12)
345# ifdef SOLVE3D
346 t3dobc(1)=iorjdim
347 t3dobc(2)=dimids( 9)
348 t3dobc(3)=dimids(14)
349 t3dobc(4)=brecdim
350 t3dobc(5)=dimids(12)
351# endif
352# endif
353# ifdef ADJUST_STFLUX
354 t4dfrc(1)=dimids( 1)
355 t4dfrc(2)=dimids( 5)
356 t4dfrc(3)=frecdim
357 t4dfrc(4)=dimids(12)
358# endif
359# ifdef ADJUST_WSTRESS
360!
361! Define dimension vectors for staggered u-momentum type variables.
362!
363 u4dfrc(1)=dimids( 2)
364 u4dfrc(2)=dimids( 6)
365 u4dfrc(3)=frecdim
366 u4dfrc(4)=dimids(12)
367# endif
368# ifdef ADJUST_WSTRESS
369!
370! Define dimension vectors for staggered v-momentum type variables.
371!
372 v4dfrc(1)=dimids( 3)
373 v4dfrc(2)=dimids( 7)
374 v4dfrc(3)=frecdim
375 v4dfrc(4)=dimids(12)
376# endif
377!
378! Initialize local information variable arrays.
379!
380 DO i=1,natt
381 DO j=1,len(vinfo(1))
382 vinfo(i)(j:j)=' '
383 END DO
384 END DO
385 DO i=1,6
386 aval(i)=0.0_r8
387 END DO
388!
389!-----------------------------------------------------------------------
390! Define additional variables. Notice that these variables have their
391! own fixed time-dimension to allow 4DVAR adjustments at other times
392! in addition to initialization time.
393!-----------------------------------------------------------------------
394
395# ifdef ADJUST_BOUNDARY
396!
397! Define free-surface open boundaries.
398!
399 IF (any(lobc(:,isfsur,ng))) THEN
400 ifield=idsbry(isfsur)
401 vinfo( 1)=vname(1,ifield)
402 vinfo( 2)=vname(2,ifield)
403 vinfo( 3)=vname(3,ifield)
404 vinfo(14)=vname(4,ifield)
405 vinfo(16)=vname(1,idtime)
406 vinfo(21)=vname(6,ifield)
407 aval(5)=real(iinfo(1,ifield,ng),r8)
408 status=def_var(ng, irpm, irp(ng)%ncid, irp(ng)%Vid(ifield), &
409 & nf_fout, 4, t2dobc, aval, vinfo, ncname, &
410 & setfillval = .false.)
411 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
412 END IF
413!
414! Define 2D U-momentum component open boundaries.
415!
416 IF (any(lobc(:,isubar,ng))) THEN
417 ifield=idsbry(isubar)
418 vinfo( 1)=vname(1,ifield)
419 vinfo( 2)=vname(2,ifield)
420 vinfo( 3)=vname(3,ifield)
421 vinfo(14)=vname(4,ifield)
422 vinfo(16)=vname(1,idtime)
423 vinfo(21)=vname(6,ifield)
424 aval(5)=real(iinfo(1,ifield,ng),r8)
425 status=def_var(ng, irpm, irp(ng)%ncid, irp(ng)%Vid(ifield), &
426 & nf_fout, 4, t2dobc, aval, vinfo, ncname, &
427 & setfillval = .false.)
428 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
429 END IF
430!
431! Define 2D V-momentum component open boundaries.
432!
433 IF (any(lobc(:,isvbar,ng))) THEN
434 ifield=idsbry(isvbar)
435 vinfo( 1)=vname(1,ifield)
436 vinfo( 2)=vname(2,ifield)
437 vinfo( 3)=vname(3,ifield)
438 vinfo(14)=vname(4,ifield)
439 vinfo(16)=vname(1,idtime)
440 vinfo(21)=vname(6,ifield)
441 aval(5)=real(iinfo(1,ifield,ng),r8)
442 status=def_var(ng, irpm, irp(ng)%ncid, irp(ng)%Vid(ifield), &
443 & nf_fout, 4, t2dobc, aval, vinfo, ncname, &
444 & setfillval = .false.)
445 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
446 END IF
447
448# ifdef SOLVE3D
449!
450! Define 3D U-momentum component open boundaries.
451!
452 IF (any(lobc(:,isuvel,ng))) THEN
453 ifield=idsbry(isuvel)
454 vinfo( 1)=vname(1,ifield)
455 vinfo( 2)=vname(2,ifield)
456 vinfo( 3)=vname(3,ifield)
457 vinfo(14)=vname(4,ifield)
458 vinfo(16)=vname(1,idtime)
459 vinfo(21)=vname(6,ifield)
460 aval(5)=real(iinfo(1,ifield,ng),r8)
461 status=def_var(ng, irpm, irp(ng)%ncid, irp(ng)%Vid(ifield), &
462 & nf_fout, 5, t3dobc, aval, vinfo, ncname, &
463 & setfillval = .false.)
464 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
465 END IF
466!
467! Define 3D V-momentum component open boundaries.
468!
469 IF (any(lobc(:,isvvel,ng))) THEN
470 ifield=idsbry(isvvel)
471 vinfo( 1)=vname(1,ifield)
472 vinfo( 2)=vname(2,ifield)
473 vinfo( 3)=vname(3,ifield)
474 vinfo(14)=vname(4,ifield)
475 vinfo(16)=vname(1,idtime)
476 vinfo(21)=vname(6,ifield)
477 aval(5)=real(iinfo(1,ifield,ng),r8)
478 status=def_var(ng, irpm, irp(ng)%ncid, irp(ng)%Vid(ifield), &
479 & nf_fout, 5, t3dobc, aval, vinfo, ncname, &
480 & setfillval = .false.)
481 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
482 END IF
483!
484! Define tracer type variables open boundaries.
485!
486 DO itrc=1,nt(ng)
487 IF (any(lobc(:,istvar(itrc),ng))) THEN
488 ifield=idsbry(istvar(itrc))
489 vinfo( 1)=vname(1,ifield)
490 vinfo( 2)=vname(2,ifield)
491 vinfo( 3)=vname(3,ifield)
492 vinfo(14)=vname(4,ifield)
493 vinfo(16)=vname(1,idtime)
494 vinfo(21)=vname(6,ifield)
495 aval(5)=real(iinfo(1,ifield,ng),r8)
496 status=def_var(ng, irpm, irp(ng)%ncid, irp(ng)%Vid(ifield), &
497 & nf_fout, 5, t3dobc, aval, vinfo, ncname, &
498 & setfillval = .false.)
499 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
500 END IF
501 END DO
502# endif
503# endif
504# ifdef ADJUST_WSTRESS
505!
506! Define surface U-momentum stress.
507!
508 IF (.not.got_var(idusms)) THEN
509 vinfo( 1)=vname(1,idusms)
510 vinfo( 2)=vname(2,idusms)
511 vinfo( 3)=vname(3,idusms)
512# if defined WRITE_WATER && defined MASKING
513 vinfo(20)='mask_u'
514# endif
515 vinfo(21)=vname(6,idusms)
516 vinfo(22)='coordinates'
517 aval(5)=real(u2dvar,r8)
518 status=def_var(ng, irpm, irp(ng)%ncid, irp(ng)%Vid(idusms), &
519 & nf_fout, 4, u4dfrc, aval, vinfo, ncname)
520 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
521 END IF
522!
523! Define surface V-momentum stress.
524!
525 IF (.not.got_var(idvsms)) THEN
526 vinfo( 1)=vname(1,idvsms)
527 vinfo( 2)=vname(2,idvsms)
528 vinfo( 3)=vname(3,idvsms)
529# if defined WRITE_WATER && defined MASKING
530 vinfo(20)='mask_v'
531# endif
532 vinfo(21)=vname(6,idvsms)
533 vinfo(22)='coordinates'
534 aval(5)=real(v2dvar,r8)
535 status=def_var(ng, irpm, irp(ng)%ncid, irp(ng)%Vid(idvsms), &
536 & nf_fout, 4, v4dfrc, aval, vinfo, ncname)
537 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
538 END IF
539# endif
540# if defined ADJUST_STFLUX && defined SOLVE3D
541!
542! Define surface tracer flux.
543!
544 DO itrc=1,nt(ng)
545 IF (.not.got_var(idtsur(itrc)).and.lstflux(itrc,ng)) THEN
546 vinfo( 1)=vname(1,idtsur(itrc))
547 vinfo( 2)=trim(vname(2,idtsur(itrc)))
548 vinfo( 3)=vname(3,idtsur(itrc))
549 IF (itrc.eq.itemp) THEN
550 vinfo(11)='upward flux, cooling'
551 vinfo(12)='downward flux, heating'
552 ELSE IF (itrc.eq.isalt) THEN
553 vinfo(11)='upward flux, freshening (net precipitation)'
554 vinfo(12)='downward flux, salting (net evaporation)'
555 END IF
556# if defined WRITE_WATER && defined MASKING
557 vinfo(20)='mask_rho'
558# endif
559 vinfo(21)=vname(6,idtsur(itrc))
560 vinfo(22)='coordinates'
561 aval(5)=real(r2dvar,r8)
562 status=def_var(ng, irpm, irp(ng)%ncid, &
563 & irp(ng)%Vid(idtsur(itrc)), nf_fout, &
564 & 4, t4dfrc, aval, vinfo, ncname)
565 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
566 END IF
567 END DO
568# endif
569!
570!-----------------------------------------------------------------------
571! Leave definition mode.
572!-----------------------------------------------------------------------
573!
574 CALL netcdf_enddef (ng, irpm, ncname, irp(ng)%ncid)
575 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
576
577 END IF define
578# endif
579!
580!=======================================================================
581! Open an existing tangent linear initial file, check its contents, and
582! prepare for appending data.
583!=======================================================================
584!
585 IF (.not.ldefirp(ng)) THEN
586 ncname=irp(ng)%name
587
588# if !(defined ADJUST_STFLUX || defined ADJUST_WSTRESS)
589!
590! Open representer model initialization file for read/write.
591!
592 CALL netcdf_open (ng, irpm, trim(ncname), 1, irp(ng)%ncid)
593 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
594 WRITE (stdout,10) trim(ncname)
595 RETURN
596 END IF
597!
598! Inquire about the dimensions and check for consistency.
599!
600 CALL netcdf_check_dim (ng, irpm, ncname, &
601 & ncid = irp(ng)%ncid)
602 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
603!
604! Inquire about the variables.
605!
606 CALL netcdf_inq_var (ng, irpm, ncname, &
607 & ncid = irp(ng)%ncid)
608 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
609!
610! Initialize logical switches.
611!
612 DO i=1,nv
613 got_var(i)=.false.
614 END DO
615# endif
616!
617! Scan variable list from input NetCDF and activate switches for
618! initialization variables. Get variable IDs.
619!
620 DO i=1,n_var
621 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
622 got_var(idtime)=.true.
623 irp(ng)%Vid(idtime)=var_id(i)
624 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idfsur))) THEN
625 got_var(idfsur)=.true.
626 irp(ng)%Vid(idfsur)=var_id(i)
627 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idubar))) THEN
628 got_var(idubar)=.true.
629 irp(ng)%Vid(idubar)=var_id(i)
630 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvbar))) THEN
631 got_var(idvbar)=.true.
632 irp(ng)%Vid(idvbar)=var_id(i)
633# ifdef ADJUST_BOUNDARY
634 ELSE IF (trim(var_name(i)).eq. &
635 & trim(vname(1,idsbry(isfsur)))) THEN
636 got_var(idsbry(isfsur))=.true.
637 irp(ng)%Vid(idsbry(isfsur))=var_id(i)
638 ELSE IF (trim(var_name(i)).eq. &
639 & trim(vname(1,idsbry(isubar)))) THEN
640 got_var(idsbry(isubar))=.true.
641 irp(ng)%Vid(idsbry(isubar))=var_id(i)
642 ELSE IF (trim(var_name(i)).eq. &
643 & trim(vname(1,idsbry(isvbar)))) THEN
644 got_var(idsbry(isvbar))=.true.
645 irp(ng)%Vid(idsbry(isvbar))=var_id(i)
646# endif
647# ifdef ADJUST_WSTRESS
648 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idusms))) THEN
649 got_var(idusms)=.true.
650 irp(ng)%Vid(idusms)=var_id(i)
651 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvsms))) THEN
652 got_var(idvsms)=.true.
653 irp(ng)%Vid(idvsms)=var_id(i)
654# endif
655# ifdef SOLVE3D
656 ELSE IF (trim(var_name(i)).eq.trim(vname(1,iduvel))) THEN
657 got_var(iduvel)=.true.
658 irp(ng)%Vid(iduvel)=var_id(i)
659 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvel))) THEN
660 got_var(idvvel)=.true.
661 irp(ng)%Vid(idvvel)=var_id(i)
662# ifdef ADJUST_BOUNDARY
663 ELSE IF (trim(var_name(i)).eq. &
664 & trim(vname(1,idsbry(isuvel)))) THEN
665 got_var(idsbry(isuvel))=.true.
666 irp(ng)%Vid(idsbry(isuvel))=var_id(i)
667 ELSE IF (trim(var_name(i)).eq. &
668 & trim(vname(1,idsbry(isvvel)))) THEN
669 got_var(idsbry(isvvel))=.true.
670 irp(ng)%Vid(idsbry(isvvel))=var_id(i)
671# endif
672# if defined BVF_MIXING || defined LMD_MIXING || \
673 defined gls_mixing || defined my25_mixing
674 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvis))) THEN
675 got_var(idvvis)=.true.
676 irp(ng)%Vid(idvvis)=var_id(i)
677 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idtdif))) THEN
678 got_var(idtdif)=.true.
679 irp(ng)%Vid(idtdif)=var_id(i)
680 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsdif))) THEN
681 got_var(idsdif)=.true.
682 irp(ng)%Vid(idsdif)=var_id(i)
683# endif
684# endif
685 END IF
686# ifdef SOLVE3D
687 DO itrc=1,nt(ng)
688 IF (trim(var_name(i)).eq.trim(vname(1,idtvar(itrc)))) THEN
689 got_var(idtvar(itrc))=.true.
690 irp(ng)%Tid(itrc)=var_id(i)
691# ifdef ADJUST_BOUNDARY
692 ELSE IF (trim(var_name(i)).eq. &
693 & trim(vname(1,idsbry(istvar(itrc))))) THEN
694 got_var(idsbry(istvar(itrc)))=.true.
695 irp(ng)%Vid(idsbry(istvar(itrc)))=var_id(i)
696# endif
697# ifdef ADJUST_STFLUX
698 ELSE IF (trim(var_name(i)).eq. &
699 & trim(vname(1,idtsur(itrc)))) THEN
700 got_var(idtsur(itrc))=.true.
701 irp(ng)%Vid(idtsur(itrc))=var_id(i)
702# endif
703 END IF
704 END DO
705# endif
706 END DO
707!
708! Check if initialization variables are available in input NetCDF
709! file.
710!
711 IF (.not.got_var(idtime)) THEN
712 IF (master) WRITE (stdout,30) trim(vname(1,idtime)), &
713 & trim(ncname)
714 exit_flag=3
715 RETURN
716 END IF
717 IF (.not.got_var(idfsur)) THEN
718 IF (master) WRITE (stdout,30) trim(vname(1,idfsur)), &
719 & trim(ncname)
720 exit_flag=3
721 RETURN
722 END IF
723 IF (.not.got_var(idubar)) THEN
724 IF (master) WRITE (stdout,30) trim(vname(1,idubar)), &
725 & trim(ncname)
726 exit_flag=3
727 RETURN
728 END IF
729 IF (.not.got_var(idvbar)) THEN
730 IF (master) WRITE (stdout,30) trim(vname(1,idvbar)), &
731 & trim(ncname)
732 exit_flag=3
733 RETURN
734 END IF
735# ifdef ADJUST_BOUNDARY
736 IF (.not.got_var(idsbry(isfsur)).and. &
737 & any(lobc(:,isfsur,ng))) THEN
738 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isfsur))), &
739 & trim(ncname)
740 exit_flag=3
741 RETURN
742 END IF
743 IF (.not.got_var(idsbry(isubar)).and. &
744 & any(lobc(:,isubar,ng))) THEN
745 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isubar))), &
746 & trim(ncname)
747 exit_flag=3
748 RETURN
749 END IF
750 IF (.not.got_var(idsbry(isvbar)).and. &
751 & any(lobc(:,isvbar,ng))) THEN
752 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isvbar))), &
753 & trim(ncname)
754 exit_flag=3
755 RETURN
756 END IF
757# endif
758# ifdef ADJUST_WSTRESS
759 IF (.not.got_var(idusms)) THEN
760 IF (master) WRITE (stdout,30) trim(vname(1,idusms)), &
761 & trim(ncname)
762 exit_flag=3
763 RETURN
764 END IF
765 IF (.not.got_var(idvsms)) THEN
766 IF (master) WRITE (stdout,30) trim(vname(1,idvsms)), &
767 & trim(ncname)
768 exit_flag=3
769 RETURN
770 END IF
771# endif
772# ifdef SOLVE3D
773 IF (.not.got_var(iduvel)) THEN
774 IF (master) WRITE (stdout,30) trim(vname(1,iduvel)), &
775 & trim(ncname)
776 exit_flag=3
777 RETURN
778 END IF
779 IF (.not.got_var(idvvel)) THEN
780 IF (master) WRITE (stdout,30) trim(vname(1,idvvel)), &
781 & trim(ncname)
782 exit_flag=3
783 RETURN
784 END IF
785# ifdef ADJUST_BOUNDARY
786 IF (.not.got_var(idsbry(isuvel)).and. &
787 & any(lobc(:,isuvel,ng))) THEN
788 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isuvel))), &
789 & trim(ncname)
790 exit_flag=3
791 RETURN
792 END IF
793 IF (.not.got_var(idsbry(isvvel)).and. &
794 & any(lobc(:,isvvel,ng))) THEN
795 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isvvel))), &
796 & trim(ncname)
797 exit_flag=3
798 RETURN
799 END IF
800# endif
801# ifdef ADJUST_BOUNDARY
802 IF (.not.got_var(idsbry(isuvel)).and. &
803 & any(lobc(:,isuvel,ng))) THEN
804 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isuvel))), &
805 & trim(ncname)
806 exit_flag=3
807 RETURN
808 END IF
809 IF (.not.got_var(idsbry(isvvel)).and. &
810 & any(lobc(:,isvvel,ng))) THEN
811 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isvvel))), &
812 & trim(ncname)
813 exit_flag=3
814 RETURN
815 END IF
816# endif
817# endif
818# ifdef SOLVE3D
819 DO itrc=1,nt(ng)
820 IF (.not.got_var(idtvar(itrc))) THEN
821 IF (master) WRITE (stdout,30) trim(vname(1,idtvar(itrc))), &
822 & trim(ncname)
823 exit_flag=3
824 RETURN
825 END IF
826# ifdef ADJUST_BOUNDARY
827 IF (.not.got_var(idsbry(istvar(itrc))).and. &
828 & any(lobc(:,istvar(itrc),ng))) THEN
829 IF (master) WRITE (stdout,30) &
830 & trim(vname(1,idsbry(istvar(itrc)))), &
831 & trim(ncname)
832 exit_flag=3
833 RETURN
834 END IF
835# endif
836# ifdef ADJUST_STFLUX
837 IF (.not.got_var(idtsur(itrc)).and.lstflux(itrc,ng)) THEN
838 IF (master) WRITE (stdout,30) trim(vname(1,idtsur(itrc))), &
839 & trim(ncname)
840 exit_flag=3
841 RETURN
842 END IF
843# endif
844 END DO
845# endif
846!
847! Set unlimited time record dimension to the appropriate value.
848!
849 irp(ng)%Rindex=rec_size
850 fcount=irp(ng)%Fcount
851 irp(ng)%Nrec(fcount)=rec_size
852 END IF
853!
854 10 FORMAT (/,' RP_DEF_INI_NF90 - unable to open initial NetCDF', &
855 & ' file: ',a)
856 20 FORMAT (/,' RP_DEF_INI_NF90 - unable to put in define mode', &
857 & ' initial NetCDF file: ',a)
858 30 FORMAT (/,' RP_DEF_INI_NF90 - unable to find variable: ',a,2x, &
859 & ' in file: ',a)
860!
861 RETURN
subroutine, public netcdf_redef(ng, model, ncname, ncid)
integer, dimension(mdims) dim_id
Definition mod_netcdf.F:158
subroutine, public netcdf_check_dim(ng, model, ncname, ncid)
Definition mod_netcdf.F:538
character(len=100), dimension(mdims) dim_name
Definition mod_netcdf.F:168
subroutine, public netcdf_open(ng, model, ncname, omode, ncid)
integer, parameter nf_fout
Definition mod_netcdf.F:188
subroutine, public netcdf_enddef(ng, model, ncname, ncid)
character(len=100), dimension(mvars) var_name
Definition mod_netcdf.F:169
integer n_dim
Definition mod_netcdf.F:151
integer, dimension(mvars) var_id
Definition mod_netcdf.F:160
integer n_var
Definition mod_netcdf.F:152
integer rec_size
Definition mod_netcdf.F:156
integer rec_id
Definition mod_netcdf.F:155
subroutine, public netcdf_inq_var(ng, model, ncname, ncid, myvarname, searchvar, varid, nvardim, nvaratt)

References mod_netcdf::dim_id, mod_netcdf::dim_name, mod_scalars::exit_flag, strings_mod::founderror(), mod_ncparam::idfsur, mod_ncparam::idsbry, mod_ncparam::idsdif, mod_ncparam::idtdif, mod_ncparam::idtime, mod_ncparam::idtsur, mod_ncparam::idtvar, mod_ncparam::idubar, mod_ncparam::idusms, mod_ncparam::iduvel, mod_ncparam::idvbar, mod_ncparam::idvsms, mod_ncparam::idvvel, mod_ncparam::idvvis, mod_ncparam::iinfo, mod_param::iobounds, mod_iounits::irp, mod_param::irpm, mod_scalars::isalt, mod_ncparam::isfsur, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::itemp, mod_scalars::ldefirp, mod_scalars::lobc, mod_scalars::lstflux, mod_parallel::master, mod_netcdf::n_dim, mod_netcdf::n_var, mod_scalars::nbrec, mod_netcdf::netcdf_check_dim(), mod_netcdf::netcdf_enddef(), mod_netcdf::netcdf_inq_var(), mod_netcdf::netcdf_open(), mod_netcdf::netcdf_redef(), mod_netcdf::nf_fout, mod_scalars::nfrec, mod_scalars::noerror, mod_param::nt, mod_ncparam::nv, mod_param::r2dvar, mod_netcdf::rec_id, mod_netcdf::rec_size, mod_iounits::sourcefile, mod_iounits::stdout, mod_param::u2dvar, mod_param::v2dvar, mod_netcdf::var_id, mod_netcdf::var_name, and mod_ncparam::vname.

Referenced by rp_def_ini().

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

◆ rp_def_ini_pio()

subroutine, private rp_def_ini_mod::rp_def_ini_pio ( integer, intent(in) ng)
private

Definition at line 867 of file rp_def_ini.F.

868!***********************************************************************
869!
871!
872! Imported variable declarations.
873!
874 integer, intent(in) :: ng
875!
876! Local variable declarations.
877!
878 logical :: Ldefine = .false.
879# ifdef ADJUST_BOUNDARY
880 logical :: got_IorJ, got_boundary, got_obc_adjust
881# endif
882# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
883 logical :: got_frc_adjust
884# endif
885 logical :: got_var(NV)
886!
887 integer, parameter :: Natt = 25
888
889 integer :: i, j, ifield, itrc, status
890 integer :: Fcount
891# ifdef ADJUST_BOUNDARY
892 integer :: IorJdim, brecdim
893# endif
894# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
895 integer :: frecdim
896# endif
897 integer :: DimIDs(nDimID)
898# ifdef ADJUST_BOUNDARY
899 integer :: t2dobc(4)
900# ifdef SOLVE3D
901 integer :: t3dobc(5)
902# endif
903# endif
904# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
905 integer :: t4dfrc(4), u4dfrc(4), v4dfrc(4)
906# endif
907!
908 real(r8) :: Aval(6)
909!
910 character (len=256) :: ncname
911 character (len=MaxLen) :: Vinfo(Natt)
912
913 character (len=*), parameter :: MyFile = &
914 & __FILE__//", rp_def_ini_pio"
915!
916 sourcefile=myfile
917
918# if defined ADJUST_BOUNDARY || \
919 defined adjust_stflux || defined adjust_wstress
920!
921!=======================================================================
922! Open existing representer model initial conditions and define new
923! variables.
924!=======================================================================
925!
926 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
927 ncname=irp(ng)%name
928!
929! Open representer model initialization file for read/write.
930!
931 CALL pio_netcdf_open (ng, irpm, ncname, 1, irp(ng)%pioFile)
932 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
933 IF (master) WRITE (stdout,10) trim(ncname)
934 RETURN
935 END IF
936!
937! Inquire about the dimensions and check for consistency.
938!
939 CALL pio_netcdf_check_dim (ng, irpm, ncname, &
940 & piofile = irp(ng)%pioFile)
941 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
942!
943! Inquire about the variables.
944!
945 CALL pio_netcdf_inq_var (ng, irpm, ncname, &
946 & piofile = irp(ng)%pioFile)
947 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
948!
949! Check if surface forcing variables have been already defined.
950!
951 DO i=1,nv
952 got_var(i)=.false.
953 END DO
954!
955 DO i=1,n_var
956# ifdef ADJUST_BOUNDARY
957 IF (trim(var_name(i)).eq. &
958 & trim(vname(1,idsbry(isfsur)))) THEN
959 got_var(idsbry(isfsur))=.true.
960 irp(ng)%Vid(idsbry(isfsur))=var_id(i)
961 ELSE IF (trim(var_name(i)).eq. &
962 & trim(vname(1,idsbry(isubar)))) THEN
963 got_var(idsbry(isubar))=.true.
964 irp(ng)%Vid(idsbry(isubar))=var_id(i)
965 ELSE IF (trim(var_name(i)).eq. &
966 & trim(vname(1,idsbry(isvbar)))) THEN
967 got_var(idsbry(isvbar))=.true.
968 irp(ng)%Vid(idsbry(isvbar))=var_id(i)
969# ifdef SOLVE3D
970 ELSE IF (trim(var_name(i)).eq. &
971 & trim(vname(1,idsbry(isuvel)))) THEN
972 got_var(idsbry(isuvel))=.true.
973 irp(ng)%Vid(idsbry(isuvel))=var_id(i)
974 ELSE IF (trim(var_name(i)).eq. &
975 & trim(vname(1,idsbry(isvvel)))) THEN
976 got_var(idsbry(isvvel))=.true.
977 irp(ng)%Vid(idsbry(isvvel))=var_id(i)
978# endif
979 END IF
980# ifdef SOLVE3D
981 DO itrc=1,nt(ng)
982 IF (trim(var_name(i)).eq. &
983 & trim(vname(1,idsbry(istvar(itrc))))) THEN
984 got_var(idsbry(istvar(itrc)))=.true.
985 irp(ng)%Vid(idsbry(istvar(itrc)))=var_id(i)
986 END IF
987 END DO
988# endif
989# endif
990# ifdef ADJUST_WSTRESS
991 IF (trim(var_name(i)).eq.trim(vname(1,idusms))) THEN
992 got_var(idusms)=.true.
993 irp(ng)%Vid(idusms)=var_id(i)
994 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvsms))) THEN
995 got_var(idvsms)=.true.
996 irp(ng)%Vid(idvsms)=var_id(i)
997 END IF
998# endif
999# if defined ADJUST_STFLUX && defined SOLVE3D
1000 DO itrc=1,nt(ng)
1001 IF (lstflux(itrc,ng)) THEN
1002 IF (trim(var_name(i)).eq.trim(vname(1,idtsur(itrc)))) THEN
1003 got_var(idtsur(itrc))=.true.
1004 irp(ng)%Vid(idtsur(itrc))=var_id(i)
1005 END IF
1006 END IF
1007 END DO
1008# endif
1009 END DO
1010
1011# ifdef ADJUST_BOUNDARY
1012 IF (.not.got_var(idsbry(isfsur))) ldefine=.true.
1013 IF (.not.got_var(idsbry(isubar))) ldefine=.true.
1014 IF (.not.got_var(idsbry(isvbar))) ldefine=.true.
1015# ifdef SOLVE3D
1016 IF (.not.got_var(idsbry(isuvel))) ldefine=.true.
1017 IF (.not.got_var(idsbry(isvvel))) ldefine=.true.
1018 DO itrc=1,nt(ng)
1019 IF (.not.got_var(idsbry(istvar(itrc)))) ldefine=.true.
1020 END DO
1021# endif
1022# endif
1023# ifdef ADJUST_WSTRESS
1024 IF (.not.got_var(idusms)) ldefine=.true.
1025 IF (.not.got_var(idvsms)) ldefine=.true.
1026# endif
1027# if defined ADJUST_STFLUX && defined SOLVE3D
1028 DO itrc=1,nt(ng)
1029 IF (lstflux(itrc,ng)) THEN
1030 IF (.not.got_var(idtsur(itrc))) ldefine=.true.
1031 END IF
1032 END DO
1033# endif
1034!
1035! Put existing file into define mode so new variables can be added.
1036!
1037 IF (ldefine) THEN
1038 CALL pio_netcdf_redef (ng, irpm, ncname, irp(ng)%pioFile)
1039 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1040 IF (master) WRITE (stdout,20) trim(ncname)
1041 RETURN
1042 END IF
1043 END IF
1044!
1045!-----------------------------------------------------------------------
1046! Define the dimensions of staggered fields.
1047!-----------------------------------------------------------------------
1048!
1049 define: IF (ldefine) THEN
1050
1051# ifdef ADJUST_BOUNDARY
1052 got_iorj=.false.
1053 got_boundary=.false.
1054 got_obc_adjust=.false.
1055# endif
1056# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
1057 got_frc_adjust=.false.
1058# endif
1059 DO i=1,n_dim
1060 SELECT CASE (trim(adjustl(dim_name(i))))
1061 CASE ('xi_rho')
1062 dimids( 1)=dim_id(i)
1063 CASE ('xi_u')
1064 dimids( 2)=dim_id(i)
1065 CASE ('xi_v')
1066 dimids( 3)=dim_id(i)
1067 CASE ('eta_rho')
1068 dimids( 5)=dim_id(i)
1069 CASE ('eta_u')
1070 dimids( 6)=dim_id(i)
1071 CASE ('eta_v')
1072 dimids( 7)=dim_id(i)
1073# ifdef SOLVE3D
1074 CASE ('s_rho')
1075 dimids( 9)=dim_id(i)
1076 CASE ('s_w')
1077 dimids(10)=dim_id(i)
1078# endif
1079# ifdef ADJUST_BOUNDARY
1080 CASE ('boundary')
1081 dimids(14)=dim_id(i)
1082 got_boundary=.true.
1083 CASE ('IorJ')
1084 iorjdim=dim_id(i)
1085 got_iorj=.true.
1086# endif
1087# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
1088 CASE ('frc_adjust')
1089 frecdim=dim_id(i)
1090 got_frc_adjust=.true.
1091# endif
1092# ifdef ADJUST_BOUNDARY
1093 CASE ('obc_adjust')
1094 brecdim=dim_id(i)
1095 got_obc_adjust=.true.
1096# endif
1097 END SELECT
1098 END DO
1099!
1100 dimids(12)=rec_id
1101# ifdef ADJUST_BOUNDARY
1102 IF (.not.got_boundary) THEN
1103 status=def_dim(ng, irpm, irp(ng)%pioFile, ncname, &
1104 & 'boundary', 4, dimids(14))
1105 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1106 END IF
1107 IF (.not.got_iorj) THEN
1108 status=def_dim(ng, irpm, irp(ng)%pioFile, ncname, &
1109 & 'IorJ', iobounds(ng)%IorJ, iorjdim)
1110 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1111 END IF
1112# endif
1113# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
1114 IF (.not.got_frc_adjust) THEN
1115 status=def_dim(ng, irpm, irp(ng)%pioFile, ncname, &
1116 & 'frc_adjust', nfrec(ng), frecdim)
1117 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1118 END IF
1119# endif
1120# ifdef ADJUST_BOUNDARY
1121 IF (.not.got_obc_adjust) THEN
1122 status=def_dim(ng, irpm, irp(ng)%pioFile, ncname, &
1123 & 'obc_adjust', nbrec(ng), brecdim)
1124 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1125 END IF
1126# endif
1127!
1128! Define dimension vectors for staggered tracer type variables.
1129!
1130# ifdef ADJUST_BOUNDARY
1131 t2dobc(1)=iorjdim
1132 t2dobc(2)=dimids(14)
1133 t2dobc(3)=brecdim
1134 t2dobc(4)=dimids(12)
1135# ifdef SOLVE3D
1136 t3dobc(1)=iorjdim
1137 t3dobc(2)=dimids( 9)
1138 t3dobc(3)=dimids(14)
1139 t3dobc(4)=brecdim
1140 t3dobc(5)=dimids(12)
1141# endif
1142# endif
1143# ifdef ADJUST_STFLUX
1144 t4dfrc(1)=dimids( 1)
1145 t4dfrc(2)=dimids( 5)
1146 t4dfrc(3)=frecdim
1147 t4dfrc(4)=dimids(12)
1148# endif
1149# ifdef ADJUST_WSTRESS
1150!
1151! Define dimension vectors for staggered u-momentum type variables.
1152!
1153 u4dfrc(1)=dimids( 2)
1154 u4dfrc(2)=dimids( 6)
1155 u4dfrc(3)=frecdim
1156 u4dfrc(4)=dimids(12)
1157# endif
1158# ifdef ADJUST_WSTRESS
1159!
1160! Define dimension vectors for staggered v-momentum type variables.
1161!
1162 v4dfrc(1)=dimids( 3)
1163 v4dfrc(2)=dimids( 7)
1164 v4dfrc(3)=frecdim
1165 v4dfrc(4)=dimids(12)
1166# endif
1167!
1168! Initialize local information variable arrays.
1169!
1170 DO i=1,natt
1171 DO j=1,len(vinfo(1))
1172 vinfo(i)(j:j)=' '
1173 END DO
1174 END DO
1175 DO i=1,6
1176 aval(i)=0.0_r8
1177 END DO
1178!
1179!-----------------------------------------------------------------------
1180! Define additional variables. Notice that these variables have their
1181! own fixed time-dimension to allow 4DVAR adjustments at other times
1182! in addition to initialization time.
1183!-----------------------------------------------------------------------
1184
1185# ifdef ADJUST_BOUNDARY
1186!
1187! Define free-surface open boundaries.
1188!
1189 IF (any(lobc(:,isfsur,ng))) THEN
1190 ifield=idsbry(isfsur)
1191 vinfo( 1)=vname(1,ifield)
1192 vinfo( 2)=vname(2,ifield)
1193 vinfo( 3)=vname(3,ifield)
1194 vinfo(14)=vname(4,ifield)
1195 vinfo(16)=vname(1,idtime)
1196 vinfo(21)=vname(6,ifield)
1197 aval(5)=real(iinfo(1,ifield,ng),r8)
1198 irp(ng)%pioVar(ifield)%dkind=pio_fout
1199 irp(ng)%pioVar(ifield)%gtype=r2dobc
1200!
1201 status=def_var(ng, irpm, irp(ng)%pioFile, &
1202 & irp(ng)%pioVar(ifield)%vd, pio_fout, &
1203 & 4, t2dobc, aval, vinfo, ncname, &
1204 & setfillval = .false.)
1205 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1206 END IF
1207!
1208! Define 2D U-momentum component open boundaries.
1209!
1210 IF (any(lobc(:,isubar,ng))) THEN
1211 ifield=idsbry(isubar)
1212 vinfo( 1)=vname(1,ifield)
1213 vinfo( 2)=vname(2,ifield)
1214 vinfo( 3)=vname(3,ifield)
1215 vinfo(14)=vname(4,ifield)
1216 vinfo(16)=vname(1,idtime)
1217 vinfo(21)=vname(6,ifield)
1218 aval(5)=real(iinfo(1,ifield,ng),r8)
1219 irp(ng)%pioVar(ifield)%dkind=pio_fout
1220 irp(ng)%pioVar(ifield)%gtype=u2dobc
1221!
1222 status=def_var(ng, irpm, irp(ng)%pioFile, &
1223 & irp(ng)%pioVar(ifield)%vd, pio_fout, &
1224 & 4, t2dobc, aval, vinfo, ncname, &
1225 & setfillval = .false.)
1226 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1227 END IF
1228!
1229! Define 2D V-momentum component open boundaries.
1230!
1231 IF (any(lobc(:,isvbar,ng))) THEN
1232 ifield=idsbry(isvbar)
1233 vinfo( 1)=vname(1,ifield)
1234 vinfo( 2)=vname(2,ifield)
1235 vinfo( 3)=vname(3,ifield)
1236 vinfo(14)=vname(4,ifield)
1237 vinfo(16)=vname(1,idtime)
1238 vinfo(21)=vname(6,ifield)
1239 aval(5)=real(iinfo(1,ifield,ng),r8)
1240 irp(ng)%pioVar(ifield)%dkind=pio_fout
1241 irp(ng)%pioVar(ifield)%gtype=v2dobc
1242!
1243 status=def_var(ng, irpm, irp(ng)%pioFile, &
1244 & irp(ng)%pioVar(ifield)%vd, pio_fout, &
1245 & 4, t2dobc, aval, vinfo, ncname, &
1246 & setfillval = .false.)
1247 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1248 END IF
1249
1250# ifdef SOLVE3D
1251!
1252! Define 3D U-momentum component open boundaries.
1253!
1254 IF (any(lobc(:,isuvel,ng))) THEN
1255 ifield=idsbry(isuvel)
1256 vinfo( 1)=vname(1,ifield)
1257 vinfo( 2)=vname(2,ifield)
1258 vinfo( 3)=vname(3,ifield)
1259 vinfo(14)=vname(4,ifield)
1260 vinfo(16)=vname(1,idtime)
1261 vinfo(21)=vname(6,ifield)
1262 aval(5)=real(iinfo(1,ifield,ng),r8)
1263 irp(ng)%pioVar(ifield)%dkind=pio_fout
1264 irp(ng)%pioVar(ifield)%gtype=u3dobc
1265!
1266 status=def_var(ng, irpm, irp(ng)%pioFile, &
1267 & irp(ng)%pioVar(ifield)%vd, pio_fout, &
1268 & 5, t3dobc, aval, vinfo, ncname, &
1269 & setfillval = .false.)
1270 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1271 END IF
1272!
1273! Define 3D V-momentum component open boundaries.
1274!
1275 IF (any(lobc(:,isvvel,ng))) THEN
1276 ifield=idsbry(isvvel)
1277 vinfo( 1)=vname(1,ifield)
1278 vinfo( 2)=vname(2,ifield)
1279 vinfo( 3)=vname(3,ifield)
1280 vinfo(14)=vname(4,ifield)
1281 vinfo(16)=vname(1,idtime)
1282 vinfo(21)=vname(6,ifield)
1283 aval(5)=real(iinfo(1,ifield,ng),r8)
1284 irp(ng)%pioVar(ifield)%dkind=pio_fout
1285 irp(ng)%pioVar(ifield)%gtype=v3dobc
1286!
1287 status=def_var(ng, irpm, irp(ng)%pioFile, &
1288 & irp(ng)%pioVar(ifield)%vd, pio_fout, &
1289 & 5, t3dobc, aval, vinfo, ncname, &
1290 & setfillval = .false.)
1291 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1292 END IF
1293!
1294! Define tracer type variables open boundaries.
1295!
1296 DO itrc=1,nt(ng)
1297 IF (any(lobc(:,istvar(itrc),ng))) THEN
1298 ifield=idsbry(istvar(itrc))
1299 vinfo( 1)=vname(1,ifield)
1300 vinfo( 2)=vname(2,ifield)
1301 vinfo( 3)=vname(3,ifield)
1302 vinfo(14)=vname(4,ifield)
1303 vinfo(16)=vname(1,idtime)
1304 vinfo(21)=vname(6,ifield)
1305 aval(5)=real(iinfo(1,ifield,ng),r8)
1306 irp(ng)%pioVar(ifield)%dkind=pio_fout
1307 irp(ng)%pioVar(ifield)%gtype=r3dobc
1308!
1309 status=def_var(ng, irpm, irp(ng)%pioFile, &
1310 & irp(ng)%pioVar(ifield)%vd, pio_fout, &
1311 & 5, t3dobc, aval, vinfo, ncname, &
1312 & setfillval = .false.)
1313 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1314 END IF
1315 END DO
1316# endif
1317# endif
1318# ifdef ADJUST_WSTRESS
1319!
1320! Define surface U-momentum stress.
1321!
1322 IF (.not.got_var(idusms)) THEN
1323 vinfo( 1)=vname(1,idusms)
1324 vinfo( 2)=vname(2,idusms)
1325 vinfo( 3)=vname(3,idusms)
1326# if defined WRITE_WATER && defined MASKING
1327 vinfo(20)='mask_u'
1328# endif
1329 vinfo(21)=vname(6,idusms)
1330 vinfo(22)='coordinates'
1331 aval(5)=real(u2dvar,r8)
1332 irp(ng)%pioVar(idusms)%dkind=pio_fout
1333 irp(ng)%pioVar(idusms)%gtype=u2dvar
1334!
1335 status=def_var(ng, irpm, irp(ng)%pioFile, &
1336 & irp(ng)%pioVar(idusms)%vd, pio_fout, &
1337 & 4, u4dfrc, aval, vinfo, ncname)
1338 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1339 END IF
1340!
1341! Define surface V-momentum stress.
1342!
1343 IF (.not.got_var(idvsms)) THEN
1344 vinfo( 1)=vname(1,idvsms)
1345 vinfo( 2)=vname(2,idvsms)
1346 vinfo( 3)=vname(3,idvsms)
1347# if defined WRITE_WATER && defined MASKING
1348 vinfo(20)='mask_v'
1349# endif
1350 vinfo(21)=vname(6,idvsms)
1351 vinfo(22)='coordinates'
1352 aval(5)=real(v2dvar,r8)
1353 irp(ng)%pioVar(idvsms)%dkind=pio_fout
1354 irp(ng)%pioVar(idvsms)%gtype=v2dvar
1355!
1356 status=def_var(ng, irpm, irp(ng)%pioFile, &
1357 & irp(ng)%pioVar(idvsms)%vd, pio_fout, &
1358 & 4, v4dfrc, aval, vinfo, ncname)
1359 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1360 END IF
1361# endif
1362# if defined ADJUST_STFLUX && defined SOLVE3D
1363!
1364! Define surface tracer flux.
1365!
1366 DO itrc=1,nt(ng)
1367 IF (.not.got_var(idtsur(itrc)).and.lstflux(itrc,ng)) THEN
1368 vinfo( 1)=vname(1,idtsur(itrc))
1369 vinfo( 2)=trim(vname(2,idtsur(itrc)))
1370 vinfo( 3)=vname(3,idtsur(itrc))
1371 IF (itrc.eq.itemp) THEN
1372 vinfo(11)='upward flux, cooling'
1373 vinfo(12)='downward flux, heating'
1374 ELSE IF (itrc.eq.isalt) THEN
1375 vinfo(11)='upward flux, freshening (net precipitation)'
1376 vinfo(12)='downward flux, salting (net evaporation)'
1377 END IF
1378# if defined WRITE_WATER && defined MASKING
1379 vinfo(20)='mask_rho'
1380# endif
1381 vinfo(21)=vname(6,idtsur(itrc))
1382 vinfo(22)='coordinates'
1383 aval(5)=real(r2dvar,r8)
1384 irp(ng)%pioVar(idtsur(itrc))%dkind=pio_fout
1385 irp(ng)%pioVar(idtsur(itrc))%gtype=r2dvar
1386!
1387 status=def_var(ng, irpm, irp(ng)%pioFile, &
1388 & irp(ng)%pioVar(idtsur(itrc))%vd, pio_fout, &
1389 & 4, t4dfrc, aval, vinfo, ncname)
1390 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1391 END IF
1392 END DO
1393# endif
1394!
1395!-----------------------------------------------------------------------
1396! Leave definition mode.
1397!-----------------------------------------------------------------------
1398!
1399 CALL pio_netcdf_enddef (ng, irpm, ncname, irp(ng)%pioFile)
1400 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1401
1402 END IF define
1403# endif
1404!
1405!=======================================================================
1406! Open an existing tangent linear initial file, check its contents, and
1407! prepare for appending data.
1408!=======================================================================
1409!
1410 IF (.not.ldefirp(ng)) THEN
1411 ncname=irp(ng)%name
1412
1413# if !(defined ADJUST_STFLUX || defined ADJUST_WSTRESS)
1414!
1415! Open representer model initialization file for read/write.
1416!
1417 CALL pio_netcdf_open (ng, irpm, ncname, 1, irp(ng)%pioFile)
1418 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
1419 WRITE (stdout,10) trim(ncname)
1420 RETURN
1421 END IF
1422!
1423! Inquire about the dimensions and check for consistency.
1424!
1425 CALL pio_netcdf_check_dim (ng, irpm, ncname, &
1426 & piofile = irp(ng)%pioFile)
1427 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1428!
1429! Inquire about the variables.
1430!
1431 CALL pio_netcdf_inq_var (ng, irpm, ncname, &
1432 & piofile = irp(ng)%pioFile)
1433 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1434!
1435! Initialize logical switches.
1436!
1437 DO i=1,nv
1438 got_var(i)=.false.
1439 END DO
1440# endif
1441!
1442! Scan variable list from input NetCDF and activate switches for
1443! initialization variables. Get variable IDs.
1444!
1445 DO i=1,n_var
1446 IF (trim(var_name(i)).eq.trim(vname(1,idtime))) THEN
1447 got_var(idtime)=.true.
1448 irp(ng)%pioVar(idtime)%vd=var_desc(i)
1449 irp(ng)%pioVar(idtime)%dkind=pio_tout
1450 irp(ng)%pioVar(idtime)%gtype=0
1451 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idfsur))) THEN
1452 got_var(idfsur)=.true.
1453 irp(ng)%pioVar(idfsur)%vd=var_desc(i)
1454 irp(ng)%pioVar(idfsur)%dkind=pio_fout
1455 irp(ng)%pioVar(idfsur)%gtype=r2dvar
1456 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idubar))) THEN
1457 got_var(idubar)=.true.
1458 irp(ng)%pioVar(idubar)%vd=var_desc(i)
1459 irp(ng)%pioVar(idubar)%dkind=pio_fout
1460 irp(ng)%pioVar(idubar)%gtype=u2dvar
1461 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvbar))) THEN
1462 got_var(idvbar)=.true.
1463 irp(ng)%pioVar(idvbar)%vd=var_desc(i)
1464 irp(ng)%pioVar(idvbar)%dkind=pio_fout
1465 irp(ng)%pioVar(idvbar)%gtype=v2dvar
1466# ifdef ADJUST_BOUNDARY
1467 ELSE IF (trim(var_name(i)).eq. &
1468 & trim(vname(1,idsbry(isfsur)))) THEN
1469 got_var(idsbry(isfsur))=.true.
1470 irp(ng)%pioVar(idsbry(isfsur))%vd=var_desc(i)
1471 irp(ng)%pioVar(idsbry(isfsur))%dkind=pio_fout
1472 irp(ng)%pioVar(idsbry(isfsur))%gtype=r2dobc
1473 ELSE IF (trim(var_name(i)).eq. &
1474 & trim(vname(1,idsbry(isubar)))) THEN
1475 got_var(idsbry(isubar))=.true.
1476 irp(ng)%pioVar(idsbry(isubar))%vd=var_desc(i)
1477 irp(ng)%pioVar(idsbry(isubar))%dkind=pio_fout
1478 irp(ng)%pioVar(idsbry(isubar))%gtype=u2dobc
1479 ELSE IF (trim(var_name(i)).eq. &
1480 & trim(vname(1,idsbry(isvbar)))) THEN
1481 got_var(idsbry(isvbar))=.true.
1482 irp(ng)%pioVar(idsbry(isvbar))%vd=var_desc(i)
1483 irp(ng)%pioVar(idsbry(isvbar))%dkind=pio_fout
1484 irp(ng)%pioVar(idsbry(isvbar))%gtype=v2dobc
1485# endif
1486# ifdef ADJUST_WSTRESS
1487 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idusms))) THEN
1488 got_var(idusms)=.true.
1489 irp(ng)%pioVar(idusms)%vd=var_desc(i)
1490 irp(ng)%pioVar(idusms)%dkind=pio_fout
1491 irp(ng)%pioVar(idusms)%gtype=u2dvar
1492 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvsms))) THEN
1493 got_var(idvsms)=.true.
1494 irp(ng)%pioVar(idvsms)%vd=var_desc(i)
1495 irp(ng)%pioVar(idvsms)%dkind=pio_fout
1496 irp(ng)%pioVar(idvsms)%gtype=v2dvar
1497# endif
1498# ifdef SOLVE3D
1499 ELSE IF (trim(var_name(i)).eq.trim(vname(1,iduvel))) THEN
1500 got_var(iduvel)=.true.
1501 irp(ng)%pioVar(iduvel)%vd=var_desc(i)
1502 irp(ng)%pioVar(iduvel)%dkind=pio_fout
1503 irp(ng)%pioVar(iduvel)%gtype=u3dvar
1504 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvel))) THEN
1505 got_var(idvvel)=.true.
1506 irp(ng)%pioVar(idvvel)%vd=var_desc(i)
1507 irp(ng)%pioVar(idvvel)%dkind=pio_fout
1508 irp(ng)%pioVar(idvvel)%gtype=v3dvar
1509# ifdef ADJUST_BOUNDARY
1510 ELSE IF (trim(var_name(i)).eq. &
1511 & trim(vname(1,idsbry(isuvel)))) THEN
1512 got_var(idsbry(isuvel))=.true.
1513 irp(ng)%pioVar(idsbry(isuvel))%vd=var_desc(i)
1514 irp(ng)%pioVar(idsbry(isuvel))%dkind=pio_fout
1515 irp(ng)%pioVar(idsbry(isuvel))%gtype=u3dobc
1516 ELSE IF (trim(var_name(i)).eq. &
1517 & trim(vname(1,idsbry(isvvel)))) THEN
1518 got_var(idsbry(isvvel))=.true.
1519 irp(ng)%pioVar(idsbry(isvvel))%vd=var_desc(i)
1520 irp(ng)%pioVar(idsbry(isvvel))%dkind=pio_fout
1521 irp(ng)%pioVar(idsbry(isvvel))%gtype=v3dobc
1522# endif
1523# if defined BVF_MIXING || defined LMD_MIXING || \
1524 defined gls_mixing || defined my25_mixing
1525 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvvis))) THEN
1526 got_var(idvvis)=.true.
1527 irp(ng)%pioVar(idvvis)%vd=var_desc(i)
1528 irp(ng)%pioVar(idvvis)%dkind=pio_fout
1529 irp(ng)%pioVar(idvvis)%gtype=w3dvar
1530 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idtdif))) THEN
1531 got_var(idtdif)=.true.
1532 irp(ng)%pioVar(idtdif)%vd=var_desc(i)
1533 irp(ng)%pioVar(idtdif)%dkind=pio_fout
1534 irp(ng)%pioVar(idtdif)%gtype=w3dvar
1535 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idsdif))) THEN
1536 got_var(idsdif)=.true.
1537 irp(ng)%pioVar(idsdif)%vd=var_desc(i)
1538 irp(ng)%pioVar(idsdif)%dkind=pio_fout
1539 irp(ng)%pioVar(idsdif)%gtype=w3dvar
1540# endif
1541# endif
1542 END IF
1543# ifdef SOLVE3D
1544 DO itrc=1,nt(ng)
1545 IF (trim(var_name(i)).eq.trim(vname(1,idtvar(itrc)))) THEN
1546 got_var(idtvar(itrc))=.true.
1547 irp(ng)%pioTrc(itrc)%vd=var_desc(i)
1548 irp(ng)%pioTrc(itrc)%dkind=pio_fout
1549 irp(ng)%pioTrc(itrc)%gtype=r3dvar
1550# ifdef ADJUST_BOUNDARY
1551 ELSE IF (trim(var_name(i)).eq. &
1552 & trim(vname(1,idsbry(istvar(itrc))))) THEN
1553 got_var(idsbry(istvar(itrc)))=.true.
1554 irp(ng)%pioVar(idsbry(istvar(itrc)))%vd=var_desc(i)
1555 irp(ng)%pioVar(idsbry(istvar(itrc)))%dkind=pio_fout
1556 irp(ng)%pioVar(idsbry(istvar(itrc)))%gtype=r3dobc
1557# endif
1558# ifdef ADJUST_STFLUX
1559 ELSE IF (trim(var_name(i)).eq. &
1560 & trim(vname(1,idtsur(itrc)))) THEN
1561 got_var(idtsur(itrc))=.true.
1562 irp(ng)%pioVar(idtsur(itrc))%vd=var_desc(i)
1563 irp(ng)%pioVar(idtsur(itrc))%dkind=pio_fout
1564 irp(ng)%pioVar(idtsur(itrc))%gtype=r2dvar
1565# endif
1566 END IF
1567 END DO
1568# endif
1569 END DO
1570!
1571! Check if initialization variables are available in input NetCDF
1572! file.
1573!
1574 IF (.not.got_var(idtime)) THEN
1575 IF (master) WRITE (stdout,30) trim(vname(1,idtime)), &
1576 & trim(ncname)
1577 exit_flag=3
1578 RETURN
1579 END IF
1580 IF (.not.got_var(idfsur)) THEN
1581 IF (master) WRITE (stdout,30) trim(vname(1,idfsur)), &
1582 & trim(ncname)
1583 exit_flag=3
1584 RETURN
1585 END IF
1586 IF (.not.got_var(idubar)) THEN
1587 IF (master) WRITE (stdout,30) trim(vname(1,idubar)), &
1588 & trim(ncname)
1589 exit_flag=3
1590 RETURN
1591 END IF
1592 IF (.not.got_var(idvbar)) THEN
1593 IF (master) WRITE (stdout,30) trim(vname(1,idvbar)), &
1594 & trim(ncname)
1595 exit_flag=3
1596 RETURN
1597 END IF
1598# ifdef ADJUST_BOUNDARY
1599 IF (.not.got_var(idsbry(isfsur)).and. &
1600 & any(lobc(:,isfsur,ng))) THEN
1601 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isfsur))), &
1602 & trim(ncname)
1603 exit_flag=3
1604 RETURN
1605 END IF
1606 IF (.not.got_var(idsbry(isubar)).and. &
1607 & any(lobc(:,isubar,ng))) THEN
1608 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isubar))), &
1609 & trim(ncname)
1610 exit_flag=3
1611 RETURN
1612 END IF
1613 IF (.not.got_var(idsbry(isvbar)).and. &
1614 & any(lobc(:,isvbar,ng))) THEN
1615 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isvbar))), &
1616 & trim(ncname)
1617 exit_flag=3
1618 RETURN
1619 END IF
1620# endif
1621# ifdef ADJUST_WSTRESS
1622 IF (.not.got_var(idusms)) THEN
1623 IF (master) WRITE (stdout,30) trim(vname(1,idusms)), &
1624 & trim(ncname)
1625 exit_flag=3
1626 RETURN
1627 END IF
1628 IF (.not.got_var(idvsms)) THEN
1629 IF (master) WRITE (stdout,30) trim(vname(1,idvsms)), &
1630 & trim(ncname)
1631 exit_flag=3
1632 RETURN
1633 END IF
1634# endif
1635# ifdef SOLVE3D
1636 IF (.not.got_var(iduvel)) THEN
1637 IF (master) WRITE (stdout,30) trim(vname(1,iduvel)), &
1638 & trim(ncname)
1639 exit_flag=3
1640 RETURN
1641 END IF
1642 IF (.not.got_var(idvvel)) THEN
1643 IF (master) WRITE (stdout,30) trim(vname(1,idvvel)), &
1644 & trim(ncname)
1645 exit_flag=3
1646 RETURN
1647 END IF
1648# ifdef ADJUST_BOUNDARY
1649 IF (.not.got_var(idsbry(isuvel)).and. &
1650 & any(lobc(:,isuvel,ng))) THEN
1651 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isuvel))), &
1652 & trim(ncname)
1653 exit_flag=3
1654 RETURN
1655 END IF
1656 IF (.not.got_var(idsbry(isvvel)).and. &
1657 & any(lobc(:,isvvel,ng))) THEN
1658 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isvvel))), &
1659 & trim(ncname)
1660 exit_flag=3
1661 RETURN
1662 END IF
1663# endif
1664# ifdef ADJUST_BOUNDARY
1665 IF (.not.got_var(idsbry(isuvel)).and. &
1666 & any(lobc(:,isuvel,ng))) THEN
1667 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isuvel))), &
1668 & trim(ncname)
1669 exit_flag=3
1670 RETURN
1671 END IF
1672 IF (.not.got_var(idsbry(isvvel)).and. &
1673 & any(lobc(:,isvvel,ng))) THEN
1674 IF (master) WRITE (stdout,30) trim(vname(1,idsbry(isvvel))), &
1675 & trim(ncname)
1676 exit_flag=3
1677 RETURN
1678 END IF
1679# endif
1680# endif
1681# ifdef SOLVE3D
1682 DO itrc=1,nt(ng)
1683 IF (.not.got_var(idtvar(itrc))) THEN
1684 IF (master) WRITE (stdout,30) trim(vname(1,idtvar(itrc))), &
1685 & trim(ncname)
1686 exit_flag=3
1687 RETURN
1688 END IF
1689# ifdef ADJUST_BOUNDARY
1690 IF (.not.got_var(idsbry(istvar(itrc))).and. &
1691 & any(lobc(:,istvar(itrc),ng))) THEN
1692 IF (master) WRITE (stdout,30) &
1693 & trim(vname(1,idsbry(istvar(itrc)))), &
1694 & trim(ncname)
1695 exit_flag=3
1696 RETURN
1697 END IF
1698# endif
1699# ifdef ADJUST_STFLUX
1700 IF (.not.got_var(idtsur(itrc)).and.lstflux(itrc,ng)) THEN
1701 IF (master) WRITE (stdout,30) trim(vname(1,idtsur(itrc))), &
1702 & trim(ncname)
1703 exit_flag=3
1704 RETURN
1705 END IF
1706# endif
1707 END DO
1708# endif
1709!
1710! Set unlimited time record dimension to the appropriate value.
1711!
1712 irp(ng)%Rindex=rec_size
1713 fcount=irp(ng)%Fcount
1714 irp(ng)%Nrec(fcount)=rec_size
1715 END IF
1716!
1717 10 FORMAT (/,' RP_DEF_INI_PIO - unable to open initial NetCDF', &
1718 & ' file: ',a)
1719 20 FORMAT (/,' RP_DEF_INI_PIO - unable to put in define mode', &
1720 & ' initial NetCDF file: ',a)
1721 30 FORMAT (/,' RP_DEF_INI_PIO - unable to find variable: ',a,2x, &
1722 & ' in file: ',a)
1723!
1724 RETURN
subroutine, public pio_netcdf_redef(ng, model, ncname, piofile)
integer, parameter pio_fout
type(var_desc_t), dimension(:), pointer var_desc
subroutine, public pio_netcdf_inq_var(ng, model, ncname, piofile, myvarname, searchvar, piovar, nvardim, nvaratt)
subroutine, public pio_netcdf_open(ng, model, ncname, omode, piofile)
subroutine, public pio_netcdf_check_dim(ng, model, ncname, piofile)
integer, parameter pio_tout
subroutine, public pio_netcdf_enddef(ng, model, ncname, piofile)

References mod_scalars::exit_flag, strings_mod::founderror(), mod_ncparam::idfsur, mod_ncparam::idsbry, mod_ncparam::idsdif, mod_ncparam::idtdif, mod_ncparam::idtime, mod_ncparam::idtsur, mod_ncparam::idtvar, mod_ncparam::idubar, mod_ncparam::idusms, mod_ncparam::iduvel, mod_ncparam::idvbar, mod_ncparam::idvsms, mod_ncparam::idvvel, mod_ncparam::idvvis, mod_ncparam::iinfo, mod_param::iobounds, mod_iounits::irp, mod_param::irpm, mod_scalars::isalt, mod_ncparam::isfsur, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::itemp, mod_scalars::ldefirp, mod_scalars::lobc, mod_scalars::lstflux, mod_parallel::master, mod_scalars::nbrec, mod_scalars::nfrec, mod_scalars::noerror, mod_param::nt, mod_ncparam::nv, mod_pio_netcdf::pio_fout, mod_pio_netcdf::pio_netcdf_check_dim(), mod_pio_netcdf::pio_netcdf_enddef(), mod_pio_netcdf::pio_netcdf_inq_var(), mod_pio_netcdf::pio_netcdf_open(), mod_pio_netcdf::pio_netcdf_redef(), mod_pio_netcdf::pio_tout, mod_param::r2dobc, mod_param::r2dvar, mod_param::r3dobc, mod_param::r3dvar, mod_iounits::sourcefile, mod_iounits::stdout, mod_param::u2dobc, mod_param::u2dvar, mod_param::u3dobc, mod_param::u3dvar, mod_param::v2dobc, mod_param::v2dvar, mod_param::v3dobc, mod_param::v3dvar, mod_pio_netcdf::var_desc, mod_ncparam::vname, and mod_param::w3dvar.

Referenced by rp_def_ini().

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