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

Annotation of /MITgcm/pkg/ecco/cost_obcss.F

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


Revision 1.3 - (hide annotations) (download)
Tue Oct 9 00:02:50 2007 UTC (16 years, 8 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.2: +6 -4 lines
add missing cvs $Header:$ or $Name:$

1 jmc 1.3 C $Header: $
2     C $Name: $
3 heimbach 1.1
4     #include "COST_CPPOPTIONS.h"
5     #ifdef ALLOW_OBCS
6     # include "OBCS_OPTIONS.h"
7     #endif
8    
9     subroutine cost_obcss(
10     I myiter,
11     I mytime,
12     I mythid
13     & )
14    
15     c ==================================================================
16     c SUBROUTINE cost_obcss
17     c ==================================================================
18     c
19     c o cost function contribution obc
20     c
21     c ==================================================================
22     c SUBROUTINE cost_obcss
23     c ==================================================================
24    
25     implicit none
26    
27     c == global variables ==
28    
29     #include "EEPARAMS.h"
30     #include "SIZE.h"
31     #include "PARAMS.h"
32     #include "GRID.h"
33     #include "DYNVARS.h"
34     #ifdef ALLOW_OBCS
35     # include "OBCS.h"
36     #endif
37    
38     #include "cal.h"
39     #include "ecco_cost.h"
40     #include "ctrl.h"
41     #include "ctrl_dummy.h"
42     #include "optim.h"
43    
44     c == routine arguments ==
45    
46     integer myiter
47     _RL mytime
48     integer mythid
49    
50     c == local variables ==
51    
52     integer bi,bj
53     integer i,j,k
54     integer itlo,ithi
55     integer jtlo,jthi
56     integer jmin,jmax
57     integer imin,imax
58     integer irec
59     integer il
60     integer iobcs
61     integer jp1
62    
63     _RL fctile
64     _RL fcthread
65     _RL dummy
66    
67 heimbach 1.2 character*(80) fnametheta
68     character*(80) fnamesalt
69     character*(80) fnameuvel
70     character*(80) fnamevvel
71 heimbach 1.1
72     logical doglobalread
73     logical ladinit
74    
75     #ifdef ECCO_VERBOSE
76     character*(MAX_LEN_MBUF) msgbuf
77     #endif
78    
79     c == external functions ==
80    
81     integer ilnblnk
82     external ilnblnk
83    
84     c == end of interface ==
85    
86     jtlo = mybylo(mythid)
87     jthi = mybyhi(mythid)
88     itlo = mybxlo(mythid)
89     ithi = mybxhi(mythid)
90     jmin = 1
91     jmax = sny
92     imin = 1
93     imax = snx
94    
95     c-- Read tiled data.
96     doglobalread = .false.
97     ladinit = .false.
98    
99     #ifdef ALLOW_OBCSS_COST_CONTRIBUTION
100    
101     jp1 = 1
102     fcthread = 0. _d 0
103    
104 heimbach 1.2 c-- Loop over records.
105     do irec = 1,nmonsrec
106 heimbach 1.1
107 heimbach 1.2 c-- temperature
108     iobcs = 1
109     c-- Read time averages and the monthly mean data.
110     il = ilnblnk( tbarfile )
111     write(fnametheta(1:80),'(2a,i10.10)')
112     & tbarfile(1:il),'.',optimcycle
113     call active_read_xyz( fnametheta, tbar, irec,
114     & doglobalread, ladinit,
115     & optimcycle, mythid,
116     & xx_tbar_mean_dummy )
117 heimbach 1.1
118 jmc 1.3 call mdsreadfieldxz( OBStFile, readBinaryPrec, 'RS',
119 heimbach 1.2 & nr, OBSt, irec, mythid)
120 heimbach 1.1
121 heimbach 1.2 do bj = jtlo,jthi
122     do bi = itlo,ithi
123     c-- Loop over the model layers
124     fctile = 0. _d 0
125     do k = 1,nr
126     c-- Compute model data misfit and cost function term for
127     c the temperature field.
128     do i = imin,imax
129     j = OB_Js(I,bi,bj)
130     if (maskS(i,j+jp1,k,bi,bj) .ne. 0.) then
131     fctile = fctile +
132     & wobcss(k,iobcs)*cosphi(i,j,bi,bj)*
133     & (tbar(i,j,k,bi,bj) - OBSt(i,k,bi,bj))*
134     & (tbar(i,j,k,bi,bj) - OBSt(i,k,bi,bj))
135     endif
136     enddo
137     enddo
138     c-- End of loop over layers.
139     fcthread = fcthread + fctile
140     objf_obcss(bi,bj) = objf_obcss(bi,bj) + fctile
141     enddo
142 heimbach 1.1 enddo
143    
144 heimbach 1.2 c-- salt
145     iobcs = 2
146     c-- Read time averages and the monthly mean data.
147     il = ilnblnk( sbarfile )
148     write(fnamesalt(1:80),'(2a,i10.10)')
149     & sbarfile(1:il),'.',optimcycle
150     call active_read_xyz( fnamesalt, sbar, irec,
151     & doglobalread, ladinit,
152     & optimcycle, mythid,
153     & xx_sbar_mean_dummy )
154    
155 jmc 1.3 call mdsreadfieldxz( OBSsFile, readBinaryPrec, 'RS',
156 heimbach 1.2 & nr, OBSs, irec, mythid)
157    
158     do bj = jtlo,jthi
159     do bi = itlo,ithi
160     c-- Loop over the model layers
161     fctile = 0. _d 0
162     do k = 1,nr
163     c-- Compute model data misfit and cost function term for
164     c the temperature field.
165     do i = imin,imax
166     j = OB_Js(I,bi,bj)
167     if (maskS(i,j+jp1,k,bi,bj) .ne. 0.) then
168     fctile = fctile +
169     & wobcss(k,iobcs)*cosphi(i,j,bi,bj)*
170     & (sbar(i,j,k,bi,bj) - OBSs(i,k,bi,bj))*
171     & (sbar(i,j,k,bi,bj) - OBSs(i,k,bi,bj))
172     endif
173     enddo
174 heimbach 1.1 enddo
175 heimbach 1.2 c-- End of loop over layers.
176     fcthread = fcthread + fctile
177     objf_obcss(bi,bj) = objf_obcss(bi,bj) + fctile
178 heimbach 1.1 enddo
179 heimbach 1.2 enddo
180 heimbach 1.1
181 heimbach 1.2 c-- uvel
182     iobcs = 3
183     c-- Read time averages and the monthly mean data.
184     il = ilnblnk( ubarfile )
185     write(fnameuvel(1:80),'(2a,i10.10)')
186     & ubarfile(1:il),'.',optimcycle
187     call active_read_xyz( fnameuvel, ubar, irec,
188     & doglobalread, ladinit,
189     & optimcycle, mythid,
190     & dummy )
191 heimbach 1.1
192 jmc 1.3 call mdsreadfieldxz( OBSuFile, readBinaryPrec, 'RS',
193 heimbach 1.2 & nr, OBSu, irec, mythid)
194 heimbach 1.1 do bj = jtlo,jthi
195     do bi = itlo,ithi
196 heimbach 1.2 c-- Loop over the model layers
197     fctile = 0. _d 0
198     do k = 1,nr
199     c-- Compute model data misfit and cost function term for
200     c the temperature field.
201     do i = imin,imax
202     j = OB_Js(I,bi,bj)
203     if (maskW(i,j,k,bi,bj) .ne. 0.) then
204     fctile = fctile +
205     & wobcss(k,iobcs)*cosphi(i,j,bi,bj)*
206     & (ubar(i,j,k,bi,bj) - OBSu(i,k,bi,bj))*
207     & (ubar(i,j,k,bi,bj) - OBSu(i,k,bi,bj))
208     endif
209     enddo
210     enddo
211     c-- End of loop over layers.
212     fcthread = fcthread + fctile
213     objf_obcss(bi,bj) = objf_obcss(bi,bj) + fctile
214     enddo
215     enddo
216    
217     c-- vvel
218     iobcs = 4
219     c-- Read time averages and the monthly mean data.
220     il = ilnblnk( vbarfile )
221     write(fnamevvel(1:80),'(2a,i10.10)')
222     & vbarfile(1:il),'.',optimcycle
223     call active_read_xyz( fnamevvel, vbar, irec,
224     & doglobalread, ladinit,
225     & optimcycle, mythid,
226     & dummy )
227 heimbach 1.1
228 jmc 1.3 call mdsreadfieldxz( OBSvFile, readBinaryPrec, 'RS',
229 heimbach 1.2 & nr, OBSv, irec, mythid)
230    
231     do bj = jtlo,jthi
232     do bi = itlo,ithi
233     c-- Loop over the model layers
234 heimbach 1.1 fctile = 0. _d 0
235 heimbach 1.2 do k = 1,nr
236     c-- Compute model data misfit and cost function term for
237     c the temperature field.
238     do i = imin,imax
239     j = OB_Js(I,bi,bj)
240     if (maskS(i,j+jp1,k,bi,bj) .ne. 0.) then
241     fctile = fctile +
242     & wobcss(k,iobcs)*cosphi(i,j,bi,bj)*
243     & (vbar(i,j,k,bi,bj) - OBSv(i,k,bi,bj))*
244     & (vbar(i,j,k,bi,bj) - OBSv(i,k,bi,bj))
245     endif
246     enddo
247 heimbach 1.1 enddo
248 heimbach 1.2 c-- End of loop over layers.
249     fcthread = fcthread + fctile
250 heimbach 1.1 objf_obcss(bi,bj) = objf_obcss(bi,bj) + fctile
251     enddo
252     enddo
253    
254     #ifdef ECCO_VERBOSE
255     c-- Print cost function for all tiles.
256     _GLOBAL_SUM_R8( fcthread , myThid )
257     write(msgbuf,'(a)') ' '
258     call print_message( msgbuf, standardmessageunit,
259     & SQUEEZE_RIGHT , mythid)
260     write(msgbuf,'(a,i8.8)')
261     & ' cost_obcss: irec = ',irec
262     call print_message( msgbuf, standardmessageunit,
263     & SQUEEZE_RIGHT , mythid)
264     write(msgbuf,'(a,a,d22.15)')
265     & ' global cost function value',
266     & ' (obcss) = ',fcthread
267     call print_message( msgbuf, standardmessageunit,
268     & SQUEEZE_RIGHT , mythid)
269     write(msgbuf,'(a)') ' '
270     call print_message( msgbuf, standardmessageunit,
271     & SQUEEZE_RIGHT , mythid)
272     #endif
273    
274     enddo
275     c-- End of loop over records.
276    
277     #endif
278    
279     return
280     end
281    
282    
283    
284    
285    
286    
287    

  ViewVC Help
Powered by ViewVC 1.1.22