FV3-JEDI
fv3jedi_fields_mod.f90
Go to the documentation of this file.
1 ! (C) Copyright 2020 UCAR
2 !
3 ! This software is licensed under the terms of the Apache Licence Version 2.0
4 ! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
5 
7 
8 ! fckit
9 use fckit_mpi_module
10 use fckit_configuration_module
11 
12 ! oops
13 use datetime_mod
14 use oops_variables_mod
15 use string_utils
16 
17 ! fv3jedi
26 
27 implicit none
28 private
29 public :: fv3jedi_fields
30 
31 ! --------------------------------------------------------------------------------------------------
32 
33 ! Fields type (base for State and Increment)
35 
36  integer :: isc, iec, jsc, jec, npx, npy, npz, nf ! Geometry convenience
37  integer :: calendar_type, date_init(6) ! FMS style datetime
38  type(fckit_mpi_comm) :: f_comm ! Communicator
39  type(fv3jedi_field), allocatable :: fields(:) ! Array of fields
40 
41  contains
42 
43  ! Methods needed by both state and increment classes
44  procedure, public :: create
45  procedure, public :: delete
46  procedure, public :: copy
47  procedure, public :: zero
48  procedure, public :: norm
49  procedure, public :: change_resol
50  procedure, public :: minmaxrms
51  procedure, public :: read
52  procedure, public :: write
53  procedure, public :: accumul
54  procedure, public :: serialize
55  procedure, public :: deserialize
56 
57  ! Public array/field accessor functions
58  procedure, public :: has_field => has_field_
59  generic, public :: get_field => get_field_return_type_pointer, &
62  procedure, public :: put_field => put_field_
63 
64  ! Private array/field accessor functions
65  procedure, private :: get_field_return_type_pointer
66  procedure, private :: get_field_return_array_pointer
67  procedure, private :: get_field_return_array_allocatable
68 
69 endtype fv3jedi_fields
70 
71 ! --------------------------------------------------------------------------------------------------
72 
73 contains
74 
75 ! --------------------------------------------------------------------------------------------------
76 
77 subroutine create(self, geom, vars, increment)
78 
79  class(fv3jedi_fields), intent(inout) :: self
80  type(fv3jedi_geom), intent(in) :: geom
81  type(oops_variables), intent(in) :: vars
82  logical, optional, intent(in) :: increment
83 
84  integer :: var, fc
85  type(field_metadata) :: fmd
86  logical :: is_increment, field_fail
87 
88  ! Allocate fields structure
89  ! -------------------------
90  self%nf = vars%nvars()
91  allocate(self%fields(self%nf))
92 
93  ! Check if this is an increment rather than state
94  ! -----------------------------------------------
95  if (.not. present(increment)) then
96  is_increment = .false.
97  else
98  is_increment = increment
99  endif
100 
101  ! Loop through and allocate actual fields
102  ! ---------------------------------------
103  fc = 0
104  do var = 1, vars%nvars()
105 
106  fmd = geom%fields%get_field(trim(vars%variable(var)))
107 
108  ! Uptick counter
109  fc=fc+1;
110 
111  self%fields(fc)%isc = geom%isc
112  self%fields(fc)%iec = geom%iec
113  self%fields(fc)%jsc = geom%jsc
114  self%fields(fc)%jec = geom%jec
115  self%fields(fc)%npz = fmd%levels
116 
117  field_fail = len(trim(fmd%field_name)) > field_clen
118  if (field_fail) call abor1_ftn("fv3jedi_fields.create: " //trim(fmd%field_name)// " too long")
119 
120  if(.not.self%fields(fc)%lalloc) then
121 
122  if (trim(fmd%stagger_loc) == 'center') then
123  allocate(self%fields(fc)%array(geom%isc:geom%iec,geom%jsc:geom%jec,1:fmd%levels))
124  elseif (trim(fmd%stagger_loc) == 'northsouth') then
125  allocate(self%fields(fc)%array(geom%isc:geom%iec,geom%jsc:geom%jec+1,1:fmd%levels))
126  elseif (trim(fmd%stagger_loc) == 'eastwest') then
127  allocate(self%fields(fc)%array(geom%isc:geom%iec+1,geom%jsc:geom%jec,1:fmd%levels))
128  endif
129 
130  endif
131 
132  self%fields(fc)%lalloc = .true.
133 
134  self%fields(fc)%short_name = trim(fmd%field_io_name)
135  if (is_increment) then
136  self%fields(fc)%long_name = "increment_of_"//trim(fmd%long_name)
137  else
138  self%fields(fc)%long_name = trim(fmd%long_name)
139  endif
140  self%fields(fc)%fv3jedi_name = trim(fmd%field_name)
141  self%fields(fc)%units = trim(fmd%units)
142  self%fields(fc)%io_file = trim(fmd%io_file)
143  self%fields(fc)%space = trim(fmd%space)
144  self%fields(fc)%staggerloc = trim(fmd%stagger_loc)
145  self%fields(fc)%tracer = fmd%tracer
146  self%fields(fc)%integerfield = trim(fmd%array_kind)=='integer'
147 
148  enddo
149 
150  ! Check field count
151  if (fc .ne. self%nf) call abor1_ftn("fv3jedi_fields_mod create: fc does not equal self%nf")
152 
153  ! Initialize all arrays to zero
154  call self%zero()
155 
156  ! Copy some geometry for convenience
157  self%isc = geom%isc
158  self%iec = geom%iec
159  self%jsc = geom%jsc
160  self%jec = geom%jec
161  self%npx = geom%npx
162  self%npy = geom%npy
163  self%npz = geom%npz
164 
165  ! Pointer to fv3jedi communicator
166  self%f_comm = geom%f_comm
167 
168  ! Check winds
169  field_fail = self%has_field('ua') .and. .not.self%has_field('va')
170  if (field_fail) call abor1_ftn("fv3jedi_fields.create: found A-Grid u but not v")
171  field_fail = .not.self%has_field('ua') .and. self%has_field('va')
172  if (field_fail) call abor1_ftn("fv3jedi_fields.create: found A-Grid v but not u")
173  field_fail = self%has_field('ud') .and. .not.self%has_field('vd')
174  if (field_fail) call abor1_ftn("fv3jedi_fields.create: found D-Grid u but not v")
175  field_fail = .not.self%has_field('ud') .and. self%has_field('vd')
176  if (field_fail) call abor1_ftn("fv3jedi_fields.create: found D-Grid v but not u")
177 
178  !Check User's choice of ozone variables.
179  field_fail = self%has_field('o3mr') .and. self%has_field('o3ppmv')
180  if (field_fail) call abor1_ftn("fv3jedi_fields.create: o3mr and o3ppmv created")
181 
182 endsubroutine create
183 
184 ! --------------------------------------------------------------------------------------------------
185 
186 subroutine delete(self)
187 
188 class(fv3jedi_fields), intent(inout) :: self
189 
190 integer :: var
191 
192 do var = 1, self%nf
193  if(self%fields(var)%lalloc) deallocate(self%fields(var)%array)
194  self%fields(var)%lalloc = .false.
195 enddo
196 deallocate(self%fields)
197 
198 end subroutine delete
199 
200 ! --------------------------------------------------------------------------------------------------
201 
202 subroutine copy(self, other)
203 
204 class(fv3jedi_fields), intent(inout) :: self
205 class(fv3jedi_fields), intent(in) :: other
206 
207 integer :: var
208 
209 call checksame(self%fields, other%fields, "fv3jedi_fields_mod.copy")
210 
211 do var = 1, self%nf
212  self%fields(var)%array = other%fields(var)%array
213 enddo
214 
215 self%calendar_type = other%calendar_type
216 self%date_init = other%date_init
217 
218 end subroutine copy
219 
220 ! --------------------------------------------------------------------------------------------------
221 
222 subroutine zero(self)
223 
224 class(fv3jedi_fields), intent(inout) :: self
225 
226 integer :: var
227 
228 do var = 1, self%nf
229  self%fields(var)%array = 0.0_kind_real
230 enddo
231 
232 endsubroutine zero
233 
234 ! --------------------------------------------------------------------------------------------------
235 
236 subroutine norm(self, normout)
237 
238 class(fv3jedi_fields), intent(inout) :: self
239 real(kind=kind_real), intent(out) :: normout
240 
241 integer :: i, j, k, ii, iisum, var
242 real(kind=kind_real) :: zz
243 
244 zz = 0.0_kind_real
245 ii = 0
246 
247 do var = 1, self%nf
248 
249  do k = 1, self%fields(var)%npz
250  do j = self%fields(var)%jsc, self%fields(var)%jec
251  do i = self%fields(var)%isc, self%fields(var)%iec
252  zz = zz + self%fields(var)%array(i,j,k)**2
253  ii = ii + 1
254  enddo
255  enddo
256  enddo
257 
258 enddo
259 
260 !Get global values
261 call self%f_comm%allreduce(zz,normout,fckit_mpi_sum())
262 call self%f_comm%allreduce(ii,iisum,fckit_mpi_sum())
263 normout = sqrt(normout/real(iisum,kind_real))
264 
265 endsubroutine norm
266 
267 ! --------------------------------------------------------------------------------------------------
268 
269 subroutine change_resol(self, geom, other, geom_other)
270 
271 implicit none
272 class(fv3jedi_fields), intent(inout) :: self
273 type(fv3jedi_geom), intent(inout) :: geom
274 class(fv3jedi_fields), intent(in) :: other
275 type(fv3jedi_geom), intent(inout) :: geom_other
276 
277 ! Interpolation
278 integer :: var
279 type(field2field_interp) :: interp
280 logical :: integer_interp = .false.
281 
282 call checksame(self%fields, other%fields, "fv3jedi_fields_mod.change_resol")
283 
284 if ((other%iec-other%isc+1)-(self%iec-self%isc+1) == 0) then
285 
286  call self%copy(other)
287 
288 else
289 
290  ! Check if integer interp needed
291  do var = 1, self%nf
292  if (other%fields(var)%integerfield) integer_interp = .true.
293  enddo
294 
295  call interp%create(geom%interp_method, integer_interp, geom_other, geom)
296  call interp%apply(self%nf, geom_other, other%fields, geom, self%fields)
297  call interp%delete()
298 
299  self%calendar_type = other%calendar_type
300  self%date_init = other%date_init
301 
302 endif
303 
304 end subroutine change_resol
305 
306 ! --------------------------------------------------------------------------------------------------
307 
308 subroutine minmaxrms(self, field_num, field_name, minmaxrmsout)
309 
310 class(fv3jedi_fields), intent(inout) :: self
311 integer, intent(in) :: field_num
312 character(len=*), intent(inout) :: field_name
313 real(kind=kind_real), intent(out) :: minmaxrmsout(3)
314 
315 integer :: isc, iec, jsc, jec, npz
316 real(kind=kind_real) :: tmp(3), gs3, gs3g
317 
318 isc = self%fields(field_num)%isc
319 iec = self%fields(field_num)%iec
320 jsc = self%fields(field_num)%jsc
321 jec = self%fields(field_num)%jec
322 npz = self%fields(field_num)%npz
323 
324 field_name = self%fields(field_num)%short_name
325 
326 ! Compute global sum over the field
327 gs3 = real((iec-isc+1)*(jec-jsc+1)*npz, kind_real)
328 call self%f_comm%allreduce(gs3,gs3g,fckit_mpi_sum())
329 
330 ! Min/Max/SumSquares
331 tmp(1) = minval(self%fields(field_num)%array(isc:iec,jsc:jec,1:npz))
332 tmp(2) = maxval(self%fields(field_num)%array(isc:iec,jsc:jec,1:npz))
333 tmp(3) = sum(self%fields(field_num)%array(isc:iec,jsc:jec,1:npz)**2)
334 
335 ! Get global min/max/sum
336 call self%f_comm%allreduce(tmp(1), minmaxrmsout(1), fckit_mpi_min())
337 call self%f_comm%allreduce(tmp(2), minmaxrmsout(2), fckit_mpi_max())
338 call self%f_comm%allreduce(tmp(3), minmaxrmsout(3), fckit_mpi_sum())
339 
340 ! SumSquares to rms
341 minmaxrmsout(3) = sqrt(minmaxrmsout(3)/gs3g)
342 
343 endsubroutine minmaxrms
344 
345 ! --------------------------------------------------------------------------------------------------
346 
347 subroutine read(self, geom, conf, vdate)
348 
349 class(fv3jedi_fields), intent(inout) :: self
350 type(fv3jedi_geom), intent(inout) :: geom
351 type(fckit_configuration), intent(in) :: conf
352 type(datetime), intent(inout) :: vdate
353 
354 type(fv3jedi_io_gfs) :: gfs
355 type(fv3jedi_io_geos) :: geos
356 
357 character(len=10) :: filetype
358 character(len=:), allocatable :: str
359 
360 
361 ! IO type
362 call conf%get_or_die("filetype",str)
363 filetype = str
364 deallocate(str)
365 
366 
367 if (trim(filetype) == 'gfs') then
368 
369  call gfs%setup_conf(conf)
370  call gfs%read_meta(geom, vdate, self%calendar_type, self%date_init)
371  call gfs%read_fields(geom, self%fields)
372 
373 elseif (trim(filetype) == 'geos') then
374 
375  call geos%setup_conf(geom, conf)
376  call geos%read_meta(geom, vdate, self%calendar_type, self%date_init, self%fields)
377  call geos%read_fields(geom, self%fields)
378  call geos%delete()
379 
380 else
381 
382  call abor1_ftn("fv3jedi_fields_mod.read: restart type not supported")
383 
384 endif
385 
386 endsubroutine read
387 
388 ! --------------------------------------------------------------------------------------------------
389 
390 subroutine write(self, geom, conf, vdate)
391 
392 class(fv3jedi_fields), intent(inout) :: self
393 type(fv3jedi_geom), intent(inout) :: geom
394 type(fckit_configuration), intent(in) :: conf
395 type(datetime), intent(inout) :: vdate
396 
397 type(fv3jedi_io_gfs) :: gfs
398 type(fv3jedi_io_geos) :: geos
399 type(fv3jedi_llgeom) :: latlon
400 
401 character(len=10) :: filetype
402 character(len=:), allocatable :: str
403 
404 ! IO type
405 call conf%get_or_die("filetype",str)
406 filetype = str
407 deallocate(str)
408 
409 if (trim(filetype) == 'gfs') then
410 
411  call gfs%setup_conf(conf)
412  call gfs%setup_date(vdate)
413  call gfs%write(geom, self%fields, vdate, self%calendar_type, self%date_init)
414 
415 elseif (trim(filetype) == 'geos') then
416 
417  call geos%setup_conf(geom, conf)
418  call geos%setup_date(vdate)
419  call geos%write(geom, self%fields, vdate)
420  call geos%delete()
421 
422 elseif (trim(filetype) == 'latlon') then
423 
424  call latlon%setup_conf(geom)
425  call latlon%setup_date(vdate)
426  call latlon%write(geom, self%fields, conf, vdate)
427 
428 else
429 
430  call abor1_ftn("fv3jedi_fields.write: restart type not supported")
431 
432 endif
433 
434 endsubroutine write
435 
436 ! --------------------------------------------------------------------------------------------------
437 
438 subroutine accumul(self, zz, rhs)
439 
440 implicit none
441 class(fv3jedi_fields), intent(inout) :: self
442 real(kind=kind_real), intent(in) :: zz
443 class(fv3jedi_fields), intent(in) :: rhs
444 
445 integer :: var
446 
447 call checksame(self%fields,rhs%fields,"fv3jedi_fields.accumul")
448 
449 do var = 1, self%nf
450  self%fields(var)%array = self%fields(var)%array + zz * rhs%fields(var)%array
451 enddo
452 
453 end subroutine accumul
454 
455 ! --------------------------------------------------------------------------------------------------
456 
457 subroutine serialize(self,vsize,vect_inc)
458 
459 implicit none
460 
461 ! Passed variables
462 class(fv3jedi_fields), intent(in) :: self
463 integer, intent(in) :: vsize
464 real(kind_real), intent(out) :: vect_inc(vsize)
465 
466 ! Local variables
467 integer :: ind, var, i, j, k
468 
469 ! Initialize
470 ind = 0
471 
472 ! Copy
473 do var = 1, self%nf
474  do k = 1,self%fields(var)%npz
475  do j = self%fields(var)%jsc,self%fields(var)%jec
476  do i = self%fields(var)%isc,self%fields(var)%iec
477  ind = ind + 1
478  vect_inc(ind) = self%fields(var)%array(i, j, k)
479  enddo
480  enddo
481  enddo
482 enddo
483 
484 end subroutine serialize
485 
486 ! --------------------------------------------------------------------------------------------------
487 
488 subroutine deserialize(self, vsize, vect_inc, index)
489 
490 implicit none
491 
492 ! Passed variables
493 class(fv3jedi_fields), intent(inout) :: self
494 integer, intent(in) :: vsize
495 real(kind_real), intent(in) :: vect_inc(vsize)
496 integer, intent(inout) :: index
497 
498 ! Local variables
499 integer :: ind, var, i, j, k
500 
501 ! Copy
502 do var = 1, self%nf
503  do k = 1,self%fields(var)%npz
504  do j = self%fields(var)%jsc,self%fields(var)%jec
505  do i = self%fields(var)%isc,self%fields(var)%iec
506  self%fields(var)%array(i, j, k) = vect_inc(index + 1)
507  index = index + 1
508  enddo
509  enddo
510  enddo
511 enddo
512 
513 end subroutine deserialize
514 
515 ! --------------------------------------------------------------------------------------------------
516 
517 logical function has_field_(self, field_name, field_index)
518 
519 class(fv3jedi_fields), intent(in) :: self
520 character(len=*), intent(in) :: field_name
521 integer, optional, intent(out) :: field_index
522 
523 has_field_ = has_field(self%fields, field_name, field_index)
524 
525 end function has_field_
526 
527 ! --------------------------------------------------------------------------------------------------
528 
529 subroutine get_field_return_type_pointer(self, field_name, field)
530 
531 class(fv3jedi_fields), target, intent(in) :: self
532 character(len=*), intent(in) :: field_name
533 type(fv3jedi_field), pointer, intent(inout) :: field
534 
535 call get_field(self%fields, field_name, field)
536 
537 endsubroutine get_field_return_type_pointer
538 
539 ! --------------------------------------------------------------------------------------------------
540 
541 subroutine get_field_return_array_pointer(self, field_name, field)
542 
543 class(fv3jedi_fields), target, intent(in) :: self
544 character(len=*), intent(in) :: field_name
545 real(kind=kind_real), pointer, intent(inout) :: field(:,:,:)
546 
547 call get_field(self%fields, field_name, field)
548 
549 endsubroutine get_field_return_array_pointer
550 
551 ! --------------------------------------------------------------------------------------------------
552 
553 subroutine get_field_return_array_allocatable(self, field_name, field)
554 
555 class(fv3jedi_fields), target, intent(in) :: self
556 character(len=*), intent(in) :: field_name
557 real(kind=kind_real), allocatable, intent(inout) :: field(:,:,:)
558 
559 call get_field(self%fields, field_name, field)
560 
562 
563 ! --------------------------------------------------------------------------------------------------
564 
565 subroutine put_field_(self, field_name, field)
566 
567 class(fv3jedi_fields), target, intent(inout) :: self
568 character(len=*), intent(in) :: field_name
569 real(kind=kind_real), allocatable, intent(in) :: field(:,:,:)
570 
571 call put_field(self%fields, field_name, field)
572 
573 endsubroutine put_field_
574 
575 ! --------------------------------------------------------------------------------------------------
576 
577 end module fv3jedi_fields_mod
fv3jedi_field_mod::checksame
subroutine, public checksame(fields1, fields2, calling_method)
Definition: fv3jedi_field_mod.f90:212
fv3jedi_fields_mod::copy
subroutine copy(self, other)
Definition: fv3jedi_fields_mod.f90:203
fv3jedi_field_mod
Definition: fv3jedi_field_mod.f90:6
fv3jedi_fields_mod::read
subroutine read(self, geom, conf, vdate)
Definition: fv3jedi_fields_mod.f90:348
fv3jedi_fields_mod::change_resol
subroutine change_resol(self, geom, other, geom_other)
Definition: fv3jedi_fields_mod.f90:270
fv3jedi_field_mod::has_field
logical function, public has_field(fields, field_name, field_index)
Definition: fv3jedi_field_mod.f90:58
fv3jedi_io_latlon_mod::fv3jedi_llgeom
Definition: fv3jedi_io_latlon_mod.f90:39
fv3jedi_fields_mod::create
subroutine create(self, geom, vars, increment)
Definition: fv3jedi_fields_mod.f90:78
fv3jedi_fields_mod::norm
subroutine norm(self, normout)
Definition: fv3jedi_fields_mod.f90:237
fv3jedi_interpolation_mod
Definition: fv3jedi_interpolation_mod.f90:6
fv3jedi_geom_mod
Fortran module handling geometry for the FV3 model.
Definition: fv3jedi_geom_mod.f90:8
fv3jedi_fields_mod::delete
subroutine delete(self)
Definition: fv3jedi_fields_mod.f90:187
fv3jedi_io_gfs_mod
Definition: fv3jedi_io_gfs_mod.f90:1
fv3jedi_fields_mod::get_field_return_array_pointer
subroutine get_field_return_array_pointer(self, field_name, field)
Definition: fv3jedi_fields_mod.f90:542
fv3jedi_field_mod::put_field
subroutine, public put_field(fields, field_name, field)
Definition: fv3jedi_field_mod.f90:175
fv3jedi_io_latlon_mod
Definition: fv3jedi_io_latlon_mod.f90:6
fields_metadata_mod::field_metadata
Definition: fields_metadata_mod.f90:28
fv3jedi_field_mod::get_field
Definition: fv3jedi_field_mod.f90:25
fields_metadata_mod
Definition: fields_metadata_mod.f90:6
fv3jedi_geom_mod::fv3jedi_geom
Fortran derived type to hold geometry data for the FV3JEDI model.
Definition: fv3jedi_geom_mod.f90:46
fv3jedi_fields_mod::has_field_
logical function has_field_(self, field_name, field_index)
Definition: fv3jedi_fields_mod.f90:518
fv3jedi_fields_mod::get_field_return_type_pointer
subroutine get_field_return_type_pointer(self, field_name, field)
Definition: fv3jedi_fields_mod.f90:530
fv3jedi_fields_mod::put_field_
subroutine put_field_(self, field_name, field)
Definition: fv3jedi_fields_mod.f90:566
fv3jedi_io_gfs_mod::fv3jedi_io_gfs
Definition: fv3jedi_io_gfs_mod.f90:35
fv3jedi_fields_mod::get_field_return_array_allocatable
subroutine get_field_return_array_allocatable(self, field_name, field)
Definition: fv3jedi_fields_mod.f90:554
fv3jedi_fields_mod::write
subroutine write(self, geom, conf, vdate)
Definition: fv3jedi_fields_mod.f90:391
fv3jedi_kinds_mod::kind_real
integer, parameter, public kind_real
Definition: fv3jedi_kinds_mod.f90:14
fv3jedi_fields_mod::serialize
subroutine serialize(self, vsize, vect_inc)
Definition: fv3jedi_fields_mod.f90:458
fv3jedi_io_geos_mod
Definition: fv3jedi_io_geos_mod.f90:6
fv3jedi_interpolation_mod::field2field_interp
Definition: fv3jedi_interpolation_mod.f90:28
fv3jedi_field_mod::fv3jedi_field
Definition: fv3jedi_field_mod.f90:36
fv3jedi_fields_mod
Definition: fv3jedi_fields_mod.f90:6
fv3jedi_fields_mod::deserialize
subroutine deserialize(self, vsize, vect_inc, index)
Definition: fv3jedi_fields_mod.f90:489
fv3jedi_io_geos_mod::fv3jedi_io_geos
Definition: fv3jedi_io_geos_mod.f90:36
fv3jedi_kinds_mod
Definition: fv3jedi_kinds_mod.f90:6
fv3jedi_fields_mod::zero
subroutine zero(self)
Definition: fv3jedi_fields_mod.f90:223
fv3jedi_fields_mod::minmaxrms
subroutine minmaxrms(self, field_num, field_name, minmaxrmsout)
Definition: fv3jedi_fields_mod.f90:309
fv3jedi_fields_mod::accumul
subroutine accumul(self, zz, rhs)
Definition: fv3jedi_fields_mod.f90:439
fv3jedi_field_mod::field_clen
integer, parameter, public field_clen
Definition: fv3jedi_field_mod.f90:31
fv3jedi_fields_mod::fv3jedi_fields
Definition: fv3jedi_fields_mod.f90:34