/[MITgcm]/MITgcm/pkg/ecco/cost_ctdsclim.F
ViewVC logotype

Contents of /MITgcm/pkg/ecco/cost_ctdsclim.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.4 - (show annotations) (download)
Mon Oct 11 16:38:53 2004 UTC (19 years, 7 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint58l_post, checkpoint57t_post, checkpoint57o_post, checkpoint58e_post, checkpoint57v_post, checkpoint58u_post, checkpoint57f_post, checkpoint57s_post, checkpoint57k_post, checkpoint57d_post, checkpoint57g_post, checkpoint57b_post, checkpoint57c_pre, checkpoint58r_post, checkpoint55j_post, checkpoint56b_post, checkpoint57i_post, checkpoint57y_post, checkpoint58g_post, checkpoint57x_post, checkpoint57m_post, checkpoint55h_post, checkpoint58n_post, checkpoint58x_post, checkpoint57g_pre, checkpoint58t_post, checkpoint58h_post, checkpoint57e_post, checkpoint58w_post, checkpoint56c_post, checkpoint57y_pre, checkpoint57f_pre, checkpoint57a_post, checkpoint58q_post, checkpoint55g_post, checkpoint58j_post, checkpoint59e, checkpoint59d, checkpoint59g, checkpoint59f, checkpoint59a, checkpoint55f_post, checkpoint59c, checkpoint59b, checkpoint59h, checkpoint57r_post, checkpoint59, checkpoint58, checkpoint57a_pre, checkpoint55i_post, checkpoint57, checkpoint56, eckpoint57e_pre, checkpoint57h_done, checkpoint58f_post, checkpoint57n_post, checkpoint58d_post, checkpoint58c_post, checkpoint57w_post, checkpoint57p_post, checkpint57u_post, checkpoint58a_post, checkpoint58i_post, checkpoint57q_post, checkpoint58o_post, checkpoint57z_post, checkpoint57c_post, checkpoint58y_post, checkpoint55e_post, checkpoint58k_post, checkpoint58v_post, checkpoint58s_post, checkpoint58p_post, checkpoint57j_post, checkpoint58b_post, checkpoint57h_pre, checkpoint58m_post, checkpoint57l_post, checkpoint57h_post, checkpoint56a_post, checkpoint55d_post
Changes since 1.3: +2 -22 lines
o ECCO specific cost function terms (up-to-date with 1x1 runs)
o ecco_cost_weights is modified to 1x1 runs
o modifs to allow observations to be read in as
  single file or yearly files

1 C $Header: /u/gcmpack/MITgcm/pkg/cost/Attic/cost_salt.F,v 1.1.2.3 2002/04/04 10:58:59 heimbach Exp $
2
3 #include "COST_CPPOPTIONS.h"
4
5
6 subroutine cost_ctdsclim(
7 I myiter,
8 I mytime,
9 I mythid
10 & )
11
12 c ==================================================================
13 c SUBROUTINE cost_Ctdsclim
14 c ==================================================================
15 c
16 c o Evaluate cost function contribution of salinity.
17 c
18 c started: Christian Eckert eckert@mit.edu 30-Jun-1999
19 c
20 c changed: Christian Eckert eckert@mit.edu 25-Feb-2000
21 c
22 c - Restructured the code in order to create a package
23 c for the MITgcmUV.
24 c
25 c changed: Patrick Heimbach heimbach@mit.edu 27-May-2000
26 c
27 c - set ladinit to .true. to initialise adsbar file
28 c
29 c ==================================================================
30 c SUBROUTINE cost_Ctdsclim
31 c ==================================================================
32
33 implicit none
34
35 c == global variables ==
36
37 #include "EEPARAMS.h"
38 #include "SIZE.h"
39 #include "GRID.h"
40 #include "DYNVARS.h"
41
42 #include "cal.h"
43 #include "ecco_cost.h"
44 #include "ctrl.h"
45 #include "ctrl_dummy.h"
46 #include "optim.h"
47
48 c == routine arguments ==
49
50 integer myiter
51 _RL mytime
52 integer mythid
53
54 c == local variables ==
55
56 _RS one_rs
57 parameter( one_rs = 1. )
58
59 integer bi,bj
60 integer i,j,k
61 integer itlo,ithi
62 integer jtlo,jthi
63 integer jmin,jmax
64 integer imin,imax
65 integer irec
66 integer levmon
67 integer levoff
68 integer ilctdsclim
69
70 _RL fctile
71 _RL fcthread
72
73 _RL cmask (1-olx:snx+olx,1-oly:sny+oly)
74 _RL spval
75
76 character*(80) fnamesalt
77
78 logical doglobalread
79 logical ladinit
80
81 character*(MAX_LEN_MBUF) msgbuf
82
83 #ifdef GENERIC_BAR_MONTH
84 integer mrec, nyears, iyear
85 #endif
86 c == external functions ==
87
88 integer ilnblnk
89 external ilnblnk
90
91 c == end of interface ==
92
93 jtlo = mybylo(mythid)
94 jthi = mybyhi(mythid)
95 itlo = mybxlo(mythid)
96 ithi = mybxhi(mythid)
97 jmin = 1
98 jmax = sny
99 imin = 1
100 imax = snx
101
102 spval = -9990.
103
104 c-- Read tiled data.
105 doglobalread = .false.
106 ladinit = .false.
107
108 #ifdef ALLOW_CTDSCLIM_COST_CONTRIBUTION
109
110 if (optimcycle .ge. 0) then
111 ilctdsclim = ilnblnk( sbarfile )
112 write(fnamesalt(1:80),'(2a,i10.10)')
113 & sbarfile(1:ilctdsclim),'.',optimcycle
114 endif
115
116 fcthread = 0. _d 0
117
118 #ifdef GENERIC_BAR_MONTH
119 c-- Loop over month
120 do irec = 1,12
121 nyears=int((nmonsrec-irec)/12)+1
122 if(nyears.gt.0) then
123 do iyear=1,nyears
124 mrec=irec+(iyear-1)*12
125 c-- Read time averages and the monthly mean data.
126 call active_read_xyz( fnamesalt, sbar, mrec,
127 & doglobalread, ladinit,
128 & optimcycle, mythid,
129 & xx_sbar_mean_dummy )
130 do bj = jtlo,jthi
131 do bi = itlo,ithi
132 do k = 1,nr
133 do j = jmin,jmax
134 do i = imin,imax
135 if(iyear.eq.1) then
136 sbar_gen(i,j,k,bi,bj) =sbar(i,j,k,bi,bj)
137 elseif(iyear.eq.nyears) then
138 sbar(i,j,k,bi,bj) =(sbar_gen(i,j,k,bi,bj)
139 $ +sbar(i,j,k,bi,bj))/float(nyears)
140 else
141 sbar_gen(i,j,k,bi,bj) =sbar_gen(i,j,k,bi,bj)
142 $ +sbar(i,j,k,bi,bj)
143 endif
144 enddo
145 enddo
146 enddo
147 enddo
148 enddo
149 enddo
150 #else
151 c-- Loop over records.
152 do irec = 1,nmonsrec
153
154 c-- Read time averages and the monthly mean data.
155 call active_read_xyz( fnamesalt, sbar, irec,
156 & doglobalread, ladinit,
157 & optimcycle, mythid,
158 & xx_sbar_mean_dummy )
159 #endif
160 c-- Determine the month to be read.
161 levoff = mod(modelstartdate(1)/100,100)
162 levmon = (irec-1) + levoff
163 levmon = mod(levmon-1,12)+1
164
165 call mdsreadfield( ctdsclimfile, cost_iprec, cost_yftype,
166 & nr, sdat, levmon, mythid)
167
168 do bj = jtlo,jthi
169 do bi = itlo,ithi
170
171 c-- Loop over the model layers
172 fctile = 0. _d 0
173 do k = 1,nr
174
175 c-- Determine the mask or weights
176 do j = jmin,jmax
177 do i = imin,imax
178 cmask(i,j) = 1. _d 0
179 if (sdat(i,j,k,bi,bj) .eq. 0.) then
180 cmask(i,j) = 0. _d 0
181 endif
182 if (sdat(i,j,k,bi,bj) .lt. spval) then
183 cmask(i,j) = 0. _d 0
184 endif
185
186 cph(
187 cph print *, 'WARNING: SPECIFIC SETUP FOR ECCO'
188 cph below statement could be replaced by following
189 cph to make it independnet of Nr:
190 cph
191 cph if ( rC(K) .GT. -1000. ) then
192 cph)
193 c set cmask=0 in areas shallower than 1000m
194 if (_hFacC(i,j,13,bi,bj) .eq. 0.) then
195 cmask(i,j) = 0. _d 0
196 endif
197
198 enddo
199 enddo
200
201 c-- Compute model data misfit and cost function term for
202 c the salinity field.
203 do j = jmin,jmax
204 do i = imin,imax
205 if (_hFacC(i,j,k,bi,bj) .ne. 0.) then
206 fctile = fctile +
207 & (wsalt2(i,j,k,bi,bj)*cosphi(i,j,bi,bj)*
208 & cmask(i,j)*
209 & (sbar(i,j,k,bi,bj) - sdat(i,j,k,bi,bj))*
210 & (sbar(i,j,k,bi,bj) - sdat(i,j,k,bi,bj)) )
211 endif
212 enddo
213 enddo
214
215 enddo
216 c-- End of loop over layers.
217
218 fcthread = fcthread + fctile
219 objf_ctdsclim(bi,bj) = objf_ctdsclim(bi,bj) + fctile
220
221 #ifdef ECCO_VERBOSE
222 c-- Print cost function for each tile in each thread.
223 write(msgbuf,'(a)') ' '
224 call print_message( msgbuf, standardmessageunit,
225 & SQUEEZE_RIGHT , mythid)
226 write(msgbuf,'(a,i8.8,1x,i3.3,1x,i3.3)')
227 & ' cost_Ctdsclim: irec,bi,bj = ',irec,bi,bj
228 call print_message( msgbuf, standardmessageunit,
229 & SQUEEZE_RIGHT , mythid)
230 write(msgbuf,'(a,d22.15)')
231 & ' cost function (salinity) = ',
232 & fctile
233 call print_message( msgbuf, standardmessageunit,
234 & SQUEEZE_RIGHT , mythid)
235 write(msgbuf,'(a)') ' '
236 call print_message( msgbuf, standardmessageunit,
237 & SQUEEZE_RIGHT , mythid)
238 #endif
239
240 enddo
241 enddo
242
243 #ifdef ECCO_VERBOSE
244 c-- Print cost function for all tiles.
245 _GLOBAL_SUM_R8( fcthread , myThid )
246 write(msgbuf,'(a)') ' '
247 call print_message( msgbuf, standardmessageunit,
248 & SQUEEZE_RIGHT , mythid)
249 write(msgbuf,'(a,i8.8)')
250 & ' cost_Ctdsclim: irec = ',irec
251 call print_message( msgbuf, standardmessageunit,
252 & SQUEEZE_RIGHT , mythid)
253 write(msgbuf,'(a,a,d22.15)')
254 & ' global cost function value',
255 & ' (salinity) = ',fcthread
256 call print_message( msgbuf, standardmessageunit,
257 & SQUEEZE_RIGHT , mythid)
258 write(msgbuf,'(a)') ' '
259 call print_message( msgbuf, standardmessageunit,
260 & SQUEEZE_RIGHT , mythid)
261 #endif
262
263 #ifdef GENERIC_BAR_MONTH
264 endif
265 #endif
266 enddo
267 c-- End of loop over records.
268
269 #else
270 c-- Do not enter the calculation of the salinity contribution to
271 c-- the final cost function.
272
273 _BEGIN_MASTER( mythid )
274 write(msgbuf,'(a)') ' '
275 call print_message( msgbuf, standardmessageunit,
276 & SQUEEZE_RIGHT , mythid)
277 write(msgbuf,'(a,a)')
278 & ' cost_Ctdsclim: no contribution of salinity field ',
279 & 'to cost function.'
280 call print_message( msgbuf, standardmessageunit,
281 & SQUEEZE_RIGHT , mythid)
282 write(msgbuf,'(a,a,i9.8)')
283 & ' cost_Ctdsclim: number of records that would have',
284 & ' been processed: ',nmonsrec
285 call print_message( msgbuf, standardmessageunit,
286 & SQUEEZE_RIGHT , mythid)
287 write(msgbuf,'(a)') ' '
288 call print_message( msgbuf, standardmessageunit,
289 & SQUEEZE_RIGHT , mythid)
290 _END_MASTER( mythid )
291 #endif
292
293 return
294 end
295

  ViewVC Help
Powered by ViewVC 1.1.22