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

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

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


Revision 1.9 - (show annotations) (download)
Tue Oct 9 00:02:50 2007 UTC (16 years, 7 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint63h, checkpoint63i, checkpoint63j, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63g, checkpoint63a, checkpoint63b, checkpoint63c, checkpoint60, checkpoint61, checkpoint62, checkpoint63, checkpoint59q, checkpoint59p, checkpoint59r, checkpoint59m, checkpoint59l, checkpoint59o, checkpoint59n, checkpoint59i, checkpoint59k, checkpoint59j, 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, checkpoint61f, checkpoint61g, checkpoint61d, checkpoint61e, checkpoint61b, checkpoint61c, checkpoint61a, checkpoint61n, checkpoint61o, checkpoint61l, checkpoint61m, checkpoint61j, checkpoint61k, checkpoint61h, checkpoint61i, checkpoint61v, checkpoint61w, checkpoint61t, checkpoint61u, checkpoint61r, checkpoint61s, checkpoint61p, checkpoint61q, checkpoint61z, checkpoint61x, checkpoint61y
Changes since 1.8: +5 -4 lines
add missing cvs $Header:$ or $Name:$

1 C $Header: $
2 C $Name: $
3
4 #include "COST_CPPOPTIONS.h"
5
6
7 subroutine cost_salt( myiter, mytime, mythid )
8
9 c ==================================================================
10 c SUBROUTINE cost_salt
11 c ==================================================================
12 c
13 c o Evaluate cost function contribution of salinity.
14 c
15 c started: Christian Eckert eckert@mit.edu 30-Jun-1999
16 c
17 c changed: Christian Eckert eckert@mit.edu 25-Feb-2000
18 c
19 c - Restructured the code in order to create a package
20 c for the MITgcmUV.
21 c
22 c changed: Patrick Heimbach heimbach@mit.edu 27-May-2000
23 c
24 c - set ladinit to .true. to initialise adsbar file
25 c
26 c ==================================================================
27 c SUBROUTINE cost_salt
28 c ==================================================================
29
30 implicit none
31
32 c == global variables ==
33
34 #include "EEPARAMS.h"
35 #include "SIZE.h"
36 #include "PARAMS.h"
37 #include "GRID.h"
38 #include "DYNVARS.h"
39
40 #include "cal.h"
41 #include "ecco_cost.h"
42 #include "ctrl.h"
43 #include "ctrl_dummy.h"
44 #include "optim.h"
45
46 c == routine arguments ==
47
48 integer myiter
49 _RL mytime
50 integer mythid
51
52 c == local variables ==
53
54 integer bi,bj
55 integer i,j,k
56 integer itlo,ithi
57 integer jtlo,jthi
58 integer jmin,jmax
59 integer imin,imax
60 integer irec, irectmp
61 integer levmon
62 integer levoff
63 integer ilsalt
64
65 _RL fctile
66 _RL fcthread
67
68 _RL cmask (1-olx:snx+olx,1-oly:sny+oly)
69 _RL spval
70 _RL spmax
71
72 character*(80) fnamesalt
73
74 logical doglobalread
75 logical ladinit
76
77 character*(MAX_LEN_MBUF) msgbuf
78
79 #ifdef GENERIC_BAR_MONTH
80 integer mrec, nyears, iyear
81 #endif
82
83 _RL diagnosfld3d(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)
84
85 c == external functions ==
86
87 integer ilnblnk
88 external ilnblnk
89
90 c == end of interface ==
91
92 jtlo = mybylo(mythid)
93 jthi = mybyhi(mythid)
94 itlo = mybxlo(mythid)
95 ithi = mybxhi(mythid)
96 jmin = 1
97 jmax = sny
98 imin = 1
99 imax = snx
100
101 spval = 25.
102 spmax = 40.
103
104 c-- Read tiled data.
105 doglobalread = .false.
106 ladinit = .false.
107
108 #ifdef ALLOW_SALT_COST_CONTRIBUTION
109
110 if (optimcycle .ge. 0) then
111 ilsalt = ilnblnk( sbarfile )
112 write(fnamesalt(1:80),'(2a,i10.10)')
113 & sbarfile(1:ilsalt),'.',optimcycle
114 endif
115
116 fcthread = 0. _d 0
117
118 #ifdef GENERIC_BAR_MONTH
119 c-- Loop over month
120 do irec = 1,min(nmonsrec,12)
121 nyears=int((nmonsrec-irec)/12)+1
122 do iyear=1,nyears
123 mrec=irec+(iyear-1)*12
124 irectmp=mrec
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 irectmp=irec
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 = (irectmp-1) + levoff
164 levmon = mod(levmon-1,12)+1
165
166 call mdsreadfield( sdatfile, 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) = cosphi(i,j,bi,bj)
180 if (sdat(i,j,k,bi,bj) .eq. 0.) then
181 cmask(i,j) = 0. _d 0
182 else if (sdat(i,j,k,bi,bj) .lt. spval) then
183 cmask(i,j) = 0. _d 0
184 else if (sdat(i,j,k,bi,bj) .gt. spmax) then
185 cmask(i,j) = 0. _d 0
186 endif
187 enddo
188 enddo
189
190 c-- Compute model data misfit and cost function term for
191 c the salinity field.
192 do j = jmin,jmax
193 do i = imin,imax
194 if ( _hFacC(i,j,k,bi,bj) .ne. 0. ) then
195 fctile = fctile +
196 & (wsaltLev(i,j,k,bi,bj)*cmask(i,j)*
197 & (sbar(i,j,k,bi,bj) - sdat(i,j,k,bi,bj))*
198 & (sbar(i,j,k,bi,bj) - sdat(i,j,k,bi,bj)) )
199 if ( wsaltLev(i,j,k,bi,bj)*cmask(i,j) .ne. 0. )
200 & num_salt(bi,bj) = num_salt(bi,bj) + 1. _d 0
201 diagnosfld3d(i,j,k,bi,bj) =
202 & (wsaltLev(i,j,k,bi,bj)*cmask(i,j)*
203 & (sbar(i,j,k,bi,bj) - sdat(i,j,k,bi,bj))*
204 & (sbar(i,j,k,bi,bj) - sdat(i,j,k,bi,bj)) )
205 else
206 diagnosfld3d(i,j,k,bi,bj) = 0.
207 endif
208 enddo
209 enddo
210
211 enddo
212 c-- End of loop over layers.
213
214 call mdswritefield( 'DiagnosCost_ClimSalt',
215 & writeBinaryPrec, globalfiles, 'RL', Nr,
216 & diagnosfld3d, irec, optimcycle, mythid )
217
218 fcthread = fcthread + fctile
219 objf_salt(bi,bj) = objf_salt(bi,bj) + fctile
220
221 enddo
222 enddo
223
224 enddo
225 c-- End of loop over records.
226
227 #endif
228
229 return
230 end
231

  ViewVC Help
Powered by ViewVC 1.1.22