ROMS
Loading...
Searching...
No Matches
ice_output.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined ICE_MODEL && defined SOLVE3D
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This module defines/writes the ice model prognostic variables into !
13! output NetCDF files. !
14! !
15!=======================================================================
16!
17 USE mod_param
18 USE mod_parallel
19 USE mod_grid
20 USE mod_ice
21 USE mod_ncparam
22 USE mod_stepping
23!
24 USE def_var_mod, ONLY : def_var
25# ifdef STATIONS
27# endif
29 USE uv_rotate_mod, ONLY : uv_rotate2d
30 USE strings_mod, ONLY : founderror
31!
32 implicit none
33!
34 PUBLIC :: ice_def_nf90
35# if defined PIO_LIB && defined DISTRIBUTE
36 PUBLIC :: ice_def_pio
37# endif
38# ifdef STATIONS
39 PUBLIC :: ice_def_station_nf90
40# if defined PIO_LIB && defined DISTRIBUTE
41 PUBLIC :: ice_def_station_pio
42# endif
43# endif
44 PUBLIC :: ice_wrt_nf90
45# if defined PIO_LIB && defined DISTRIBUTE
46 PUBLIC :: ice_wrt_pio
47# endif
48# ifdef STATIONS
49 PUBLIC :: ice_wrt_station_nf90
50# if defined PIO_LIB && defined DISTRIBUTE
51 PUBLIC :: ice_wrt_station_pio
52# endif
53# endif
54!
55 CONTAINS
56!
57!***********************************************************************
58 SUBROUTINE ice_def_nf90 (ng, model, ldef, VarOut, S, &
59 & t2dgrd, u2dgrd, v2dgrd)
60!***********************************************************************
61!
62 USE mod_netcdf
63!
64! Imported variable declarations.
65!
66 logical, intent(in) :: ldef, VarOut(NV,Ngrids)
67!
68 integer, intent(in) :: ng, model
69 integer, intent(in), optional :: t2dgrd(:), u2dgrd(:), v2dgrd(:)
70!
71 TYPE(T_IO), intent(inout) :: S(Ngrids)
72!
73! Local variable declarations.
74!
75 logical :: got_var(NV)
76 logical :: LdefVar
77!
78 integer, parameter :: Natt = 25
79
80 integer :: i, ifield, j, nf, nvd3, status, vtype
81
82 integer :: icegrd(3)
83!
84
85 real(r8) :: Aval(6)
86!
87 character (len=120) :: Vinfo(Natt)
88 character (len=256) :: ncname
89!
90 character (len=*), parameter :: MyFile = &
91 & __FILE__//", ice_def_nf90"
92!
93 sourcefile=myfile
94!
95!-----------------------------------------------------------------------
96! Define prognostic ice variables into specified output NetCDF.
97!-----------------------------------------------------------------------
98!
99 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
100 ncname=s(ng)%name
101!
102! Set number of dimensions for output variables.
103!
104# if defined WRITE_WATER && defined MASKING
105 nvd3=2
106# else
107 nvd3=3
108# endif
109!
110! Initialize local information variable arrays.
111!
112 DO i=1,natt
113 DO j=1,len(vinfo(1))
114 vinfo(i)(j:j)=' '
115 END DO
116 END DO
117 DO i=1,6
118 aval(i)=0.0_r8
119 END DO
120!
121!-----------------------------------------------------------------------
122! Define ice model state prognostic variables.
123!-----------------------------------------------------------------------
124!
125 define : IF (ldef) THEN
126
127 DO nf=1,nices
128 IF (isice(nf).gt.0) THEN
129 ifield=isice(nf)
130 IF (s(ng)%ncid.eq.rst(ng)%ncid) THEN
131 ldefvar=.true.
132 vtype=nf_frst
133 ELSE
134 ldefvar=varout(ifield,ng)
135 vtype=nf_fout
136 END IF
137 IF (ldefvar) THEN
138 vinfo( 1)=vname(1,ifield)
139 vinfo( 2)=vname(2,ifield)
140 vinfo( 3)=vname(3,ifield)
141 vinfo(14)=vname(4,ifield)
142 vinfo(16)=vname(1,idtime)
143 vinfo(21)=vname(6,ifield)
144 vinfo(22)='coordinates'
145 aval(5)=real(iinfo(1,ifield,ng),r8)
146!
147 SELECT CASE (nf)
148 CASE (isuice)
149 icegrd(1:3)=u2dgrd(1:3)
150# if defined WRITE_WATER && defined MASKING
151 vinfo(20)='mask_u'
152# endif
153 CASE (isvice)
154 icegrd(1:3)=v2dgrd(1:3)
155# if defined WRITE_WATER && defined MASKING
156 vinfo(20)='mask_v'
157# endif
158 CASE DEFAULT
159 icegrd(1:3)=t2dgrd(1:3)
160# if defined WRITE_WATER && defined MASKING
161 vinfo(20)='mask_rho'
162# endif
163 END SELECT
164!
165 status=def_var(ng, model, s(ng)%ncid, &
166 & s(ng)%Vid(ifield), &
167 & vtype, nvd3, icegrd, aval, vinfo, ncname)
168 IF (founderror(exit_flag, noerror, &
169 & __line__, myfile)) RETURN
170 END IF
171 END IF
172 END DO
173!
174! Define ice EASTward velocity at RHO-points.
175!
176 IF (varout(iduier,ng)) THEN
177 vinfo( 1)=vname(1,iduier)
178 vinfo( 2)=vname(2,iduier)
179 vinfo( 3)=vname(3,iduier)
180 vinfo(14)=vname(4,iduier)
181 vinfo(16)=vname(1,idtime)
182# if defined WRITE_WATER && defined MASKING
183 vinfo(20)='mask_rho'
184# endif
185 vinfo(21)=vname(6,iduier)
186 vinfo(22)='coordinates'
187 aval(5)=real(iinfo(1,iduier,ng),r8)
188 status=def_var(ng, model, s(ng)%ncid, &
189 & s(ng)%Vid(iduier), &
190 & vtype, nvd3, t2dgrd, aval, vinfo, ncname)
191 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
192 END IF
193!
194! Define ice NORTHward velocity at RHO-points.
195!
196 IF (varout(idvinr,ng)) THEN
197 vinfo( 1)=vname(1,idvinr)
198 vinfo( 2)=vname(2,idvinr)
199 vinfo( 3)=vname(3,idvinr)
200 vinfo(14)=vname(4,idvinr)
201 vinfo(16)=vname(1,idtime)
202# if defined WRITE_WATER && defined MASKING
203 vinfo(20)='mask_rho'
204# endif
205 vinfo(21)=vname(6,idvinr)
206 vinfo(22)='coordinates'
207 aval(5)=real(iinfo(1,idvinr,ng),r8)
208 status=def_var(ng, model, s(ng)%ncid, &
209 & s(ng)%Vid(idvinr), &
210 & vtype, nvd3, t2dgrd, aval, vinfo, ncname)
211 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
212 END IF
213!
214!-----------------------------------------------------------------------
215! Define ice model internal variables.
216!-----------------------------------------------------------------------
217!
218 DO nf=1,nicef
219 IF (ifice(nf).gt.0) THEN
220 ifield=ifice(nf)
221 IF (varout(ifield,ng)) THEN
222 vinfo( 1)=vname(1,ifield)
223 vinfo( 2)=vname(2,ifield)
224 vinfo( 3)=vname(3,ifield)
225 vinfo(14)=vname(4,ifield)
226 vinfo(16)=vname(1,idtime)
227# if defined WRITE_WATER && defined MASKING
228 vinfo(20)='mask_rho'
229# endif
230 vinfo(21)=vname(6,ifield)
231 vinfo(22)='coordinates'
232 aval(5)=real(iinfo(1,ifield,ng),r8)
233!
234 SELECT CASE (nf)
235 CASE (icwdiv)
236 vinfo(11)='increase ice thickness'
237 vinfo(12)='decrease ice concentration'
238 CASE (icw_ai)
239 vinfo(11)='freezing'
240 vinfo(12)='melting'
241 CASE (icw_ao, icw_io)
242 vinfo(11)='melting'
243 vinfo(12)='freezing'
244 END SELECT
245!
246 status=def_var(ng, model, s(ng)%ncid, &
247 & s(ng)%Vid(ifield), &
248 & vtype, nvd3, t2dgrd, aval, vinfo, ncname)
249 IF (founderror(exit_flag, noerror, &
250 & __line__, myfile)) RETURN
251 END IF
252 END IF
253 END DO
254 END IF define
255!
256!-----------------------------------------------------------------------
257! Otherwise, check existing output file and prepare for appending
258! data.
259!-----------------------------------------------------------------------
260!
261 query : IF (.not.ldef) THEN
262!
263! Initialize locallogical switches.
264!
265 DO i=1,nv
266 got_var(i)=.false.
267 END DO
268!
269! Scan variable list from input NetCDF and activate switches for
270! sea-ice variables. Get variable IDs.
271!
272 DO i=1,n_var
273 DO nf=1,nices
274 IF (isice(nf).gt.0) THEN
275 ifield=isice(nf)
276 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
277 got_var(ifield)=.true.
278 s(ng)%Vid(ifield)=var_id(i)
279 cycle
280 END IF
281 END IF
282 END DO
283!
284 IF (trim(var_name(i)).eq.trim(vname(1,iduier))) THEN
285 got_var(iduier)=.true.
286 s(ng)%Vid(iduier)=var_id(i)
287 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvinr))) THEN
288 got_var(idvinr)=.true.
289 s(ng)%Vid(idvinr)=var_id(i)
290 END IF
291!
292 DO nf=1,nicef
293 IF (ifice(nf).gt.0) THEN
294 ifield=ifice(nf)
295 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
296 got_var(ifield)=.true.
297 s(ng)%Vid(ifield)=var_id(i)
298 cycle
299 END IF
300 END IF
301 END DO
302 END DO
303!
304! Check if output variables are available in input NetCDF file.
305!
306 DO nf=1,nices
307 IF (isice(nf).gt.0) THEN
308 ifield=isice(nf)
309 IF (.not.got_var(ifield).and.varout(ifield,ng)) THEN
310 IF (master) WRITE (stdout,10) trim(vname(1,ifield)), &
311 & trim(ncname)
312 exit_flag=3
313 RETURN
314 END IF
315 END IF
316 END DO
317!
318 IF (.not.got_var(iduier).and.varout(iduier,ng)) THEN
319 IF (master) WRITE (stdout,10) trim(vname(1,iduier)), &
320 & trim(ncname)
321 exit_flag=3
322 RETURN
323 END IF
324 IF (.not.got_var(idvinr).and.varout(idvinr,ng)) THEN
325 IF (master) WRITE (stdout,10) trim(vname(1,idvinr)), &
326 & trim(ncname)
327 exit_flag=3
328 RETURN
329 END IF
330!
331 DO nf=1,nicef
332 IF (ifice(nf).gt.0) THEN
333 ifield=ifice(nf)
334 IF (.not.got_var(ifield).and.varout(ifield,ng)) THEN
335 IF (master) WRITE (stdout,10) trim(vname(1,ifield)), &
336 & trim(ncname)
337 exit_flag=3
338 RETURN
339 END IF
340 END IF
341 END DO
342 END IF query
343!
344 10 FORMAT (/,' ICE_DEF_NF90 - unable to find variable: ',a,2x, &
345 & ' in output NetCDF file: ',a)
346 RETURN
347 END SUBROUTINE ice_def_nf90
348!
349!***********************************************************************
350 SUBROUTINE ice_wrt_nf90 (ng, model, tile, &
351 & LBi, UBi, LBj, UBj, &
352 & VarOut, S)
353!***********************************************************************
354!
355 USE mod_netcdf
356!
357! Imported variable declarations.
358!
359 logical, intent(in) :: VarOut(NV,Ngrids)
360!
361 integer, intent(in) :: ng, model, tile
362 integer, intent(in) :: LBi, UBi, LBj, UBj
363!
364 TYPE(T_IO), intent(inout) :: S(Ngrids)
365!
366! Local variable declarations.
367!
368 logical :: LwrtVar
369!
370 integer :: ifield, ifld
371 integer :: gfactor, gtype, status
372!
373 real(dp) :: scale
374!
375 real(r8), pointer :: iceField(:,:)
376 real(r8), pointer :: iceMask(:,:)
377!
378 real(r8), allocatable :: Ur2d(:,:)
379 real(r8), allocatable :: Vr2d(:,:)
380!
381 character (len=*), parameter :: MyFile = &
382 & __FILE__//", ice_wrt_nf90"
383!
384 sourcefile=myfile
385!
386!-----------------------------------------------------------------------
387! Write out ice prognostic variables into specified output NetCDF file.
388!-----------------------------------------------------------------------
389!
390 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
391!
392! Set grid type factor to write full (gfactor=1) fields or water
393! points (gfactor=-1) fields only.
394!
395# if defined WRITE_WATER && defined MASKING
396 gfactor=-1
397# else
398 gfactor=1
399# endif
400!
401!-----------------------------------------------------------------------
402! write out ice ice model state prognostic variables.
403!-----------------------------------------------------------------------
404!
405 DO ifld=1,nices
406 IF (isice(ifld).gt.0) THEN
407 ifield=isice(ifld)
408 IF (s(ng)%ncid.eq.rst(ng)%ncid) THEN
409 lwrtvar=.true.
410 ELSE
411 lwrtvar=varout(ifield,ng)
412 END IF
413!
414 IF (lwrtvar) THEN
415 IF ((model.eq.inlm).and. &
416 ((s(ng)%ncid.eq.his(ng)%ncid).or. &
417 & (s(ng)%ncid.eq.qck(ng)%ncid).or. &
418 & (s(ng)%ncid.eq.rst(ng)%ncid))) THEN
419 icefield => ice(ng) % Si(lbi:ubi,lbj:ubj,iuout,ifld)
420# ifdef AVERAGES
421 ELSE IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
422 icefield => ice_savg(ifld,ng) % var(lbi:ubi,lbj:ubj)
423# endif
424 END IF
425!
426 SELECT CASE (ifld)
427 CASE (isuice)
428 gtype=gfactor*u2dvar
429# ifdef MASKING
430 icemask => grid(ng) % umask_full
431# endif
432 CASE (isvice)
433 gtype=gfactor*v2dvar
434# ifdef MASKING
435 icemask => grid(ng) % vmask_full
436# endif
437 CASE DEFAULT
438 gtype=gfactor*r2dvar
439# ifdef MASKING
440 icemask => grid(ng) % rmask_full
441# endif
442 END SELECT
443 scale=1.0_dp
444!
445 status=nf_fwrite2d(ng, model, s(ng)%ncid, ifield, &
446 & s(ng)%Vid(ifield), &
447 & s(ng)%Rindex, gtype, &
448 & lbi, ubi, lbj, ubj, scale, &
449# ifdef MASKING
450 & icemask, &
451# endif
452 & icefield)
453 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
454 IF (master) WRITE (stdout,10) trim(vname(1,ifield)), &
455 & s(ng)%Rindex, &
456 & trim(s(ng)%name)
457 exit_flag=3
458 ioerror=status
459 RETURN
460 END IF
461 END IF
462 END IF
463 END DO
464!
465! Write ice EASTward and NORTHward velocity component (m/s) at
466! RHO-points.
467!
468 IF (varout(iduier,ng).and.varout(idvinr,ng)) THEN
469 IF (.not.allocated(ur2d)) THEN
470 allocate (ur2d(lbi:ubi,lbj:ubj))
471 ur2d(lbi:ubi,lbj:ubj)=0.0_r8
472 END IF
473 IF (.not.allocated(vr2d)) THEN
474 allocate (vr2d(lbi:ubi,lbj:ubj))
475 vr2d(lbi:ubi,lbj:ubj)=0.0_r8
476 END IF
477 CALL uv_rotate2d (ng, tile, .false., .true., &
478 & lbi, ubi, lbj, ubj, &
479 & grid(ng) % CosAngler, &
480 & grid(ng) % SinAngler, &
481# ifdef MASKING
482 & grid(ng) % rmask_full, &
483# endif
484 & ice(ng) % Si(:,:,iuout,isuice), &
485 & ice(ng) % Si(:,:,iuout,isvice), &
486 & ur2d, vr2d)
487!
488 scale=1.0_dp
489 gtype=gfactor*r2dvar
490 status=nf_fwrite2d(ng, model, s(ng)%ncid, iduier, &
491 & s(ng)%Vid(iduier), &
492 & s(ng)%Rindex, gtype, &
493 & lbi, ubi, lbj, ubj, scale, &
494# ifdef MASKING
495 & grid(ng) % rmask_full, &
496# endif
497 & ur2d)
498 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
499 IF (master) WRITE (stdout,10) trim(vname(1,iduier)), &
500 & s(ng)%Rindex, trim(s(ng)%name)
501 exit_flag=3
502 ioerror=status
503 RETURN
504 END IF
505!
506 status=nf_fwrite2d(ng, model, s(ng)%ncid, idvinr, &
507 & s(ng)%Vid(idvinr), &
508 & s(ng)%Rindex, gtype, &
509 & lbi, ubi, lbj, ubj, scale, &
510# ifdef MASKING
511 & grid(ng) % rmask_full, &
512# endif
513 & vr2d)
514 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
515 IF (master) WRITE (stdout,10) trim(vname(1,idvinr)), &
516 & s(ng)%Rindex, trim(s(ng)%name)
517 exit_flag=3
518 ioerror=status
519 RETURN
520 END IF
521 deallocate (ur2d)
522 deallocate (vr2d)
523 END IF
524!
525!-----------------------------------------------------------------------
526! write out ice model internal variables.
527!-----------------------------------------------------------------------
528!
529 DO ifld=1,nicef
530 IF (ifice(ifld).gt.0) THEN
531 ifield=ifice(ifld)
532 IF (s(ng)%ncid.eq.rst(ng)%ncid) THEN
533 lwrtvar=.false.
534 ELSE
535 lwrtvar=varout(ifield,ng)
536 END IF
537!
538 IF (lwrtvar) THEN
539 IF ((model.eq.inlm).and. &
540 ((s(ng)%ncid.eq.his(ng)%ncid).or. &
541 & (s(ng)%ncid.eq.qck(ng)%ncid).or. &
542 & (s(ng)%ncid.eq.rst(ng)%ncid))) THEN
543 icefield => ice(ng) % Fi(lbi:ubi,lbj:ubj,ifld)
544# ifdef AVERAGES
545 ELSE IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
546 icefield => ice_favg(ifld,ng) % var(lbi:ubi,lbj:ubj)
547# endif
548 END IF
549!
550 scale=1.0_dp
551 gtype=gfactor*r2dvar
552 status=nf_fwrite2d(ng, model, s(ng)%ncid, ifield, &
553 & s(ng)%Vid(ifield), &
554 & s(ng)%Rindex, gtype, &
555 & lbi, ubi, lbj, ubj, scale, &
556# ifdef MASKING
557 & grid(ng) % rmask_full, &
558# endif
559 & icefield)
560 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
561 IF (master) WRITE (stdout,10) trim(vname(1,ifield)), &
562 & s(ng)%Rindex, &
563 & trim(s(ng)%name)
564 exit_flag=3
565 ioerror=status
566 RETURN
567 END IF
568 END IF
569 END IF
570 END DO
571!
572 10 FORMAT (/," ICE_WRT_NF90 - error while writing variable '",a, &
573 & "', time record = ",i0,/,11x,'into file: ',a)
574!
575 RETURN
576 END SUBROUTINE ice_wrt_nf90
577
578# ifdef STATIONS
579!
580!***********************************************************************
581 SUBROUTINE ice_def_station_nf90 (ng, model, ldef, VarOut, S, &
582 & pgrd)
583!***********************************************************************
584!
585 USE mod_netcdf
586!
587! Imported variable declarations.
588!
589 logical, intent(in) :: ldef, VarOut(NV,Ngrids)
590!
591 integer, intent(in) :: ng, model
592 integer, intent(in), optional :: pgrd(:)
593!
594 TYPE(T_IO), intent(inout) :: S(Ngrids)
595!
596! Local variable declarations.
597!
598 logical :: got_var(NV)
599!
600 integer, parameter :: Natt = 25
601
602 integer :: i, ifield, j, nf, status
603!
604 real(r8) :: Aval(6)
605!
606 character (len=120) :: Vinfo(Natt)
607 character (len=256) :: ncname
608!
609 character (len=*), parameter :: MyFile = &
610 & __FILE__//", ice_def_stations_nf90"
611!
612 sourcefile=myfile
613!
614!-----------------------------------------------------------------------
615! Define prognostic ice variables into specified output NetCDF.
616!-----------------------------------------------------------------------
617!
618 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
619 ncname=s(ng)%name
620!
621! Initialize local information variable arrays.
622!
623 DO i=1,natt
624 DO j=1,len(vinfo(1))
625 vinfo(i)(j:j)=' '
626 END DO
627 END DO
628 DO i=1,6
629 aval(i)=0.0_r8
630 END DO
631!
632!-----------------------------------------------------------------------
633! Define ice model state prognostic variables.
634!-----------------------------------------------------------------------
635!
636 define : IF (ldef) THEN
637
638 DO nf=1,nices
639 IF (isice(nf).gt.0) THEN
640 ifield=isice(nf)
641 IF (varout(ifield,ng)) THEN
642 vinfo( 1)=vname(1,ifield)
643 vinfo( 2)=vname(2,ifield)
644 vinfo( 3)=vname(3,ifield)
645 vinfo(14)=vname(4,ifield)
646 vinfo(16)=vname(1,idtime)
647!
648 status=def_var(ng, model, s(ng)%ncid, &
649 & s(ng)%Vid(ifield), &
650 & nf_fout, 2, pgrd, aval, vinfo, ncname, &
651 & setfillval = .true., &
652 & setparaccess = .true.)
653 IF (founderror(exit_flag, noerror, &
654 & __line__, myfile)) RETURN
655 END IF
656 END IF
657 END DO
658!
659!-----------------------------------------------------------------------
660! Define ice model internal variables.
661!-----------------------------------------------------------------------
662!
663 DO nf=1,nicef
664 IF (ifice(nf).gt.0) THEN
665 ifield=ifice(nf)
666 IF (varout(ifield,ng)) THEN
667 vinfo( 1)=vname(1,ifield)
668 vinfo( 2)=vname(2,ifield)
669 vinfo( 3)=vname(3,ifield)
670 vinfo(14)=vname(4,ifield)
671 vinfo(16)=vname(1,idtime)
672!
673 status=def_var(ng, model, s(ng)%ncid, &
674 & s(ng)%Vid(ifield), &
675 & nf_fout, 2, pgrd, aval, vinfo, ncname, &
676 & setfillval = .true., &
677 & setparaccess = .true.)
678 IF (founderror(exit_flag, noerror, &
679 & __line__, myfile)) RETURN
680 END IF
681 END IF
682 END DO
683 END IF define
684!
685!-----------------------------------------------------------------------
686! Otherwise, check existing output file and prepare for appending
687! data.
688!-----------------------------------------------------------------------
689!
690 query : IF (.not.ldef) THEN
691!
692! Initialize locallogical switches.
693!
694 DO i=1,nv
695 got_var(i)=.false.
696 END DO
697!
698! Scan variable list from input NetCDF and activate switches for
699! sea-ice variables. Get variable IDs.
700!
701 DO i=1,n_var
702 DO nf=1,nices
703 IF (isice(nf).gt.0) THEN
704 ifield=isice(nf)
705 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
706 got_var(ifield)=.true.
707 s(ng)%Vid(ifield)=var_id(i)
708 cycle
709 END IF
710 END IF
711 END DO
712!
713 DO nf=1,nicef
714 IF (ifice(nf).gt.0) THEN
715 ifield=ifice(nf)
716 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
717 got_var(ifield)=.true.
718 s(ng)%Vid(ifield)=var_id(i)
719 cycle
720 END IF
721 END IF
722 END DO
723 END DO
724!
725! Check if output variables are available in input NetCDF file.
726!
727 DO nf=1,nices
728 IF (isice(nf).gt.0) THEN
729 ifield=isice(nf)
730 IF (.not.got_var(ifield).and.varout(ifield,ng)) THEN
731 IF (master) WRITE (stdout,10) trim(vname(1,ifield)), &
732 & trim(ncname)
733 exit_flag=3
734 RETURN
735 END IF
736 END IF
737 END DO
738!
739 DO nf=1,nicef
740 IF (ifice(nf).gt.0) THEN
741 ifield=ifice(nf)
742 IF (.not.got_var(ifield).and.varout(ifield,ng)) THEN
743 IF (master) WRITE (stdout,10) trim(vname(1,ifield)), &
744 & trim(ncname)
745 exit_flag=3
746 RETURN
747 END IF
748 END IF
749 END DO
750 END IF query
751!
752 10 FORMAT (/,' ICE_DEF_STATION_NF90 - unable to find variable: ', &
753 & a,2x,' in output NetCDF file: ',a)
754 RETURN
755 END SUBROUTINE ice_def_station_nf90
756!
757!***********************************************************************
758 SUBROUTINE ice_wrt_station_nf90 (ng, model, tile, &
759 & LBi, UBi, LBj, UBj, &
760 & VarOut, S)
761!***********************************************************************
762!
763 USE mod_netcdf
764!
765! Imported variable declarations.
766!
767 logical, intent(in) :: VarOut(NV,Ngrids)
768!
769 integer, intent(in) :: ng, model, tile
770 integer, intent(in) :: LBi, UBi, LBj, UBj
771!
772 TYPE(T_IO), intent(inout) :: S(Ngrids)
773!
774! Local variable declarations.
775!
776 logical :: Cgrid
777!
778 integer :: i, ifield, ifld
779!
780 real(dp) :: scale
781!
782 real(r8), dimension(Nstation(ng)) :: Xpos, Ypos, psta
783!
784 character (len=*), parameter :: MyFile = &
785 & __FILE__//", ice_wrt_station_nf90"
786!
787 sourcefile=myfile
788!
789!-----------------------------------------------------------------------
790! Write out ice prognostic variables into specified output NetCDF file.
791!-----------------------------------------------------------------------
792!
793 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
794!
795! Set switch to extract station data at native C-grid position (TRUE)
796! or at RHO-points (FALSE).
797!
798# ifdef STATIONS_CGRID
799 cgrid=.true.
800# else
801 cgrid=.false.
802# endif
803!
804! Set positions for generic extraction routine.
805!
806 DO i=1,nstation(ng)
807 xpos(i)=scalars(ng)%SposX(i)
808 ypos(i)=scalars(ng)%SposY(i)
809 END DO
810!
811!-----------------------------------------------------------------------
812! write out ice ice model state prognostic variables.
813!-----------------------------------------------------------------------
814!
815 DO ifld=1,nices
816 IF (isice(ifld).gt.0) THEN
817 ifield=isice(ifld)
818 IF (varout(ifield,ng)) THEN
819 scale=1.0_dp
820 CALL extract_sta2d (ng, model, cgrid, ifield, r2dvar, &
821 & lbi, ubi, lbj, ubj, &
822 & scale, ice(ng)%Si(:,:,iuout,ifld), &
823 & nstation(ng), xpos, ypos, psta)
824!
825 CALL netcdf_put_fvar (ng, model, s(ng)%name, &
826 & trim(vname(1,ifield)), psta, &
827 & (/1,s(ng)%Rindex/), &
828 & (/nstation(ng),1/), &
829 & ncid = s(ng)%ncid, &
830 & varid = s(ng)%Vid(ifield))
831 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
832 END IF
833 END IF
834 END DO
835!
836!-----------------------------------------------------------------------
837! write out ice model internal variables.
838!-----------------------------------------------------------------------
839!
840 DO ifld=1,nicef
841 IF (ifice(ifld).gt.0) THEN
842 ifield=ifice(ifld)
843 IF (varout(ifield,ng)) THEN
844 scale=1.0_dp
845 CALL extract_sta2d (ng, model, cgrid, ifield, r2dvar, &
846 & lbi, ubi, lbj, ubj, &
847 & scale, ice(ng)%Fi(:,:,ifld), &
848 & nstation(ng), xpos, ypos, psta)
849!
850 CALL netcdf_put_fvar (ng, model, s(ng)%name, &
851 & trim(vname(1,ifield)), psta, &
852 & (/1,s(ng)%Rindex/), &
853 & (/nstation(ng),1/), &
854 & ncid = s(ng)%ncid, &
855 & varid = s(ng)%Vid(ifield))
856 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
857 END IF
858 END IF
859 END DO
860!
861 10 FORMAT (/,' ICE_WRT_STATION_NF90 - error while writing ', &
862 & "variable '",a,"', time record = ",i0,/,11x, &
863 & 'into file: ',a)
864!
865 RETURN
866 END SUBROUTINE ice_wrt_station_nf90
867# endif
868
869# if defined PIO_LIB && defined DISTRIBUTE
870!
871!***********************************************************************
872 SUBROUTINE ice_def_pio (ng, model, ldef, VarOut, S, &
873 & t2dgrd, u2dgrd, v2dgrd)
874!***********************************************************************
875!
877!
878! Imported variable declarations.
879!
880 logical, intent(in) :: ldef, VarOut(NV,Ngrids)
881!
882 integer, intent(in) :: ng, model
883 integer, intent(in), optional :: t2dgrd(:), u2dgrd(:), v2dgrd(:)
884!
885 TYPE(T_IO), intent(inout) :: S(Ngrids)
886!
887! Local variable declarations.
888!
889 logical :: LdefVar
890!
891 integer, parameter :: Natt = 25
892
893 integer :: i, ifield, j, nf, nvd3, status, vtype
894
895 integer :: icegrd(3)
896!
897 real(r8) :: Aval(6)
898!
899 character (len=120) :: Vinfo(Natt)
900 character (len=256) :: ncname
901!
902 character (len=*), parameter :: MyFile = &
903 & __FILE__//", ice_def_pio"
904!
905 sourcefile=myfile
906!
907!-----------------------------------------------------------------------
908! Define ice model state prognostic variables.
909!-----------------------------------------------------------------------
910!
911 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
912 ncname=s(ng)%name
913!!
914 define : IF (ldef) THEN
915!
916! Set number of dimensions for output variables.
917!
918# if defined WRITE_WATER && defined MASKING
919 nvd3=2
920# else
921 nvd3=3
922# endif
923!
924! Initialize local information variable arrays.
925!
926 DO i=1,natt
927 DO j=1,len(vinfo(1))
928 vinfo(i)(j:j)=' '
929 END DO
930 END DO
931 DO i=1,6
932 aval(i)=0.0_r8
933 END DO
934!
935! Define ice model state prognostic variables.
936!
937 DO nf=1,nices
938 IF (isice(nf).gt.0) THEN
939 ifield=isice(nf)
940 IF (s(ng)%ncid.eq.rst(ng)%ncid) THEN
941 ldefvar=.true.
942 vtype=pio_frst
943 ELSE
944 ldefvar=varout(ifield,ng)
945 vtype=pio_fout
946 END IF
947
948 IF (ldefvar) THEN
949 vinfo( 1)=vname(1,ifield)
950 vinfo( 2)=vname(2,ifield)
951 vinfo( 3)=vname(3,ifield)
952 vinfo(14)=vname(4,ifield)
953 vinfo(16)=vname(1,idtime)
954 vinfo(21)=vname(6,ifield)
955 vinfo(22)='coordinates'
956 aval(5)=real(iinfo(1,ifield,ng),r8)
957 s(ng)%pioVar(ifield)%dkind=pio_fout
958!
959 SELECT CASE (nf)
960 CASE (isuice)
961 icegrd(1:3)=u2dgrd(1:3)
962# if defined WRITE_WATER && defined MASKING
963 vinfo(20)='mask_u'
964# endif
965 s(ng)%pioVar(ifield)%gtype=u2dvar
966 CASE (isvice)
967 icegrd(1:3)=v2dgrd(1:3)
968# if defined WRITE_WATER && defined MASKING
969 vinfo(20)='mask_v'
970# endif
971 s(ng)%pioVar(ifield)%gtype=v2dvar
972 CASE DEFAULT
973 icegrd(1:3)=t2dgrd(1:3)
974# if defined WRITE_WATER && defined MASKING
975 vinfo(20)='mask_rho'
976# endif
977 s(ng)%pioVar(ifield)%gtype=r2dvar
978 END SELECT
979!
980 status=def_var(ng, model, s(ng)%pioFile, &
981 & s(ng)%pioVar(iduice)%vd, &
982 & vtype, nvd3, icegrd, &
983 & aval, vinfo, ncname)
984 IF (founderror(exit_flag, noerror, &
985 & __line__, myfile)) RETURN
986 END IF
987 END IF
988 END DO
989!
990! Define ice EASTward velocity at RHO-points.
991!
992 IF (varout(iduier,ng)) THEN
993 vinfo( 1)=vname(1,iduier)
994 vinfo( 2)=vname(2,iduier)
995 vinfo( 3)=vname(3,iduier)
996 vinfo(14)=vname(4,iduier)
997 vinfo(16)=vname(1,idtime)
998# if defined WRITE_WATER && defined MASKING
999 vinfo(20)='mask_rho'
1000# endif
1001 vinfo(21)=vname(6,iduier)
1002 vinfo(22)='coordinates'
1003 aval(5)=real(iinfo(1,iduier,ng),r8)
1004 s(ng)%pioVar(iduier)%dkind=pio_fout
1005 s(ng)%pioVar(iduier)%gtype=r2dvar
1006!
1007 status=def_var(ng, model, s(ng)%pioFile, &
1008 & s(ng)%pioVar(iduier)%vd, &
1009 & vtype, nvd3, t2dgrd, aval, vinfo, ncname)
1010 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1011 END IF
1012!
1013! Define ice NORTHward velocity at RHO-points.
1014!
1015 IF (varout(idvinr,ng)) THEN
1016 vinfo( 1)=vname(1,idvinr)
1017 vinfo( 2)=vname(2,idvinr)
1018 vinfo( 3)=vname(3,idvinr)
1019 vinfo(14)=vname(4,idvinr)
1020 vinfo(16)=vname(1,idtime)
1021# if defined WRITE_WATER && defined MASKING
1022 vinfo(20)='mask_rho'
1023# endif
1024 vinfo(21)=vname(6,idvinr)
1025 vinfo(22)='coordinates'
1026 aval(5)=real(iinfo(1,idvinr,ng),r8)
1027 s(ng)%pioVar(idvinr)%dkind=pio_fout
1028 s(ng)%pioVar(idvinr)%gtype=r2dvar
1029!
1030 status=def_var(ng, model, s(ng)%pioFile, &
1031 & s(ng)%pioVar(idvinr)%vd, &
1032 & vtype, nvd3, t2dgrd, aval, vinfo, ncname)
1033 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1034 END IF
1035!
1036!-----------------------------------------------------------------------
1037! Define ice model internal variables.
1038!-----------------------------------------------------------------------
1039!
1040 DO nf=1,nicef
1041 IF (ifice(nf).gt.0) THEN
1042 ifield=ifice(nf)
1043 IF (varout(ifield,ng)) THEN
1044 vinfo( 1)=vname(1,ifield)
1045 vinfo( 2)=vname(2,ifield)
1046 vinfo( 3)=vname(3,ifield)
1047 vinfo(14)=vname(4,ifield)
1048 vinfo(16)=vname(1,idtime)
1049# if defined WRITE_WATER && defined MASKING
1050 vinfo(20)='mask_rho'
1051# endif
1052 vinfo(21)=vname(6,ifield)
1053 vinfo(22)='coordinates'
1054 aval(5)=real(iinfo(1,ifield,ng),r8)
1055 s(ng)%pioVar(ifield)%dkind=pio_fout
1056 s(ng)%pioVar(ifield)%gtype=r2dvar
1057!
1058 SELECT CASE (nf)
1059 CASE (icwdiv)
1060 vinfo(11)='increase ice thickness'
1061 vinfo(12)='decrease ice concentration'
1062 CASE (icw_ai)
1063 vinfo(11)='freezing'
1064 vinfo(12)='melting'
1065 CASE (icw_ao, icw_io)
1066 vinfo(11)='melting'
1067 vinfo(12)='freezing'
1068 END SELECT
1069!
1070 status=def_var(ng, model, s(ng)%pioFile, &
1071 & s(ng)%pioVar(idaice)%vd, &
1072 & vtype, nvd3, t2dgrd, &
1073 & aval, vinfo, ncname)
1074 IF (founderror(exit_flag, noerror, &
1075 & __line__, myfile)) RETURN
1076 END IF
1077 END IF
1078 END DO
1079 END IF define
1080!
1081!-----------------------------------------------------------------------
1082! Otherwise, check existing output file and prepare for appending
1083! data.
1084!-----------------------------------------------------------------------
1085!
1086 query : IF (.not.ldef) THEN
1087!
1088! Initialize locallogical switches.
1089!
1090 DO i=1,nv
1091 got_var(i)=.false.
1092 END DO
1093!
1094! Scan variable list from input NetCDF and activate switches for
1095! sea-ice variables. Get variable IDs.
1096!
1097 IF (s(ng)%ncid.eq.rst(ng)%ncid) THEN
1098 vtype=pio_frst
1099 ELSE
1100 vtype=pio_fout
1101 END IF
1102!
1103 DO i=1,n_var
1104 DO nf=1,nices
1105 IF (isice(nf).gt.0) THEN
1106 ifield=isice(nf)
1107 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
1108 got_var(ifield)=.true.
1109 s(ng)%pioVar(ifield)%vd=var_desc(i)
1110 s(ng)%pioVar(ifield)%dkind=vtype
1111!
1112 SELECT CASE (nf)
1113 CASE (isuice)
1114 s(ng)%pioVar(ifield)%gtype=u2dvar
1115 CASE (isvice)
1116 s(ng)%pioVar(ifield)%gtype=v2dvar
1117 CASE DEFAULT
1118 s(ng)%pioVar(ifield)%gtype=r2dvar
1119 END SELECT
1120!
1121 cycle
1122 END IF
1123 END IF
1124 END DO
1125!
1126 IF (trim(var_name(i)).eq.trim(vname(1,iduier))) THEN
1127 got_var(iduier)=.true.
1128 s(ng)%pioVar(iduier)%vd=var_desc(i)
1129 s(ng)%pioVar(iduier)%dkind=pio_fout
1130 s(ng)%pioVar(iduier)%gtype=r2dvar
1131 ELSE IF (trim(var_name(i)).eq.trim(vname(1,idvinr))) THEN
1132 got_var(idvinr)=.true.
1133 s(ng)%pioVar(idvinr)%vd=var_desc(i)
1134 s(ng)%pioVar(idvinr)%dkind=vtype
1135 s(ng)%pioVar(idvinr)%gtype=r2dvar
1136 END IF
1137!
1138 DO nf=1,nicef
1139 IF (ifice(nf).gt.0) THEN
1140 ifield=ifice(nf)
1141 IF (trim(var_name(i)).eq.trim(vname(1,ifield))) THEN
1142 got_var(ifield)=.true.
1143 s(ng)%pioVar(ifield)%vd=var_desc(i)
1144 s(ng)%pioVar(ifield)%dkind=vtype
1145 s(ng)%pioVar(ifield)%gtype=r2dvar
1146 cycle
1147 END IF
1148 END IF
1149 END DO
1150 END DO
1151!
1152! Check if output variables are available in input NetCDF file.
1153!
1154 DO nf=1,nices
1155 IF (isice(nf).gt.0) THEN
1156 ifield=isice(nf)
1157 IF (.not.got_var(ifield).and.varout(ifield,ng)) THEN
1158 IF (master) WRITE (stdout,10) trim(vname(1,ifield)), &
1159 & trim(ncname)
1160 exit_flag=3
1161 RETURN
1162 END IF
1163 END IF
1164 END DO
1165!
1166 IF (.not.got_var(iduier).and.varout(iduier,ng)) THEN
1167 IF (master) WRITE (stdout,10) trim(vname(1,iduier)), &
1168 & trim(ncname)
1169 exit_flag=3
1170 RETURN
1171 END IF
1172 IF (.not.got_var(idvinr).and.varout(idvinr,ng)) THEN
1173 IF (master) WRITE (stdout,10) trim(vname(1,idvinr)), &
1174 & trim(ncname)
1175 exit_flag=3
1176 RETURN
1177 END IF
1178!
1179 DO nf=1,nicef
1180 IF (ifice(nf).gt.0) THEN
1181 ifield=ifice(nf)
1182 IF (.not.got_var(ifield).and.varout(ifield,ng)) THEN
1183 IF (master) WRITE (stdout,10) trim(vname(1,ifield)), &
1184 & trim(ncname)
1185 exit_flag=3
1186 RETURN
1187 END IF
1188 END IF
1189 END DO
1190 END IF query
1191!
1192 10 FORMAT (/,' ICE_DEF_PIO - unable to find variable: ',a,2x, &
1193 & ' in output NetCDF file: ',a)
1194!
1195 RETURN
1196 END SUBROUTINE ice_def_pio
1197!
1198!***********************************************************************
1199 SUBROUTINE ice_wrt_pio (ng, model, tile, &
1200 & LBi, UBi, LBj, UBj, &
1201 & VarOut, S)
1202!***********************************************************************
1203!
1204 USE mod_pio_netcdf
1205!
1206! Imported variable declarations.
1207!
1208 logical, intent(in) :: VarOut(NV,Ngrids)
1209!
1210 integer, intent(in) :: ng, model, tile
1211 integer, intent(in) :: LBi, UBi, LBj, UBj
1212!
1213 TYPE(T_IO), intent(inout) :: S(Ngrids)
1214!
1215! Local variable declarations.
1216!
1217 logical :: LwrtVar
1218!
1219 integer :: ifield, ifld
1220 integer :: gfactor, gtype, status
1221!
1222 real(dp) :: scale
1223!
1224 real(r8), pointer :: iceField(:,:)
1225 real(r8), pointer :: iceMask(:,:)
1226!
1227 real(r8), allocatable :: Ur2d(:,:)
1228 real(r8), allocatable :: Vr2d(:,:)
1229!
1230 character (len=*), parameter :: MyFile = &
1231 & __FILE__//", ice_wrt_pio"
1232!
1233 TYPE (IO_desc_t), pointer :: ioDesc
1234!
1235 sourcefile=myfile
1236!
1237!-----------------------------------------------------------------------
1238! Write out ice prognostic variables into specified output NetCDF file.
1239!-----------------------------------------------------------------------
1240!
1241 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
1242!
1243! Set grid type factor to write full (gfactor=1) fields or water
1244! points (gfactor=-1) fields only.
1245!
1246# if defined WRITE_WATER && defined MASKING
1247 gfactor=-1
1248# else
1249 gfactor=1
1250# endif
1251!
1252!-----------------------------------------------------------------------
1253! write out oce ice model state prognostic variables.
1254!-----------------------------------------------------------------------
1255!
1256 DO ifld=1,nices
1257 IF (isice(ifld).gt.0) THEN
1258 ifield=isice(ifld)
1259 IF (s(ng)%ncid.eq.rst(ng)%ncid) THEN
1260 lwrtvar=.true.
1261 ELSE
1262 lwrtvar=varout(ifield,ng)
1263 END IF
1264!
1265 IF (lwrtvar) THEN
1266 IF ((model.eq.inlm).and. &
1267 ((s(ng)%ncid.eq.his(ng)%ncid).or. &
1268 & (s(ng)%ncid.eq.qck(ng)%ncid).or. &
1269 & (s(ng)%ncid.eq.rst(ng)%ncid))) THEN
1270 icefield => ice(ng) % Si(lbi:ubi,lbj:ubj,iuout,ifld)
1271# ifdef AVERAGES
1272 ELSE IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1273 icefield => ice_savg(ifld,ng) % var(lbi:ubi,lbj:ubj)
1274# endif
1275 END IF
1276!
1277 SELECT CASE (ifld)
1278 CASE (isuice)
1279 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1280 iodesc => iodesc_dp_u2dvar(ng)
1281 ELSE
1282 iodesc => iodesc_sp_u2dvar(ng)
1283 END IF
1284# ifdef MASKING
1285 icemask => grid(ng) % umask_full
1286# endif
1287 CASE (isvice)
1288 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1289 iodesc => iodesc_dp_v2dvar(ng)
1290 ELSE
1291 iodesc => iodesc_sp_v2dvar(ng)
1292 END IF
1293# ifdef MASKING
1294 icemask => grid(ng) % vmask_full
1295# endif
1296 CASE DEFAULT
1297 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1298 iodesc => iodesc_dp_r2dvar(ng)
1299 ELSE
1300 iodesc => iodesc_sp_r2dvar(ng)
1301 END IF
1302# ifdef MASKING
1303 icemask => grid(ng) % rmask_full
1304# endif
1305 END SELECT
1306 scale=1.0_dp
1307!
1308 status=nf_fwrite2d(ng, model, s(ng)%pioFile, ifield, &
1309 & s(ng)%pioVar(ifield), s(ng)%Rindex, &
1310 & iodesc, &
1311 & lbi, ubi, lbj, ubj, scale, &
1312# ifdef MASKING
1313 & icemask, &
1314# endif
1315 & icefield)
1316 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
1317 IF (master) WRITE (stdout,10) trim(vname(1,ifield)), &
1318 & s(ng)%Rindex, &
1319 & trim(s(ng)%name)
1320 exit_flag=3
1321 ioerror=status
1322 RETURN
1323 END IF
1324 END IF
1325 END IF
1326 END DO
1327!
1328! Write ice EASTward and NORTHward velocity component (m/s) at
1329! RHO-points.
1330!
1331 IF (varout(iduier,ng).and.varout(idvinr,ng)) THEN
1332 IF (.not.allocated(ur2d)) THEN
1333 allocate (ur2d(lbi:ubi,lbj:ubj))
1334 ur2d(lbi:ubi,lbj:ubj)=0.0_r8
1335 END IF
1336 IF (.not.allocated(vr2d)) THEN
1337 allocate (vr2d(lbi:ubi,lbj:ubj))
1338 vr2d(lbi:ubi,lbj:ubj)=0.0_r8
1339 END IF
1340 CALL uv_rotate2d (ng, tile, .false., .true., &
1341 & lbi, ubi, lbj, ubj, &
1342 & grid(ng) % CosAngler, &
1343 & grid(ng) % SinAngler, &
1344# ifdef MASKING
1345 & grid(ng) % rmask_full, &
1346# endif
1347 & ice(ng) % Si(:,:,iuout,isuice), &
1348 & ice(ng) % Si(:,:,iuout,isvice), &
1349 & ur2d, vr2d)
1350!
1351 scale=1.0_dp
1352 IF (s(ng)%pioVar(iduier)%dkind.eq.pio_double) THEN
1353 iodesc => iodesc_dp_r2dvar(ng)
1354 ELSE
1355 iodesc => iodesc_sp_r2dvar(ng)
1356 END IF
1357 status=nf_fwrite2d(ng, model, s(ng)%pioFile, iduier, &
1358 & s(ng)%pioVAR(iduier), s(ng)%Rindex, &
1359 & iodesc, &
1360 & lbi, ubi, lbj, ubj, scale, &
1361# ifdef MASKING
1362 & grid(ng) % rmask_full, &
1363# endif
1364 & ur2d)
1365 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
1366 IF (master) WRITE (stdout,10) trim(vname(1,iduier)), &
1367 & s(ng)%Rindex, trim(s(ng)%name)
1368 exit_flag=3
1369 ioerror=status
1370 RETURN
1371 END IF
1372!
1373 IF (s(ng)%pioVar(idvinr)%dkind.eq.pio_double) THEN
1374 iodesc => iodesc_dp_r2dvar(ng)
1375 ELSE
1376 iodesc => iodesc_sp_r2dvar(ng)
1377 END IF
1378 status=nf_fwrite2d(ng, model, s(ng)%pioFile, idvinr, &
1379 & s(ng)%pioVar(idvinr), s(ng)%Rindex, &
1380 & iodesc, &
1381 & lbi, ubi, lbj, ubj, scale, &
1382# ifdef MASKING
1383 & grid(ng) % rmask_full, &
1384# endif
1385 & vr2d)
1386 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
1387 IF (master) WRITE (stdout,10) trim(vname(1,idvinr)), &
1388 & s(ng)%Rindex, trim(s(ng)%name)
1389 exit_flag=3
1390 ioerror=status
1391 RETURN
1392 END IF
1393 deallocate (ur2d)
1394 deallocate (vr2d)
1395 END IF
1396!
1397!-----------------------------------------------------------------------
1398! write out ice model internal variables.
1399!-----------------------------------------------------------------------
1400!
1401 DO ifld=1,nicef
1402 IF (ifice(ifld).gt.0) THEN
1403 ifield=ifice(ifld)
1404 IF (s(ng)%ncid.eq.rst(ng)%ncid) THEN
1405 lwrtvar=.false.
1406 ELSE
1407 lwrtvar=varout(ifield,ng)
1408 END IF
1409!
1410 IF (lwrtvar) THEN
1411 IF ((model.eq.inlm).and. &
1412 ((s(ng)%ncid.eq.his(ng)%ncid).or. &
1413 & (s(ng)%ncid.eq.qck(ng)%ncid).or. &
1414 & (s(ng)%ncid.eq.rst(ng)%ncid))) THEN
1415 icefield => ice(ng) % Fi(lbi:ubi,lbj:ubj,ifld)
1416# ifdef AVERAGES
1417 ELSE IF (s(ng)%ncid.eq.avg(ng)%ncid) THEN
1418 icefield => ice_favg(ifld,ng) % var(lbi:ubi,lbj:ubj)
1419# endif
1420 END IF
1421!
1422 scale=1.0_dp
1423 IF (s(ng)%pioVar(ifield)%dkind.eq.pio_double) THEN
1424 iodesc => iodesc_dp_r2dvar(ng)
1425 ELSE
1426 iodesc => iodesc_sp_r2dvar(ng)
1427 END IF
1428 status=nf_fwrite2d(ng, model, s(ng)%pioFile, ifield, &
1429 & s(ng)%pioVar(ifield), s(ng)%Rindex, &
1430 & iodesc, &
1431 & lbi, ubi, lbj, ubj, scale, &
1432# ifdef MASKING
1433 & grid(ng) % rmask_full, &
1434# endif
1435 & icefield)
1436 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
1437 IF (master) WRITE (stdout,10) trim(vname(1,ifield)), &
1438 & s(ng)%Rindex, &
1439 & trim(s(ng)%name)
1440 exit_flag=3
1441 ioerror=status
1442 RETURN
1443 END IF
1444 END IF
1445 END IF
1446 END DO
1447!
1448 10 FORMAT (/," ICE_WRT - error while writing variable '",a, &
1449 & "', time record = ",i0,/,11x,'into file: ',a)
1450!
1451 RETURN
1452 END SUBROUTINE ice_wrt_pio
1453# endif
1454!
1455#endif
1456 END MODULE ice_output_mod
subroutine extract_sta2d(ng, model, cgrid, ifield, gtype, lbi, ubi, lbj, ubj, ascl, a, npos, xpos, ypos, apos)
Definition extract_sta.F:65
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer, parameter isvice
Definition ice_mod.h:147
integer idvinr
Definition ice_mod.h:119
integer, parameter icwdiv
Definition ice_mod.h:189
integer, dimension(nicef) ifice
Definition ice_mod.h:169
integer, parameter icw_ao
Definition ice_mod.h:191
type(t_ice), dimension(:), allocatable ice
Definition ice_mod.h:283
integer, parameter nicef
Definition ice_mod.h:167
integer, dimension(nices) isice
Definition ice_mod.h:135
integer idaice
Definition ice_mod.h:96
integer, parameter isuice
Definition ice_mod.h:146
type(t_ice_avg), dimension(:,:), allocatable ice_savg
Definition ice_mod.h:320
integer iduice
Definition ice_mod.h:114
type(t_ice_avg), dimension(:,:), allocatable ice_favg
Definition ice_mod.h:319
integer, parameter nices
Definition ice_mod.h:130
integer, parameter icw_io
Definition ice_mod.h:193
integer, parameter icw_ai
Definition ice_mod.h:190
integer iduier
Definition ice_mod.h:116
character(len=maxlen), dimension(6, 0:nv) vname
integer idtime
integer, dimension(:,:,:), allocatable iinfo
integer, parameter nf_fout
Definition mod_netcdf.F:188
character(len=100), dimension(mvars) var_name
Definition mod_netcdf.F:169
integer, dimension(mvars) var_id
Definition mod_netcdf.F:160
integer n_var
Definition mod_netcdf.F:152
integer, parameter nf_frst
Definition mod_netcdf.F:193
logical master
integer, parameter inlm
Definition mod_param.F:662
integer, dimension(:), allocatable nstation
Definition mod_param.F:557
integer, parameter u2dvar
Definition mod_param.F:718
integer, parameter r2dvar
Definition mod_param.F:717
integer, parameter v2dvar
Definition mod_param.F:719
type(io_desc_t), dimension(:), pointer iodesc_sp_u2dvar
integer, parameter pio_fout
type(var_desc_t), dimension(:), pointer var_desc
type(io_desc_t), dimension(:), pointer iodesc_sp_r2dvar
integer, parameter pio_frst
type(io_desc_t), dimension(:), pointer iodesc_dp_u2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_v2dvar
type(io_desc_t), dimension(:), pointer iodesc_dp_r2dvar
type(io_desc_t), dimension(:), pointer iodesc_sp_v2dvar
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine, public uv_rotate2d(ng, tile, add, lboundary, lbi, ubi, lbj, ubj, cosangler, sinangler, rmask_full, uinp, vinp, uout, vout)
Definition uv_rotate.F:35