ROMS
Loading...
Searching...
No Matches
obs_k2z.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#ifdef SOLVE3D
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine converts observations vertical fractional coordinate !
13! (Zobs => obs_grid) to depth in meters (obs_depths). The depths are !
14! negative downwards. This needed by the Ensemble Kalman (EnKF) for !
15! localization. !
16! !
17! On Input: !
18! !
19! ng Nested grid number. !
20! Imin Global I-coordinate lower bound of RHO-points. !
21! Imax Global I-coordinate upper bound of RHO-points. !
22! Jmin Global J-coordinate lower bound of RHO-points. !
23! Jmax Global J-coordinate upper bound of RHO-points. !
24! LBi I-dimension Lower bound. !
25! UBi I-dimension Upper bound. !
26! LBj J-dimension Lower bound. !
27! UBj J-dimension Upper bound. !
28! LBk K-dimension Lower bound. !
29! UBk K-dimension Upper bound. !
30! Xmin Global minimum fractional I-coordinate to consider. !
31! Xmax Global maximum fractional I-coordinate to consider. !
32! Ymin Global minimum fractional J-coordinate to consider. !
33! Ymax Global maximum fractional J-coordinate to consider. !
34! Mobs Number of observations. !
35! Xobs Observations X-locations (fractional coordinates). !
36! Yobs Observations Y-locations (fractional coordinates). !
37! Zobs Observations Z-locations (fractional coordinates or !
38! or actual meters). !
39! obs_scale Observation screenning flag. !
40! z Model grid depths of W-points (meters, 3D array). !
41! !
42! On Output: !
43! !
44! obs_depths Observations depth (meters, negative downwards. !
45! !
46! The interpolation weights matrix, Hmat(1:8), is as follows: !
47! !
48! 8____________7 !
49! /. /| (i2,j2,k2) !
50! / . / | !
51! 5/___________/6 | !
52! | . | | !
53! | . | | Grid Cell !
54! | 4.........|..|3 !
55! | . | / !
56! |. | / !
57! (i1,j1,k1) |___________|/ !
58! 1 2 !
59! !
60! Notice that the indices i2 and j2 are reset when observations are !
61! located exactly at the eastern and/or northern boundaries. This is !
62! needed to avoid out-of-range array computations. !
63! !
64! All the observations are assumed to in fractional coordinates with !
65! respect to RHO-points: !
66! !
67! !
68! M r..u..r..u..r..u..r..u..r..u..r..u..r..u..r..u..r..u..r !
69! : : !
70! Mm+.5 v p++v++p++v++p++v++p++v++p++v++p++v++p++v++p++v++p v !
71! : + | | | | | | | + : !
72! Mm r u r u r u r u r u r u r u r u r u r !
73! : + | | | | | | | + : !
74! Mm-.5 v p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p v !
75! : + | | | | | | | + : !
76! r u r u r u r u r u r u r u r u r u r !
77! : + | | | | | | | + : !
78! v p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p v !
79! : + | | | | | | | + : !
80! r u r u r u r u r u r u r u r u r u r !
81! : + | | | | | | | + : !
82! 2.5 v p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p v !
83! : + | | | | | | | + : !
84! 2.0 r u r u r u r u r u r u r u r u r u r !
85! : + | | | | | | | + : !
86! 1.5 v p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p v !
87! : + | | | | | | | + : !
88! 1.0 r u r u r u r u r u r u r u r u r u r !
89! : + | | | | | | | + : !
90! 0.5 v p++v++p++v++p++v++p++v++p++v++p++v++p++v++p++v++p v !
91! : : !
92! 0.0 r..u..r..u..r..u..r..u..r..u..r..u..r..u..r..u..r..u..r !
93! !
94! 0.5 1.5 2.5 Lm-.5 Lm+.5 !
95! !
96! 0.0 1.0 2.0 Lm L !
97! !
98!=======================================================================
99!
100 USE mod_kinds
101!
102 implicit none
103!
104 CONTAINS
105!
106!***********************************************************************
107 SUBROUTINE obs_k2z (ng, Imin, Imax, Jmin, Jmax, &
108 & LBi, UBi, LBj, UBj, LBk, UBk, &
109 & Xmin, Xmax, Ymin, Ymax, &
110 & Mobs, Xobs, Yobs, Zobs, obs_scale, &
111 & z, obs_depths)
112!***********************************************************************
113!
114 USE mod_param
115!
116! Imported variable declarations.
117!
118 integer, intent(in) :: ng, Imin, Imax, Jmin, Jmax
119 integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk
120 integer, intent(in) :: Mobs
121!
122 real(r8), intent(in) :: Xmin, Xmax, Ymin, Ymax
123!
124# ifdef ASSUMED_SHAPE
125 real(r8), intent(in) :: obs_scale(:)
126 real(r8), intent(in) :: Xobs(:)
127 real(r8), intent(in) :: Yobs(:)
128 real(r8), intent(in) :: Zobs(:)
129 real(r8), intent(in) :: z(LBi:,LBj:,LBk:)
130
131 real(r8), intent(out) :: obs_depths(:)
132# else
133 real(r8), intent(in) :: obs_scale(Mobs)
134 real(r8), intent(in) :: Xobs(Mobs)
135 real(r8), intent(in) :: Yobs(Mobs)
136 real(r8), intent(in) :: Zobs(Mobs)
137 real(r8), intent(in) :: z(LBi:UBi,LBj:UBj,LBk:UBk)
138
139 real(r8), intent(out) :: obs_depths(Mobs)
140# endif
141!
142! Local variable declarations.
143!
144 logical :: Linterpolate
145
146 integer :: i, ic, iobs, i1, i2, j1, j2, k, k1, k2
147
148 real(r8) :: Zbot, Ztop, dz, p1, p2, q1, q2, r1, r2
149 real(r8) :: w11, w12, w21, w22
150
151 real(r8), dimension(8) :: Hmat
152!
153!-----------------------------------------------------------------------
154! Interpolate vertical fractional coordinate to depths.
155!-----------------------------------------------------------------------
156!
157 DO iobs=1,mobs
158 IF (((xmin.le.xobs(iobs)).and.(xobs(iobs).lt.xmax)).and. &
159 & ((ymin.le.yobs(iobs)).and.(yobs(iobs).lt.ymax))) THEN
160 i1=int(xobs(iobs))
161 j1=int(yobs(iobs))
162 i2=i1+1
163 j2=j1+1
164 IF (i2.gt.imax) THEN
165 i2=i1 ! Observation at the eastern boundary
166 END IF
167 IF (j2.gt.jmax) THEN
168 j2=j1 ! Observation at the northern boundary
169 END IF
170 p2=real(i2-i1,r8)*(xobs(iobs)-real(i1,r8))
171 q2=real(j2-j1,r8)*(yobs(iobs)-real(j1,r8))
172 p1=1.0_r8-p2
173 q1=1.0_r8-q2
174 w11=p1*q1
175 w21=p2*q1
176 w22=p2*q2
177 w12=p1*q2
178 IF (zobs(iobs).gt.0.0_r8) THEN
179 IF (abs(real(ubk,r8)-zobs(iobs)).lt.1.0e-8_r8) THEN
180 linterpolate=.false. ! surface observation
181 obs_depths(iobs)=0.0_r8
182 ELSE
183 linterpolate=.true. ! fractional level
184 k1=max(lbk,int(zobs(iobs)-0.5_r8)) ! W-point
185 k2=min(k1+1,ubk)
186 r2=real(k2-k1,r8)*(zobs(iobs)-real(k1,r8))
187 r1=1.0_r8-r2
188 END IF
189 ELSE
190 linterpolate=.false. ! already depths in meters
191 obs_depths(iobs)=zobs(iobs)
192 END IF
193 IF (linterpolate) THEN
194 IF ((r1+r2).gt.0.0_r8) THEN
195 hmat(1)=w11*r1
196 hmat(2)=w21*r1
197 hmat(3)=w22*r1
198 hmat(4)=w12*r1
199 hmat(5)=w11*r2
200 hmat(6)=w21*r2
201 hmat(7)=w22*r2
202 hmat(8)=w12*r2
203 obs_depths(iobs)=hmat(1)*z(i1,j1,k1)+ &
204 & hmat(2)*z(i2,j1,k1)+ &
205 & hmat(3)*z(i2,j2,k1)+ &
206 & hmat(4)*z(i1,j2,k1)+ &
207 & hmat(5)*z(i1,j1,k2)+ &
208 & hmat(6)*z(i2,j1,k2)+ &
209 & hmat(7)*z(i2,j2,k2)+ &
210 & hmat(8)*z(i1,j2,k2)
211 END IF
212 END IF
213 END IF
214 END DO
215
216 RETURN
217 END SUBROUTINE obs_k2z
218#endif
219 END MODULE obs_k2z_mod
subroutine obs_k2z(ng, imin, imax, jmin, jmax, lbi, ubi, lbj, ubj, lbk, ubk, xmin, xmax, ymin, ymax, mobs, xobs, yobs, zobs, obs_scale, z, obs_depths)
Definition obs_k2z.F:112