11 use atlas_module,
only: atlas_field
12 use fckit_configuration_module,
only: fckit_configuration
13 use fckit_log_module,
only: fckit_log
41 real(kind_real),
allocatable :: values(:,:)
45 #define LISTED_TYPE qg_gom
48 #include "oops/util/linkedList_i.f"
58 #include "oops/util/linkedList_c.f"
66 type(
qg_gom),
intent(inout) :: self
67 integer,
intent(in) :: nobs
77 self%ix = 0; self%iq = 0; self%iu = 0; self%iv = 0
78 do ivar = 1, vars%nvars()
79 if (vars%variable(ivar) ==
'x')
then
83 if (vars%variable(ivar) ==
'q')
then
87 if (vars%variable(ivar) ==
'u')
then
91 if (vars%variable(ivar) ==
'v')
then
98 allocate(self%values(self%nv,self%nobs))
109 type(
qg_gom),
intent(inout) :: self
112 self%lalloc = .false.
122 type(
qg_gom),
intent(inout) :: self
125 if (self%lalloc)
then
126 deallocate(self%values)
137 type(
qg_gom),
intent(inout) :: self
138 type(
qg_gom),
intent(in) :: other
141 self%nobs = other%nobs
147 self%used = other%used
149 if (.not.self%lalloc)
then
150 allocate(self%values(self%nv,self%nobs))
155 self%values = other%values
165 type(
qg_gom),
intent(inout) :: self
178 type(
qg_gom),
intent(inout) :: self
181 self%values = abs(self%values)
191 type(
qg_gom),
intent(inout) :: self
194 call normal_distribution(self%values,0.0_kind_real,1.0_kind_real)
204 type(
qg_gom),
intent(inout) :: self
205 real(kind_real),
intent(in) :: zz
213 self%values(jv,jo) = zz*self%values(jv,jo)
225 type(
qg_gom),
intent(inout) :: self
226 type(
qg_gom),
intent(in) :: other
234 self%values(jv,jo) = self%values(jv,jo)+other%values(jv,jo)
246 type(
qg_gom),
intent(inout) :: self
247 type(
qg_gom),
intent(in) :: other
255 self%values(jv,jo) = self%values(jv,jo)-other%values(jv,jo)
267 type(
qg_gom),
intent(inout) :: self
268 type(
qg_gom),
intent(in) :: other
276 self%values(jv,jo) = self%values(jv,jo)*other%values(jv,jo)
288 type(
qg_gom),
intent(inout) :: self
289 type(
qg_gom),
intent(in) :: other
292 real(kind_real) :: tol
301 if (abs(other%values(jvar,jloc))>tol)
then
302 self%values(jvar,jloc) = self%values(jvar,jloc)/other%values(jvar,jloc)
304 self%values(jvar,jloc) = 0.0
317 type(
qg_gom),
intent(inout) :: self
318 real(kind_real),
intent(inout) :: rms
329 rms = rms+self%values(jv,jo)**2
334 rms = sqrt(rms/(self%nobs*self%nv))
344 type(
qg_gom),
intent(in) :: gom1
345 type(
qg_gom),
intent(in) :: gom2
346 real(kind_real),
intent(inout) :: prod
352 if ((gom1%nv/=gom2%nv).or.(gom1%nobs/=gom2%nobs))
call abor1_ftn(
'qg_gom_dotprod: inconsistent GOM sizes')
360 prod = prod+gom1%values(jv,jo)*gom2%values(jv,jo)
372 type(
qg_gom),
intent(inout) :: self
373 integer,
intent(inout) :: kobs
374 real(kind_real),
intent(inout) :: scaling
375 real(kind_real),
intent(inout) :: pmin
376 real(kind_real),
intent(inout) :: pmax
377 real(kind_real),
intent(inout) :: prms
380 real(kind_real) :: expo,scalinginv,svalues(self%nv,self%nobs)
381 character(len=1024) :: msg
384 if (maxval(abs(self%values))>0.0)
then
385 expo = aint(log(maxval(abs(self%values)))/log(10.0_kind_real))-1
390 scalinginv = 1.0/scaling
393 svalues = self%values*scalinginv
397 pmin = minval(svalues)
398 pmax = maxval(svalues)
399 prms = sqrt(sum(svalues**2)/real(self%nobs*self%nv,kind_real))
409 type(
qg_gom),
intent(inout) :: self
410 real(kind_real),
intent(inout) :: mxval
411 integer,
intent(inout) :: iloc
412 integer,
intent(inout) :: ivar
418 mxval = maxval(self%values)
419 mxloc = maxloc(self%values)
433 type(
qg_gom),
intent(inout) :: self
434 type(fckit_configuration),
intent(in) :: f_conf
437 integer :: ncid,nobs_id,nv_id,values_id
438 character(len=1024) :: filename
439 character(len=:),
allocatable :: str
442 if (self%lalloc)
call abor1_ftn(
'qg_gom_read_file: gom alredy allocated')
445 call f_conf%get_or_die(
"filename",str)
447 call fckit_log%info(
'qg_gom_read_file: reading '//trim(filename))
450 call ncerr(nf90_open(trim(filename)//
'.nc',nf90_nowrite,ncid))
453 call ncerr(nf90_inq_dimid(ncid,
'nobs',nobs_id))
454 call ncerr(nf90_inq_dimid(ncid,
'nv',nv_id))
457 call ncerr(nf90_inquire_dimension(ncid,nobs_id,len=self%nobs))
458 call ncerr(nf90_inquire_dimension(ncid,nv_id,len=self%nv))
461 call ncerr(nf90_get_att(ncid,nf90_global,
'used',self%used))
462 call ncerr(nf90_get_att(ncid,nf90_global,
'ix',self%ix))
463 call ncerr(nf90_get_att(ncid,nf90_global,
'iq',self%iq))
464 call ncerr(nf90_get_att(ncid,nf90_global,
'iu',self%iu))
465 call ncerr(nf90_get_att(ncid,nf90_global,
'iv',self%iv))
466 call ncerr(nf90_get_att(ncid,nf90_global,
'nv',self%nv))
469 allocate(self%values(self%nv,self%nobs))
476 call ncerr(nf90_inq_varid(ncid,
'values',values_id))
479 call ncerr(nf90_get_var(ncid,values_id,self%values))
482 call ncerr(nf90_close(ncid))
492 type(
qg_gom),
intent(inout) :: self
493 type(fckit_configuration),
intent(in) :: f_conf
496 integer :: ncid,nobs_id,nv_id,values_id
497 character(len=1024) :: filename
498 character(len=:),
allocatable :: str
501 if (.not.self%lalloc)
call abor1_ftn(
'qg_gom_write_file: gom not allocated')
504 call f_conf%get_or_die(
"filename",str)
506 call fckit_log%info(
'qg_gom_write_file: writing '//trim(filename))
509 call ncerr(nf90_create(trim(filename)//
'.nc',or(nf90_clobber,nf90_64bit_offset),ncid))
512 call ncerr(nf90_def_dim(ncid,
'nobs',self%nobs,nobs_id))
513 call ncerr(nf90_def_dim(ncid,
'nv',self%nv,nv_id))
516 call ncerr(nf90_put_att(ncid,nf90_global,
'used',self%used))
517 call ncerr(nf90_put_att(ncid,nf90_global,
'ix',self%ix))
518 call ncerr(nf90_put_att(ncid,nf90_global,
'iq',self%iq))
519 call ncerr(nf90_put_att(ncid,nf90_global,
'iu',self%iu))
520 call ncerr(nf90_put_att(ncid,nf90_global,
'iv',self%iv))
521 call ncerr(nf90_put_att(ncid,nf90_global,
'nv',self%nv))
524 call ncerr(nf90_def_var(ncid,
'values',nf90_double,(/nv_id,nobs_id/),values_id))
527 call ncerr(nf90_enddef(ncid))
530 call ncerr(nf90_put_var(ncid,values_id,self%values))
533 call ncerr(nf90_close(ncid))
543 type(
qg_gom),
intent(inout) :: self
544 type(
qg_locs),
intent(inout) :: locs
545 type(fckit_configuration),
intent(in) :: f_conf
549 real(kind_real) :: x,y
550 character(len=30) :: ic
551 character(len=:),
allocatable :: str
552 real(kind_real),
pointer :: lonlat(:,:), z(:)
553 type(atlas_field) :: lonlat_field, z_field
556 lonlat_field = locs%lonlat()
557 call lonlat_field%data(lonlat)
559 z_field = locs%altitude()
563 if (.not. self%lalloc)
call abor1_ftn(
'qg_gom_analytic init: gom not allocated')
566 call f_conf%get_or_die(
"analytic_init",str)
568 call fckit_log%info(
'qg_gom_analytic_init: ic = '//trim(ic))
569 do iloc=1,locs%nlocs()
570 select case (trim(ic))
571 case (
'baroclinic-instability')
580 case (
'large-vortices')
585 if (self%ix>0)
call large_vortices(x,y,z(iloc),
'x',self%values(self%ix,iloc))
586 if (self%iq>0)
call large_vortices(x,y,z(iloc),
'q',self%values(self%iq,iloc))
587 if (self%iu>0)
call large_vortices(x,y,z(iloc),
'u',self%values(self%iu,iloc))
588 if (self%iv>0)
call large_vortices(x,y,z(iloc),
'v',self%values(self%iv,iloc))
590 call abor1_ftn(
'qg_gom_analytic_init: unknown initialization')
594 call lonlat_field%final()