ocnice_prep 1.14.0
Loading...
Searching...
No Matches
utils_mod.F90
Go to the documentation of this file.
1
7module utils_mod
8
9 use netcdf
10 use arrays_mod, only : eta
11 use init_mod, only : debug, logunit, vardefs, fsrc
12
13 implicit none
14
15 private
16
17 interface getfield
18 module procedure getfield2d
19 module procedure getfield3d
20 end interface getfield
21
22 interface packarrays
23 module procedure packarrays2d
24 module procedure packarrays3d
25 end interface packarrays
26
27 interface remap
28 module procedure remap1d
29 module procedure remap2d
30 module procedure remap3d
31 end interface remap
32
33 interface getvecpair
34 module procedure getvecpair2d
35 module procedure getvecpair3d
36 end interface getvecpair
37
38 interface dumpnc
39 module procedure dumpnc1d
40 module procedure dumpnc2d
41 module procedure dumpnc3d
42 module procedure dumpnc3dk
43 end interface dumpnc
44
45 public :: getfield
46 public :: packarrays
47 public :: remap
48 public :: dumpnc
49 public :: nf90_err
50 public :: zero_out_land_ice
51 public :: zero_out_phantom_ice
52
53contains
66 subroutine packarrays2d(filesrc, wgtsdir, cosrot, sinrot, vars, dims, nflds, fields)
67
68 character(len=*), intent(in) :: filesrc,wgtsdir
69 real(kind=8), intent(in) :: cosrot(:),sinrot(:)
70 type(vardefs), intent(in) :: vars(:)
71 integer, intent(in) :: dims(:)
72 integer, intent(in) :: nflds
73 real(kind=8), intent(out) :: fields(:,:)
74
75 ! local variables
76 integer :: n, nn
77 real(kind=8), allocatable :: vecpair(:,:)
78 character(len=20) :: subname = 'packarrays2d'
79 !----------------------------------------------------------------------------
80
81 if (debug)write(logunit,'(a)')'enter '//trim(subname)
82
83 fields=0.0
84 ! obtain vector pairs
85 do n = 1,nflds
86 if (trim(vars(n)%var_grid) == 'Cu' .or. trim(vars(n)%var_grid) == 'Bu_x') then
87 allocate(vecpair(dims(1)*dims(2),2)); vecpair = 0.0
88 call getvecpair(trim(filesrc), trim(wgtsdir), cosrot, sinrot, &
89 trim(vars(n)%var_name), trim(vars(n)%var_grid(1:2)), &
90 trim(vars(n)%var_pair), trim(vars(n)%var_pair_grid(1:2)), &
91 dims=(/dims(1),dims(2)/), vecpair=vecpair)
92 end if
93 end do
94
95 ! create packed array
96 nn = 0
97 do n = 1,nflds
98 if (len_trim(vars(n)%var_pair) == 0) then
99 nn = nn + 1
100 call getfield(trim(filesrc), trim(vars(n)%var_name), dims=(/dims(1),dims(2)/), &
101 field=fields(nn,:))
102 else ! fill with vector pairs
103 nn = nn+1
104 ! ocn vectors
105 if (trim(vars(n)%var_grid) == 'Cu')fields(nn,:) = vecpair(:,1)
106 if (trim(vars(n)%var_grid) == 'Cv')fields(nn,:) = vecpair(:,2)
107 ! ice vectors
108 if (trim(vars(n)%var_grid) == 'Bu_x')fields(nn,:) = vecpair(:,1)
109 if (trim(vars(n)%var_grid) == 'Bu_y')fields(nn,:) = vecpair(:,2)
110 end if
111 end do
112
113 if (debug)write(logunit,'(a)')'exit '//trim(subname)
114 end subroutine packarrays2d
115
128 subroutine packarrays3d(filesrc, wgtsdir, cosrot, sinrot, vars, dims, nflds, fields)
129
130 character(len=*), intent(in) :: filesrc,wgtsdir
131 real(kind=8), intent(in) :: cosrot(:),sinrot(:)
132 type(vardefs), intent(in) :: vars(:)
133 integer, intent(in) :: dims(:)
134 integer, intent(in) :: nflds
135 real(kind=8), intent(out) :: fields(:,:,:)
136
137 ! local variables
138 integer :: n, nn
139 real(kind=8), allocatable :: vecpair(:,:,:)
140 character(len=20) :: subname = 'packarrays3d'
141 !----------------------------------------------------------------------------
142
143 if (debug)write(logunit,'(a)')'enter '//trim(subname)
144
145 fields=0.0
146 ! obtain vector pairs
147 do n = 1,nflds
148 if (trim(vars(n)%var_grid) == 'Cu') then
149 allocate(vecpair(dims(3),dims(1)*dims(2),2)); vecpair = 0.0
150 call getvecpair(trim(filesrc), trim(wgtsdir), cosrot, sinrot, &
151 trim(vars(n)%var_name), trim(vars(n)%var_grid), &
152 trim(vars(n)%var_pair), trim(vars(n)%var_pair_grid), &
153 dims=(/dims(1),dims(2),dims(3)/), vecpair=vecpair)
154 end if
155 end do
156
157 ! create packed array
158 nn = 0
159 do n = 1,nflds
160 if (len_trim(vars(n)%var_pair) == 0) then
161 nn = nn + 1
162 if (trim(vars(n)%var_name) .eq. 'eta') then
163 fields(nn,:,:) = eta(:,:)
164 else
165 call getfield(trim(filesrc), trim(vars(n)%var_name), dims=(/dims(1),dims(2),dims(3)/), &
166 field=fields(nn,:,:))
167 end if
168 else ! fill with vector pairs
169 nn = nn+1
170 if (trim(vars(n)%var_grid) == 'Cu')fields(nn,:,:) = vecpair(:,:,1)
171 if (trim(vars(n)%var_grid) == 'Cv')fields(nn,:,:) = vecpair(:,:,2)
172 end if
173 end do
174
175 if (debug)write(logunit,'(a)')'exit '//trim(subname)
176 end subroutine packarrays3d
177
192 subroutine getvecpair2d(fname, wdir, cosrot, sinrot, vname1, vgrid1, vname2, vgrid2, dims, vecpair)
193
194 character(len=*), intent(in) :: fname
195 character(len=*), intent(in) :: wdir
196 real(kind=8), intent(in) :: cosrot(:), sinrot(:)
197 character(len=*), intent(in) :: vname1, vgrid1, vname2, vgrid2
198 integer, intent(in) :: dims(:)
199 real(kind=8), intent(out) :: vecpair(:,:)
200
201 ! local variables
202 integer :: ii
203 real(kind=8), dimension(dims(1)*dims(2)) :: urot, vrot
204 character(len=240) :: wgtsfile
205 character(len=20) :: subname = 'getvecpair2d'
206 !----------------------------------------------------------------------------
207
208 if (debug)write(logunit,'(a)')'enter '//trim(subname)
209
210 wgtsfile = trim(wdir)//'tripole.'//trim(fsrc)//'.'//vgrid1//'.to.Ct.bilinear.nc'
211 call getfield(fname, vname1, dims=dims, field=vecpair(:,1), wgts=trim(wgtsfile))
212 if (debug)write(logunit,'(a)')'wgtsfile for 2d vector '//trim(vname1)//' '//trim(wgtsfile)
213 wgtsfile = trim(wdir)//'tripole.'//trim(fsrc)//'.'//vgrid2//'.to.Ct.bilinear.nc'
214 call getfield(fname, vname2, dims=dims, field=vecpair(:,2), wgts=trim(wgtsfile))
215 if (debug)write(logunit,'(a)')'wgtsfile for 2d vector '//trim(vname2)//' '//trim(wgtsfile)
216
217 urot = 0.0; vrot = 0.0
218 do ii = 1,dims(1)*dims(2)
219 urot(ii) = vecpair(ii,1)*cosrot(ii) + vecpair(ii,2)*sinrot(ii)
220 vrot(ii) = vecpair(ii,2)*cosrot(ii) - vecpair(ii,1)*sinrot(ii)
221 end do
222 vecpair(:,1) = urot(:)
223 vecpair(:,2) = vrot(:)
224
225 if (debug) write(logunit,'(a)')'exit '//trim(subname)
226 end subroutine getvecpair2d
227
242 subroutine getvecpair3d(fname, wdir, cosrot, sinrot, vname1, vgrid1, vname2, vgrid2, dims, vecpair)
243
244 character(len=*), intent(in) :: fname
245 character(len=*), intent(in) :: wdir
246 real(kind=8), intent(in) :: cosrot(:), sinrot(:)
247 character(len=*), intent(in) :: vname1, vgrid1, vname2, vgrid2
248 integer, intent(in) :: dims(:)
249 real(kind=8), intent(out) :: vecpair(:,:,:)
250
251 ! local variables
252 integer :: ii,k
253 real(kind=8), dimension(dims(1)*dims(2)) :: urot, vrot
254 character(len=240) :: wgtsfile
255 character(len=20) :: subname = 'getvecpair3d'
256 !----------------------------------------------------------------------------
257
258 if (debug)write(logunit,'(a)')'enter '//trim(subname)
259
260 wgtsfile = trim(wdir)//'tripole.'//trim(fsrc)//'.'//vgrid1//'.to.Ct.bilinear.nc'
261 call getfield(fname, vname1, dims=dims, field=vecpair(:,:,1), wgts=trim(wgtsfile))
262 wgtsfile = trim(wdir)//'tripole.'//trim(fsrc)//'.'//vgrid2//'.to.Ct.bilinear.nc'
263 call getfield(fname, vname2, dims=dims, field=vecpair(:,:,2), wgts=trim(wgtsfile))
264
265 do k = 1,dims(3)
266 urot = 0.0; vrot = 0.0
267 do ii = 1,dims(1)*dims(2)
268 urot(ii)= vecpair(k,ii,1)*cosrot(ii) + vecpair(k,ii,2)*sinrot(ii)
269 vrot(ii)= vecpair(k,ii,2)*cosrot(ii) - vecpair(k,ii,1)*sinrot(ii)
270 end do
271 vecpair(k,:,1) = urot(:)
272 vecpair(k,:,2) = vrot(:)
273 end do
274
275 if (debug) write(logunit,'(a)')'exit '//trim(subname)
276 end subroutine getvecpair3d
277
287 subroutine getfield2d(fname, vname, dims, field, wgts)
288
289 character(len=*), intent(in) :: fname, vname
290 integer, intent(in) :: dims(:)
291 real(kind=8), intent(out) :: field(:)
292 character(len=*), optional, intent(in) :: wgts
293
294 ! local variable
295 integer :: ncid, varid
296 real(kind=8), allocatable :: a2d(:,:)
297 real(kind=8), allocatable :: atmp(:)
298 character(len=20) :: subname = 'getfield2d'
299 !----------------------------------------------------------------------------
300
301 if (debug)write(logunit,'(a)')'enter '//trim(subname)//' variable '//vname
302
303 allocate(a2d(dims(1),dims(2))); a2d = 0.0
304 allocate(atmp(dims(1)*dims(2))); atmp = 0.0
305
306 call nf90_err(nf90_open(fname, nf90_nowrite, ncid), 'nf90_open: '//fname)
307 call nf90_err(nf90_inq_varid(ncid, vname, varid), &
308 'get variable ID: '//vname)
309 call nf90_err(nf90_get_var(ncid, varid, a2d), &
310 'get variable: '//vname)
311 call nf90_err(nf90_close(ncid), 'close: '//fname)
312
313 atmp(:) = reshape(a2d, (/dims(1)*dims(2)/))
314 if(present(wgts)) then
315 call remap(trim(wgts), src_field=atmp, dst_field=field)
316 else
317 field = atmp
318 end if
319
320 if (debug) write(logunit,'(a)')'exit '//trim(subname)//' variable '//vname
321 end subroutine getfield2d
322
332 subroutine getfield3d(fname, vname, dims, field, wgts)
333
334 character(len=*), intent(in) :: fname, vname
335 integer, intent(in) :: dims(:)
336 real(kind=8), intent(out) :: field(:,:)
337 character(len=*), optional, intent(in) :: wgts
338
339 ! local variable
340 integer :: k, ncid, varid
341 real(kind=8), allocatable :: a3d(:,:,:)
342 real(kind=8), allocatable :: atmp(:,:)
343 character(len=20) :: subname = 'getfield3d'
344 !----------------------------------------------------------------------------
345
346 if (debug)write(logunit,'(a)')'enter '//trim(subname)//' variable '//vname
347
348 allocate(a3d(dims(1),dims(2),dims(3))); a3d = 0.0
349 allocate(atmp(dims(3),dims(1)*dims(2))); atmp = 0.0
350
351 call nf90_err(nf90_open(fname, nf90_nowrite, ncid), 'nf90_open: '//fname)
352 call nf90_err(nf90_inq_varid(ncid, vname, varid), &
353 'get variable ID: '//vname)
354 call nf90_err(nf90_get_var(ncid, varid, a3d), &
355 'get variable: '//vname)
356 call nf90_err(nf90_close(ncid), 'close: '//fname)
357
358 do k = 1,dims(3)
359 atmp(k,:) = reshape(a3d(1:dims(1),1:dims(2),k), (/dims(1)*dims(2)/))
360 end do
361 if(present(wgts)) then
362 call remap(trim(wgts), dim2=dims(3), src_field=atmp, dst_field=field)
363 else
364 field = atmp
365 end if
366
367 if (debug) write(logunit,'(a)')'exit '//trim(subname)//' variable '//vname
368 end subroutine getfield3d
369
377 subroutine remap1d(fname, src_field, dst_field)
378
379 character(len=*), intent(in) :: fname
380 real(kind=8), intent(in) :: src_field(:)
381 real(kind=8), intent(out) :: dst_field(:)
382
383 ! local variables
384 integer :: ncid, id
385 integer :: i,ii,jj
386 integer :: n_a, n_b, n_s
387 integer(kind=4), allocatable, dimension(:) :: col, row
388 real(kind=8), allocatable, dimension(:) :: s
389 character(len=20) :: subname = 'remap1d'
390 !----------------------------------------------------------------------------
391
392 if (debug)write(logunit,'(a)')'enter '//trim(subname)
393
394 ! retrieve the weights
395 call nf90_err(nf90_open(trim(fname), nf90_nowrite, ncid), 'open: '//fname)
396 call nf90_err(nf90_inq_dimid(ncid, 'n_s', id), 'get dimension Id: n_s')
397 call nf90_err(nf90_inquire_dimension(ncid, id, len=n_s), 'get dimension: n_s' )
398 call nf90_err(nf90_inq_dimid(ncid, 'n_a', id), 'get dimension Id: n_a')
399 call nf90_err(nf90_inquire_dimension(ncid, id, len=n_a), 'get dimension: n_a' )
400 call nf90_err(nf90_inq_dimid(ncid, 'n_b', id), 'get dimension Id: n_b')
401 call nf90_err(nf90_inquire_dimension(ncid, id, len=n_b), 'get dimension: n_b' )
402
403 allocate(col(1:n_s)); col = 0
404 allocate(row(1:n_s)); row = 0
405 allocate( s(1:n_s)); s = 0.0
406
407 call nf90_err(nf90_inq_varid(ncid, 'col', id),'get variable Id: col')
408 call nf90_err(nf90_get_var(ncid, id, col),'get variable: col')
409 call nf90_err(nf90_inq_varid(ncid, 'row', id),'get variable Id: row')
410 call nf90_err(nf90_get_var(ncid, id, row),'get variable: row')
411 call nf90_err(nf90_inq_varid(ncid, 'S', id),'get variable Id: S')
412 call nf90_err(nf90_get_var(ncid, id, s),'get variable: S')
413 call nf90_err(nf90_close(ncid), 'close: '//fname)
414
415 dst_field = 0.0
416 do i = 1,n_s
417 ii = row(i); jj = col(i)
418 dst_field(ii) = dst_field(ii) + s(i)*src_field(jj)
419 enddo
420
421 if (debug) write(logunit,'(a)')'exit '//trim(subname)
422 end subroutine remap1d
423
432 subroutine remap2d(fname, dim2, src_field, dst_field)
433
434 character(len=*), intent(in) :: fname
435 integer, intent(in) :: dim2
436 real(kind=8), intent(in) :: src_field(:,:)
437 real(kind=8), intent(out) :: dst_field(:,:)
438
439 ! local variables
440 integer :: ncid, id
441 integer :: i,ii,jj
442 integer :: n_a, n_b, n_s
443 integer(kind=4), allocatable, dimension(:) :: col, row
444 real(kind=8), allocatable, dimension(:) :: s
445 character(len=20) :: subname = 'remap2d'
446 !----------------------------------------------------------------------------
447
448 if (debug)write(logunit,'(a)')'enter '//trim(subname)//' weights = '//trim(fname)
449
450 ! retrieve the weights
451 call nf90_err(nf90_open(trim(fname), nf90_nowrite, ncid), 'open: '//fname)
452 call nf90_err(nf90_inq_dimid(ncid, 'n_s', id), 'get dimension Id: n_s')
453 call nf90_err(nf90_inquire_dimension(ncid, id, len=n_s), 'get dimension: n_s')
454 call nf90_err(nf90_inq_dimid(ncid, 'n_a', id), 'get dimension Id: n_a')
455 call nf90_err(nf90_inquire_dimension(ncid, id, len=n_a), 'get dimension: n_a')
456 call nf90_err(nf90_inq_dimid(ncid, 'n_b', id), 'get dimension Id: n_b')
457 call nf90_err(nf90_inquire_dimension(ncid, id, len=n_b), 'get dimension: n_b')
458
459 allocate(col(1:n_s)); col = 0
460 allocate(row(1:n_s)); row = 0
461 allocate( s(1:n_s)); s = 0.0
462
463 call nf90_err(nf90_inq_varid(ncid, 'col', id),'get variable Id: col')
464 call nf90_err(nf90_get_var(ncid, id, col),'get variable: col')
465 call nf90_err(nf90_inq_varid(ncid, 'row', id),'get variable Id: row')
466 call nf90_err(nf90_get_var(ncid, id, row),'get variable: row')
467 call nf90_err(nf90_inq_varid(ncid, 'S', id),'get variable Id: S')
468 call nf90_err(nf90_get_var(ncid, id, s),'get variable: S')
469 call nf90_err(nf90_close(ncid), 'close: '//fname)
470
471 dst_field = 0.0
472 do i = 1,n_s
473 ii = row(i); jj = col(i)
474 dst_field(:,ii) = dst_field(:,ii) + s(i)*src_field(:,jj)
475 enddo
476
477 if (debug) write(logunit,'(a)')'exit '//trim(subname)
478 end subroutine remap2d
479
489 subroutine remap3d(fname, nk, nflds, src_field, dst_field)
490
491 character(len=*), intent(in) :: fname
492 integer, intent(in) :: nk, nflds
493 real(kind=8), intent(in) :: src_field(:,:,:)
494 real(kind=8), intent(out) :: dst_field(:,:,:)
495
496 ! local variables
497 integer :: ncid, id
498 integer :: i,ii,jj
499 integer :: n_a, n_b, n_s
500 integer(kind=4), allocatable, dimension(:) :: col, row
501 real(kind=8), allocatable, dimension(:) :: s
502 character(len=20) :: subname = 'remap3d'
503 !----------------------------------------------------------------------------
504
505 if (debug)write(logunit,'(a)')'enter '//trim(subname)//' weights = '//trim(fname)
506
507 ! retrieve the weights
508 call nf90_err(nf90_open(trim(fname), nf90_nowrite, ncid), 'open: '//fname)
509 call nf90_err(nf90_inq_dimid(ncid, 'n_s', id), 'get dimension Id: n_s')
510 call nf90_err(nf90_inquire_dimension(ncid, id, len=n_s), 'get dimension: n_s')
511 call nf90_err(nf90_inq_dimid(ncid, 'n_a', id), 'get dimension Id: n_a')
512 call nf90_err(nf90_inquire_dimension(ncid, id, len=n_a), 'get dimension: n_a')
513 call nf90_err(nf90_inq_dimid(ncid, 'n_b', id), 'get dimension Id: n_b')
514 call nf90_err(nf90_inquire_dimension(ncid, id, len=n_b), 'get dimension: n_b')
515
516 allocate(col(1:n_s)); col = 0
517 allocate(row(1:n_s)); row = 0
518 allocate( s(1:n_s)); s = 0.0
519
520 call nf90_err(nf90_inq_varid(ncid, 'col', id),'get variable Id: col')
521 call nf90_err(nf90_get_var(ncid, id, col),'get variable: col')
522 call nf90_err(nf90_inq_varid(ncid, 'row', id),'get variable Id: row')
523 call nf90_err(nf90_get_var(ncid, id, row),'get variable: row')
524 call nf90_err(nf90_inq_varid(ncid, 'S', id),'get variable Id: S')
525 call nf90_err(nf90_get_var(ncid, id, s),'get variable: S')
526 call nf90_err(nf90_close(ncid), 'close: '//fname)
527
528 dst_field = 0.0
529 do i = 1,n_s
530 ii = row(i); jj = col(i)
531 dst_field(:,:,ii) = dst_field(:,:,ii) + s(i)*src_field(:,:,jj)
532 enddo
533
534 if (debug) write(logunit,'(a)')'exit '//trim(subname)
535 end subroutine remap3d
536
546 subroutine dumpnc2d(fname, vname, dims, nflds, field)
547
548 character(len=*), intent(in) :: fname, vname
549 integer, intent(in) :: dims(:)
550 integer, intent(in) :: nflds
551 real(kind=8), intent(in) :: field(:,:)
552
553 ! local variable
554 integer :: n, ncid, varid, idimid, jdimid, fdimid
555 real(kind=8), allocatable :: a3d(:,:,:)
556 character(len=20) :: subname = 'dumpnc2d'
557 !----------------------------------------------------------------------------
558
559 if (debug)write(logunit,'(a)')'enter '//trim(subname)//' variable '//vname
560 allocate(a3d(dims(1),dims(2),nflds)); a3d = 0.0
561
562 call nf90_err(nf90_create(trim(fname), nf90_64bit_offset, ncid), 'create: '//fname)
563 call nf90_err(nf90_def_dim(ncid, 'nx', dims(1), idimid), 'define dimension: nx')
564 call nf90_err(nf90_def_dim(ncid, 'ny', dims(2), jdimid), 'define dimension: ny')
565 call nf90_err(nf90_def_dim(ncid, 'nf', nflds, fdimid), 'define dimension: nf')
566 call nf90_err(nf90_def_var(ncid, vname, nf90_double, (/idimid,jdimid,fdimid/), varid), &
567 'define variable: '//vname)
568 call nf90_err(nf90_enddef(ncid), 'nf90_enddef: '//fname)
569
570 do n = 1,nflds
571 a3d(:,:,n) = reshape(field(n,1:dims(1)*dims(2)), (/dims(1),dims(2)/))
572 end do
573 call nf90_err(nf90_put_var(ncid, varid, a3d), 'put variable: '//vname)
574 call nf90_err(nf90_close(ncid), 'close: '//fname)
575
576 if (debug)write(logunit,'(a)')'exit '//trim(subname)//' variable '//vname
577 end subroutine dumpnc2d
578
589 subroutine dumpnc3d(fname, vname, dims, nk, nflds, field)
590
591 character(len=*), intent(in) :: fname, vname
592 integer, intent(in) :: dims(:)
593 integer, intent(in) :: nk, nflds
594 real(kind=8), intent(in) :: field(:,:,:)
595
596 ! local variable
597 integer :: n, k, ncid, varid, idimid, jdimid, kdimid, fdimid
598 real(kind=8), allocatable :: a4d(:,:,:,:)
599 character(len=20) :: subname = 'dumpnc3d'
600 !----------------------------------------------------------------------------
601
602 if (debug)write(logunit,'(a)')'enter '//trim(subname)//' variable '//vname
603 allocate(a4d(dims(1),dims(2),dims(3),nflds)); a4d = 0.0
604
605 call nf90_err(nf90_create(trim(fname), nf90_64bit_offset, ncid), 'create: '//fname)
606 call nf90_err(nf90_def_dim(ncid, 'nx', dims(1), idimid), 'define dimension: nx')
607 call nf90_err(nf90_def_dim(ncid, 'ny', dims(2), jdimid), 'define dimension: ny')
608 call nf90_err(nf90_def_dim(ncid, 'nk', dims(3), kdimid), 'define dimension: nk')
609 call nf90_err(nf90_def_dim(ncid, 'nf', nflds, fdimid), 'define dimension: nf')
610 call nf90_err(nf90_def_var(ncid, vname, nf90_double, (/idimid,jdimid,kdimid,fdimid/), varid), &
611 'define variable: '//vname)
612 call nf90_err(nf90_enddef(ncid), 'nf90_enddef: '//fname)
613
614 do n = 1,nflds
615 do k = 1,nk
616 a4d(:,:,k,n) = reshape(field(n,k,1:dims(1)*dims(2)), (/dims(1),dims(2)/))
617 end do
618 end do
619 call nf90_err(nf90_put_var(ncid, varid, a4d), 'put variable: '//vname)
620 call nf90_err(nf90_close(ncid), 'close: '//fname)
621
622 if (debug)write(logunit,'(a)')'exit '//trim(subname)//' variable '//vname
623 end subroutine dumpnc3d
624
633 subroutine dumpnc3dk(fname, vname, dims, field)
634
635 character(len=*), intent(in) :: fname, vname
636 integer, intent(in) :: dims(:)
637 real(kind=8), intent(in) :: field(:,:)
638
639 ! local variable
640 integer :: k, ncid, varid, idimid, jdimid, kdimid
641 real(kind=8), allocatable :: a3d(:,:,:)
642 character(len=20) :: subname = 'dumpnc3dk'
643 !----------------------------------------------------------------------------
644
645 if (debug)write(logunit,'(a)')'enter '//trim(subname)//' variable '//vname
646 allocate(a3d(dims(1),dims(2),dims(3))); a3d = 0.0
647
648 call nf90_err(nf90_create(trim(fname), nf90_64bit_offset, ncid), 'nf90_create: '//fname)
649 call nf90_err(nf90_def_dim(ncid, 'nx', dims(1), idimid), 'define dimension: nx')
650 call nf90_err(nf90_def_dim(ncid, 'ny', dims(2), jdimid), 'define dimension: ny')
651 call nf90_err(nf90_def_dim(ncid, 'nk', dims(3), kdimid), 'define dimension: nk')
652 call nf90_err(nf90_def_var(ncid, vname, nf90_double, (/idimid,jdimid,kdimid/), varid), &
653 'define variable: '//vname)
654 call nf90_err(nf90_enddef(ncid), 'nf90_enddef: '//fname)
655
656 do k = 1,dims(3)
657 a3d(:,:,k) = reshape(field(k,1:dims(1)*dims(2)), (/dims(1),dims(2)/))
658 end do
659 call nf90_err(nf90_put_var(ncid, varid, a3d), 'put variable: '//vname)
660 call nf90_err(nf90_close(ncid), 'close: '//fname)
661
662 if (debug)write(logunit,'(a)')'exit '//trim(subname)//' variable '//vname
663
664 end subroutine dumpnc3dk
665
674 subroutine dumpnc1d(fname, vname, dims, field)
675
676 character(len=*), intent(in) :: fname, vname
677 integer, intent(in) :: dims(:)
678 real(kind=8), intent(in) :: field(:)
679
680 ! local variable
681 integer :: ncid, varid, idimid, jdimid
682 real(kind=8), allocatable :: a2d(:,:)
683 character(len=20) :: subname = 'dumpnc1d'
684 !----------------------------------------------------------------------------
685
686 if (debug)write(logunit,'(a)')'enter '//trim(subname)//' variable '//vname
687 allocate(a2d(dims(1),dims(2))); a2d = 0.0
688
689 call nf90_err(nf90_create(trim(fname), nf90_64bit_offset, ncid), 'nf90_create: '//fname)
690 call nf90_err(nf90_def_dim(ncid, 'nx', dims(1), idimid), 'define dimension: nx')
691 call nf90_err(nf90_def_dim(ncid, 'ny', dims(2), jdimid), 'define dimension: ny')
692 call nf90_err(nf90_def_var(ncid, vname, nf90_double, (/idimid,jdimid/), varid), &
693 'define variable: '//vname)
694 call nf90_err(nf90_enddef(ncid), 'nf90_enddef: '//fname)
695
696 a2d(:,:) = reshape(field(1:dims(1)*dims(2)), (/dims(1),dims(2)/))
697 call nf90_err(nf90_put_var(ncid, varid, a2d), 'put variable: '//vname)
698 call nf90_err(nf90_close(ncid), 'close: '//fname)
699
700 if (debug)write(logunit,'(a)')'exit '//trim(subname)//' variable '//vname
701
702 end subroutine dumpnc1d
703
711 subroutine zero_out_land_ice(mask, fout, icnt)
712
713 integer(kind=4), intent(in) :: mask(:)
714 real(kind=8), intent(inout) :: fout(:,:)
715 integer, intent(out) :: icnt
716
717 !local variables
718 integer :: ij
719 !----------------------------------------------------------------------------
720
721 icnt = 0
722 do ij = 1,size(fout,2)
723 if ( mask(ij) .eq. 0 .and. sum(fout(:,ij)) .ne. 0.0) then
724 icnt = icnt + 1
725 fout(:,ij) = 0.0
726 end if
727 end do
728 end subroutine zero_out_land_ice
729
738 subroutine zero_out_phantom_ice(fin, fout, icnt)
739
740 real(kind=8), intent(in) :: fin(:,:)
741 real(kind=8), intent(inout) :: fout(:,:)
742 integer, intent(out) :: icnt
743
744 !local variables
745 integer :: ij
746 !----------------------------------------------------------------------------
747
748 icnt = 0
749 do ij = 1,size(fout,2)
750 if (sum(fin(:,ij)) .eq. 0.0 .and. sum(fout(:,ij)) .ne. 0.0 ) then
751 icnt = icnt + 1
752 fout(:,ij) = 0.0
753 end if
754 end do
755 end subroutine zero_out_phantom_ice
756
763 subroutine nf90_err(ierr, string)
764
765 integer , intent(in) :: ierr
766 character(len=*), intent(in) :: string
767 !----------------------------------------------------------------------------
768
769 if (ierr /= nf90_noerr) then
770 write(0, '(a)') 'FATAL ERROR: ' // trim(string)// ' : ' // trim(nf90_strerror(ierr))
771 ! This fails on WCOSS2 with Intel 19 compiler. See
772 ! https://community.intel.com/t5/Intel-Fortran-Compiler/STOP-and-ERROR-STOP-with-variable-stop-codes/m-p/1182521#M149254
773 ! When WCOSS2 moves to Intel 2020+, uncomment the next line and remove stop 99
774 !stop ierr
775 stop 99
776 end if
777 end subroutine nf90_err
778
779end module utils_mod