/[MITgcm]/MITgcm_contrib/mlosch/optim_m1qn3/optim_readdata.F
ViewVC logotype

Contents of /MITgcm_contrib/mlosch/optim_m1qn3/optim_readdata.F

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


Revision 1.6 - (show annotations) (download)
Thu May 3 11:26:05 2018 UTC (2 years ago) by mlosch
Branch: MAIN
CVS Tags: HEAD
Changes since 1.5: +85 -104 lines
spring cleaning

- adjust some debugging output
- reduce amount of output
- code cleaning (mainly indentation) for better readability

1 C $Header: /u/gcmpack/MITgcm_contrib/mlosch/optim_m1qn3/optim_readdata.F,v 1.5 2016/05/09 09:37:16 mlosch Exp $
2 C $Name: $
3
4 C ECCO_CPPOPTIONS used to affect maxcvars and defined ALLOW_OBCS?_CONTROL
5 C#include "ECCO_CPPOPTIONS.h"
6 C now:
7 C CTRL_OPTIONS affects maxcvars and may define ALLOW_OBCS?_CONTROL
8 #include "CTRL_OPTIONS.h"
9
10 subroutine optim_readdata(
11 I nn,
12 I dfile,
13 I lheaderonly,
14 O ff,
15 O vv
16 & )
17
18 c ==================================================================
19 c SUBROUTINE optim_readdata
20 c ==================================================================
21 c
22 c o Read the data written by the MITgcmUV state estimation setup and
23 c join them to one vector that is subsequently used by the minimi-
24 c zation algorithm "lsopt". Depending on the specified file name
25 c either the control vector or the gradient vector can be read.
26 c
27 c *dfile* should be the radix of the file: ecco_ctrl or ecco_cost
28 c
29 c started: Christian Eckert eckert@mit.edu 12-Apr-2000
30 c
31 c changed: Patrick Heimbach heimbach@mit.edu 19-Jun-2000
32 c - finished, revised and debugged
33 c
34 c ==================================================================
35 c SUBROUTINE optim_readdata
36 c ==================================================================
37
38 implicit none
39
40 c == global variables ==
41
42 #include "EEPARAMS.h"
43 #include "SIZE.h"
44 #if (defined (ALLOW_GENARR2D_CONTROL) || defined (ALLOW_GENARR3D_CONTROL) || defined (ALLOW_GENTIM2D_CONTROL))
45 # include "CTRL_SIZE.h"
46 #endif
47 #include "ctrl.h"
48 #include "optim.h"
49
50 c == routine arguments ==
51
52 integer nn
53 _RL ff
54
55 #ifdef DYNAMIC
56 _RL vv(nn)
57 #else
58 integer nmax
59 parameter( nmax = MAX_INDEPEND )
60 _RL vv(nmax)
61 #endif
62
63 character*(9) dfile
64 logical lheaderonly
65
66 c == local variables ==
67
68 integer bi,bj
69 integer biG,bjG
70 integer i,j
71 integer ii,k
72 integer icvar
73 integer icvrec
74 integer icvcomp
75 integer icvoffset
76 integer nopt
77 integer funit
78
79 integer cbuffindex
80 real*4 cbuff( sNx*nSx*nPx*sNy*nSy*nPy )
81
82 character*(128) fname
83 character*(17) prefix
84 parameter ( prefix = " OPTIM_READDATA: " )
85
86 c integer filei
87 c integer filej
88 c integer filek
89 c integer fileiG
90 c integer filejG
91 c integer filensx
92 c integer filensy
93 integer filenopt
94 _RL fileff
95
96 cgg(
97 _RL gg
98 integer igg
99 integer iobcs
100 cgg)
101
102 c == end of interface ==
103
104 c-- The reference i/o unit.
105 funit = 20
106
107 c-- Next optimization cycle.
108 nopt = optimcycle
109
110 if ( dfile .eq. ctrlname ) then
111 print*
112 print*,' OPTIM_READDATA: Reading control vector'
113 print*,' for optimization cycle: ',nopt
114 print*
115 else if ( dfile .eq. costname ) then
116 print*
117 print*,' OPTIM_READDATA: Reading cost function and'//
118 & ' gradient of cost function'
119 print*,' for optimization cycle: ',nopt
120 print*
121 else
122 print*
123 print*,' OPTIM_READDATA: subroutine called by a false *dfile*'
124 print*,' argument. *dfile* = ',dfile
125 print*
126 stop ' ... stopped in OPTIM_READDATA.'
127 endif
128
129 c-- Read the data.
130
131 bjG = 1 + (myygloballo - 1)/sny
132 biG = 1 + (myxgloballo - 1)/snx
133
134 c-- Generate file name and open the file.
135 write(fname(1:128),'(4a,i4.4)')
136 & dfile,'_',yctrlid(1:10),'.opt', nopt
137 open( funit, file = fname,
138 & status = 'old',
139 & form = 'unformatted',
140 & access = 'sequential' )
141 print*, prefix, 'opened file ', fname
142
143 c-- Read the header.
144 read( funit ) nvartype
145 read( funit ) nvarlength
146 read( funit ) yctrlid
147 read( funit ) filenopt
148 read( funit ) fileff
149 read( funit ) fileiG
150 read( funit ) filejG
151 read( funit ) filensx
152 read( funit ) filensy
153
154 read( funit ) (nWetcGlobal(k), k=1,nr)
155 read( funit ) (nWetsGlobal(k), k=1,nr)
156 read( funit ) (nWetwGlobal(k), k=1,nr)
157 #ifdef ALLOW_CTRL_WETV
158 read( funit ) (nWetvGlobal(k), k=1,nr)
159 #endif
160 #ifdef ALLOW_SHIFWFLX_CONTROL
161 read(funit) (nWetiGlobal(k), k=1,nr)
162 c read(funit) nWetiGlobal(1)
163 #endif
164
165 cgg( Add OBCS Mask information into the header section for optimization.
166 #ifdef ALLOW_OBCSN_CONTROL
167 read( funit ) ((nWetobcsnGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
168 #endif
169 #ifdef ALLOW_OBCSS_CONTROL
170 read( funit ) ((nWetobcssGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
171 #endif
172 #ifdef ALLOW_OBCSW_CONTROL
173 read( funit ) ((nWetobcswGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
174 #endif
175 #ifdef ALLOW_OBCSE_CONTROL
176 read( funit ) ((nWetobcseGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
177 #endif
178 cgg)
179 read( funit ) (ncvarindex(i), i=1,maxcvars)
180 read( funit ) (ncvarrecs(i), i=1,maxcvars)
181 read( funit ) (ncvarxmax(i), i=1,maxcvars)
182 read( funit ) (ncvarymax(i), i=1,maxcvars)
183 read( funit ) (ncvarnrmax(i), i=1,maxcvars)
184 read( funit ) (ncvargrd(i), i=1,maxcvars)
185 read( funit )
186
187 print *, prefix, 'nvartype ', nvartype
188 print *, prefix, 'nvarlength ', nvarlength
189 print *, prefix, 'yctrlid ', yctrlid
190 print *, prefix, 'filenopt ', filenopt
191 print *, prefix, 'fileff ', fileff
192 print *, prefix, 'fileiG ', fileiG
193 print *, prefix, 'filejG ', filejG
194 print *, prefix, 'filensx ', filensx
195 print *, prefix, 'filensy ', filensy
196
197 if (lheaderonly) then
198 print *, prefix, 'nWetcGlobal ', (nWetcGlobal(k), k=1,nr)
199 print *, prefix, 'nWetsGlobal ', (nWetsGlobal(k), k=1,nr)
200 print *, prefix, 'nWetwGlobal ', (nWetwGlobal(k), k=1,nr)
201 print *, prefix, 'nWetvGlobal ', (nWetvGlobal(k), k=1,nr)
202 #ifdef ALLOW_SHIFWFLX_CONTROL
203 print *, prefix, 'nWetiGlobal ', (nWetiGlobal(k), k=1,nr)
204 #endif
205 #ifdef ALLOW_OBCSN_CONTROL
206 do iobcs=1,nobcs
207 print *, prefix, 'nWetobcsnGlo (iobcs=', iobcs,')',
208 & (nWetobcsnGlo(k,iobcs), k=1,nr)
209 enddo
210 #endif
211 #ifdef ALLOW_OBCSS_CONTROL
212 do iobcs=1,nobcs
213 print *, prefix, 'nWetobcssGlo (iobcs=', iobcs,')',
214 & (nWetobcssGlo(k,iobcs), k=1,nr)
215 enddo
216 #endif
217 #ifdef ALLOW_OBCSW_CONTROL
218 do iobcs=1,nobcs
219 print *, prefix, 'nWetobcswGlo (iobcs=', iobcs,')',
220 & (nWetobcswGlo(k,iobcs), k=1,nr)
221 enddo
222 #endif
223 #ifdef ALLOW_OBCSE_CONTROL
224 do iobcs=1,nobcs
225 print *, prefix, 'nWetobcseGlo (iobcs=', iobcs,')',
226 & (nWetobcseGlo(k,iobcs), k=1,nr)
227 enddo
228 #endif
229 print *, prefix, 'ncvarindex ', (ncvarindex(i), i=1,maxcvars)
230 print *, prefix, 'ncvarrecs ', (ncvarrecs(i), i=1,maxcvars)
231 print *, prefix, 'ncvarxmax ', (ncvarxmax(i), i=1,maxcvars)
232 print *, prefix, 'ncvarymax ', (ncvarymax(i), i=1,maxcvars)
233 print *, prefix, 'ncvarnrmax ', (ncvarnrmax(i), i=1,maxcvars)
234 print *, prefix, 'ncvargrd ', (ncvargrd(i), i=1,maxcvars)
235 end if
236 c-- Check the header information for consistency.
237
238 cph if ( filenopt .ne. nopt ) then
239 cph print*
240 cph print*,' READ_HEADER: Input data belong to the wrong'
241 cph print*,' optimization cycle.'
242 cph print*,' optimization cycle = ',nopt
243 cph print*,' input optim cycle = ',filenopt
244 cph print*
245 cph stop ' ... stopped in READ_HEADER.'
246 cph endif
247
248 if ( (fileiG .ne. biG) .or. (filejG .ne. bjG) ) then
249 print*
250 print*, prefix, 'Tile indices of loop and data do not match.'
251 print*,' loop x/y component = ',
252 & biG,bjG
253 print*,' data x/y component = ',
254 & fileiG,filejG
255 print*
256 stop ' ... stopped in OPTIM_READDATA.'
257 endif
258
259 if ( (filensx .ne. nsx) .or. (filensy .ne. nsy) ) then
260 print*
261 print*, prefix, ' Numbers of tiles do not match.'
262 print*,' parameter x/y no. of tiles = ',
263 & bi,bj
264 print*,' data x/y no. of tiles = ',
265 & filensx,filensy
266 print*
267 stop ' ... stopped in OPTIM_READDATA.'
268 endif
269
270 ce Add some more checks. ...
271
272 if (.NOT. lheaderonly) then
273 c-- Read the data.
274 icvoffset = 0
275 do icvar = 1,maxcvars
276 if ( ncvarindex(icvar) .ne. -1 ) then
277 do icvrec = 1,ncvarrecs(icvar)
278 cph do bj = 1,nsy
279 cph do bi = 1,nsx
280 read( funit ) ncvarindex(icvar)
281 read( funit ) filej
282 read( funit ) filei
283 do k = 1,ncvarnrmax(icvar)
284 cbuffindex = 0
285 if (ncvargrd(icvar) .eq. 'c') then
286 cbuffindex = nWetcGlobal(k)
287 else if (ncvargrd(icvar) .eq. 's') then
288 cbuffindex = nWetsGlobal(k)
289 else if (ncvargrd(icvar) .eq. 'w') then
290 cbuffindex = nWetwGlobal(k)
291 else if (ncvargrd(icvar) .eq. 'v') then
292 cbuffindex = nWetvGlobal(k)
293 #ifdef ALLOW_SHIFWFLX_CONTROL
294 else if (ncvargrd(icvar) .eq. 'i') then
295 cbuffindex = nWetiGlobal(k)
296 #endif
297 cgg( O.B. points have the grid mask "m".
298 else if (ncvargrd(icvar) .eq. 'm') then
299 cgg From "icvrec", calculate what iobcs must be.
300 gg = (icvrec-1)/nobcs
301 igg = int(gg)
302 iobcs= icvrec - igg*nobcs
303 #ifdef ALLOW_OBCSN_CONTROL
304 if (icvar .eq. 11) cbuffindex = nWetobcsnGlo(k,iobcs)
305 #endif
306 #ifdef ALLOW_OBCSS_CONTROL
307 if (icvar .eq. 12) cbuffindex = nWetobcssGlo(k,iobcs)
308 #endif
309 #ifdef ALLOW_OBCSW_CONTROL
310 if (icvar .eq. 13) cbuffindex = nWetobcswGlo(k,iobcs)
311 #endif
312 #ifdef ALLOW_OBCSE_CONTROL
313 if (icvar .eq. 14) cbuffindex = nWetobcseGlo(k,iobcs)
314 #endif
315 cgg)
316 endif
317 if ( icvoffset + cbuffindex .gt. nvarlength ) then
318 print*
319 print *, ' ERROR:'
320 print *, ' There are at least ', icvoffset+cbuffindex,
321 & ' records in '//fname(1:28)//'.'
322 print *, ' This is more than expected from nvarlength =',
323 & nvarlength, '.'
324 print *, ' Something is wrong in the computation of '//
325 & 'the wet points or'
326 print *, ' in computing the number of records in '//
327 & 'some variable(s).'
328 print *, ' ... stopped in OPTIM_READDATA.'
329 stop ' ... stopped in OPTIM_READDATA.'
330 endif
331 if (cbuffindex .gt. 0) then
332 read( funit ) cbuffindex
333 read( funit ) filek
334 read( funit ) (cbuff(ii), ii=1,cbuffindex)
335 do icvcomp = 1,cbuffindex
336 vv(icvoffset+icvcomp) = cbuff(icvcomp)
337 c If you want to optimize with respect to just O.B. T and S
338 c uncomment the next two lines.
339 c if (iobcs .eq. 3) vv(icvoffset+icvcomp)=0.
340 c if (iobcs .eq. 4) vv(icvoffset+icvcomp)=0.
341 enddo
342 icvoffset = icvoffset + cbuffindex
343 endif
344 enddo
345 cph enddo
346 cph enddo
347 enddo
348 endif
349 enddo
350
351 else
352
353 c-- Assign the number of control variables.
354 nn = nvarlength
355
356 endif
357
358 close( funit )
359
360 c-- Assign the cost function value in case we read the cost file.
361
362 if ( dfile .eq. ctrlname ) then
363 ff = 0. d 0
364 else if ( dfile .eq. costname ) then
365 ff = fileff
366 endif
367 c-- Always return the cost function value if lheaderonly
368 if ( lheaderonly) ff = fileff
369
370 print *, prefix, 'end of optim_readdata'
371 print *, ' '
372
373 return
374 end

  ViewVC Help
Powered by ViewVC 1.1.22