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

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

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


Revision 1.7 - (show annotations) (download)
Fri Sep 2 16:02:11 2005 UTC (18 years, 9 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint57t_post, checkpoint57v_post, checkpoint57s_post, checkpoint57r_post, checkpint57u_post
Changes since 1.6: +11 -92 lines
Clean up. Allow flag GENERIC_BAR_MONTH to be used for less than a year.

1 C $Header: /u/gcmpack/MITgcm/pkg/ecco/cost_theta.F,v 1.6 2005/04/29 10:31:56 heimbach Exp $
2
3 #include "COST_CPPOPTIONS.h"
4
5
6 subroutine cost_theta( myiter, mytime, mythid )
7
8 c ==================================================================
9 c SUBROUTINE cost_theta
10 c ==================================================================
11 c
12 c o Evaluate cost function contribution of temperature.
13 c
14 c started: Christian Eckert eckert@mit.edu 30-Jun-1999
15 c
16 c changed: Christian Eckert eckert@mit.edu 25-Feb-2000
17 c
18 c - Restructured the code in order to create a package
19 c for the MITgcmUV.
20 c
21 c changed: Patrick Heimbach heimbach@mit.edu 27-May-2000
22 c
23 c - set ladinit to .true. to initialise adtbar file
24 c
25 c ==================================================================
26 c SUBROUTINE cost_theta
27 c ==================================================================
28
29 implicit none
30
31 c == global variables ==
32
33 #include "EEPARAMS.h"
34 #include "SIZE.h"
35 #include "GRID.h"
36 #include "DYNVARS.h"
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, irectmp
59 integer levmon
60 integer levoff
61 integer iltheta
62
63 _RL fctile
64 _RL fcthread
65
66 _RL cmask (1-olx:snx+olx,1-oly:sny+oly)
67 _RL spval
68 _RL spmax
69
70 character*(80) fnametheta
71
72 logical doglobalread
73 logical ladinit
74
75 character*(MAX_LEN_MBUF) msgbuf
76 #ifdef GENERIC_BAR_MONTH
77 integer mrec, nyears, iyear
78 #endif
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 spval = -1.8
96 spmax = 40.
97
98 c-- Read tiled data.
99 doglobalread = .false.
100 ladinit = .false.
101
102 #ifdef ALLOW_THETA_COST_CONTRIBUTION
103
104 if (optimcycle .ge. 0) then
105 iltheta = ilnblnk( tbarfile )
106 write(fnametheta(1:80),'(2a,i10.10)')
107 & tbarfile(1:iltheta),'.',optimcycle
108 endif
109
110 fcthread = 0. _d 0
111
112 #ifdef GENERIC_BAR_MONTH
113 c-- Loop over month
114 do irec = 1,min(nmonsrec,12)
115 nyears=int((nmonsrec-irec)/12)+1
116 do iyear=1,nyears
117 mrec=irec+(iyear-1)*12
118 irectmp=mrec
119 c-- Read time averages and the monthly mean data.
120 call active_read_xyz( fnametheta, tbar, mrec,
121 & doglobalread, ladinit,
122 & optimcycle, mythid,
123 & xx_tbar_mean_dummy )
124 do bj = jtlo,jthi
125 do bi = itlo,ithi
126 do k = 1,nr
127 do j = jmin,jmax
128 do i = imin,imax
129 if(iyear.eq.1) then
130 tbar_gen(i,j,k,bi,bj) =tbar(i,j,k,bi,bj)
131 elseif(iyear.eq.nyears) then
132 tbar(i,j,k,bi,bj) =(tbar_gen(i,j,k,bi,bj)
133 $ +tbar(i,j,k,bi,bj))/float(nyears)
134 else
135 tbar_gen(i,j,k,bi,bj) =tbar_gen(i,j,k,bi,bj)
136 $ +tbar(i,j,k,bi,bj)
137 endif
138 enddo
139 enddo
140 enddo
141 enddo
142 enddo
143 enddo
144 #else
145 c-- Loop over records.
146 do irec = 1,nmonsrec
147
148 irectmp = irec
149 c-- Read time averages and the monthly mean data.
150 call active_read_xyz( fnametheta, tbar, irec,
151 & doglobalread, ladinit,
152 & optimcycle, mythid,
153 & xx_tbar_mean_dummy )
154 #endif
155 c-- Determine the month to be read.
156 levoff = mod(modelstartdate(1)/100,100)
157 levmon = (irectmp-1) + levoff
158 levmon = mod(levmon-1,12)+1
159
160 call mdsreadfield( tdatfile, cost_iprec, cost_yftype,
161 & nr, tdat, levmon, mythid)
162
163 do bj = jtlo,jthi
164 do bi = itlo,ithi
165
166 c-- Loop over the model layers
167 fctile = 0. _d 0
168 do k = 1,nr
169
170 c-- Determine the mask on weights
171 do j = jmin,jmax
172 do i = imin,imax
173 cmask(i,j) = cosphi(i,j,bi,bj)
174 if (tdat(i,j,k,bi,bj) .eq. 0.) then
175 cmask(i,j) = 0. _d 0
176 else if (tdat(i,j,k,bi,bj) .lt. spval) then
177 cmask(i,j) = 0. _d 0
178 else if (tdat(i,j,k,bi,bj) .gt. spmax) then
179 cmask(i,j) = 0. _d 0
180 endif
181 enddo
182 enddo
183
184 c-- Compute model data misfit and cost function term for
185 c the temperature field.
186 do j = jmin,jmax
187 do i = imin,imax
188 if ( _hFacC(i,j,k,bi,bj) .ne. 0. ) then
189 fctile = fctile +
190 & (wtheta(k,bi,bj)*cmask(i,j)*
191 & (tbar(i,j,k,bi,bj) - tdat(i,j,k,bi,bj))*
192 & (tbar(i,j,k,bi,bj) - tdat(i,j,k,bi,bj)) )
193 if ( wtheta(k,bi,bj)*cmask(i,j) .ne. 0. )
194 & num_temp(bi,bj) = num_temp(bi,bj) + 1. _d 0
195 endif
196 enddo
197 enddo
198
199 enddo
200 c-- End of loop over layers.
201
202 fcthread = fcthread + fctile
203 objf_temp(bi,bj) = objf_temp(bi,bj) + fctile
204
205 enddo
206 enddo
207
208 enddo
209 c-- End of loop over records.
210
211 #endif
212
213 return
214 end
215

  ViewVC Help
Powered by ViewVC 1.1.22