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

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

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


Revision 1.1 - (hide annotations) (download)
Thu Nov 6 22:10:07 2003 UTC (20 years, 7 months ago) by heimbach
Branch: MAIN
CVS Tags: branch-netcdf, checkpoint52b_pre, checkpoint52, checkpoint52a_post, checkpoint52b_post, checkpoint52c_post, ecco_c52_e35, checkpoint52a_pre, checkpoint51u_post
o merging from ecco-branch
o pkg/ecco now containes ecco-specific part of cost function
o top level routines the_main_loop, forward_step
  supersede those in model/src/
  previous input data.cost now in data.ecco
  (new namelist ecco_cost_nml)

1 heimbach 1.1 C $Header: /u/gcmpack/MITgcm/pkg/cost/Attic/cost_ctdt.F,v 1.1.2.3 2002/04/04 10:58:58 heimbach Exp $
2    
3     #include "COST_CPPOPTIONS.h"
4    
5    
6     subroutine cost_CTDT(
7     I myiter,
8     I mytime,
9     I mythid
10     & )
11    
12     c ==================================================================
13     c SUBROUTINE cost_CTDT
14     c ==================================================================
15     c
16     c o Evaluate cost function contribution of CTD temperature data.
17     c
18     c started: Elisabeth Remy eremy@ucsd.edu 30-Aug-2000
19     c
20     c
21     c ==================================================================
22     c SUBROUTINE cost_CTDT
23     c ==================================================================
24    
25     implicit none
26    
27     c == global variables ==
28    
29     #include "EEPARAMS.h"
30     #include "SIZE.h"
31     #include "GRID.h"
32     #include "DYNVARS.h"
33    
34     #include "cal.h"
35     #include "ecco_cost.h"
36     #include "ctrl.h"
37     #include "ctrl_dummy.h"
38     #include "optim.h"
39    
40     c == routine arguments ==
41    
42     integer myiter
43     _RL mytime
44     integer mythid
45    
46     c == local variables ==
47    
48     integer bi,bj
49     integer i,j,k
50     integer itlo,ithi
51     integer jtlo,jthi
52     integer jmin,jmax
53     integer imin,imax
54     integer nrec
55     integer irec
56     integer ilu
57    
58     _RL fctile_ctdt
59     _RL fcthread_ctdt
60     _RL www (1-olx:snx+olx,1-oly:sny+oly)
61     _RL wtmp (1-olx:snx+olx,1-oly:sny+oly)
62     _RL tmpobs (1-olx:snx+olx,1-oly:sny+oly)
63     _RL tmpbar (1-olx:snx+olx,1-oly:sny+oly)
64     _RL cmask (1-olx:snx+olx,1-oly:sny+oly)
65     _RL spval
66    
67     character*(80) fnametheta
68    
69     logical doglobalread
70     logical ladinit
71    
72     character*(MAX_LEN_MBUF) msgbuf
73    
74     c == external functions ==
75    
76     integer ilnblnk
77     external ilnblnk
78    
79     c == end of interface ==
80    
81     jtlo = mybylo(mythid)
82     jthi = mybyhi(mythid)
83     itlo = mybxlo(mythid)
84     ithi = mybxhi(mythid)
85     jmin = 1
86     jmax = sny
87     imin = 1
88     imax = snx
89    
90     spval = -9990.
91    
92     c-- Read state record from global file.
93     doglobalread = .false.
94     ladinit = .false.
95    
96     #ifdef ALLOW_CTDT_COST_CONTRIBUTION
97    
98     #ifdef ECCO_VERBOSE
99     write(msgbuf,'(a)') ' '
100     call print_message( msgbuf, standardmessageunit,
101     & SQUEEZE_RIGHT , mythid)
102     write(msgbuf,'(a,i8.8)')
103     & ' cost_CTDT: number of records to process =', nmonsrec
104     call print_message( msgbuf, standardmessageunit,
105     & SQUEEZE_RIGHT , mythid)
106     write(msgbuf,'(a)') ' '
107     call print_message( msgbuf, standardmessageunit,
108     & SQUEEZE_RIGHT , mythid)
109     #endif
110    
111     if (optimcycle .ge. 0) then
112     ilu=ilnblnk( tbarfile )
113     write(fnametheta(1:80),'(2a,i10.10)') tbarfile(1:ilu),'.',
114     & optimcycle
115     else
116     print*
117     print*,' cost_CTDT: optimcycle has a wrong value.'
118     print*,' optimcycle = ',optimcycle
119     print*
120     stop ' ... stopped in cost_CTDT.'
121     endif
122    
123     fcthread_ctdt = 0. _d 0
124    
125     c-- Loop over records.
126     do irec = 1,nmonsrec
127    
128     c-- Read time averages and the monthly mean data.
129     call active_read_xyz( fnametheta, tbar, irec,
130     & doglobalread, ladinit,
131     & optimcycle, mythid, xx_tbar_mean_dummy )
132    
133     call mdsreadfield( ctdtfile, cost_iprec, 'RL', nr, ctdtobs,
134     & irec, mythid)
135    
136     c-- Loop over this thread's tiles.
137     do bj = jtlo,jthi
138     do bi = itlo,ithi
139     c-- Loop over the model layers
140    
141     fctile_ctdt = 0. _d 0
142    
143     do k = 1,nr
144    
145     c-- Determine the weights to be used.
146     do j = jmin,jmax
147     do i = imin,imax
148     cmask(i,j) = 1. _d 0
149     if (ctdtobs(i,j,k,bi,bj) .eq. 0.) then
150     cmask(i,j) = 0. _d 0
151     endif
152    
153     if (ctdtobs(i,j,k,bi,bj) .lt. spval) then
154     cmask(i,j) = 0. _d 0
155     endif
156    
157     cph(
158     cph print *, 'WARNING: SPECIFIC SETUP FOR ECCO'
159     cph below statement could be replaced by following
160     cph to make it independnet of Nr:
161     cph
162     cph if ( rC(K) .GT. -1000. ) then
163     cph)
164     c set cmask=0 in areas shallower than 1000m
165     if (_hFacC(i,j,13,bi,bj) .eq. 0.) then
166     cmask(i,j) = 0. _d 0
167     endif
168    
169     if (_hFacC(i,j,k,bi,bj) .eq. 0.) then
170     cmask(i,j) = 0. _d 0
171     endif
172    
173     enddo
174     enddo
175    
176     do j = jmin,jmax
177     do i = imin,imax
178     www(i,j) = cosphi(i,j,bi,bj)*cmask(i,j)
179     tmpobs(i,j) = ctdtobs(i,j,k,bi,bj)
180     tmpbar(i,j) = tbar(i,j,k,bi,bj)
181     c wtmp(i,j) = wtheta(k,bi,bj)
182     wtmp(i,j) = wtheta2(i,j,k,bi,bj)
183     enddo
184     enddo
185    
186     do j = jmin,jmax
187     do i = imin,imax
188     c-- The array ctdtobs contains CTD temperature.
189     fctile_ctdt = fctile_ctdt +
190     & (wtmp(i,j)*www(i,j))*
191     & (tmpbar(i,j)-tmpobs(i,j))*
192     & (tmpbar(i,j)-tmpobs(i,j))
193    
194     enddo
195     enddo
196     enddo
197     c-- End of loop over layers.
198    
199     fcthread_ctdt = fcthread_ctdt + fctile_ctdt
200     objf_ctdt(bi,bj) = objf_ctdt(bi,bj) + fctile_ctdt
201    
202     #ifdef ECCO_VERBOSE
203     write(msgbuf,'(a)') ' '
204     call print_message( msgbuf, standardmessageunit,
205     & SQUEEZE_RIGHT , mythid)
206     write(msgbuf,'(a,i8.8,1x,i3.3,1x,i3.3)')
207     & ' COST_CTDT: irec,bi,bj = ',irec,bi,bj
208     call print_message( msgbuf, standardmessageunit,
209     & SQUEEZE_RIGHT , mythid)
210     write(msgbuf,'(a,d22.15)')
211     & ' COST_CTDT: cost function = ', fctile_ctdt
212     call print_message( msgbuf, standardmessageunit,
213     & SQUEEZE_RIGHT , mythid)
214     write(msgbuf,'(a)') ' '
215     call print_message( msgbuf, standardmessageunit,
216     & SQUEEZE_RIGHT , mythid)
217     #endif
218    
219     enddo
220     enddo
221    
222     #ifdef ECCO_VERBOSE
223     c-- Print cost function for all tiles.
224     _GLOBAL_SUM_R8( fcthread_ctdt , myThid )
225     write(msgbuf,'(a)') ' '
226     call print_message( msgbuf, standardmessageunit,
227     & SQUEEZE_RIGHT , mythid)
228     write(msgbuf,'(a,i8.8)')
229     & ' cost_CTDT: irec = ',irec
230     call print_message( msgbuf, standardmessageunit,
231     & SQUEEZE_RIGHT , mythid)
232     write(msgbuf,'(a,a,d22.15)')
233     & ' global cost function value',
234     & ' ( CTD temp. ) = ',fcthread_ctdt
235     call print_message( msgbuf, standardmessageunit,
236     & SQUEEZE_RIGHT , mythid)
237     write(msgbuf,'(a)') ' '
238     call print_message( msgbuf, standardmessageunit,
239     & SQUEEZE_RIGHT , mythid)
240     #endif
241    
242     enddo
243     c-- End of second loop over records.
244    
245     #else
246     c-- Do not enter the calculation of the CTD temperature contribution
247     c-- to the final cost function.
248    
249     fctile_ctdt = 0. _d 0
250     fcthread_ctdt = 0. _d 0
251    
252     crg
253     nrec = 1
254     crg
255    
256     _BEGIN_MASTER( mythid )
257     write(msgbuf,'(a)') ' '
258     call print_message( msgbuf, standardmessageunit,
259     & SQUEEZE_RIGHT , mythid)
260     write(msgbuf,'(a,a)')
261     & ' cost_CTDT: no contribution of CTD temperature ',
262     & ' to cost function.'
263     call print_message( msgbuf, standardmessageunit,
264     & SQUEEZE_RIGHT , mythid)
265     write(msgbuf,'(a,a,i9.8)')
266     & ' cost_CTDT: number of records that would have',
267     & ' been processed: ',nrec
268     call print_message( msgbuf, standardmessageunit,
269     & SQUEEZE_RIGHT , mythid)
270     write(msgbuf,'(a)') ' '
271     call print_message( msgbuf, standardmessageunit,
272     & SQUEEZE_RIGHT , mythid)
273     _END_MASTER( mythid )
274     #endif
275    
276     return
277     end

  ViewVC Help
Powered by ViewVC 1.1.22