1 |
gforget |
1.12 |
C $Header: /u/gcmpack/MITgcm/pkg/ecco/cost_gencost_sshv4.F,v 1.11 2013/02/05 14:51:09 gforget Exp $ |
2 |
gforget |
1.1 |
C $Name: $ |
3 |
|
|
|
4 |
jmc |
1.8 |
#include "ECCO_OPTIONS.h" |
5 |
gforget |
1.1 |
|
6 |
|
|
|
7 |
|
|
subroutine cost_gencost_sshv4( |
8 |
|
|
I myiter, |
9 |
|
|
I mytime, |
10 |
|
|
I mythid |
11 |
|
|
& ) |
12 |
|
|
|
13 |
|
|
c ================================================================== |
14 |
|
|
c SUBROUTINE cost_gencost_sshv4 |
15 |
|
|
c ================================================================== |
16 |
|
|
c |
17 |
|
|
c o Evaluate cost function contributions of sea surface height. |
18 |
|
|
c |
19 |
|
|
c started: Gael Forget, Oct-2009 |
20 |
|
|
c |
21 |
|
|
c working assumption for the time mean dynamic topography (MDT) constraint: |
22 |
|
|
c the various SLA data sets (tp, ers, gfo) have been consistenty cross-referenced, |
23 |
|
|
c as done in the RADS data sets. We do not need to know the reference dynamic |
24 |
|
|
c topography (or SSH/Geoid). Simply it is assumed that the biases |
25 |
|
|
c between instruments have been taken care of. This is only a matter |
26 |
|
|
c for the MDT constraint, not for SLA constraints (see below). |
27 |
|
|
c |
28 |
|
|
cgf 1) there are a few hardcoded numbers that will eventually be put in common |
29 |
|
|
cgf blocks/namelists |
30 |
|
|
cgf 2) there are a several refinements that should be considered, such as |
31 |
|
|
cgf modulating weights with respect to numbers of samples |
32 |
|
|
c |
33 |
|
|
c ================================================================== |
34 |
|
|
c SUBROUTINE cost_gencost_sshv4 |
35 |
|
|
c ================================================================== |
36 |
|
|
|
37 |
|
|
implicit none |
38 |
|
|
|
39 |
|
|
c == global variables == |
40 |
|
|
|
41 |
|
|
#include "EEPARAMS.h" |
42 |
|
|
#include "SIZE.h" |
43 |
|
|
#include "PARAMS.h" |
44 |
|
|
#include "GRID.h" |
45 |
|
|
|
46 |
|
|
#include "ecco_cost.h" |
47 |
heimbach |
1.7 |
#include "CTRL_SIZE.h" |
48 |
gforget |
1.1 |
#include "ctrl.h" |
49 |
|
|
#include "ctrl_dummy.h" |
50 |
|
|
#include "optim.h" |
51 |
|
|
#include "DYNVARS.h" |
52 |
|
|
#ifdef ALLOW_PROFILES |
53 |
|
|
#include "profiles.h" |
54 |
|
|
#endif |
55 |
gforget |
1.4 |
#include "cal.h" |
56 |
gforget |
1.1 |
|
57 |
|
|
c == routine arguments == |
58 |
|
|
|
59 |
|
|
integer myiter |
60 |
|
|
_RL mytime |
61 |
|
|
integer mythid |
62 |
|
|
|
63 |
|
|
#ifdef ALLOW_SSH_COST_CONTRIBUTION |
64 |
|
|
#ifdef ALLOW_GENCOST_CONTRIBUTION |
65 |
|
|
c == local variables == |
66 |
|
|
|
67 |
|
|
integer bi,bj |
68 |
|
|
integer i,j,k |
69 |
|
|
integer itlo,ithi |
70 |
|
|
integer jtlo,jthi |
71 |
|
|
integer jmin,jmax |
72 |
|
|
integer imin,imax |
73 |
|
|
integer irec,jrec,krec |
74 |
|
|
integer ilps |
75 |
|
|
integer gwunit |
76 |
|
|
|
77 |
|
|
logical doglobalread |
78 |
|
|
logical ladinit |
79 |
|
|
|
80 |
|
|
c mapping to gencost |
81 |
|
|
integer igen_mdt, igen_lsc |
82 |
jmc |
1.8 |
integer igen_tp, igen_ers, igen_gfo |
83 |
gforget |
1.1 |
|
84 |
|
|
_RL offset,fac |
85 |
|
|
_RL offset_sum |
86 |
|
|
_RL psmean ( 1-olx:snx+olx, 1-oly:sny+oly, nsx, nsy ) |
87 |
|
|
_RL diagnosfld ( 1-olx:snx+olx, 1-oly:sny+oly, nsx, nsy ) |
88 |
|
|
|
89 |
|
|
c for PART 1: re-reference MDT (tpmean) to the inferred SLA reference field |
90 |
|
|
_RL mean_slaobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
91 |
|
|
_RL mean_slaobs_NUM(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
92 |
|
|
|
93 |
|
|
c for PART 2: compute time mean differences over the model period |
94 |
|
|
_RL mean_slaobs2(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
95 |
|
|
_RL mean_psMssh_all(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
96 |
|
|
_RL mean_psMssh_all_NUM(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
97 |
|
|
_RL mean_psMssh_all_MSK(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
98 |
|
|
|
99 |
|
|
_RL mean_psMtpobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
100 |
|
|
_RL mean_psMtpobs_NUM(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
101 |
|
|
_RL mean_psMtpobs_MSK(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
102 |
|
|
|
103 |
|
|
_RL mean_psMersobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
104 |
|
|
_RL mean_psMersobs_NUM(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
105 |
|
|
_RL mean_psMersobs_MSK(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
106 |
|
|
|
107 |
|
|
_RL mean_psMgfoobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
108 |
|
|
_RL mean_psMgfoobs_NUM(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
109 |
|
|
_RL mean_psMgfoobs_MSK(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
110 |
|
|
|
111 |
|
|
c for PART 4/5: compute smooth/raw anomalies |
112 |
|
|
_RL anom_psMslaobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
113 |
|
|
_RL anom_slaobs (1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
114 |
|
|
_RL anom_psMslaobs_NUM (1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
115 |
|
|
|
116 |
|
|
_RL anom_psMtpobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
117 |
|
|
_RL anom_psMtpobs_NUM (1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
118 |
|
|
_RL anom_tpobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
119 |
|
|
|
120 |
|
|
_RL anom_psMersobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
121 |
|
|
_RL anom_psMersobs_NUM(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
122 |
|
|
_RL anom_ersobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
123 |
|
|
|
124 |
|
|
_RL anom_psMgfoobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
125 |
|
|
_RL anom_psMgfoobs_NUM(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
126 |
|
|
_RL anom_gfoobs(1-olx:snx+olx,1-oly:sny+oly,nsx,nsy) |
127 |
|
|
|
128 |
|
|
integer tpmean_y0,tpmean_y1,year,day |
129 |
|
|
integer num_var |
130 |
|
|
|
131 |
|
|
_RL junk,junkweight |
132 |
|
|
|
133 |
|
|
integer ndaysave |
134 |
|
|
_RL ndaysaveRL |
135 |
|
|
|
136 |
|
|
character*(80) fname |
137 |
|
|
character*(80) fname4test |
138 |
|
|
character*(MAX_LEN_MBUF) msgbuf |
139 |
|
|
|
140 |
gforget |
1.4 |
LOGICAL doReference |
141 |
|
|
|
142 |
gforget |
1.1 |
c == external functions == |
143 |
|
|
|
144 |
|
|
integer ilnblnk |
145 |
|
|
external ilnblnk |
146 |
|
|
|
147 |
|
|
c == end of interface == |
148 |
|
|
|
149 |
|
|
jtlo = mybylo(mythid) |
150 |
|
|
jthi = mybyhi(mythid) |
151 |
|
|
itlo = mybxlo(mythid) |
152 |
|
|
ithi = mybxhi(mythid) |
153 |
|
|
jmin = 1 |
154 |
|
|
jmax = sny |
155 |
|
|
imin = 1 |
156 |
|
|
imax = snx |
157 |
|
|
|
158 |
|
|
c-- detect the relevant gencost indices |
159 |
|
|
igen_mdt=0 |
160 |
|
|
igen_tp =0 |
161 |
|
|
igen_ers=0 |
162 |
|
|
igen_gfo=0 |
163 |
|
|
igen_lsc=0 |
164 |
|
|
do k=1,NGENCOST |
165 |
|
|
if (gencost_name(k).EQ.'sshv4-mdt') igen_mdt=k |
166 |
|
|
if (gencost_name(k).EQ.'sshv4-tp') igen_tp=k |
167 |
|
|
if (gencost_name(k).EQ.'sshv4-ers') igen_ers=k |
168 |
|
|
if (gencost_name(k).EQ.'sshv4-gfo') igen_gfo=k |
169 |
|
|
if (gencost_name(k).EQ.'sshv4-lsc') igen_lsc=k |
170 |
|
|
enddo |
171 |
|
|
|
172 |
|
|
c-- First, read tiled data. |
173 |
|
|
doglobalread = .false. |
174 |
|
|
ladinit = .false. |
175 |
|
|
|
176 |
|
|
write(fname(1:80),'(80a)') ' ' |
177 |
|
|
ilps=ilnblnk( psbarfile ) |
178 |
|
|
write(fname(1:80),'(2a,i10.10)') |
179 |
|
|
& psbarfile(1:ilps),'.',optimcycle |
180 |
|
|
|
181 |
|
|
cgf ======================================================= |
182 |
|
|
cgf PART 1: |
183 |
|
|
cgf x Get the MDT (tpmean) ... compute the sample mean |
184 |
|
|
cgf (mean_slaobs) of the SLA data (i.e. RADS for tp, ers, and gfo |
185 |
|
|
cgf together) over the time interval of the MDT ... subtract |
186 |
|
|
cgf mean_slaobs from tpmean. |
187 |
|
|
cgf x At this point, tpmean is the inferred SLA reference field. |
188 |
|
|
cgf x From there on, tpmean+sla will be directly comparable to |
189 |
|
|
cgf the model SSH (psbar). |
190 |
|
|
cgf ======================================================= |
191 |
|
|
|
192 |
|
|
c-- Read mean field and mask |
193 |
|
|
call cost_ReadTopexMean( mythid ) |
194 |
|
|
|
195 |
|
|
c-- Compute mean_slaobs: sample mean SLA over the time period of tpmean. |
196 |
|
|
|
197 |
|
|
c pavlis and ecco/rio |
198 |
|
|
tpmean_y0=1993 |
199 |
|
|
tpmean_y1=2004 |
200 |
|
|
c maximenko |
201 |
|
|
c tpmean_y0=1992 |
202 |
|
|
c tpmean_y1=2002 |
203 |
|
|
c rio |
204 |
|
|
c tpmean_y0=1993 |
205 |
|
|
c tpmean_y1=1999 |
206 |
|
|
|
207 |
gforget |
1.4 |
doReference=.FALSE. |
208 |
gforget |
1.12 |
if ((modelstartdate(1).GT.1992*10000).AND. |
209 |
|
|
& (modelstartdate(1).LT.2011*10000).AND. |
210 |
gforget |
1.4 |
& (ndaysrec.GE.365)) doReference=.TRUE. |
211 |
|
|
|
212 |
|
|
write(msgbuf,'(a,l)') ' sshv4:re-reference MDT=',doReference |
213 |
|
|
call print_message( msgbuf, standardmessageunit, |
214 |
|
|
& SQUEEZE_RIGHT , mythid) |
215 |
|
|
|
216 |
gforget |
1.1 |
do bj = jtlo,jthi |
217 |
|
|
do bi = itlo,ithi |
218 |
|
|
do j = jmin,jmax |
219 |
|
|
do i = imin,imax |
220 |
|
|
mean_slaobs(i,j,bi,bj) = 0. _d 0 |
221 |
|
|
mean_slaobs_NUM(i,j,bi,bj) = 0. _d 0 |
222 |
|
|
enddo |
223 |
|
|
enddo |
224 |
|
|
enddo |
225 |
|
|
enddo |
226 |
|
|
|
227 |
|
|
do year=tpmean_y0,tpmean_y1 |
228 |
|
|
do day=1,366 |
229 |
|
|
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
230 |
gforget |
1.5 |
call cost_sla_read_yd( topexfile, topexstartdate, |
231 |
gforget |
1.1 |
& tpobs, tpmask, |
232 |
|
|
& year, day, mythid ) |
233 |
|
|
#endif |
234 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
235 |
gforget |
1.5 |
call cost_sla_read_yd( ersfile, ersstartdate, |
236 |
gforget |
1.1 |
& ersobs, ersmask, |
237 |
|
|
& year, day, mythid ) |
238 |
|
|
#endif |
239 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
240 |
gforget |
1.5 |
call cost_sla_read_yd( gfofile, gfostartdate, |
241 |
gforget |
1.1 |
& gfoobs, gfomask, |
242 |
|
|
& year, day, mythid ) |
243 |
|
|
#endif |
244 |
|
|
do bj = jtlo,jthi |
245 |
|
|
do bi = itlo,ithi |
246 |
|
|
do j = jmin,jmax |
247 |
|
|
do i = imin,imax |
248 |
|
|
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
249 |
|
|
if ( tpmask(i,j,bi,bj)*tpmeanmask(i,j,bi,bj)* |
250 |
|
|
& gencost_weight(i,j,bi,bj,igen_tp) .NE. 0. ) then |
251 |
|
|
mean_slaobs(i,j,bi,bj)= mean_slaobs(i,j,bi,bj)+ |
252 |
|
|
& tpobs(i,j,bi,bj) |
253 |
|
|
mean_slaobs_NUM(i,j,bi,bj)= mean_slaobs_NUM(i,j,bi,bj)+1. _d 0 |
254 |
|
|
endif |
255 |
|
|
#endif |
256 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
257 |
|
|
if ( ersmask(i,j,bi,bj)*tpmeanmask(i,j,bi,bj)* |
258 |
|
|
& gencost_weight(i,j,bi,bj,igen_ers) .NE. 0. ) then |
259 |
|
|
mean_slaobs(i,j,bi,bj)= mean_slaobs(i,j,bi,bj)+ |
260 |
|
|
& ersobs(i,j,bi,bj) |
261 |
|
|
mean_slaobs_NUM(i,j,bi,bj)= mean_slaobs_NUM(i,j,bi,bj)+1. _d 0 |
262 |
|
|
endif |
263 |
|
|
#endif |
264 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
265 |
|
|
if ( gfomask(i,j,bi,bj)*tpmeanmask(i,j,bi,bj)* |
266 |
|
|
& gencost_weight(i,j,bi,bj,igen_gfo) .NE. 0. ) then |
267 |
|
|
mean_slaobs(i,j,bi,bj)= mean_slaobs(i,j,bi,bj)+ |
268 |
|
|
& gfoobs(i,j,bi,bj) |
269 |
|
|
mean_slaobs_NUM(i,j,bi,bj)= mean_slaobs_NUM(i,j,bi,bj)+1. _d 0 |
270 |
|
|
endif |
271 |
|
|
#endif |
272 |
|
|
enddo |
273 |
|
|
enddo |
274 |
|
|
enddo |
275 |
|
|
enddo |
276 |
|
|
|
277 |
|
|
enddo !do day=1,366 |
278 |
|
|
enddo !do year=tpmean_y0,tpmean_y1 |
279 |
|
|
|
280 |
|
|
do bj = jtlo,jthi |
281 |
|
|
do bi = itlo,ithi |
282 |
|
|
do j = jmin,jmax |
283 |
|
|
do i = imin,imax |
284 |
|
|
if ( ( mean_slaobs_NUM(i,j,bi,bj) .NE. 0. ).AND. |
285 |
|
|
& ( maskc(i,j,1,bi,bj) .NE. 0. ).AND. |
286 |
gforget |
1.11 |
#ifndef ALLOW_HIGHLAT_ALTIMETRY |
287 |
gforget |
1.1 |
& ( abs(YC(i,j,bi,bj)) .LE. 66. ).AND. |
288 |
gforget |
1.11 |
#endif |
289 |
gforget |
1.1 |
& ( tpmeanmask(i,j,bi,bj) .NE. 0. ) ) then |
290 |
|
|
mean_slaobs(i,j,bi,bj) = mean_slaobs(i,j,bi,bj) / |
291 |
|
|
& mean_slaobs_NUM(i,j,bi,bj) |
292 |
|
|
else |
293 |
|
|
mean_slaobs(i,j,bi,bj) = 0. _d 0 |
294 |
|
|
endif |
295 |
|
|
enddo |
296 |
|
|
enddo |
297 |
|
|
enddo |
298 |
|
|
enddo |
299 |
|
|
|
300 |
|
|
|
301 |
|
|
c-- smooth mean_slaobs: |
302 |
|
|
|
303 |
|
|
write(fname4test(1:80),'(1a)') 'sla2mdt_raw' |
304 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
305 |
|
|
& 1,mean_slaobs,1,1,mythid) |
306 |
|
|
|
307 |
gforget |
1.9 |
#ifdef ALLOW_SMOOTH |
308 |
|
|
if ( useSMOOTH ) |
309 |
|
|
& call smooth_hetero2d(mean_slaobs,maskc, |
310 |
gforget |
1.1 |
& gencost_scalefile(igen_mdt),300,mythid) |
311 |
gforget |
1.9 |
#endif |
312 |
gforget |
1.1 |
|
313 |
|
|
write(fname4test(1:80),'(1a)') 'sla2mdt_smooth' |
314 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
315 |
|
|
& 1,mean_slaobs,1,1,mythid) |
316 |
|
|
|
317 |
|
|
c-- re-reference tpmean: |
318 |
|
|
do bj = jtlo,jthi |
319 |
|
|
do bi = itlo,ithi |
320 |
|
|
do j = jmin,jmax |
321 |
|
|
do i = imin,imax |
322 |
|
|
if ( ( tpmeanmask(i,j,bi,bj) .NE. 0. ).AND. |
323 |
gforget |
1.4 |
& ( maskc(i,j,1,bi,bj) .NE. 0. ).AND. |
324 |
|
|
& ( doReference ) ) then |
325 |
gforget |
1.1 |
tpmean(i,j,bi,bj) = tpmean(i,j,bi,bj) |
326 |
|
|
& -mean_slaobs(i,j,bi,bj) |
327 |
|
|
endif |
328 |
|
|
enddo |
329 |
|
|
enddo |
330 |
|
|
enddo |
331 |
|
|
enddo |
332 |
|
|
|
333 |
|
|
|
334 |
|
|
cgf ======================================================= |
335 |
|
|
cgf PART 2: compute sample means of psbar-slaobs over the |
336 |
|
|
cgf period that is covered by the model (i.e. psbar). |
337 |
|
|
cgf x for all SLA data sets together: mean_psMssh_all, mean_psMssh_all_MSK, |
338 |
|
|
cgf and offset will be used in PART 3 (MDT cost term). |
339 |
|
|
cgf x for each SLA data individually. mean_psMtpobs, mean_psMtpobs_MS, etc. |
340 |
|
|
cgf will be used in PART 4&5 (SLA cost terms). |
341 |
|
|
cgf ======================================================= |
342 |
|
|
|
343 |
|
|
do bj = jtlo,jthi |
344 |
|
|
do bi = itlo,ithi |
345 |
|
|
do j = jmin,jmax |
346 |
|
|
do i = imin,imax |
347 |
|
|
|
348 |
|
|
psmean(i,j,bi,bj) = 0. _d 0 |
349 |
|
|
mean_psMtpobs(i,j,bi,bj) = 0. _d 0 |
350 |
|
|
mean_psMersobs(i,j,bi,bj) = 0. _d 0 |
351 |
|
|
mean_psMgfoobs(i,j,bi,bj) = 0. _d 0 |
352 |
|
|
mean_psMssh_all(i,j,bi,bj) = 0. _d 0 |
353 |
|
|
mean_slaobs2(i,j,bi,bj) = 0. _d 0 |
354 |
|
|
|
355 |
|
|
mean_psMtpobs_NUM(i,j,bi,bj) = 0. _d 0 |
356 |
|
|
mean_psMersobs_NUM(i,j,bi,bj) = 0. _d 0 |
357 |
|
|
mean_psMgfoobs_NUM(i,j,bi,bj) = 0. _d 0 |
358 |
|
|
mean_psMssh_all_NUM(i,j,bi,bj) = 0. _d 0 |
359 |
|
|
|
360 |
|
|
mean_psMtpobs_MSK(i,j,bi,bj) = 0. _d 0 |
361 |
|
|
mean_psMersobs_MSK(i,j,bi,bj) = 0. _d 0 |
362 |
|
|
mean_psMgfoobs_MSK(i,j,bi,bj) = 0. _d 0 |
363 |
|
|
|
364 |
|
|
enddo |
365 |
|
|
enddo |
366 |
|
|
enddo |
367 |
|
|
enddo |
368 |
|
|
offset = 0. _d 0 |
369 |
|
|
offset_sum = 0. _d 0 |
370 |
|
|
|
371 |
|
|
|
372 |
|
|
do irec = 1, ndaysrec |
373 |
|
|
|
374 |
|
|
call active_read_xy( fname, psbar, irec, doglobalread, |
375 |
|
|
& ladinit, optimcycle, mythid, |
376 |
|
|
& xx_psbar_mean_dummy ) |
377 |
|
|
|
378 |
gforget |
1.3 |
#ifndef ALLOW_PSBAR_MEAN |
379 |
|
|
CALL REMOVE_MEAN_RL( 1, psbar, maskInC, maskInC, rA, drF, |
380 |
|
|
& 'psbar', myTime, myThid ) |
381 |
|
|
#endif |
382 |
|
|
|
383 |
gforget |
1.1 |
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
384 |
|
|
call cost_sla_read( topexfile, topexstartdate, topexperiod, |
385 |
|
|
& topexintercept, topexslope, |
386 |
|
|
& tpobs, tpmask, |
387 |
|
|
& irec, mythid ) |
388 |
|
|
#endif |
389 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
390 |
|
|
call cost_sla_read( ersfile, ersstartdate, ersperiod, |
391 |
|
|
& ersintercept, ersslope, |
392 |
|
|
& ersobs, ersmask, |
393 |
|
|
& irec, mythid ) |
394 |
|
|
#endif |
395 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
396 |
|
|
call cost_sla_read( gfofile, gfostartdate, gfoperiod, |
397 |
|
|
& gfointercept, gfoslope, |
398 |
|
|
& gfoobs, gfomask, |
399 |
|
|
& irec, mythid ) |
400 |
|
|
#endif |
401 |
|
|
|
402 |
|
|
do bj = jtlo,jthi |
403 |
|
|
do bi = itlo,ithi |
404 |
|
|
do j = jmin,jmax |
405 |
|
|
do i = imin,imax |
406 |
|
|
psmean(i,j,bi,bj) = psmean(i,j,bi,bj) + |
407 |
|
|
& psbar(i,j,bi,bj) / float(ndaysrec) |
408 |
|
|
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
409 |
|
|
if ( tpmask(i,j,bi,bj)* |
410 |
|
|
& gencost_weight(i,j,bi,bj,igen_tp) .NE. 0. ) then |
411 |
|
|
mean_slaobs2(i,j,bi,bj)= |
412 |
|
|
& mean_slaobs2(i,j,bi,bj)+tpobs(i,j,bi,bj) |
413 |
|
|
mean_psMtpobs(i,j,bi,bj) = |
414 |
|
|
& mean_psMtpobs(i,j,bi,bj) + |
415 |
|
|
& psbar(i,j,bi,bj)-tpobs(i,j,bi,bj) |
416 |
|
|
mean_psMtpobs_NUM(i,j,bi,bj) = |
417 |
|
|
& mean_psMtpobs_NUM(i,j,bi,bj) + 1. _d 0 |
418 |
|
|
endif |
419 |
|
|
#endif |
420 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
421 |
|
|
if ( ersmask(i,j,bi,bj)* |
422 |
|
|
& gencost_weight(i,j,bi,bj,igen_ers) .NE. 0. ) then |
423 |
|
|
mean_slaobs2(i,j,bi,bj)= |
424 |
|
|
& mean_slaobs2(i,j,bi,bj)+ersobs(i,j,bi,bj) |
425 |
|
|
mean_psMersobs(i,j,bi,bj) = |
426 |
|
|
& mean_psMersobs(i,j,bi,bj) + |
427 |
|
|
& psbar(i,j,bi,bj)-ersobs(i,j,bi,bj) |
428 |
|
|
mean_psMersobs_NUM(i,j,bi,bj) = |
429 |
|
|
& mean_psMersobs_NUM(i,j,bi,bj) + 1. _d 0 |
430 |
|
|
endif |
431 |
|
|
#endif |
432 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
433 |
|
|
if ( gfomask(i,j,bi,bj)* |
434 |
|
|
& gencost_weight(i,j,bi,bj,igen_gfo) .NE. 0. ) then |
435 |
|
|
mean_slaobs2(i,j,bi,bj)= |
436 |
|
|
& mean_slaobs2(i,j,bi,bj)+gfoobs(i,j,bi,bj) |
437 |
|
|
mean_psMgfoobs(i,j,bi,bj) = |
438 |
|
|
& mean_psMgfoobs(i,j,bi,bj) + |
439 |
|
|
& psbar(i,j,bi,bj)-gfoobs(i,j,bi,bj) |
440 |
|
|
mean_psMgfoobs_NUM(i,j,bi,bj) = |
441 |
|
|
& mean_psMgfoobs_NUM(i,j,bi,bj) + 1. _d 0 |
442 |
|
|
endif |
443 |
|
|
#endif |
444 |
|
|
enddo |
445 |
|
|
enddo |
446 |
|
|
enddo |
447 |
|
|
enddo |
448 |
|
|
|
449 |
|
|
c-- END loop over records for the first time. |
450 |
|
|
enddo |
451 |
|
|
|
452 |
|
|
do bj = jtlo,jthi |
453 |
|
|
do bi = itlo,ithi |
454 |
|
|
do j = jmin,jmax |
455 |
|
|
do i = imin,imax |
456 |
|
|
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
457 |
gforget |
1.11 |
if ( ( mean_psMtpobs_NUM(i,j,bi,bj) .NE. 0. ) |
458 |
|
|
#ifndef ALLOW_HIGHLAT_ALTIMETRY |
459 |
|
|
& .AND.( abs(YC(i,j,bi,bj)) .LE. 66. ) |
460 |
|
|
#endif |
461 |
|
|
& ) then |
462 |
gforget |
1.1 |
mean_psMssh_all(i,j,bi,bj) = |
463 |
|
|
& mean_psMssh_all(i,j,bi,bj) + |
464 |
|
|
& mean_psMtpobs(i,j,bi,bj) |
465 |
|
|
mean_psMssh_all_NUM(i,j,bi,bj) = |
466 |
|
|
& mean_psMssh_all_NUM(i,j,bi,bj) + |
467 |
|
|
& mean_psMtpobs_NUM(i,j,bi,bj) |
468 |
|
|
mean_psMtpobs(i,j,bi,bj) = |
469 |
|
|
& mean_psMtpobs(i,j,bi,bj) / |
470 |
|
|
& mean_psMtpobs_NUM(i,j,bi,bj) |
471 |
|
|
mean_psMtpobs_MSK(i,j,bi,bj) = 1. _d 0 |
472 |
|
|
endif |
473 |
|
|
#endif |
474 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
475 |
gforget |
1.11 |
if ( ( mean_psMersobs_NUM(i,j,bi,bj) .NE. 0. ) |
476 |
|
|
#ifndef ALLOW_HIGHLAT_ALTIMETRY |
477 |
|
|
& .AND.( abs(YC(i,j,bi,bj)) .LE. 66. ) |
478 |
|
|
#endif |
479 |
|
|
& ) then |
480 |
gforget |
1.1 |
mean_psMssh_all(i,j,bi,bj) = |
481 |
|
|
& mean_psMssh_all(i,j,bi,bj) + |
482 |
|
|
& mean_psMersobs(i,j,bi,bj) |
483 |
|
|
mean_psMssh_all_NUM(i,j,bi,bj) = |
484 |
|
|
& mean_psMssh_all_NUM(i,j,bi,bj) + |
485 |
|
|
& mean_psMersobs_NUM(i,j,bi,bj) |
486 |
|
|
mean_psMersobs(i,j,bi,bj) = |
487 |
|
|
& mean_psMersobs(i,j,bi,bj) / |
488 |
|
|
& mean_psMersobs_NUM(i,j,bi,bj) |
489 |
|
|
mean_psMersobs_MSK(i,j,bi,bj) = 1. _d 0 |
490 |
|
|
endif |
491 |
|
|
#endif |
492 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
493 |
gforget |
1.11 |
if ( ( mean_psMgfoobs_NUM(i,j,bi,bj) .NE. 0. ) |
494 |
|
|
#ifndef ALLOW_HIGHLAT_ALTIMETRY |
495 |
|
|
& .AND.( abs(YC(i,j,bi,bj)) .LE. 66. ) |
496 |
|
|
#endif |
497 |
|
|
& ) then |
498 |
gforget |
1.1 |
mean_psMssh_all(i,j,bi,bj) = |
499 |
|
|
& mean_psMssh_all(i,j,bi,bj) + |
500 |
|
|
& mean_psMgfoobs(i,j,bi,bj) |
501 |
|
|
mean_psMssh_all_NUM(i,j,bi,bj) = |
502 |
|
|
& mean_psMssh_all_NUM(i,j,bi,bj) + |
503 |
|
|
& mean_psMgfoobs_NUM(i,j,bi,bj) |
504 |
|
|
mean_psMgfoobs(i,j,bi,bj) = |
505 |
|
|
& mean_psMgfoobs(i,j,bi,bj) / |
506 |
|
|
& mean_psMgfoobs_NUM(i,j,bi,bj) |
507 |
|
|
mean_psMgfoobs_MSK(i,j,bi,bj) = 1. _d 0 |
508 |
|
|
endif |
509 |
|
|
#endif |
510 |
|
|
if ( ( mean_psMssh_all_NUM(i,j,bi,bj) .NE. 0. ).AND. |
511 |
|
|
& ( maskc(i,j,1,bi,bj) .NE. 0. ) .AND. |
512 |
gforget |
1.11 |
#ifndef ALLOW_HIGHLAT_ALTIMETRY |
513 |
gforget |
1.1 |
& ( abs(YC(i,j,bi,bj)) .LE. 66. ).AND. |
514 |
gforget |
1.11 |
#endif |
515 |
gforget |
1.4 |
& ( tpmeanmask(i,j,bi,bj) .NE. 0. ).AND. |
516 |
|
|
& ( doReference ) ) then |
517 |
gforget |
1.1 |
mean_slaobs2(i,j,bi,bj) = |
518 |
|
|
& mean_slaobs2(i,j,bi,bj) / |
519 |
|
|
& mean_psMssh_all_NUM(i,j,bi,bj) |
520 |
|
|
mean_psMssh_all(i,j,bi,bj) = |
521 |
|
|
& mean_psMssh_all(i,j,bi,bj) / |
522 |
|
|
& mean_psMssh_all_NUM(i,j,bi,bj)-tpmean(i,j,bi,bj) |
523 |
|
|
mean_psMssh_all_MSK(i,j,bi,bj) = 1. _d 0 |
524 |
|
|
offset=offset+mean_psMssh_all(i,j,bi,bj)* |
525 |
|
|
& mean_psMssh_all_NUM(i,j,bi,bj) |
526 |
|
|
offset_sum=offset_sum+mean_psMssh_all_NUM(i,j,bi,bj) |
527 |
gforget |
1.4 |
elseif ( ( tpmeanmask(i,j,bi,bj) .NE. 0. ) .AND. |
528 |
|
|
& ( maskc(i,j,1,bi,bj) .NE. 0. ).AND. |
529 |
|
|
& ( .NOT.doReference ) ) then |
530 |
|
|
mean_slaobs2(i,j,bi,bj) = 0.d0 |
531 |
|
|
mean_psMssh_all(i,j,bi,bj) = 0. _d 0 |
532 |
|
|
mean_psMssh_all_MSK(i,j,bi,bj) = 1. _d 0 |
533 |
|
|
offset=offset+RA(i,j,bi,bj) |
534 |
|
|
& *( psmean(i,j,bi,bj) -tpmean(i,j,bi,bj) ) |
535 |
|
|
offset_sum=offset_sum+RA(i,j,bi,bj) |
536 |
gforget |
1.1 |
else |
537 |
|
|
mean_slaobs2(i,j,bi,bj) = 0.d0 |
538 |
|
|
mean_psMssh_all(i,j,bi,bj) = 0. _d 0 |
539 |
|
|
mean_psMssh_all_MSK(i,j,bi,bj) = 0. _d 0 |
540 |
|
|
endif |
541 |
|
|
enddo |
542 |
|
|
enddo |
543 |
|
|
enddo |
544 |
|
|
enddo |
545 |
|
|
|
546 |
|
|
c-- Do a global summation. |
547 |
|
|
_GLOBAL_SUM_RL( offset , mythid ) |
548 |
|
|
_GLOBAL_SUM_RL( offset_sum , mythid ) |
549 |
|
|
|
550 |
gforget |
1.4 |
write(msgbuf,'(a,2d22.15)') ' sshv4:offset=',offset,offset_sum |
551 |
|
|
call print_message( msgbuf, standardmessageunit, |
552 |
|
|
& SQUEEZE_RIGHT , mythid) |
553 |
|
|
|
554 |
gforget |
1.1 |
if (offset_sum .eq. 0.0) then |
555 |
|
|
_BEGIN_MASTER( mythid ) |
556 |
gforget |
1.5 |
write(msgbuf,'(a)') ' sshv4: offset_sum = zero!' |
557 |
gforget |
1.1 |
call print_message( msgbuf, standardmessageunit, |
558 |
|
|
& SQUEEZE_RIGHT , mythid) |
559 |
|
|
_END_MASTER( mythid ) |
560 |
|
|
c stop ' ... stopped in cost_ssh.' |
561 |
|
|
else |
562 |
|
|
_BEGIN_MASTER( mythid ) |
563 |
|
|
write(msgbuf,'(a,d22.15)') |
564 |
gforget |
1.5 |
& ' sshv4: offset_sum = ',offset_sum |
565 |
gforget |
1.1 |
call print_message( msgbuf, standardmessageunit, |
566 |
|
|
& SQUEEZE_RIGHT , mythid) |
567 |
|
|
_END_MASTER( mythid ) |
568 |
|
|
endif |
569 |
|
|
|
570 |
|
|
c-- Compute (average) offset |
571 |
|
|
offset = offset / offset_sum |
572 |
|
|
|
573 |
|
|
c-- subtract offset from mean_psMssh_all and psmean |
574 |
|
|
do bj = jtlo,jthi |
575 |
|
|
do bi = itlo,ithi |
576 |
|
|
do j = jmin,jmax |
577 |
|
|
do i = imin,imax |
578 |
|
|
|
579 |
|
|
if ( (mean_psMssh_all_MSK(i,j,bi,bj) .NE. 0.) .AND. |
580 |
gforget |
1.11 |
& ( maskc(i,j,1,bi,bj) .NE. 0. ) |
581 |
|
|
#ifndef ALLOW_HIGHLAT_ALTIMETRY |
582 |
|
|
& .AND.( abs(YC(i,j,bi,bj)) .LE. 66. ) |
583 |
|
|
#endif |
584 |
|
|
& ) then |
585 |
gforget |
1.1 |
c use the re-referencing approach |
586 |
|
|
mean_psMssh_all(i,j,bi,bj) = |
587 |
|
|
& mean_psMssh_all(i,j,bi,bj) - offset |
588 |
|
|
elseif ( ( tpmeanmask(i,j,bi,bj) .NE. 0. ) .AND. |
589 |
gforget |
1.11 |
& ( maskc(i,j,1,bi,bj) .NE. 0. ) |
590 |
|
|
#ifndef ALLOW_HIGHLAT_ALTIMETRY |
591 |
|
|
& .AND.( abs(YC(i,j,bi,bj)) .LE. 66. ) |
592 |
|
|
#endif |
593 |
|
|
& ) then |
594 |
gforget |
1.1 |
c use the simpler approach |
595 |
|
|
mean_psMssh_all(i,j,bi,bj) = |
596 |
|
|
& psmean(i,j,bi,bj) -tpmean(i,j,bi,bj) - offset |
597 |
|
|
else |
598 |
|
|
mean_psMssh_all(i,j,bi,bj) = 0. _d 0 |
599 |
|
|
endif |
600 |
|
|
|
601 |
gforget |
1.4 |
c use the simpler approach |
602 |
gforget |
1.12 |
if ( .NOT.doReference ) then |
603 |
|
|
if ( ( tpmeanmask(i,j,bi,bj) .NE. 0. ) .AND. |
604 |
gforget |
1.4 |
& ( maskc(i,j,1,bi,bj) .NE. 0. ) ) then |
605 |
|
|
mean_psMssh_all(i,j,bi,bj) = |
606 |
|
|
& psmean(i,j,bi,bj) -tpmean(i,j,bi,bj) - offset |
607 |
|
|
else |
608 |
|
|
mean_psMssh_all(i,j,bi,bj) = 0. _d 0 |
609 |
|
|
endif |
610 |
gforget |
1.12 |
endif |
611 |
gforget |
1.4 |
|
612 |
gforget |
1.1 |
if ( maskc(i,j,1,bi,bj) .NE. 0. ) |
613 |
|
|
& psmean(i,j,bi,bj)=psmean(i,j,bi,bj)-offset |
614 |
gforget |
1.12 |
|
615 |
gforget |
1.1 |
enddo |
616 |
|
|
enddo |
617 |
|
|
enddo |
618 |
|
|
enddo |
619 |
|
|
|
620 |
|
|
c-- smooth mean_psMssh_all |
621 |
|
|
write(fname4test(1:80),'(1a)') 'mdtdiff_raw' |
622 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
623 |
|
|
& 1,mean_psMssh_all,1,1,mythid) |
624 |
|
|
|
625 |
gforget |
1.9 |
#ifdef ALLOW_SMOOTH |
626 |
|
|
if ( useSMOOTH ) |
627 |
|
|
& call smooth_hetero2d(mean_psMssh_all,maskc, |
628 |
gforget |
1.1 |
& gencost_scalefile(igen_mdt),300,mythid) |
629 |
gforget |
1.9 |
#endif |
630 |
gforget |
1.1 |
|
631 |
|
|
write(fname4test(1:80),'(1a)') 'mdtdiff_smooth' |
632 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
633 |
|
|
& 1,mean_psMssh_all,1,1,mythid) |
634 |
|
|
|
635 |
|
|
write(fname4test(1:80),'(1a)') 'sla2model_raw' |
636 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
637 |
|
|
& 1,mean_slaobs2,1,1,mythid) |
638 |
|
|
|
639 |
gforget |
1.9 |
#ifdef ALLOW_SMOOTH |
640 |
|
|
if ( useSMOOTH ) |
641 |
|
|
& call smooth_hetero2d(mean_slaobs2,maskc, |
642 |
gforget |
1.1 |
& gencost_scalefile(igen_mdt),300,mythid) |
643 |
gforget |
1.9 |
#endif |
644 |
gforget |
1.1 |
|
645 |
|
|
write(fname4test(1:80),'(1a)') 'sla2model_smooth' |
646 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
647 |
|
|
& 1,mean_slaobs2,1,1,mythid) |
648 |
|
|
|
649 |
|
|
cgf at this point: |
650 |
|
|
cgf 1) mean_psMssh_all is the sample mean <psbar-slaobs-tpmean-offset>, |
651 |
|
|
cgf to which a smoothing filter has been applied. |
652 |
|
|
cgf 2) mean_psMtpobs is the (unsmoothed) sample mean <psbar-tpobs>. |
653 |
|
|
cgf And similarly for ers and gfo, each treated separately. |
654 |
|
|
|
655 |
|
|
#ifdef ALLOW_PROFILES |
656 |
gforget |
1.6 |
if ( usePROFILES ) then |
657 |
gforget |
1.1 |
do bj = jtlo,jthi |
658 |
|
|
do bi = itlo,ithi |
659 |
|
|
do j = 1,sny |
660 |
|
|
do i = 1,snx |
661 |
|
|
prof_etan_mean(i,j,bi,bj)=psmean(i,j,bi,bj) |
662 |
|
|
enddo |
663 |
|
|
enddo |
664 |
|
|
enddo |
665 |
|
|
enddo |
666 |
|
|
_EXCH_XY_RL( prof_etan_mean, mythid ) |
667 |
gforget |
1.6 |
endif |
668 |
gforget |
1.1 |
#endif |
669 |
|
|
|
670 |
|
|
|
671 |
|
|
cgf ======================================================= |
672 |
|
|
cgf PART 3: compute MDT cost term |
673 |
|
|
cgf ======================================================= |
674 |
|
|
|
675 |
|
|
|
676 |
|
|
#ifdef ALLOW_SSH_MEAN_COST_CONTRIBUTION |
677 |
|
|
|
678 |
|
|
do bj = jtlo,jthi |
679 |
|
|
do bi = itlo,ithi |
680 |
|
|
do j = jmin,jmax |
681 |
|
|
do i = imin,imax |
682 |
|
|
junk = mean_psMssh_all(i,j,bi,bj) |
683 |
|
|
junkweight = gencost_weight(i,j,bi,bj,igen_mdt) |
684 |
|
|
& *tpmeanmask(i,j,bi,bj) |
685 |
jmc |
1.8 |
objf_gencost(igen_mdt,bi,bj) = objf_gencost(igen_mdt,bi,bj) |
686 |
gforget |
1.1 |
& + junk*junk*junkweight |
687 |
jmc |
1.8 |
if ( junkweight .ne. 0. ) num_gencost(igen_mdt,bi,bj) = |
688 |
gforget |
1.1 |
& num_gencost(igen_mdt,bi,bj) + 1. _d 0 |
689 |
|
|
diagnosfld(i,j,bi,bj) = junk*junk*junkweight |
690 |
|
|
enddo |
691 |
|
|
enddo |
692 |
|
|
enddo |
693 |
|
|
enddo |
694 |
|
|
|
695 |
|
|
CALL WRITE_FLD_XY_RL( 'DiagnosSSHmean', ' ', diagnosfld, |
696 |
|
|
& optimcycle, mythid ) |
697 |
|
|
|
698 |
|
|
#endif /* ALLOW_SSH_MEAN_COST_CONTRIBUTION */ |
699 |
|
|
|
700 |
|
|
|
701 |
|
|
|
702 |
|
|
cgf ======================================================= |
703 |
|
|
cgf PART 4: compute smooth SLA cost term |
704 |
|
|
cgf ======================================================= |
705 |
|
|
|
706 |
|
|
|
707 |
|
|
ndaysave=35 |
708 |
|
|
ndaysaveRL=ndaysave |
709 |
|
|
|
710 |
gforget |
1.2 |
do irec = 1, ndaysrec-ndaysave+1, 7 |
711 |
gforget |
1.1 |
|
712 |
|
|
do bj = jtlo,jthi |
713 |
|
|
do bi = itlo,ithi |
714 |
|
|
do j = jmin,jmax |
715 |
|
|
do i = imin,imax |
716 |
|
|
anom_psMslaobs(i,j,bi,bj) = 0. _d 0 |
717 |
|
|
anom_slaobs(i,j,bi,bj) = 0. _d 0 |
718 |
|
|
anom_psMtpobs(i,j,bi,bj) = 0. _d 0 |
719 |
|
|
anom_psMersobs(i,j,bi,bj) = 0. _d 0 |
720 |
|
|
anom_psMgfoobs(i,j,bi,bj) = 0. _d 0 |
721 |
|
|
anom_psMslaobs_NUM(i,j,bi,bj) = 0. _d 0 |
722 |
|
|
anom_psMtpobs_NUM(i,j,bi,bj) = 0. _d 0 |
723 |
|
|
anom_psMersobs_NUM(i,j,bi,bj) = 0. _d 0 |
724 |
|
|
anom_psMgfoobs_NUM(i,j,bi,bj) = 0. _d 0 |
725 |
|
|
enddo |
726 |
|
|
enddo |
727 |
|
|
enddo |
728 |
|
|
enddo |
729 |
|
|
|
730 |
|
|
c PART 4.1: compute running sample average over ndaysave |
731 |
|
|
c ------------------------------------------------------ |
732 |
|
|
|
733 |
|
|
do jrec=1,ndaysave |
734 |
|
|
|
735 |
|
|
krec=irec+jrec-1 |
736 |
|
|
|
737 |
|
|
call active_read_xy( fname, psbar, krec, doglobalread, |
738 |
|
|
& ladinit, optimcycle, mythid, |
739 |
|
|
& xx_psbar_mean_dummy ) |
740 |
|
|
|
741 |
gforget |
1.3 |
#ifndef ALLOW_PSBAR_MEAN |
742 |
|
|
CALL REMOVE_MEAN_RL( 1, psbar, maskInC, maskInC, rA, drF, |
743 |
|
|
& 'psbar', myTime, myThid ) |
744 |
|
|
#endif |
745 |
|
|
|
746 |
gforget |
1.1 |
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
747 |
|
|
call cost_sla_read( topexfile, topexstartdate, topexperiod, |
748 |
|
|
& topexintercept, topexslope, |
749 |
|
|
& tpobs, tpmask, |
750 |
|
|
& krec, mythid ) |
751 |
|
|
#endif |
752 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
753 |
|
|
call cost_sla_read( ersfile, ersstartdate, ersperiod, |
754 |
|
|
& ersintercept, ersslope, |
755 |
|
|
& ersobs, ersmask, |
756 |
|
|
& krec, mythid ) |
757 |
|
|
#endif |
758 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
759 |
|
|
call cost_sla_read( gfofile, gfostartdate, gfoperiod, |
760 |
|
|
& gfointercept, gfoslope, |
761 |
|
|
& gfoobs, gfomask, |
762 |
|
|
& krec, mythid ) |
763 |
|
|
#endif |
764 |
|
|
|
765 |
|
|
do bj = jtlo,jthi |
766 |
|
|
do bi = itlo,ithi |
767 |
|
|
do j = jmin,jmax |
768 |
|
|
do i = imin,imax |
769 |
|
|
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
770 |
|
|
if ( tpmask(i,j,bi,bj)*mean_psMtpobs_MSK(i,j,bi,bj) |
771 |
|
|
& .NE.0. ) then |
772 |
|
|
anom_psMtpobs(i,j,bi,bj)= anom_psMtpobs(i,j,bi,bj)+ |
773 |
|
|
& psbar(i,j,bi,bj)-tpobs(i,j,bi,bj) |
774 |
|
|
& -mean_psMtpobs(i,j,bi,bj) |
775 |
|
|
anom_psMslaobs(i,j,bi,bj)= anom_psMslaobs(i,j,bi,bj)+ |
776 |
|
|
& psbar(i,j,bi,bj)-tpobs(i,j,bi,bj) |
777 |
|
|
& -mean_psMtpobs(i,j,bi,bj) |
778 |
|
|
anom_slaobs(i,j,bi,bj)= anom_slaobs(i,j,bi,bj)+ |
779 |
|
|
& tpobs(i,j,bi,bj) |
780 |
|
|
anom_psMtpobs_NUM(i,j,bi,bj)= |
781 |
|
|
& anom_psMtpobs_NUM(i,j,bi,bj)+1. _d 0 |
782 |
|
|
anom_psMslaobs_NUM(i,j,bi,bj)= |
783 |
|
|
& anom_psMslaobs_NUM(i,j,bi,bj)+1. _d 0 |
784 |
|
|
endif |
785 |
|
|
#endif |
786 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
787 |
|
|
if ( ersmask(i,j,bi,bj)*mean_psMersobs_MSK(i,j,bi,bj) |
788 |
|
|
& .NE.0. ) then |
789 |
|
|
anom_psMersobs(i,j,bi,bj)= anom_psMersobs(i,j,bi,bj)+ |
790 |
|
|
& psbar(i,j,bi,bj)-ersobs(i,j,bi,bj) |
791 |
|
|
& -mean_psMersobs(i,j,bi,bj) |
792 |
|
|
anom_psMersobs_NUM(i,j,bi,bj)= |
793 |
|
|
& anom_psMersobs_NUM(i,j,bi,bj)+1. _d 0 |
794 |
|
|
anom_psMslaobs(i,j,bi,bj)= anom_psMslaobs(i,j,bi,bj)+ |
795 |
|
|
& psbar(i,j,bi,bj)-ersobs(i,j,bi,bj) |
796 |
|
|
& -mean_psMersobs(i,j,bi,bj) |
797 |
|
|
anom_slaobs(i,j,bi,bj)= anom_slaobs(i,j,bi,bj)+ |
798 |
|
|
& ersobs(i,j,bi,bj) |
799 |
|
|
anom_psMslaobs_NUM(i,j,bi,bj)= |
800 |
|
|
& anom_psMslaobs_NUM(i,j,bi,bj)+1. _d 0 |
801 |
|
|
endif |
802 |
|
|
#endif |
803 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
804 |
|
|
if ( gfomask(i,j,bi,bj)*mean_psMgfoobs_MSK(i,j,bi,bj) |
805 |
|
|
& .NE.0. ) then |
806 |
|
|
anom_psMgfoobs(i,j,bi,bj)= anom_psMgfoobs(i,j,bi,bj)+ |
807 |
|
|
& psbar(i,j,bi,bj)-gfoobs(i,j,bi,bj) |
808 |
|
|
& -mean_psMgfoobs(i,j,bi,bj) |
809 |
|
|
anom_psMgfoobs_NUM(i,j,bi,bj)= |
810 |
|
|
& anom_psMgfoobs_NUM(i,j,bi,bj)+1. _d 0 |
811 |
|
|
anom_psMslaobs(i,j,bi,bj)= anom_psMslaobs(i,j,bi,bj)+ |
812 |
|
|
& psbar(i,j,bi,bj)-gfoobs(i,j,bi,bj) |
813 |
|
|
& -mean_psMgfoobs(i,j,bi,bj) |
814 |
|
|
anom_slaobs(i,j,bi,bj)= anom_slaobs(i,j,bi,bj)+ |
815 |
|
|
& gfoobs(i,j,bi,bj) |
816 |
|
|
anom_psMslaobs_NUM(i,j,bi,bj)= |
817 |
|
|
& anom_psMslaobs_NUM(i,j,bi,bj)+1. _d 0 |
818 |
|
|
endif |
819 |
|
|
#endif |
820 |
|
|
enddo |
821 |
|
|
enddo |
822 |
|
|
enddo |
823 |
|
|
enddo |
824 |
|
|
|
825 |
|
|
enddo !do jrec=1,ndaysave |
826 |
|
|
|
827 |
|
|
do bj = jtlo,jthi |
828 |
|
|
do bi = itlo,ithi |
829 |
|
|
do j = jmin,jmax |
830 |
|
|
do i = imin,imax |
831 |
|
|
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
832 |
|
|
if ( anom_psMtpobs_NUM(i,j,bi,bj) .NE. 0. ) then |
833 |
|
|
anom_psMtpobs(i,j,bi,bj) = |
834 |
|
|
& anom_psMtpobs(i,j,bi,bj) / |
835 |
|
|
& anom_psMtpobs_NUM(i,j,bi,bj) |
836 |
|
|
endif |
837 |
|
|
#endif |
838 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
839 |
|
|
if ( anom_psMersobs_NUM(i,j,bi,bj) .NE. 0. ) then |
840 |
|
|
anom_psMersobs(i,j,bi,bj) = |
841 |
|
|
& anom_psMersobs(i,j,bi,bj) / |
842 |
|
|
& anom_psMersobs_NUM(i,j,bi,bj) |
843 |
|
|
endif |
844 |
|
|
#endif |
845 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
846 |
|
|
if ( anom_psMgfoobs_NUM(i,j,bi,bj) .NE. 0. ) then |
847 |
|
|
anom_psMgfoobs(i,j,bi,bj) = |
848 |
|
|
& anom_psMgfoobs(i,j,bi,bj) / |
849 |
|
|
& anom_psMgfoobs_NUM(i,j,bi,bj) |
850 |
|
|
endif |
851 |
|
|
#endif |
852 |
|
|
if ( ( anom_psMslaobs_NUM(i,j,bi,bj) .NE. 0. ).AND. |
853 |
gforget |
1.11 |
& ( maskc(i,j,1,bi,bj) .NE. 0. ) |
854 |
|
|
#ifndef ALLOW_HIGHLAT_ALTIMETRY |
855 |
|
|
& .AND.( abs(YC(i,j,bi,bj)) .LE. 66. ) |
856 |
|
|
#endif |
857 |
|
|
& ) then |
858 |
gforget |
1.1 |
anom_psMslaobs(i,j,bi,bj) = |
859 |
|
|
& anom_psMslaobs(i,j,bi,bj) / |
860 |
|
|
& anom_psMslaobs_NUM(i,j,bi,bj) |
861 |
|
|
anom_slaobs(i,j,bi,bj) = |
862 |
|
|
& anom_slaobs(i,j,bi,bj) / |
863 |
|
|
& anom_psMslaobs_NUM(i,j,bi,bj) |
864 |
|
|
else |
865 |
|
|
anom_psMslaobs(i,j,bi,bj) = 0. _d 0 |
866 |
|
|
anom_slaobs(i,j,bi,bj) = 0. _d 0 |
867 |
|
|
endif |
868 |
|
|
enddo |
869 |
|
|
enddo |
870 |
|
|
enddo |
871 |
|
|
enddo |
872 |
|
|
|
873 |
|
|
c PART 4.2: smooth anom_psMslaobs in space |
874 |
|
|
c ---------------------------------------- |
875 |
|
|
|
876 |
|
|
#ifdef ALLOW_GENCOST_SSHV4_OUTPUT |
877 |
|
|
write(fname4test(1:80),'(1a)') 'sladiff_raw' |
878 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
879 |
|
|
& 1,anom_psMslaobs,irec,1,mythid) |
880 |
|
|
|
881 |
|
|
write(fname4test(1:80),'(1a)') 'slaobs_raw' |
882 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
883 |
|
|
& 1,anom_slaobs,irec,1,mythid) |
884 |
|
|
#endif |
885 |
|
|
|
886 |
gforget |
1.9 |
#ifdef ALLOW_SMOOTH |
887 |
|
|
if ( useSMOOTH ) |
888 |
|
|
& call smooth_hetero2d(anom_psMslaobs,maskc, |
889 |
gforget |
1.1 |
& gencost_scalefile(igen_lsc),300,mythid) |
890 |
gforget |
1.9 |
#endif |
891 |
gforget |
1.1 |
|
892 |
|
|
#ifdef ALLOW_GENCOST_SSHV4_OUTPUT |
893 |
gforget |
1.9 |
#ifdef ALLOW_SMOOTH |
894 |
|
|
if ( useSMOOTH ) |
895 |
|
|
& call smooth_hetero2d(anom_slaobs,maskc, |
896 |
gforget |
1.1 |
& gencost_scalefile(igen_lsc),300,mythid) |
897 |
gforget |
1.9 |
#endif |
898 |
gforget |
1.1 |
|
899 |
|
|
write(fname4test(1:80),'(1a)') 'sladiff_smooth' |
900 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
901 |
|
|
& 1,anom_psMslaobs,irec,1,mythid) |
902 |
|
|
|
903 |
|
|
write(fname4test(1:80),'(1a)') 'slaobs_smooth' |
904 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
905 |
|
|
& 1,anom_slaobs,irec,1,mythid) |
906 |
|
|
#endif |
907 |
|
|
|
908 |
|
|
c PART 4.3: compute cost function term |
909 |
|
|
c ------------------------------------ |
910 |
|
|
|
911 |
|
|
do bj = jtlo,jthi |
912 |
|
|
do bi = itlo,ithi |
913 |
|
|
do j = jmin,jmax |
914 |
|
|
do i = imin,imax |
915 |
|
|
# if (defined (ALLOW_SSH_GFOANOM_COST_CONTRIBUTION) || \ |
916 |
|
|
defined (ALLOW_SSH_TPANOM_COST_CONTRIBUTION) || \ |
917 |
|
|
defined (ALLOW_SSH_ERSANOM_COST_CONTRIBUTION)) |
918 |
|
|
junk = anom_psMslaobs(i,j,bi,bj) |
919 |
|
|
objf_gencost(igen_lsc,bi,bj) = objf_gencost(igen_lsc,bi,bj) |
920 |
|
|
& + junk*junk*gencost_weight(i,j,bi,bj,igen_lsc) |
921 |
|
|
& *maskc(i,j,1,bi,bj)/ndaysaveRL |
922 |
|
|
if ( (gencost_weight(i,j,bi,bj,igen_lsc).GT.0.).AND. |
923 |
|
|
& (anom_psMslaobs_NUM(i,j,bi,bj).GT.0.).AND. |
924 |
|
|
& (maskc(i,j,1,bi,bj) .ne. 0.) ) |
925 |
jmc |
1.8 |
& num_gencost(igen_lsc,bi,bj) = |
926 |
gforget |
1.1 |
& num_gencost(igen_lsc,bi,bj) + 1. _d 0 /ndaysaveRL |
927 |
|
|
#endif |
928 |
|
|
enddo |
929 |
|
|
enddo |
930 |
|
|
enddo |
931 |
|
|
enddo |
932 |
|
|
|
933 |
|
|
enddo |
934 |
|
|
|
935 |
|
|
|
936 |
|
|
cgf ======================================================= |
937 |
|
|
cgf PART 5: compute raw SLA cost term |
938 |
|
|
cgf ======================================================= |
939 |
|
|
|
940 |
|
|
|
941 |
|
|
do irec = 1, ndaysrec |
942 |
|
|
|
943 |
|
|
call active_read_xy( fname, psbar, irec, doglobalread, |
944 |
|
|
& ladinit, optimcycle, mythid, |
945 |
|
|
& xx_psbar_mean_dummy ) |
946 |
|
|
|
947 |
gforget |
1.3 |
#ifndef ALLOW_PSBAR_MEAN |
948 |
|
|
CALL REMOVE_MEAN_RL( 1, psbar, maskInC, maskInC, rA, drF, |
949 |
|
|
& 'psbar', myTime, myThid ) |
950 |
|
|
#endif |
951 |
|
|
|
952 |
gforget |
1.1 |
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
953 |
|
|
call cost_readtopex( irec, mythid ) |
954 |
|
|
#endif |
955 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
956 |
|
|
call cost_readers( irec, mythid ) |
957 |
|
|
#endif |
958 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
959 |
|
|
call cost_readgfo( irec, mythid ) |
960 |
|
|
#endif |
961 |
|
|
|
962 |
|
|
do bj = jtlo,jthi |
963 |
|
|
do bi = itlo,ithi |
964 |
|
|
do j = jmin,jmax |
965 |
|
|
do i = imin,imax |
966 |
|
|
anom_psMtpobs(i,j,bi,bj) = 0. _d 0 |
967 |
|
|
anom_psMersobs(i,j,bi,bj) = 0. _d 0 |
968 |
|
|
anom_psMgfoobs(i,j,bi,bj) = 0. _d 0 |
969 |
|
|
anom_tpobs(i,j,bi,bj) = 0. _d 0 |
970 |
|
|
anom_ersobs(i,j,bi,bj) = 0. _d 0 |
971 |
|
|
anom_gfoobs(i,j,bi,bj) = 0. _d 0 |
972 |
|
|
enddo |
973 |
|
|
enddo |
974 |
|
|
enddo |
975 |
|
|
enddo |
976 |
|
|
|
977 |
|
|
do bj = jtlo,jthi |
978 |
|
|
do bi = itlo,ithi |
979 |
|
|
do j = jmin,jmax |
980 |
|
|
do i = imin,imax |
981 |
|
|
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
982 |
|
|
if ( tpmask(i,j,bi,bj)*mean_psMtpobs_MSK(i,j,bi,bj).NE.0. ) |
983 |
|
|
& then |
984 |
|
|
anom_psMtpobs(i,j,bi,bj)= |
985 |
|
|
& psbar(i,j,bi,bj) - tpobs(i,j,bi,bj) |
986 |
|
|
& - mean_psMtpobs(i,j,bi,bj) |
987 |
|
|
anom_tpobs(i,j,bi,bj)=tpobs(i,j,bi,bj) |
988 |
|
|
endif |
989 |
|
|
#endif |
990 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
991 |
|
|
if ( ersmask(i,j,bi,bj)*mean_psMersobs_MSK(i,j,bi,bj).NE.0. ) |
992 |
|
|
& then |
993 |
|
|
anom_psMersobs(i,j,bi,bj)= |
994 |
|
|
& psbar(i,j,bi,bj) - ersobs(i,j,bi,bj) |
995 |
|
|
& - mean_psMersobs(i,j,bi,bj) |
996 |
|
|
anom_ersobs(i,j,bi,bj)=ersobs(i,j,bi,bj) |
997 |
|
|
endif |
998 |
|
|
#endif |
999 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
1000 |
|
|
if ( gfomask(i,j,bi,bj)*mean_psMgfoobs_MSK(i,j,bi,bj).NE.0. ) |
1001 |
|
|
& then |
1002 |
|
|
anom_psMgfoobs(i,j,bi,bj)= |
1003 |
|
|
& psbar(i,j,bi,bj) - gfoobs(i,j,bi,bj) |
1004 |
|
|
& - mean_psMgfoobs(i,j,bi,bj) |
1005 |
|
|
anom_gfoobs(i,j,bi,bj)=gfoobs(i,j,bi,bj) |
1006 |
|
|
endif |
1007 |
|
|
#endif |
1008 |
|
|
enddo |
1009 |
|
|
enddo |
1010 |
|
|
enddo |
1011 |
|
|
enddo |
1012 |
|
|
|
1013 |
|
|
#ifdef ALLOW_GENCOST_SSHV4_OUTPUT |
1014 |
|
|
write(fname4test(1:80),'(1a)') 'sladiff_tp_raw' |
1015 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
1016 |
|
|
& 1,anom_psMtpobs,irec,1,mythid) |
1017 |
|
|
write(fname4test(1:80),'(1a)') 'sladiff_ers_raw' |
1018 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
1019 |
|
|
& 1,anom_psMersobs,irec,1,mythid) |
1020 |
|
|
write(fname4test(1:80),'(1a)') 'sladiff_gfo_raw' |
1021 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
1022 |
|
|
& 1,anom_psMgfoobs,irec,1,mythid) |
1023 |
|
|
|
1024 |
|
|
write(fname4test(1:80),'(1a)') 'slaobs_tp_raw' |
1025 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
1026 |
|
|
& 1,anom_tpobs,irec,1,mythid) |
1027 |
|
|
write(fname4test(1:80),'(1a)') 'slaobs_ers_raw' |
1028 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
1029 |
|
|
& 1,anom_ersobs,irec,1,mythid) |
1030 |
|
|
write(fname4test(1:80),'(1a)') 'slaobs_gfo_raw' |
1031 |
|
|
call mdswritefield(fname4test,32,.false.,'RL', |
1032 |
|
|
& 1,anom_gfoobs,irec,1,mythid) |
1033 |
jmc |
1.8 |
#endif |
1034 |
gforget |
1.1 |
|
1035 |
|
|
do bj = jtlo,jthi |
1036 |
|
|
do bi = itlo,ithi |
1037 |
|
|
do j = jmin,jmax |
1038 |
|
|
do i = imin,imax |
1039 |
|
|
#ifdef ALLOW_SSH_TPANOM_COST_CONTRIBUTION |
1040 |
|
|
c-- The array psobs contains SSH anomalies. |
1041 |
|
|
junkweight = mean_psMtpobs_MSK(i,j,bi,bj) |
1042 |
|
|
& *gencost_weight(i,j,bi,bj,igen_tp) |
1043 |
|
|
& *tpmask(i,j,bi,bj) |
1044 |
|
|
& *cosphi(i,j,bi,bj) |
1045 |
|
|
junk = anom_psMtpobs(i,j,bi,bj) |
1046 |
jmc |
1.8 |
objf_gencost(igen_tp,bi,bj) = |
1047 |
gforget |
1.1 |
& objf_gencost(igen_tp,bi,bj)+junk*junk*junkweight |
1048 |
|
|
if ( junkweight .ne. 0. ) |
1049 |
|
|
& num_gencost(igen_tp,bi,bj) = |
1050 |
|
|
& num_gencost(igen_tp,bi,bj) + 1. _d 0 |
1051 |
|
|
#endif |
1052 |
|
|
#ifdef ALLOW_SSH_ERSANOM_COST_CONTRIBUTION |
1053 |
|
|
c-- The array ersobs contains SSH anomalies. |
1054 |
|
|
junkweight = mean_psMersobs_MSK(i,j,bi,bj) |
1055 |
|
|
& *gencost_weight(i,j,bi,bj,igen_ers) |
1056 |
|
|
& *ersmask(i,j,bi,bj) |
1057 |
|
|
& *cosphi(i,j,bi,bj) |
1058 |
|
|
junk = anom_psMersobs(i,j,bi,bj) |
1059 |
|
|
objf_gencost(igen_ers,bi,bj) = |
1060 |
|
|
& objf_gencost(igen_ers,bi,bj)+junk*junk*junkweight |
1061 |
|
|
if ( junkweight .ne. 0. ) |
1062 |
|
|
& num_gencost(igen_ers,bi,bj) = |
1063 |
|
|
& num_gencost(igen_ers,bi,bj) + 1. _d 0 |
1064 |
|
|
#endif |
1065 |
|
|
#ifdef ALLOW_SSH_GFOANOM_COST_CONTRIBUTION |
1066 |
|
|
c-- The array gfoobs contains SSH anomalies. |
1067 |
|
|
junkweight = mean_psMgfoobs_MSK(i,j,bi,bj) |
1068 |
|
|
& *gencost_weight(i,j,bi,bj,igen_gfo) |
1069 |
|
|
& *gfomask(i,j,bi,bj) |
1070 |
|
|
& *cosphi(i,j,bi,bj) |
1071 |
|
|
junk = anom_psMgfoobs(i,j,bi,bj) |
1072 |
|
|
objf_gencost(igen_gfo,bi,bj) = |
1073 |
|
|
& objf_gencost(igen_gfo,bi,bj)+junk*junk*junkweight |
1074 |
|
|
if ( junkweight .ne. 0. ) |
1075 |
|
|
& num_gencost(igen_gfo,bi,bj) = |
1076 |
|
|
& num_gencost(igen_gfo,bi,bj) + 1. _d 0 |
1077 |
|
|
#endif |
1078 |
|
|
enddo |
1079 |
|
|
enddo |
1080 |
|
|
enddo |
1081 |
|
|
enddo |
1082 |
|
|
|
1083 |
|
|
enddo |
1084 |
|
|
|
1085 |
|
|
|
1086 |
|
|
#endif /* ifdef ALLOW_GENCOST_CONTRIBUTION */ |
1087 |
|
|
#endif /* ifdef ALLOW_SSH_COST_CONTRIBUTION */ |
1088 |
|
|
|
1089 |
|
|
end |