1 |
C $Header: /u/gcmpack/MITgcm/pkg/ecco/cost_gencost_seaicev4.F,v 1.19 2015/11/08 03:12:56 atn Exp $ |
2 |
C $Name: $ |
3 |
|
4 |
#include "ECCO_OPTIONS.h" |
5 |
|
6 |
subroutine cost_gencost_seaicev4(mythid) |
7 |
|
8 |
c ================================================================== |
9 |
c SUBROUTINE cost_gencost_seaicev4 |
10 |
c ================================================================== |
11 |
c |
12 |
c o Evaluate cost function contributions of ice concentration. |
13 |
c |
14 |
c ================================================================== |
15 |
c SUBROUTINE cost_gencost_seaicev4 |
16 |
c ================================================================== |
17 |
|
18 |
implicit none |
19 |
|
20 |
c == global variables == |
21 |
|
22 |
#include "EEPARAMS.h" |
23 |
#include "SIZE.h" |
24 |
#include "PARAMS.h" |
25 |
#include "GRID.h" |
26 |
#ifdef ALLOW_CAL |
27 |
# include "cal.h" |
28 |
#endif |
29 |
#ifdef ALLOW_ECCO |
30 |
# include "ecco.h" |
31 |
#endif |
32 |
#ifdef ALLOW_SEAICE |
33 |
# include "SEAICE_SIZE.h" |
34 |
# include "SEAICE_COST.h" |
35 |
# include "SEAICE_PARAMS.h" |
36 |
#endif |
37 |
|
38 |
c == routine arguments == |
39 |
integer mythid |
40 |
|
41 |
#ifdef ALLOW_SEAICE |
42 |
#ifdef ALLOW_GENCOST_CONTRIBUTION |
43 |
|
44 |
c == local variables == |
45 |
|
46 |
integer nnzsiv4, nnzbar |
47 |
parameter (nnzsiv4 = 1 , nnzbar = 1) |
48 |
integer nrecloc |
49 |
integer localstartdate(4) |
50 |
|
51 |
catn changing names to make more self-explanatory |
52 |
c old:sst -> model has deficiency in iceconc -> new:deconc |
53 |
c old:heff -> model has excess of iceconc -> new:exconc |
54 |
|
55 |
_RL areabar (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
56 |
_RL deconcbar (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
57 |
_RL exconcbar (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
58 |
_RL localweight (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
59 |
_RL dummy |
60 |
_RL localperiod |
61 |
_RL spminloc |
62 |
_RL spmaxloc |
63 |
_RL spzeroloc |
64 |
|
65 |
character*(MAX_LEN_FNAM) areabarfile |
66 |
character*(MAX_LEN_FNAM) deconcbarfile |
67 |
character*(MAX_LEN_FNAM) exconcbarfile |
68 |
character*(MAX_LEN_FNAM) localobsfile |
69 |
|
70 |
integer igen_conc, igen_deconc, igen_exconc |
71 |
|
72 |
integer bi,bj |
73 |
integer i,j,k |
74 |
integer itlo,ithi |
75 |
integer jtlo,jthi |
76 |
integer jmin,jmax |
77 |
integer imin,imax |
78 |
integer irec, jrec |
79 |
integer il, k2 |
80 |
integer localrec |
81 |
integer obsrec |
82 |
logical dosumsq, dovarwei |
83 |
|
84 |
integer preproc_i(NGENPPROC) |
85 |
_RL preproc_r(NGENPPROC) |
86 |
character*(MAX_LEN_FNAM) preproc(NGENPPROC) |
87 |
character*(MAX_LEN_FNAM) preproc_c(NGENPPROC) |
88 |
|
89 |
_RL localmask (1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy) |
90 |
|
91 |
_RL localobs (1-olx:snx+olx,1-oly:sny+oly,nnzsiv4,nsx,nsy) |
92 |
_RL localtmp (1-olx:snx+olx,1-oly:sny+oly,nnzsiv4,nsx,nsy) |
93 |
_RL localdif (1-olx:snx+olx,1-oly:sny+oly,nnzsiv4,nsx,nsy) |
94 |
_RL difmask (1-olx:snx+olx,1-oly:sny+oly,nnzsiv4,nsx,nsy) |
95 |
_RL difmask1 (1-olx:snx+olx,1-oly:sny+oly,nnzsiv4,nsx,nsy) |
96 |
|
97 |
character*(128) fname0, fname0w, fname1 |
98 |
|
99 |
character*(MAX_LEN_FNAM) localobswfile |
100 |
logical exst |
101 |
|
102 |
c == external functions == |
103 |
|
104 |
integer ilnblnk |
105 |
external ilnblnk |
106 |
|
107 |
c == end of interface == |
108 |
|
109 |
jtlo = mybylo(mythid) |
110 |
jthi = mybyhi(mythid) |
111 |
itlo = mybxlo(mythid) |
112 |
ithi = mybxhi(mythid) |
113 |
jmin = 1 |
114 |
jmax = sny |
115 |
imin = 1 |
116 |
imax = snx |
117 |
|
118 |
c=============== PART 0: initilization =================== |
119 |
|
120 |
c-- detect the relevant gencost indices |
121 |
igen_conc=0 |
122 |
igen_deconc=0 |
123 |
igen_exconc=0 |
124 |
do k=1,NGENCOST |
125 |
if (gencost_name(k).EQ.'siv4-conc') igen_conc=k |
126 |
if (gencost_name(k).EQ.'siv4-deconc') igen_deconc=k |
127 |
if (gencost_name(k).EQ.'siv4-exconc') igen_exconc=k |
128 |
enddo |
129 |
|
130 |
c-- Dependency: |
131 |
c A) igen_conc can exist on its own |
132 |
c B) igen_deconc needs igen_conc |
133 |
c C) igen_exconc needs both igen_conc and igen_deconc |
134 |
if (igen_conc.NE.0) then |
135 |
|
136 |
c-- initialize objf and num: |
137 |
do bj = jtlo,jthi |
138 |
do bi = itlo,ithi |
139 |
objf_gencost(bi,bj,igen_conc) = 0. _d 0 |
140 |
num_gencost(bi,bj,igen_conc) = 0. _d 0 |
141 |
if(igen_deconc.ne.0) then |
142 |
objf_gencost(bi,bj,igen_deconc) = 0. _d 0 |
143 |
num_gencost(bi,bj,igen_deconc) = 0. _d 0 |
144 |
endif |
145 |
if(igen_exconc.ne.0) then |
146 |
objf_gencost(bi,bj,igen_exconc) = 0. _d 0 |
147 |
num_gencost(bi,bj,igen_exconc) = 0. _d 0 |
148 |
endif |
149 |
enddo |
150 |
enddo |
151 |
|
152 |
c-- Initialise local variables. |
153 |
nrecloc=0 |
154 |
localperiod=0. |
155 |
|
156 |
areabarfile=gencost_barfile(igen_conc) |
157 |
if(igen_deconc.ne.0) deconcbarfile=gencost_barfile(igen_deconc) |
158 |
if(igen_exconc.ne.0) exconcbarfile=gencost_barfile(igen_exconc) |
159 |
|
160 |
localobsfile=gencost_datafile(igen_conc) |
161 |
localobswfile=gencost_errfile(igen_conc) |
162 |
dummy=gencost_dummy(igen_conc) |
163 |
localstartdate(1)=modelstartdate(1) |
164 |
localstartdate(2)=modelstartdate(2) |
165 |
localstartdate(3)=modelstartdate(3) |
166 |
localstartdate(4)=modelstartdate(4) |
167 |
spminloc=gencost_spmin(igen_conc) |
168 |
spmaxloc=gencost_spmax(igen_conc) |
169 |
spzeroloc=gencost_spzero(igen_conc) |
170 |
|
171 |
localperiod=gencost_period(igen_conc) |
172 |
nrecloc=gencost_nrec(igen_conc) |
173 |
|
174 |
c-- flag to add cost: true=(obs-mod)*(obs-mod)*weight |
175 |
dosumsq=.TRUE. |
176 |
dovarwei=.FALSE. |
177 |
do k2 = 1, NGENPPROC |
178 |
preproc(k2)=gencost_preproc(k2,igen_conc) |
179 |
preproc_i(k2)=gencost_preproc_i(k2,igen_conc) |
180 |
preproc_c(k2)=gencost_preproc_c(k2,igen_conc) |
181 |
preproc_r(k2)=gencost_preproc_r(k2,igen_conc) |
182 |
if (preproc(k2).EQ.'variaweight') dovarwei=.TRUE. |
183 |
if (preproc(k2).EQ.'nosumsq') dosumsq=.FALSE. |
184 |
enddo |
185 |
|
186 |
c-- initialize arrays, copy maskC to localmask |
187 |
call ecco_zero(localobs,1,spzeroloc,myThid) |
188 |
call ecco_zero(localweight,1,zeroRL,myThid) |
189 |
call ecco_zero(localmask,Nr,zeroRL,myThid) |
190 |
call ecco_cprsrl(maskC,Nr,localmask,Nr,myThid) |
191 |
|
192 |
c=============== PART 1: main loop =================== |
193 |
if ( .NOT. ( localobsfile.EQ.' ' ) ) then |
194 |
|
195 |
c-- Loop over records for the second time. |
196 |
do irec = 1, nrecloc |
197 |
|
198 |
c==================================================== |
199 |
c--------- PART 1.1 read weights -------------------- |
200 |
c==================================================== |
201 |
exst=.FALSE. |
202 |
jrec=1 |
203 |
if( dovarwei ) jrec = irec |
204 |
call cost_gencal(areabarfile,gencost_errfile(igen_conc), |
205 |
& jrec, localstartdate, localperiod, fname1, |
206 |
& fname0w, localrec, obsrec, exst, mythid) |
207 |
call ecco_zero(localweight,nnzsiv4,zeroRL,myThid) |
208 |
#ifdef SEAICECOST_JPL |
209 |
fname0w=gencost_errfile(igen_conc) |
210 |
call ecco_readwei(fname0w,localweight,localrec,nnzsiv4,mythid) |
211 |
call ecco_readwei(gencost_errfile(igen_deconc), |
212 |
& gencost_weight(1,1,1,1,igen_deconc),localrec,nnzsiv4,mythid) |
213 |
call ecco_readwei(gencost_errfile(igen_exconc), |
214 |
& gencost_weight(1,1,1,1,igen_exconc),localrec,nnzsiv4,mythid) |
215 |
#else |
216 |
if ( (localrec. GT. 0).AND.(obsrec .GT. 0).AND.(exst) ) then |
217 |
call ecco_readwei(fname0w,localweight,localrec,nnzsiv4,mythid) |
218 |
else |
219 |
WRITE(standardMessageUnit,'(A)') |
220 |
& 'siv4cost WARNING: ALL WEIGHTS ZEROS! NO CONTRIBUTION' |
221 |
endif |
222 |
#endif |
223 |
|
224 |
c==================================================== |
225 |
c--------- PART 1.2 read barfiles ------------------ |
226 |
c==================================================== |
227 |
c-- set all bars to zeros: |
228 |
call ecco_zero(areabar,nnzbar,zeroRL,myThid) |
229 |
call ecco_zero(deconcbar,nnzbar,zeroRL,myThid) |
230 |
call ecco_zero(exconcbar,nnzbar,zeroRL,myThid) |
231 |
|
232 |
c--1.2.A sea-ice concentration barfile |
233 |
exst=.FALSE. |
234 |
call cost_gencal(areabarfile,gencost_datafile(igen_conc), |
235 |
& irec,localstartdate,localperiod,fname1, |
236 |
& fname0,localrec,obsrec,exst,mythid) |
237 |
call cost_genread(fname1,areabar,localtmp,irec,nnzbar, |
238 |
& nrecloc,preproc,preproc_c,preproc_i,preproc_r, |
239 |
& dummy,mythid) |
240 |
|
241 |
c--1.2.B sst as proxy for deconc barfile, needs igen_conc |
242 |
if(igen_deconc.ne.0) then |
243 |
exst=.FALSE. |
244 |
call cost_gencal(deconcbarfile,gencost_datafile(igen_conc), |
245 |
& irec,localstartdate,localperiod,fname1, |
246 |
& fname0,localrec,obsrec,exst,mythid) |
247 |
call cost_genread(fname1,deconcbar,localtmp,irec,nnzbar, |
248 |
& nrecloc,preproc,preproc_c,preproc_i,preproc_r, |
249 |
& dummy,mythid) |
250 |
endif |
251 |
|
252 |
c--1.2.C heff as proxy for exconc barfile, need igen_conc and igen_exconc |
253 |
if(igen_deconc.ne.0 .and. igen_exconc.ne.0) then |
254 |
exst=.FALSE. |
255 |
call cost_gencal(exconcbarfile,gencost_datafile(igen_conc), |
256 |
& irec,localstartdate,localperiod,fname1, |
257 |
& fname0,localrec,obsrec,exst,mythid) |
258 |
call cost_genread(fname1,exconcbar,localtmp,irec,nnzbar, |
259 |
& nrecloc,preproc,preproc_c,preproc_i,preproc_r, |
260 |
& dummy,mythid) |
261 |
endif |
262 |
|
263 |
c==================================================== |
264 |
c--------- PART 1.3 read data -------------------- |
265 |
c==================================================== |
266 |
c-- initialize to spzerloc = -9999. |
267 |
call ecco_zero(localobs,nnzsiv4,spzeroloc,mythid) |
268 |
if ( (localrec .GT. 0).AND.(obsrec .GT. 0).AND.(exst) ) then |
269 |
call mdsreadfield( fname0, cost_iprec, cost_yftype, nnzsiv4, |
270 |
& localobs, localrec, mythid ) |
271 |
else |
272 |
il=ilnblnk( fname0 ) |
273 |
WRITE(standardMessageUnit,'(2A)') |
274 |
& 'siv4cost WARNING: DATA MISING! NO CONTRIBUTION, ', |
275 |
& fname0(1:il) |
276 |
endif |
277 |
|
278 |
c==================================================== |
279 |
c--------- PART 1.4 Cost calculation ------------- |
280 |
c==================================================== |
281 |
c compute obs minus bar (localdif) and mask (difmask) |
282 |
call ecco_zero(localdif,nnzsiv4,zeroRL,mythid) |
283 |
call ecco_zero(difmask,nnzsiv4,zeroRL,mythid) |
284 |
call ecco_diffmsk( |
285 |
I areabar, nnzbar, localobs, nnzsiv4, localmask, |
286 |
I spminloc, spmaxloc, spzeroloc, |
287 |
O localdif, difmask, |
288 |
I myThid ) |
289 |
|
290 |
c---1.4.A area term: |
291 |
call ecco_addcost( |
292 |
I localdif,localweight,difmask,nnzsiv4,dosumsq, |
293 |
O objf_gencost(1,1,igen_conc),num_gencost(1,1,igen_conc), |
294 |
I mythid) |
295 |
|
296 |
c---1.4.B defficient ice term: (old: sst term, new: deconc) |
297 |
c Add ice: model_A==0 but obs_A > 0, calc enthalpy E: |
298 |
if(igen_deconc.ne.0) then |
299 |
call ecco_zero(difmask1,nnzsiv4,zeroRL,mythid) |
300 |
call ecco_zero(localdif,nnzsiv4,zeroRL,mythid) |
301 |
call ecco_zero(localtmp,nnzsiv4,zeroRL,mythid) |
302 |
|
303 |
call get_exconc_deconc( |
304 |
I localobs,nnzsiv4,areabar,exconcbar,deconcbar,nnzbar, |
305 |
I difmask,'de', |
306 |
O localdif,difmask1,localtmp, |
307 |
I myThid ) |
308 |
#ifdef SEAICECOST_JPL |
309 |
call ecco_cp(gencost_weight(1,1,1,1,igen_deconc),nnzsiv4, |
310 |
O localtmp,nnzsiv4,myThid) |
311 |
#endif |
312 |
call ecco_addcost( |
313 |
I localdif,localtmp,difmask1,nnzsiv4,dosumsq, |
314 |
O objf_gencost(1,1,igen_deconc),num_gencost(1,1,igen_deconc), |
315 |
I mythid) |
316 |
endif |
317 |
|
318 |
c---1.4.C excessive ice term: (old: heff and sst term, new: exconc) |
319 |
c Removing ice: model_A > 0 but obs_A==0, calc enthalpy E: |
320 |
if(igen_deconc.ne.0 .and. igen_exconc.ne.0) then |
321 |
call ecco_zero(difmask1,nnzsiv4,zeroRL,mythid) |
322 |
call ecco_zero(localdif,nnzsiv4,zeroRL,mythid) |
323 |
call ecco_zero(localtmp,nnzsiv4,zeroRL,mythid) |
324 |
|
325 |
call get_exconc_deconc( |
326 |
I localobs,nnzsiv4,areabar,exconcbar,deconcbar,nnzbar, |
327 |
I difmask,'ex', |
328 |
O localdif,difmask1,localtmp, |
329 |
I myThid ) |
330 |
#ifdef SEAICECOST_JPL |
331 |
call ecco_cp(gencost_weight(1,1,1,1,igen_exconc),nnzsiv4, |
332 |
O localtmp,nnzsiv4,myThid) |
333 |
#endif |
334 |
call ecco_addcost( |
335 |
I localdif,localtmp,difmask1,nnzsiv4,dosumsq, |
336 |
O objf_gencost(1,1,igen_exconc),num_gencost(1,1,igen_exconc), |
337 |
I mythid) |
338 |
endif |
339 |
|
340 |
enddo |
341 |
|
342 |
endif !if ( .NOT. ( localobsfile.EQ.' ' ) ) then |
343 |
endif !if (igen_conc.NE.0) |
344 |
|
345 |
#endif /* ALLOW_GENCOST_CONTRIBUTION */ |
346 |
#endif /* ALLOW_SEAICE */ |
347 |
|
348 |
RETURN |
349 |
end |
350 |
|
351 |
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| |
352 |
|
353 |
#include "ECCO_OPTIONS.h" |
354 |
|
355 |
subroutine get_exconc_deconc( |
356 |
I localobs,nnzobs,concbar,exconcbar,deconcbar,nnzbar, |
357 |
I localmask,flag_exconc_deconc, |
358 |
O localfld,localfldmsk,localfldweight, |
359 |
I myThid ) |
360 |
|
361 |
C !DESCRIPTION: \bv |
362 |
c Routine to calculate Enthalpy for the case of |
363 |
c defficient/excessive model seaice |
364 |
C \ev |
365 |
|
366 |
C !USES: |
367 |
implicit none |
368 |
|
369 |
c == global variables == |
370 |
#include "EEPARAMS.h" |
371 |
#include "SIZE.h" |
372 |
#include "PARAMS.h" |
373 |
#include "GRID.h" |
374 |
#ifdef ALLOW_SEAICE |
375 |
# include "SEAICE_SIZE.h" |
376 |
# include "SEAICE_COST.h" |
377 |
# include "SEAICE_PARAMS.h" |
378 |
#endif |
379 |
#ifdef ALLOW_ECCO |
380 |
# include "ecco.h" |
381 |
#endif |
382 |
|
383 |
c == routine arguments == |
384 |
|
385 |
integer mythid, nnzbar, nnzobs |
386 |
_RL localmask (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
387 |
|
388 |
_RL localobs (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
389 |
_RL concbar (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
390 |
_RL deconcbar (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
391 |
_RL exconcbar (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
392 |
_RL localfld (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
393 |
_RL localfldmsk (1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
394 |
_RL localfldweight(1-olx:snx+olx,1-oly:sny+oly,1,nsx,nsy) |
395 |
|
396 |
character*2 flag_exconc_deconc |
397 |
|
398 |
#ifdef ALLOW_GENCOST_CONTRIBUTION |
399 |
#ifdef ALLOW_SEAICE |
400 |
|
401 |
c == local variables == |
402 |
integer bi,bj |
403 |
integer itlo,ithi |
404 |
integer jtlo,jthi |
405 |
integer jmin,jmax |
406 |
integer imin,imax |
407 |
integer i,j,k |
408 |
|
409 |
C- jmc: params SEAICE_freeze has been retired; set it locally until someone |
410 |
C who knows what this cost-cointribution means fix it. |
411 |
C- atn: also adding 1 normalizing factor same order of magnitude as |
412 |
C rhoNil*HeatCapacity_cp*dz = SEAICE_rhoice*SEAICE_lhFusion*heff |
413 |
C = 1e3*1e3*1e1=1e7 |
414 |
C- atn: lastly, define 2 cutoff values for cost to be read in from data.seaice |
415 |
C and initialized in seaice_readparms: SEAICE_cutoff_[area,heff] |
416 |
C Reason: some iceconc data set have "bogus" mask with area>0 |
417 |
C at climatological max locations -> not real data. So either need |
418 |
C to clean up the data or take SEAICE_cutoff_area>=0.15 for example. |
419 |
C Might need to migrate into pkg/ecco instead of in pkg/seaice. |
420 |
_RL SEAICE_freeze, epsilonTemp, epsilonHEFF |
421 |
_RL localnorm, localnormsq |
422 |
_RL const1,const2 |
423 |
CEOP |
424 |
|
425 |
SEAICE_freeze = -1.96 _d 0 |
426 |
epsilonTemp = 0.0001 _d 0 |
427 |
#ifdef SEAICECOST_JPL |
428 |
epsilonHEFF = 0.3 _d 0 |
429 |
#else |
430 |
epsilonHEFF = 0.01 _d 0 |
431 |
#endif |
432 |
localnorm = 1. _d -07 |
433 |
localnormsq=localnorm*localnorm |
434 |
|
435 |
const1=HeatCapacity_Cp*rhoNil*drF(1) |
436 |
const2=SEAICE_lhFusion*SEAICE_rhoIce |
437 |
|
438 |
jtlo = mybylo(mythid) |
439 |
jthi = mybyhi(mythid) |
440 |
itlo = mybxlo(mythid) |
441 |
ithi = mybxhi(mythid) |
442 |
jmin = 1-oly |
443 |
jmax = sny+oly |
444 |
imin = 1-olx |
445 |
imax = snx+olx |
446 |
|
447 |
c intialize |
448 |
call ecco_zero(localfld,nnzobs,zeroRL,myThid) |
449 |
call ecco_zero(localfldmsk,nnzobs,zeroRL,myThid) |
450 |
call ecco_zero(localfldweight,nnzobs,zeroRL,myThid) |
451 |
|
452 |
c----------------------DECONC------------------------------- |
453 |
catn-- old: sst term, new: deconc |
454 |
c needs localconcbar and localsstbar |
455 |
c Add ice: model_A==0 but obs_A > 0, calc enthalpy E: |
456 |
c E_current = (deconcbar(i,j,k,bi,bj)-Tfreeze) |
457 |
c *HeatCapacity_Cp*rhoNil*drF(1) |
458 |
c HEFF_target = epsilon_HEFF [m] |
459 |
c E_target = -(HEFF_target*SEAICE_lhFusion*SEAICE_rhoIce) |
460 |
c cost=(Model-data)^2 |
461 |
if(flag_exconc_deconc.EQ.'de') then |
462 |
do bj = jtlo,jthi |
463 |
do bi = itlo,ithi |
464 |
do k = 1,nnzobs |
465 |
do j = jmin,jmax |
466 |
do i = imin,imax |
467 |
|
468 |
if ( (concbar(i,j,k,bi,bj) .LE. 0.).AND. |
469 |
& (localobs(i,j,k,bi,bj) .GT. 0.) ) then |
470 |
|
471 |
localfldmsk(i,j,k,bi,bj) = localmask(i,j,k,bi,bj) |
472 |
|
473 |
#ifndef SEAICECOST_JPL |
474 |
localfldweight(i,j,k,bi,bj) = localnormsq |
475 |
#endif |
476 |
localfld(i,j,k,bi,bj) = |
477 |
& (deconcbar(i,j,k,bi,bj)-SEAICE_freeze)*const1 |
478 |
& - (-1. _d 0 *epsilonHEFF*const2) |
479 |
endif |
480 |
enddo |
481 |
enddo |
482 |
enddo |
483 |
enddo |
484 |
enddo |
485 |
endif |
486 |
|
487 |
c----------------------EXCONC------------------------------- |
488 |
catn-- old: heff and sst term, new: exconc |
489 |
c needs localconcbar, localsstbar, and localheffbar |
490 |
c Removing ice: model_A > 0 but obs_A==0, calc enthalpy E: |
491 |
c E_current = [(deconcbar-SEAICE_freeze)*HeatCapacity_Cp*rhoNil*drF(1) |
492 |
c - (exconcbar * SEAICE_lhFusion * SEAICE_rhoIce) |
493 |
c - (HSNOW * SEAICE_lhFusion * SEAICE_rhoSnow)] |
494 |
c E_target = (epsilonTemp) * HeatCapacity_Cp * rhoNil * drF(1) |
495 |
c cost(Model-data)^2 |
496 |
|
497 |
if(flag_exconc_deconc.EQ.'ex') then |
498 |
do bj = jtlo,jthi |
499 |
do bi = itlo,ithi |
500 |
do k = 1,nnzobs |
501 |
do j = jmin,jmax |
502 |
do i = imin,imax |
503 |
|
504 |
if ((localobs(i,j,k,bi,bj) .LE. SEAICE_cutoff_area).AND. |
505 |
& (exconcbar(i,j,k,bi,bj) .GT. SEAICE_cutoff_heff)) then |
506 |
|
507 |
localfldmsk(i,j,k,bi,bj) = localmask(i,j,k,bi,bj) |
508 |
#ifndef SEAICECOST_JPL |
509 |
localfldweight(i,j,k,bi,bj) = localnormsq |
510 |
#endif |
511 |
localfld(i,j,k,bi,bj) = |
512 |
& ( (deconcbar(i,j,k,bi,bj)-SEAICE_freeze)*const1 |
513 |
#ifdef SEAICECOST_JPL |
514 |
& - max(exconcbar(i,j,k,bi,bj),epsilonHEFF)* |
515 |
#else |
516 |
& - exconcbar(i,j,k,bi,bj)* |
517 |
#endif |
518 |
& const2 ) - (epsilonTemp*const1) |
519 |
endif |
520 |
enddo |
521 |
enddo |
522 |
enddo |
523 |
enddo |
524 |
enddo |
525 |
endif |
526 |
|
527 |
#endif /* ALLOW_GENCOST_CONTRIBUTION */ |
528 |
#endif /* ALLOW_SEAICE */ |
529 |
RETURN |
530 |
END |
531 |
|
532 |
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| |