c ECCO_CPPOPTIONS used to affect maxcvars
c and def ALLOW_OBCSN_CONTROL etc (OBCS masks etc)
c#include ECCO_CPPOPTIONS.h
c CTRL_OPTIONS affects maxcvars and may def
c ALLOW_OBCSN_CONTROL etc (OBCS masks etc)
#include "CTRL_OPTIONS.h"
subroutine OPTIM_READDATA(
I nn,
I dfile,
I lheaderonly,
O ff,
O vv
& )
c ==================================================================
c SUBROUTINE optim_readdata
c ==================================================================
c
c o Read the data written by the MITgcmUV state estimation setup and
c join them to one vector that is subsequently used by the minimi-
c zation algorithm "lsopt". Depending on the specified file name
c either the control vector or the gradient vector can be read.
c
c *dfile* should be the radix of the file: ecco_ctrl or ecco_cost
c
c started: Christian Eckert eckert@mit.edu 12-Apr-2000
c
c changed: Patrick Heimbach heimbach@mit.edu 19-Jun-2000
c - finished, revised and debugged
c
c ==================================================================
c SUBROUTINE optim_readdata
c ==================================================================
implicit none
c == global variables ==
#include "EEPARAMS.h"
#include "SIZE.h"
#include "ctrl.h"
#include "optim.h"
#include "minimization.h"
c == routine arguments ==
integer nn
_RL ff
#if defined (DYNAMIC)
_RL vv(nn)
#elif defined (USE_POINTER) (MAX_INDEPEND == 0)
_RL vv
pointer (pvv,vv(1))
#else
integer nmax
parameter( nmax = MAX_INDEPEND )
_RL vv(nmax)
#endif
character*(9) dfile
logical lheaderonly
c == local variables ==
integer bi,bj
integer biG,bjG
integer i,j
integer ii,k
integer icvar
integer icvrec
integer icvcomp
integer icvoffset
integer nopt
integer funit
integer cbuffindex
real*4 cbuff( sNx*nSx*nPx*sNy*nSy*nPy )
character*(128) fname
c integer filei
c integer filej
c integer filek
c integer fileiG
c integer filejG
c integer filensx
c integer filensy
integer filenopt
_RL fileff
cgg(
_RL gg
integer igg
integer iobcs
cgg)
c == end of interface ==
print *, 'pathei-lsopt in optim_readdata'
c-- The reference i/o unit.
funit = 20
c-- Next optimization cycle.
nopt = optimcycle
if ( dfile .eq. ctrlname ) then
print*
print*,' OPTIM_READDATA: Reading control vector'
print*,' for optimization cycle: ',nopt
print*
else if ( dfile .eq. costname ) then
print*
print*,' OPTIM_READDATA: Reading cost function and'
print*,' gradient of cost function'
print*,' for optimization cycle: ',nopt
print*
else
print*
print*,' OPTIM_READDATA: subroutine called by a false *dfile*'
print*,' argument. *dfile* = ',dfile
print*
stop ' ... stopped in OPTIM_READDATA.'
endif
c-- Read the data.
bjG = 1 + (myygloballo - 1)/sny
biG = 1 + (myxgloballo - 1)/snx
c-- Generate file name and open the file.
write(fname(1:128),'(4a,i4.4)')
& dfile,'_',yctrlid(1:10),'.opt', nopt
open( funit, file = fname,
& status = 'old',
& form = 'unformatted',
& access = 'sequential' )
print*, 'opened file ', fname
c-- Read the header.
read( funit ) nvartype
read( funit ) nvarlength
read( funit ) yctrlid
read( funit ) filenopt
read( funit ) fileff
read( funit ) fileiG
read( funit ) filejG
read( funit ) filensx
read( funit ) filensy
read( funit ) (nWetcGlobal(k), k=1,nr)
read( funit ) (nWetsGlobal(k), k=1,nr)
read( funit ) (nWetwGlobal(k), k=1,nr)
#ifdef ALLOW_CTRL_WETV
read( funit ) (nWetvGlobal(k), k=1,nr)
#endif
#ifdef ALLOW_SHIFWFLX_CONTROL
read(funit) (nWetiGlobal(k), k=1,nr)
c read(funit) nWetiGlobal(1)
#endif
cgg( Add OBCS Mask information into the header section for optimization.
#ifdef ALLOW_OBCSN_CONTROL
read( funit ) ((nWetobcsnGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
#endif
#ifdef ALLOW_OBCSS_CONTROL
read( funit ) ((nWetobcssGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
#endif
#ifdef ALLOW_OBCSW_CONTROL
read( funit ) ((nWetobcswGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
#endif
#ifdef ALLOW_OBCSE_CONTROL
read( funit ) ((nWetobcseGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
#endif
cgg)
read( funit ) (ncvarindex(i), i=1,maxcvars)
read( funit ) (ncvarrecs(i), i=1,maxcvars)
read( funit ) (ncvarxmax(i), i=1,maxcvars)
read( funit ) (ncvarymax(i), i=1,maxcvars)
read( funit ) (ncvarnrmax(i), i=1,maxcvars)
read( funit ) (ncvargrd(i), i=1,maxcvars)
read( funit )
cph(
cph if (lheaderonly) then
print *, 'pathei: nvartype ', nvartype
print *, 'pathei: nvarlength ', nvarlength
print *, 'pathei: yctrlid ', yctrlid
print *, 'pathei: filenopt ', filenopt
print *, 'pathei: fileff ', fileff
print *, 'pathei: fileiG ', fileiG
print *, 'pathei: filejG ', filejG
print *, 'pathei: filensx ', filensx
print *, 'pathei: filensy ', filensy
print *, 'pathei: nWetcGlobal ',
& (nWetcGlobal(k), k=1,nr)
print *, 'pathei: nWetsGlobal ',
& (nWetsGlobal(k), k=1,nr)
print *, 'pathei: nWetwGlobal ',
& (nWetwGlobal(k), k=1,nr)
print *, 'pathei: nWetvGlobal ',
& (nWetvGlobal(k), k=1,nr)
#ifdef ALLOW_SHIFWFLX_CONTROL
print *, 'pathei: nWetiGlobal ',
& (nWetiGlobal(k), k=1,nr)
#endif
#ifdef ALLOW_OBCSN_CONTROL
do iobcs=1,nobcs
print *, 'pathei: nWetobcsnGlo (iobcs=', iobcs,')',
& (nWetobcsnGlo(k,iobcs), k=1,nr)
enddo
#endif
#ifdef ALLOW_OBCSS_CONTROL
do iobcs=1,nobcs
print *, 'pathei: nWetobcssGlo (iobcs=', iobcs,')',
& (nWetobcssGlo(k,iobcs), k=1,nr)
enddo
#endif
#ifdef ALLOW_OBCSW_CONTROL
do iobcs=1,nobcs
print *, 'pathei: nWetobcswGlo (iobcs=', iobcs,')',
& (nWetobcswGlo(k,iobcs), k=1,nr)
enddo
#endif
#ifdef ALLOW_OBCSE_CONTROL
do iobcs=1,nobcs
print *, 'pathei: nWetobcseGlo (iobcs=', iobcs,')',
& (nWetobcseGlo(k,iobcs), k=1,nr)
enddo
#endif
print *, 'pathei: ncvarindex ',
& (ncvarindex(i), i=1,maxcvars)
print *, 'pathei: ncvarrecs ',
& (ncvarrecs(i), i=1,maxcvars)
print *, 'pathei: ncvarxmax ',
& (ncvarxmax(i), i=1,maxcvars)
print *, 'pathei: ncvarymax ',
& (ncvarymax(i), i=1,maxcvars)
print *, 'pathei: ncvarnrmax ',
& (ncvarnrmax(i), i=1,maxcvars)
print *, 'pathei: ncvargrd ',
& (ncvargrd(i), i=1,maxcvars)
cph end if
cph)
c-- Check the header information for consistency.
cph if ( filenopt .ne. nopt ) then
cph print*
cph print*,' READ_HEADER: Input data belong to the wrong'
cph print*,' optimization cycle.'
cph print*,' optimization cycle = ',nopt
cph print*,' input optim cycle = ',filenopt
cph print*
cph stop ' ... stopped in READ_HEADER.'
cph endif
if ( (fileiG .ne. biG) .or. (filejG .ne. bjG) ) then
print*
print*,' READ_HEADER: Tile indices of loop and data '
print*,' do not match.'
print*,' loop x/y component = ',
& biG,bjG
print*,' data x/y component = ',
& fileiG,filejG
print*
stop ' ... stopped in READ_HEADER.'
endif
if ( (filensx .ne. nsx) .or. (filensy .ne. nsy) ) then
print*
print*,' READ_HEADER: Numbers of tiles do not match.'
print*,' parameter x/y no. of tiles = ',
& bi,bj
print*,' data x/y no. of tiles = ',
& filensx,filensy
print*
stop ' ... stopped in READ_HEADER.'
endif
ce Add some more checks. ...
if (.NOT. lheaderonly) then
c-- Read the data.
icvoffset = 0
do icvar = 1,maxcvars
if ( ncvarindex(icvar) .ne. -1 ) then
do icvrec = 1,ncvarrecs(icvar)
cph do bj = 1,nsy
cph do bi = 1,nsx
read( funit ) ncvarindex(icvar)
read( funit ) filej
read( funit ) filei
do k = 1,ncvarnrmax(icvar)
cbuffindex = 0
if (ncvargrd(icvar) .eq. 'c') then
cbuffindex = nWetcGlobal(k)
else if (ncvargrd(icvar) .eq. 's') then
cbuffindex = nWetsGlobal(k)
else if (ncvargrd(icvar) .eq. 'w') then
cbuffindex = nWetwGlobal(k)
else if (ncvargrd(icvar) .eq. 'v') then
cbuffindex = nWetvGlobal(k)
#ifdef ALLOW_SHIFWFLX_CONTROL
else if (ncvargrd(icvar) .eq. 'i') then
cbuffindex = nWetiGlobal(k)
#endif
cgg( O.B. points have the grid mask "m".
else if (ncvargrd(icvar) .eq. 'm') then
cgg From "icvrec", calculate what iobcs must be.
gg = (icvrec-1)/nobcs
igg = int(gg)
iobcs= icvrec - igg*nobcs
#ifdef ALLOW_OBCSN_CONTROL
if (icvar .eq. 11) then
cbuffindex = nWetobcsnGlo(k,iobcs)
endif
#endif
#ifdef ALLOW_OBCSS_CONTROL
if (icvar .eq. 12) then
cbuffindex = nWetobcssGlo(k,iobcs)
endif
#endif
#ifdef ALLOW_OBCSW_CONTROL
if (icvar .eq. 13) then
cbuffindex = nWetobcswGlo(k,iobcs)
endif
#endif
#ifdef ALLOW_OBCSE_CONTROL
if (icvar .eq. 14) then
cbuffindex = nWetobcseGlo(k,iobcs)
endif
#endif
cgg)
endif
if ( icvoffset + cbuffindex .gt. nvarlength ) then
print*
print *, ' ERROR:'
print *, ' There are at least ', icvoffset+cbuffindex,
& ' records in '//fname(1:28)//'.'
print *, ' This is more than expected from nvarlength =',
& nvarlength, '.'
print *, ' Something is wrong in the computation of '//
& 'the wet points or'
print *, ' in computing the number of records in '//
& 'some variable(s).'
print *, ' ... stopped in OPTIM_READDATA.'
stop ' ... stopped in OPTIM_READDATA.'
endif
if (cbuffindex .gt. 0) then
read( funit ) cbuffindex
read( funit ) filek
read( funit ) (cbuff(ii), ii=1,cbuffindex)
do icvcomp = 1,cbuffindex
vv(icvoffset+icvcomp) = cbuff(icvcomp)
c If you want to optimize with respect to just O.B. T and S
c uncomment the next two lines.
c if (iobcs .eq. 3) vv(icvoffset+icvcomp)=0.
c if (iobcs .eq. 4) vv(icvoffset+icvcomp)=0.
enddo
icvoffset = icvoffset + cbuffindex
endif
enddo
cph enddo
cph enddo
enddo
endif
enddo
else
c-- Assign the number of control variables.
nn = nvarlength
endif
close( funit )
c-- Assign the cost function value in case we read the cost file.
if ( dfile .eq. ctrlname ) then
ff = 0. d 0
else if ( dfile .eq. costname ) then
ff = fileff
endif
return
end