ROMS
Loading...
Searching...
No Matches
obs_depth.F
Go to the documentation of this file.
1#include "cppdefs.h"
2#if (defined FOUR_DVAR || defined VERIFICATION) && \
3 defined observations && \
4 defined distribute && defined solve3d
5 SUBROUTINE obs_depth (ng, tile, model)
6!
7!git $Id$
8!================================================== Hernan G. Arango ===
9! Copyright (c) 2002-2025 The ROMS Group !
10! Licensed under a MIT/X style license !
11! See License_ROMS.md !
12!=======================================================================
13! !
14! This routine proccess the observation vertical position, Zobs. If !
15! the input position is in meters (Zobs < 0) instead of fractional !
16! level (Zobs > 0; 1 <= Zobs <= N), zero out its Zobs value in all !
17! nodes in group but itself to faciliate exchange between all tiles !
18! with "mp_collect", after model extraction. !
19! !
20!-----------------------------------------------------------------------
21!
22 USE mod_param
23 USE mod_fourdvar
24 USE mod_ncparam
25 USE mod_scalars
26!
27 implicit none
28!
29! Imported variable declarations.
30!
31 integer, intent(in) :: ng, tile, model
32!
33! Local variable declarations.
34!
35 logical :: r_bound, u_bound, v_bound
36
37 integer :: Mstr, Mend, iobs, itrc
38
39 real(r8) :: IniVal = 0.0_r8
40
41# include "set_bounds.h"
42!
43!-----------------------------------------------------------------------
44! If distributed-memory and Zobs < 0, zero-out Zobs values in all
45! nodes but itself to facilitate exchages between tiles latter.
46!-----------------------------------------------------------------------
47!
48! Set starting index of obervation vectors for reading. In weak
49! constraint, the entire observation data is loaded. Otherwise,
50! only the observation for the current time window are loaded
51! and started from vector index one.
52!
53# ifdef WEAK_CONSTRAINT
54 mstr=nstrobs(ng)
55 mend=nendobs(ng)
56# else
57 mstr=1
58 mend=nobs(ng)
59# endif
60!
61! Set observations depth (Zobs) to zero if not bounded in the current
62! parallel tile to facilitate global reduction in "mp_collect".
63!
64 DO iobs=mstr,mend
65 r_bound=((rxmin(ng) .le.xobs(iobs)).and. &
66 & (xobs(iobs).lt.rxmax(ng))).and. &
67 & ((rymin(ng) .le.yobs(iobs)).and. &
68 & (yobs(iobs).lt.rymax(ng)))
69 u_bound=((uxmin(ng) .le.xobs(iobs)).and. &
70 & (xobs(iobs).lt.uxmax(ng))).and. &
71 & ((uymin(ng) .le.yobs(iobs)).and. &
72 & (yobs(iobs).lt.uymax(ng)))
73 v_bound=((vxmin(ng) .le.xobs(iobs)).and. &
74 & (xobs(iobs).lt.vxmax(ng))).and. &
75 & ((vymin(ng) .le.yobs(iobs)).and. &
76 & (yobs(iobs).lt.vymax(ng)))
77 IF (.not.(r_bound.or.u_bound.or.v_bound)) THEN
78 zobs(iobs)=inival ! not bouded
79 END IF
80 END DO
81!
82 RETURN
83 END SUBROUTINE obs_depth
84#else
85 SUBROUTINE obs_depth
86 RETURN
87 END SUBROUTINE obs_depth
88#endif
integer, dimension(:), allocatable nobs
real(r8), dimension(:), allocatable zobs
real(r8), dimension(:), allocatable xobs
real(r8), dimension(:), allocatable yobs
integer, dimension(:), allocatable nstrobs
integer, dimension(:), allocatable nendobs
real(r8), dimension(:), allocatable rymin
real(r8), dimension(:), allocatable vymin
real(r8), dimension(:), allocatable rymax
real(r8), dimension(:), allocatable uymin
real(r8), dimension(:), allocatable vymax
real(r8), dimension(:), allocatable uxmin
real(r8), dimension(:), allocatable uxmax
real(r8), dimension(:), allocatable rxmax
real(r8), dimension(:), allocatable uymax
real(r8), dimension(:), allocatable vxmin
real(r8), dimension(:), allocatable vxmax
real(r8), dimension(:), allocatable rxmin
subroutine obs_depth(ng, tile, model)
Definition obs_depth.F:6