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