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

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

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


Revision 1.7 - (show annotations) (download)
Wed Jun 17 15:11:32 2009 UTC (14 years, 10 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint63p, checkpoint63q, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63h, checkpoint63i, checkpoint63j, checkpoint63k, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63g, checkpoint63a, checkpoint63b, checkpoint63c, checkpoint62, checkpoint63, checkpoint62c, checkpoint62b, checkpoint62a, checkpoint62g, checkpoint62f, checkpoint62e, checkpoint62d, checkpoint62k, checkpoint62j, checkpoint62i, checkpoint62h, checkpoint62o, checkpoint62n, checkpoint62m, checkpoint62l, checkpoint62s, checkpoint62r, checkpoint62q, checkpoint62p, checkpoint62w, checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62z, checkpoint62y, checkpoint62x, checkpoint61v, checkpoint61w, checkpoint61t, checkpoint61u, checkpoint61r, checkpoint61s, checkpoint61q, checkpoint61z, checkpoint61x, checkpoint61y
Changes since 1.6: +1 -13 lines
Remove depth limitations for in-situ costs.

1 C $Header: /u/gcmpack/MITgcm/pkg/ecco/cost_ctdtclim.F,v 1.6 2009/04/28 18:13:27 jmc Exp $
2 C $Name: $
3
4 #include "COST_CPPOPTIONS.h"
5
6
7 subroutine cost_Ctdtclim(
8 I myiter,
9 I mytime,
10 I mythid
11 & )
12
13 c ==================================================================
14 c SUBROUTINE cost_Ctdtclim
15 c ==================================================================
16 c
17 c o Evaluate cost function contribution of temperature.
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 adtbar file
29 c
30 c ==================================================================
31 c SUBROUTINE cost_Ctdtclim
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 ilctdtclim
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) fnametheta
78
79 logical doglobalread
80 logical ladinit
81
82 character*(MAX_LEN_MBUF) msgbuf
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 = -1.8
103
104 c-- Read tiled data.
105 doglobalread = .false.
106 ladinit = .false.
107
108 #ifdef ALLOW_CTDTCLIM_COST_CONTRIBUTION
109
110 if (optimcycle .ge. 0) then
111 ilctdtclim = ilnblnk( tbarfile )
112 write(fnametheta(1:80),'(2a,i10.10)')
113 & tbarfile(1:ilctdtclim),'.',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( fnametheta, tbar, mrec,
127 & doglobalread, ladinit,
128 & optimcycle, mythid,
129 & xx_tbar_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 tbar_gen(i,j,k,bi,bj) =tbar(i,j,k,bi,bj)
137 elseif(iyear.eq.nyears) then
138 tbar(i,j,k,bi,bj) =(tbar_gen(i,j,k,bi,bj)
139 $ +tbar(i,j,k,bi,bj))/float(nyears)
140 else
141 tbar_gen(i,j,k,bi,bj) =tbar_gen(i,j,k,bi,bj)
142 $ +tbar(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( fnametheta, tbar, irec,
156 & doglobalread, ladinit,
157 & optimcycle, mythid,
158 & xx_tbar_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( ctdtclimfile, cost_iprec, cost_yftype,
166 & nr, tdat, 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 on weights
176 do j = jmin,jmax
177 do i = imin,imax
178 cmask(i,j) = 1. _d 0
179 if (tdat(i,j,k,bi,bj) .eq. 0.) then
180 cmask(i,j) = 0. _d 0
181 endif
182
183 if (tdat(i,j,k,bi,bj) .lt. spval) then
184 cmask(i,j) = 0. _d 0
185 endif
186 enddo
187 enddo
188
189 c-- Compute model data misfit and cost function term for
190 c the temperature field.
191 do j = jmin,jmax
192 do i = imin,imax
193 if (_hFacC(i,j,k,bi,bj) .ne. 0.) then
194 fctile = fctile +
195 & (wtheta2(i,j,k,bi,bj)*cosphi(i,j,bi,bj)*
196 & cmask(i,j)*
197 & (tbar(i,j,k,bi,bj) - tdat(i,j,k,bi,bj))*
198 & (tbar(i,j,k,bi,bj) - tdat(i,j,k,bi,bj)) )
199 endif
200 enddo
201 enddo
202
203 enddo
204 c-- End of loop over layers.
205
206 fcthread = fcthread + fctile
207 objf_ctdtclim(bi,bj) = objf_ctdtclim(bi,bj) + fctile
208
209 #ifdef ECCO_VERBOSE
210 c-- Print cost function for each tile in each thread.
211 write(msgbuf,'(a)') ' '
212 call print_message( msgbuf, standardmessageunit,
213 & SQUEEZE_RIGHT , mythid)
214 write(msgbuf,'(a,i8.8,1x,i3.3,1x,i3.3)')
215 & ' cost_Ctdtclim: irec,bi,bj = ',irec,bi,bj
216 call print_message( msgbuf, standardmessageunit,
217 & SQUEEZE_RIGHT , mythid)
218 write(msgbuf,'(a,d22.15)')
219 & ' cost function (temperature) = ',
220 & fctile
221 call print_message( msgbuf, standardmessageunit,
222 & SQUEEZE_RIGHT , mythid)
223 write(msgbuf,'(a)') ' '
224 call print_message( msgbuf, standardmessageunit,
225 & SQUEEZE_RIGHT , mythid)
226 #endif
227
228 enddo
229 enddo
230
231 #ifdef ECCO_VERBOSE
232 c-- Print cost function for all tiles.
233 _GLOBAL_SUM_RL( fcthread , myThid )
234 write(msgbuf,'(a)') ' '
235 call print_message( msgbuf, standardmessageunit,
236 & SQUEEZE_RIGHT , mythid)
237 write(msgbuf,'(a,i8.8)')
238 & ' cost_Ctdtclim: irec = ',irec
239 call print_message( msgbuf, standardmessageunit,
240 & SQUEEZE_RIGHT , mythid)
241 write(msgbuf,'(a,a,d22.15)')
242 & ' global cost function value',
243 & ' (temperature) = ',fcthread
244 call print_message( msgbuf, standardmessageunit,
245 & SQUEEZE_RIGHT , mythid)
246 write(msgbuf,'(a)') ' '
247 call print_message( msgbuf, standardmessageunit,
248 & SQUEEZE_RIGHT , mythid)
249 #endif
250
251 #ifdef GENERIC_BAR_MONTH
252 endif
253 #endif
254 enddo
255 c-- End of loop over records.
256
257 #else
258 c-- Do not enter the calculation of the temperature contribution to
259 c-- the final cost function.
260
261 #ifdef ECCO_VERBOSE
262 _BEGIN_MASTER( mythid )
263 write(msgbuf,'(a)') ' '
264 call print_message( msgbuf, standardmessageunit,
265 & SQUEEZE_RIGHT , mythid)
266 write(msgbuf,'(a,a)')
267 & ' cost_Ctdtclim: no contribution of temperature field ',
268 & 'to cost function.'
269 call print_message( msgbuf, standardmessageunit,
270 & SQUEEZE_RIGHT , mythid)
271 write(msgbuf,'(a,a,i9.8)')
272 & ' cost_Ctdtclim: number of records that would have',
273 & ' been processed: ',nmonsrec
274 call print_message( msgbuf, standardmessageunit,
275 & SQUEEZE_RIGHT , mythid)
276 write(msgbuf,'(a)') ' '
277 call print_message( msgbuf, standardmessageunit,
278 & SQUEEZE_RIGHT , mythid)
279 _END_MASTER( mythid )
280 #endif
281
282 #endif
283
284 return
285 end
286

  ViewVC Help
Powered by ViewVC 1.1.22