/[MITgcm]/MITgcm/optim/optim_readdata.F
ViewVC logotype

Diff of /MITgcm/optim/optim_readdata.F

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

revision 1.1 by heimbach, Tue Feb 5 20:34:35 2002 UTC revision 1.1.2.1 by heimbach, Tue Feb 5 20:34:35 2002 UTC
# Line 0  Line 1 
1    
2          subroutine optim_readdata(
3         I                      nn,
4         I                      dfile,
5         I                      lheaderonly,
6         O                      ff,
7         O                      vv
8         &                    )
9    
10    c     ==================================================================
11    c     SUBROUTINE optim_readdata
12    c     ==================================================================
13    c
14    c     o Read the data written by the MITgcmUV state estimation setup and
15    c       join them to one vector that is subsequently used by the minimi-
16    c       zation algorithm "lsopt". Depending on the specified file name
17    c       either the control vector or the gradient vector can be read.
18    c
19    c       *dfile* should be the radix of the file: ecco_ctrl or ecco_cost
20    c
21    c     started: Christian Eckert eckert@mit.edu 12-Apr-2000
22    c
23    c     changed:  Patrick Heimbach heimbach@mit.edu 19-Jun-2000
24    c               - finished, revised and debugged
25    c
26    c     ==================================================================
27    c     SUBROUTINE optim_readdata
28    c     ==================================================================
29    
30          implicit none
31    
32    c     == global variables ==
33    
34    #include "EEPARAMS.h"
35    #include "SIZE.h"
36    
37    #include "ecco.h"
38    #include "ctrl.h"
39    #include "optim.h"
40    #include "minimization.h"
41    
42    c     == routine arguments ==
43    
44          integer nn
45          _RL     ff
46    
47    #if defined (DYNAMIC)
48          _RL     vv(nn)
49    #elif defined (USE_POINTER) || (MAX_INDEPEND == 0)
50          _RL            vv
51          pointer (pvv,vv(1))
52    #else
53          integer nmax
54          parameter( nmax = MAX_INDEPEND )
55          _RL   vv(nmax)
56    #endif
57    
58          character*(9) dfile
59          logical lheaderonly
60    
61    c     == local variables ==
62    
63          integer bi,bj
64          integer biG,bjG
65          integer i,j
66          integer ii,k
67          integer icvar
68          integer icvrec
69          integer icvcomp
70          integer icvoffset
71          integer nopt
72          integer funit
73    
74          integer cbuffindex
75          _RL     cbuff( sNx*nSx*nPx*sNy*nSy*nPy )
76    
77          character*(128) fname
78    
79          integer         filei
80          integer         filej
81          integer         filek
82          integer         filenopt
83          integer         fileig
84          integer         filejg
85          integer         filensx
86          integer         filensy
87          _RL             fileff
88    
89    c     == end of interface ==
90    
91    c--   The reference i/o unit.
92          funit = 20
93    
94    c--   Next optimization cycle.
95          nopt = optimcycle
96    
97          if      ( dfile .eq. ctrlname ) then
98            print*
99            print*,' OPTIM_READDATA: Reading control vector'
100            print*,'            for optimization cycle: ',nopt
101            print*
102          else if ( dfile .eq. costname ) then
103            print*
104            print*,' OPTIM_READDATA: Reading cost function and'
105            print*,'            gradient of cost function'
106            print*,'            for optimization cycle: ',nopt
107            print*
108          else
109            print*
110            print*,' OPTIM_READDATA: subroutine called by a false *dfile*'
111            print*,'            argument. *dfile* = ',dfile
112            print*
113            stop   '  ...  stopped in OPTIM_READDATA.'
114          endif
115    
116    c--   Read the data.
117    
118          bjG = 1 + (myygloballo - 1)/sny
119          biG = 1 + (myxgloballo - 1)/snx
120    
121    c--   Generate file name and open the file.
122          write(fname(1:128),'(4a,i4.4)')
123         &     dfile,'_',expId(1:10),'.opt', nopt
124          open( funit, file   = fname,
125         &     status = 'old',
126         &     form   = 'unformatted',
127         &     access = 'sequential'   )
128          print*, 'opened file ', fname
129    
130    c--   Read the header.
131          read( funit ) nvartype
132          read( funit ) nvarlength
133          read( funit ) expId
134          read( funit ) filenopt
135          read( funit ) fileff
136          read( funit ) fileiG
137          read( funit ) filejG
138          read( funit ) filensx
139          read( funit ) filensy
140    
141    cph(
142          print *,'ph-opt 1 ', nvartype, nvarlength, filensx, filensy
143    cph)
144    
145          read( funit ) (((nWetcTile(i,j,k), i=1,nsx), j=1,nsy),
146         &     k=1,nr)
147          read( funit ) (((nWetsTile(i,j,k), i=1,nsx), j=1,nsy),
148         &     k=1,nr)
149          read( funit ) (((nWetwTile(i,j,k), i=1,nsx), j=1,nsy),
150         &     k=1,nr)
151          read( funit ) (ncvarindex(i), i=1,maxcvars)
152          read( funit ) (ncvarrecs(i),  i=1,maxcvars)
153          read( funit ) (ncvarxmax(i),  i=1,maxcvars)
154          read( funit ) (ncvarymax(i),  i=1,maxcvars)
155          read( funit ) (ncvarnrmax(i), i=1,maxcvars)
156          read( funit ) (ncvargrd(i),   i=1,maxcvars)
157          read( funit )
158    
159    cph(
160    cph      if (lheaderonly) then
161             print *, 'pathei: nvartype ', nvartype
162             print *, 'pathei: nvarlength ', nvarlength
163             print *, 'pathei: expId ', expId
164             print *, 'pathei: filenopt ', filenopt
165             print *, 'pathei: fileff ', fileff
166             print *, 'pathei: fileiG ', fileiG
167             print *, 'pathei: filejG ', filejG
168             print *, 'pathei: filensx ', filensx
169             print *, 'pathei: filensy ', filensy
170            
171             print *, 'pathei: nWetcTile ',
172         &        (((nWetcTile(i,j,k), i=1,nsx), j=1,nsy), k=1,nr)
173             print *, 'pathei: nWetsTile ',
174         &        (((nWetsTile(i,j,k), i=1,nsx), j=1,nsy), k=1,nr)
175             print *, 'pathei: nWetwTile ',
176         &        (((nWetwTile(i,j,k), i=1,nsx), j=1,nsy), k=1,nr)
177             print *, 'pathei: ncvarindex ',
178         &        (ncvarindex(i), i=1,maxcvars)
179             print *, 'pathei: ncvarrecs ',
180         &        (ncvarrecs(i),  i=1,maxcvars)
181             print *, 'pathei: ncvarxmax ',
182         &        (ncvarxmax(i),  i=1,maxcvars)
183             print *, 'pathei: ncvarymax ',
184         &        (ncvarymax(i),  i=1,maxcvars)
185             print *, 'pathei: ncvarnrmax ',
186         &        (ncvarnrmax(i), i=1,maxcvars)
187             print *, 'pathei: ncvargrd ',
188         &        (ncvargrd(i),   i=1,maxcvars)
189    cph      end if
190    cph)
191    c--   Check the header information for consistency.
192    
193    cph      if ( filenopt .ne. nopt ) then
194    cph         print*
195    cph         print*,' READ_HEADER: Input data belong to the wrong'
196    cph         print*,'              optimization cycle.'
197    cph         print*,'              optimization cycle = ',nopt
198    cph         print*,'              input optim  cycle = ',filenopt
199    cph         print*
200    cph         stop   ' ... stopped in READ_HEADER.'
201    cph      endif
202          
203          if ( (fileiG .ne. biG) .or. (filejG .ne. bjG) ) then
204             print*
205             print*,' READ_HEADER: Tile indices of loop and data '
206             print*,'              do not match.'
207             print*,'              loop x/y component = ',
208         &        biG,bjG
209             print*,'              data x/y component = ',
210         &        fileiG,filejG
211             print*
212             stop   ' ... stopped in READ_HEADER.'
213          endif
214          
215          if ( (filensx .ne. nsx) .or. (filensy .ne. nsy) ) then
216             print*
217             print*,' READ_HEADER: Numbers of tiles do not match.'
218             print*,'              parameter x/y no. of tiles = ',
219         &        bi,bj
220             print*,'              data      x/y no. of tiles = ',
221         &        filensx,filensy
222             print*
223             stop   ' ... stopped in READ_HEADER.'
224          endif
225    
226    ce    Add some more checks. ...
227    
228          if (.NOT. lheaderonly) then
229    c--   Read the data.
230             icvoffset = 0
231             do icvar = 1,maxcvars
232                if ( ncvarindex(icvar) .ne. -1 ) then
233                   do icvrec = 1,ncvarrecs(icvar)
234                      do bj = 1,nsy
235                         do bi = 1,nsx
236                            read( funit ) ncvarindex(icvar)
237                            read( funit ) filej
238                            read( funit ) filei
239                            do k = 1,ncvarnrmax(icvar)
240                               cbuffindex = 0
241                               if (ncvargrd(icvar) .eq. 'c') then
242                                  cbuffindex = nwetctile(bi,bj,k)
243                               else if (ncvargrd(icvar) .eq. 's') then
244                                  cbuffindex = nwetstile(bi,bj,k)
245                               else if (ncvargrd(icvar) .eq. 'w') then
246                                  cbuffindex = nwetwtile(bi,bj,k)
247                               endif
248                               if (cbuffindex .gt. 0) then
249                                  read( funit ) cbuffindex
250                                  read( funit ) filek
251                                  read( funit ) (cbuff(ii), ii=1,cbuffindex)
252                                  do icvcomp = 1,cbuffindex
253                                     vv(icvoffset+icvcomp) = cbuff(icvcomp)
254                                  enddo
255                                  icvoffset = icvoffset + cbuffindex
256                               endif
257                            enddo
258                         enddo
259                      enddo
260                   enddo
261                endif
262             enddo
263    
264          else
265    
266    c--   Assign the number of control variables.
267             nn = nvarlength
268            
269          endif
270    
271          close( funit )
272    
273    c--   Assign the cost function value in case we read the cost file.
274    
275          if      ( dfile .eq. ctrlname ) then
276            ff = 0. d 0
277          else if ( dfile .eq. costname ) then
278            ff = fileff
279          endif
280    
281          return
282          end

Legend:
Removed from v.1.1  
changed lines
  Added in v.1.1.2.1

  ViewVC Help
Powered by ViewVC 1.1.22