/[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.5 - (show annotations) (download)
Tue Oct 9 00:02:50 2007 UTC (16 years, 7 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint60, checkpoint61, checkpoint59q, checkpoint59p, checkpoint59r, checkpoint59m, checkpoint59l, checkpoint59o, checkpoint59n, checkpoint59i, checkpoint59k, checkpoint59j, checkpoint61f, checkpoint61g, checkpoint61d, checkpoint61e, checkpoint61b, checkpoint61c, checkpoint61a, checkpoint61l, checkpoint61j, checkpoint61k, checkpoint61h, checkpoint61i
Changes since 1.4: +7 -6 lines
add missing cvs $Header:$ or $Name:$

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

  ViewVC Help
Powered by ViewVC 1.1.22