Go to the documentation of this file.
29 type(c_ptr),
intent(in) :: odb
30 integer(c_int),
intent(in) :: nn
31 real(c_double),
intent(inout) :: zz
38 real(kind_real),
allocatable :: values(:,:)
41 #define LISTED_TYPE qg_obsvec
44 #include "oops/util/linkedList_i.f"
54 #include "oops/util/linkedList_c.f"
63 integer,
intent(in) :: nlev
64 integer,
intent(in) :: nobs
71 if (
allocated(self%values))
deallocate(self%values)
74 allocate(self%values(self%nlev,self%nobs))
77 self%values = 0.0_kind_real
91 self%nlev = other%nlev
92 self%nobs = other%nobs
95 allocate(self%values(self%nlev,self%nobs))
108 deallocate(self%values)
121 if ((other%nlev/=self%nlev).or.(other%nobs/=self%nobs))
then
123 deallocate(self%values)
126 self%nlev = other%nlev
127 self%nobs = other%nobs
130 allocate(self%values(self%nlev,self%nobs))
134 self%values = other%values
146 integer,
intent(in) :: idx(:)
151 if ((other%nlev/=self%nlev).or.(
size(idx)/=self%nobs))
then
153 deallocate(self%values)
156 self%nlev = other%nlev
157 self%nobs =
size(idx)
160 allocate(self%values(self%nlev,self%nobs))
165 self%values(:,i) = other%values(:,idx(i)+1)
190 real(kind_real),
intent(in) :: zz
193 self%values = zz*self%values
207 self%values = self%values+other%values
221 self%values = self%values-other%values
235 self%values = self%values*other%values
249 self%values = self%values/other%values
260 real(kind_real),
intent(in) :: zz
264 self%values = self%values+zz*other%values
277 self%values = 1.0/self%values
287 type(c_ptr),
intent(in) :: c_odb
294 nval = self%nobs*self%nlev
309 real(kind_real),
intent(inout) :: zz
315 if ((obsvec1%nobs/=obsvec2%nobs).or.(obsvec1%nlev/=obsvec2%nlev))
call abor1_ftn(
'qg_obsvec_dotprod: inconsistent sizes')
321 do jobs=1,obsvec1%nobs
322 do jlev=1,obsvec1%nlev
323 zz = zz+obsvec1%values(jlev,jobs)*obsvec2%values(jlev,jobs)
336 real(kind_real),
intent(inout) :: scaling
337 real(kind_real),
intent(inout) :: zmin
338 real(kind_real),
intent(inout) :: zmax
339 real(kind_real),
intent(inout) :: zavg
343 real(kind_real) :: expo
345 if (self%nobs>0.and.self%nlev>0)
then
347 if (.not.
allocated(self%values))
call abor1_ftn(
'qg_obsvec_stats: obs vector not allocated')
357 zmin = min(self%values(jlev,jobs),zmin)
358 zmax = max(self%values(jlev,jobs),zmax)
359 zavg = zavg+self%values(jlev,jobs)
364 zavg = zavg/real(self%nlev*self%nobs,kind_real)
373 if (abs(zavg)>0.0)
then
374 expo = aint(log(abs(zavg))/log(10.0_kind_real))
392 integer,
intent(inout) :: kobs
395 kobs = self%nobs*self%nlev
407 integer,
intent(in) :: iob
408 real(kind_real),
intent(out) :: val
412 i1 = iob / self%nobs + 1
413 i2 = iob - self%nobs*(i1-1) + 1
416 if (i1>self%nlev .or. i2>self%nobs)
call abor1_ftn (
'qg_obsvec_getat: index is out of bounds')
418 val = self%values(i1,i2)
subroutine, public qg_obsvec_copy(self, other)
Copy observation vector.
subroutine, public qg_obsvec_mul_scal(self, zz)
Multiply observation vector with a scalar.
subroutine, public qg_obsvec_dotprod(obsvec1, obsvec2, zz)
Compute dot product between observation vectors.
subroutine, public qg_obsvec_axpy(self, zz, other)
Apply axpy on observation vector.
subroutine, public qg_obsvec_nobs(self, kobs)
Get observation vector size.
subroutine, public qg_obsvec_sub(self, other)
Subtract observation vector.
subroutine, public qg_obsvec_zero(self)
Set observation vector to zero.
subroutine, public qg_obsvec_delete(self)
Delete observation vector.
subroutine, public qg_obsvec_clone(self, other)
Clone observation vector.
subroutine, public qg_obsvec_mul(self, other)
Multiply observation vector.
type(registry_t), public qg_obsvec_registry
Linked list interface - defines registry_t type.
subroutine, public qg_obsvec_random(c_odb, self)
Generate random observation vector.
subroutine, public qg_obsvec_copy_local(self, other, idx)
Copy a local subset of the observation vector.
subroutine, public qg_obsvec_setup(self, nlev, nobs)
Linked list implementation.
subroutine, public qg_obsvec_div(self, other)
Divide observation vector.
subroutine, public qg_obsvec_stats(self, scaling, zmin, zmax, zavg)
Compute observation vector statistics.
subroutine, public qg_obsvec_getat(self, iob, val)
Get value from observation vector at location (iob)
subroutine, public qg_obsvec_add(self, other)
Add observation vector.
subroutine, public qg_obsvec_invert(self)
Invert observation vector.