C $Header: /u/gcmpack/MITgcm/pkg/fizhi/update_ocean_exports.F,v 1.33 2010/09/22 22:21:41 jmc Exp $
C $Name: $
#include "FIZHI_OPTIONS.h"
subroutine UPDATE_OCEAN_EXPORTS (myTime, myIter, myThid)
c----------------------------------------------------------------------
c Subroutine update_ocean_exports - 'Wrapper' routine to update
c the fields related to the ocean surface that are needed
c by fizhi (sst and sea ice extent).
c
c Call: getsst (Return the current sst field-read dataset if needed)
c getsice (Return the current sea ice field-read data if needed)
c-----------------------------------------------------------------------
implicit none
#include "SIZE.h"
#include "GRID.h"
#include "fizhi_ocean_coms.h"
#include "EEPARAMS.h"
#include "chronos.h"
#ifdef ALLOW_EXCH2
#include "W2_EXCH2_SIZE.h"
#include "W2_EXCH2_TOPOLOGY.h"
#endif /* ALLOW_EXCH2 */
integer myIter, myThid
_RL myTime
INTEGER xySize
#if defined(ALLOW_EXCH2)
PARAMETER ( xySize = W2_ioBufferSize )
#else
PARAMETER ( xySize = Nx*Ny )
#endif
integer i, j, bi, bj, bislot, bjslot
integer im1, im2, jm1, jm2, idim1, idim2, jdim1, jdim2
integer xsize, ysize
_RL sstmin
parameter ( sstmin = 273.16 )
_RL sst1 (xySize), sst2 (xySize)
_RL sice1(xySize), sice2(xySize)
c _RL sst1(xsize,ysize),sst2(xsize,ysize)
c _RL sice1(xsize,ysize),sice2(xsize,ysize)
integer nymd1sst(nSx,nSy),nymd2sst(nSx,nSy)
integer nymd1sice(nSx,nSy),nymd2sice(nSx,nSy)
integer nhms1sst(nSx,nSy),nhms2sst(nSx,nSy)
integer nhms1sice(nSx,nSy),nhms2sice(nSx,nSy)
integer sstdates(370,nSx,nSy),sicedates(370,nSx,nSy)
integer ssttimes(370,nSx,nSy),sicetimes(370,nSx,nSy)
logical first(nSx,nSy)
integer nSxnSy
parameter(nSxnSy = nSx*nSy)
data first/nSxnSy*.true./
save nymd1sst,nymd2sst,nymd1sice,nymd2sice
save nhms1sst,nhms2sst,nhms1sice,nhms2sice
save sst1, sst2, sice1, sice2
save sstdates, sicedates
save ssttimes, sicetimes
#if defined(ALLOW_EXCH2)
xsize = exch2_global_Nx
ysize = exch2_global_Ny
#else
xsize = Nx
ysize = Ny
#endif
idim1 = 1-OLx
idim2 = sNx+OLx
jdim1 = 1-OLy
jdim2 = sNy+OLy
im1 = 1
im2 = sNx
jm1 = 1
jm2 = sNy
C***********************************************************************
DO BJ = myByLo(myThid),myByHi(myThid)
DO BI = myBxLo(myThid),myBxHi(myThid)
#if defined(ALLOW_EXCH2)
bislot = exch2_txglobalo(W2_myTileList(bi,bj))-1
bjslot = exch2_tyglobalo(W2_myTileList(bi,bj))-1
#else
bislot = myXGlobalLo-1+(bi-1)*sNx
bjslot = myYGlobalLo-1+(bj-1)*sNy
#endif
call GETSST(ksst,sstclim,idim1,idim2,jdim1,jdim2,im1,im2,
. jm1,jm2,nSx,nSy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,
. sst1,sst2,first(bi,bj),nymd1sst(bi,bj),nymd2sst(bi,bj),
. nhms1sst(bi,bj),nhms2sst(bi,bj),sstdates(1,bi,bj),
. ssttimes(1,bi,bj),sst,myThid)
call GETSICE(kice,siceclim,idim1,idim2,jdim1,jdim2,im1,im2,
. jm1,jm2,nSx,nSy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,
. sice1,sice2,first(bi,bj),nymd1sice(bi,bj),nymd2sice(bi,bj),
. nhms1sice(bi,bj),nhms2sice(bi,bj),sicedates(1,bi,bj),
. sicetimes(1,bi,bj),sice,myThid)
c Check for Minimum Open-Water SST
c --------------------------------
do j=jm1,jm2
do i=im1,im2
if(sice(i,j,bi,bj).eq.0.0 .and. sst(i,j,bi,bj).lt.sstmin)
. sst(i,j,bi,bj) = sstmin
enddo
enddo
ENDDO
ENDDO
_EXCH_XY_RL(sst,myThid)
_EXCH_XY_RL(sice,myThid)
return
end
subroutine GETSICE(iunit,clim,idim1,idim2,jdim1,jdim2,im1,im2,
. jm1,jm2,nSumx,nSumy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,
. sicebc1,sicebc2,first,nymdbc1,nymdbc2,nhmsbc1,nhmsbc2,
. nymdbc,nhmsbc,sice,mythid)
C***********************************************************************
C
C!ROUTINE: GETSICE
C!DESCRIPTION: GETSICE returns the sea ice depth.
C! This routine is adaptable for any frequency
C! data upto a daily frequency.
C! note: for diurnal data ndmax should be increased.
C
C!INPUT PARAMETERS:
C! iunit Unit number assigned to the sice data file
C! idim1 Start dimension in x-direction
C! idim2 End dimension in x-direction
C! jdim1 Start dimension in y-direction
C! jdim2 End dimension in y-direction
C! im1 Begin of x-direction span for filling sice
C! im2 End of x-direction span for filling sice
C! jm1 Begin of y-direction span for filling sice
C! jm2 End of y-direction span for filling sice
C! nSumx Number of processors in x-direction (local processor)
C! nSumy Number of processors in y-direction (local processor)
C! xsize Number of processors in x-direction (global)
C! ysize Number of processors in y-direction (global)
C! bi Processor number in x-direction (local to processor)
C! bj Processor number in y-direction (local to processor)
C! bislot Processor number in x-direction (global)
C! bjslot Processor number in y-direction (global)
C! nymd YYMMDD of the current model timestep
C! nhms HHMMSS of the model time
C
C!OUTPUT PARAMETERS:
C! sice(idim1:idim2,jdim1:jdim2,nSumx,nSumy) Sea ice depth in meters
C
C!ROUTINES CALLED:
C
C! bcdata Reads the data for a given unit number
C! bcheader Reads the header info for a given unit number
C! interp_time Returns weights for linear interpolation
C
C--------------------------------------------------------------------------
implicit none
#include "SIZE.h"
integer iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSumx,nSumy
integer xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,mythid
logical clim
_RL sicebc1(xsize,ysize)
_RL sicebc2(xsize,ysize)
_RL sice(idim1:idim2,jdim1:jdim2,nSumx,nSumy)
integer nhmsbc1,nhmsbc2,nymdbc1,nymdbc2
logical first
C Maximum number of dates in one year for the data
integer ndmax
parameter (ndmax = 370)
integer nymdbc(ndmax),nhmsbc(ndmax)
character*8 cname
character*80 cdscrip
character*22 sicedata
_RL fac1, fac2, lat0, lon0, timebc1, timebc2, timemod, undef
logical found, error
integer i,j,n,nn,iyear,iyearbc,nd,ndby3,imbc,jmbc
integer ndatebc,nrec
integer nymdmod
C--------- Variable Initialization ---------------------------------
data error /.false./
c save header info
save imbc,jmbc,lat0,lon0,ndatebc,undef
c this only works for between 1950-2050
if (nymd .lt. 500101) then
nymdmod = 20000000 + nymd
else if (nymd .le. 991231) then
nymdmod = 19000000 + nymd
else
nymdmod = nymd
endif
iyear = nymdmod/10000
if(clim) then
if(xsize.eq.192)sicedata='sice19232.weekly.clim'
if(xsize.eq.612)sicedata='sice612102.weekly.clim'
else
if(xsize.eq.192)
. WRITE(sicedata,'(A,I4)')'sice19232.weekly.y',iyear
if(xsize.eq.612)
. WRITE(sicedata,'(A,I4)')'sice612102.weekly.y',iyear
endif
c initialize so that first time through they have values for the check
c these values make the iyear .ne. iyearbc true anyways for
c for the first time so first isnt checked below.
if (first) then
nymdbc(2) = 0
nymdbc1 = 0
nymdbc2 = 0
nhmsbc1 = 0
nhmsbc2 = 0
first = .false.
endif
C---------- Read in Header file ----------------------------------
iyearbc = nymdbc(2)/10000
if( iyear.ne.iyearbc ) then
close(iunit)
open (iunit,file=sicedata,form='unformatted',access='direct',
. recl=xsize*ysize*4)
nrec = 1
call BCHEADER (iunit, ndmax, nrec,
. cname, cdscrip, imbc, jmbc, lat0, lon0,
. ndatebc, nymdbc, nhmsbc, undef, error)
C--------- Check data for Compatibility ------------------------------
C Check for correct data in boundary condition file
if (.not.error .and. cname.ne.'SICE') then
write(6,*)'Wrong data in SICE boundary condition file => ',cname
error = .true.
endif
C Check Horizontal Resolution
if(.not.error.and.imbc*jmbc.ne.xsize*ysize)then
write(6,*) ' B.C. Resolution DOES NOT match Model Resolution!'
write(6,*) ' B.C. Resolution: ',imbc*jmbc
write(6,*) 'Model Resolution: ',xsize*ysize
error = .true.
endif
C Check Year
iyearbc = nymdbc(2)/10000
if (.not.error .and. iyear.ne.iyearbc .and. iyearbc.ne.0) then
write(6,*)' B.C. Year DOES NOT match REQUESTED Year!'
write(6,*)' B.C. Year: ', iyearbc
write(6,*)'Requested Year: ', iyear
error = .true.
endif
if (.not.error) then
C if climatology, fill dates for data with current model year
if (iyearbc.eq.0) then
write(6,*)
write(6,*) 'Climatological Dataset is being used.'
write(6,*) 'Current model year to be used to fill Header Dates'
do n = 2, ndatebc-1
nymdbc(n) = nymdbc(n) +(nymdmod/10000)*10000
enddo
C For the first date subtract 1 year from the current model NYMD
n = 1
nymdbc(n) = nymdbc(n) +(nymdmod/10000-1)*10000
C For the last date add 1 year to the current model NYMD
n = ndatebc
nymdbc(n) = nymdbc(n) +(nymdmod/10000+1)*10000
endif
C Write out header info
_BEGIN_MASTER( myThid )
write(6,*) ' Updated boundary condition data'
write(6,*) ' ---------------------------------'
write(6,*) ' Variable: ',cname
write(6,*) ' Description: ',cdscrip
write(6,*) ' Resolution: x= ',imbc,' y= ',jmbc,
. ' Undefined value = ',undef
write(6,*) ' Starting latitude = ',lat0,' longitude =',lon0
write(6,*) ' Data valid at these times: '
ndby3 = ndatebc/3
do n = 1, ndby3*3,3
write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=n,n+2)
1000 format(3(2x,i3,':',i8,2x,i8))
enddo
write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=ndby3*3+1,ndatebc)
_END_MASTER( myThid )
endif
endif
C---------- Read sice data if necessary -------------------------------
found = .false.
nd = 2
c If model time is not within the times of saved sice data
c from previous call to getsice then read new data
timemod = dfloat(nymdmod) + dfloat(nhms) /1000000
timebc1 = dfloat(nymdbc1) + dfloat(nhmsbc1)/1000000
timebc2 = dfloat(nymdbc2) + dfloat(nhmsbc2)/1000000
if (timemod .lt. timebc1 .or. timemod .ge. timebc2) then
do while (.not.found .and. nd .le. ndatebc)
timebc2 = dfloat(nymdbc(nd)) + dfloat(nhmsbc(nd))/1000000
if (timebc2 .gt. timemod) then
nymdbc1 = nymdbc(nd-1)
nymdbc2 = nymdbc(nd)
nhmsbc1 = nhmsbc(nd-1)
nhmsbc2 = nhmsbc(nd)
call BCDATA (iunit,imbc,jmbc,nd,nd+1,sicebc1,sicebc2)
found = .true.
else
nd = nd + 1
endif
enddo
c Otherwise the data from the last time in getsice surrounds the
c current model time.
else
found = .true.
endif
if (.not.found) then
print *, 'STOP: Could not find SICE dates for model time.'
call MY_FINALIZE
call MY_EXIT (101)
endif
C---------- Interpolate sice data ------------------------------------
call INTERP_TIME(nymdmod,nhms,nymdbc1,nhmsbc1,nymdbc2,nhmsbc2,
. fac1,fac2)
do j = jm1,jm2
do i = im1,im2
sice(i,j,bi,bj) = sicebc1(i+bislot,j+bjslot)*fac1
. + sicebc2(i+bislot,j+bjslot)*fac2
c average to 0 or 1
c -----------------
if (sice(i,j,bi,bj) .ge. 0.5) then
sice(i,j,bi,bj) = 1.
else
sice(i,j,bi,bj) = 0.
endif
enddo
enddo
C---------- Fill sice with depth of ice ------------------------------------
do j = jm1,jm2
do i = im1,im2
if (sice(i,j,bi,bj) .eq. 1.) then
sice(i,j,bi,bj) = 3.
endif
enddo
enddo
C---------------------------------------------------------------------------
return
end
subroutine GETSST(iunit,clim,idim1,idim2,jdim1,jdim2,im1,im2,
. jm1,jm2,nSumx,nSumy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,
. sstbc1,sstbc2,first,nymdbc1,nymdbc2,nhmsbc1,nhmsbc2,
. nymdbc,nhmsbc,sst,mythid)
C***********************************************************************
C
C!ROUTINE: GETSST
C!DESCRIPTION: GETSST gets the SST data.
C! This routine is adaptable for any frequency
C! data upto a daily frequency.
C! note: for diurnal data ndmax should be increased.
C
C!INPUT PARAMETERS:
C! iunit Unit number assigned to the sice data file
C! idim1 Start dimension in x-direction
C! idim2 End dimension in x-direction
C! jdim1 Start dimension in y-direction
C! jdim2 End dimension in y-direction
C! im1 Begin of x-direction span for filling sst
C! im2 End of x-direction span for filling sst
C! jm1 Begin of y-direction span for filling sst
C! jm2 End of y-direction span for filling sst
C! nSumx Number of processors in x-direction (local processor)
C! nSumy Number of processors in y-direction (local processor)
C! xsize x-dimension of global array
C! ysize y-dimension of global array
C! bi Processor number in x-direction (local to processor)
C! bj Processor number in y-direction (local to processor)
C! bislot Slot number into global array in x-direction (global)
C! bjslot Slot number into global array in y-direction (global)
C! nymd YYMMDD of the current model timestep
C! nhms HHMMSS of the model time
C
C!OUTPUT PARAMETERS:
C! sst(idim1:idim2,jdim1:jdim2,nSumx,nSumy) Sea surface temperature (K)
C
C!ROUTINES CALLED:
C
C! bcdata Reads the data for a given unit number
C! bcheader Reads the header info for a given unit number
C! interp_time Returns weights for linear interpolation
C
C--------------------------------------------------------------------------
implicit none
#include "SIZE.h"
integer iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSumx,nSumy
integer xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,mythid
logical clim
_RL sstbc1(xsize,ysize)
_RL sstbc2(xsize,ysize)
_RL sst(idim1:idim2,jdim1:jdim2,nSumx,nSumy)
integer nhmsbc1,nhmsbc2,nymdbc1,nymdbc2
logical first
C Maximum number of dates in one year for the data
integer ndmax
parameter (ndmax = 370)
integer nymdbc(ndmax),nhmsbc(ndmax)
character*8 cname
character*80 cdscrip
character*21 sstdata
_RL fac1, fac2, lat0, lon0, timebc1, timebc2, timemod, undef
logical found, error
integer i,j,n,nn,iyear,iyearbc,nd,ndby3,imbc,jmbc
integer ndatebc,nrec
integer nymdmod
C--------- Variable Initialization ---------------------------------
data error /.false./
c save header info
save imbc,jmbc,lat0,lon0,ndatebc,undef
c this only works for between 1950-2050
if (nymd .lt. 500101) then
nymdmod = 20000000 + nymd
else if (nymd .le. 991231) then
nymdmod = 19000000 + nymd
else
nymdmod = nymd
endif
iyear = nymdmod/10000
if(clim) then
if(xsize.eq.192)sstdata='sst19232.weekly.clim'
if(xsize.eq.612)sstdata='sst612102.weekly.clim'
else
if(xsize.eq.192)
. WRITE(sstdata,'(A,I4)')'sst19232.weekly.y',iyear
if(xsize.eq.612)
. WRITE(sstdata,'(A,I4)')'sst612102.weekly.y',iyear
endif
c initialize so that first time through they have values for the check
c these vaules make the iyear .ne. iyearbc true anyways for
c for the first time so first isnt checked below.
if (first) then
nymdbc(2) = 0
nymdbc1 = 0
nymdbc2 = 0
nhmsbc1 = 0
nhmsbc2 = 0
first = .false.
endif
C---------- Read in Header file ----------------------------------
iyearbc = nymdbc(2)/10000
if( iyear.ne.iyearbc ) then
close(iunit)
open (iunit,file=sstdata,form='unformatted',access='direct',
. recl=xsize*ysize*4)
nrec = 1
call BCHEADER (iunit, ndmax, nrec,
. cname, cdscrip, imbc, jmbc, lat0, lon0,
. ndatebc, nymdbc, nhmsbc, undef, error)
C--------- Check data for Compatibility
C Check for correct data in boundary condition file
if (.not.error .and. cname.ne.'SST') then
write(6,*)'Wrong data in SST boundary condition file => ',cname
error = .true.
endif
C Check Horizontal Resolution
if(.not.error.and.imbc*jmbc.ne.xsize*ysize)then
write(6,*) ' B.C. Resolution DOES NOT match Model Resolution!'
write(6,*) ' B.C. Resolution: ',imbc*jmbc
write(6,*) 'Model Resolution: ',xsize*ysize
error = .true.
endif
C Check Year
iyearbc = nymdbc(2)/10000
if (.not.error .and. iyear.ne.iyearbc .and. iyearbc.ne.0) then
write(6,*)' B.C. Year DOES NOT match REQUESTED Year!'
write(6,*)' B.C. Year: ', iyearbc
write(6,*)'Requested Year: ', iyear
error = .true.
endif
if (.not.error) then
C if climatology, fill dates for data with current model year
if (iyearbc.eq.0) then
write(6,*)
write(6,*)'Climatological Dataset is being used.'
write(6,*)'Current model year is used to fill Header Dates'
do n = 2, ndatebc-1
nymdbc(n) = nymdbc(n) +(nymdmod/10000)*10000
enddo
C For the first date subtract 1 year from the current model NYMD
n = 1
nymdbc(n) = nymdbc(n) +(nymdmod/10000-1)*10000
C For the last date add 1 year to the current model NYMD
n = ndatebc
nymdbc(n) = nymdbc(n) +(nymdmod/10000+1)*10000
endif
C Write out header info
_BEGIN_MASTER( myThid )
write(6,*) ' Updated boundary condition data'
write(6,*) ' ---------------------------------'
write(6,*) ' Variable: ',cname
write(6,*) ' Description: ',cdscrip
write(6,*) ' Resolution: x= ',imbc,' y= ',jmbc,
. ' Undefined value = ',undef
write(6,*) ' Starting latitude = ',lat0,' longitude =',lon0
write(6,*) ' Data valid at these times: '
ndby3 = ndatebc/3
do n = 1, ndby3*3,3
write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=n,n+2)
1000 format(3(2x,i3,':',i8,2x,i8))
enddo
write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=ndby3*3+1,ndatebc)
_END_MASTER( myThid )
endif
if( error ) call MY_EXIT (101)
endif
C---------- Read SST data if necessary -------------------------------
found = .false.
nd = 2
c If model time is not within the times of saved sst data
c from previous call to getsst then read new data
timemod = dfloat(nymdmod) + dfloat(nhms) /1000000
timebc1 = dfloat(nymdbc1) + dfloat(nhmsbc1)/1000000
timebc2 = dfloat(nymdbc2) + dfloat(nhmsbc2)/1000000
if (timemod .lt. timebc1 .or. timemod .ge. timebc2) then
do while (.not.found .and. nd .le. ndatebc)
timebc2 = dfloat(nymdbc(nd)) + dfloat(nhmsbc(nd))/1000000
if (timebc2 .gt. timemod) then
nymdbc1 = nymdbc(nd-1)
nymdbc2 = nymdbc(nd)
nhmsbc1 = nhmsbc(nd-1)
nhmsbc2 = nhmsbc(nd)
call BCDATA (iunit,imbc,jmbc,nd,nd+1,sstbc1,sstbc2)
found = .true.
else
nd = nd + 1
endif
enddo
c Otherwise the data from the last time in getsst surrounds the
c current model time.
else
found = .true.
endif
if (.not.found) then
print *, 'STOP: Could not find SST dates for model time.'
call MY_FINALIZE
call MY_EXIT (101)
endif
C---------- Interpolate SST data ------------------------------------
call INTERP_TIME(nymdmod,nhms,nymdbc1,nhmsbc1,nymdbc2,nhmsbc2,
. fac1,fac2)
do j = jm1,jm2
do i = im1,im2
sst(i,j,bi,bj) = sstbc1(i+bislot,j+bjslot)*fac1
. + sstbc2(i+bislot,j+bjslot)*fac2
enddo
enddo
return
end
subroutine BCDATA (iunit,im,jm,nrec1,nrec2,field1,field2)
C************************************************************************
C
C!ROUTINE: BCDATA
C!DESCRIPTION: BCDATA reads the data from the file assigned to the
C! passed unit number and returns data from the two times
C! surrounding the current model time. The two record
C! postitions are not assumed to be next to each other.
C
C!INPUT PARAMETERS:
C! im number of x points
C! im number of x points
C! nrec1 record number of the time before the model time
C! nrec2 record number of the time after the model time
C
C!OUTPUT PARAMETERS:
C! field1(im,jm) data field before the model time
C! field2(im,jm) data field after the model time
C
C--------------------------------------------------------------------------
implicit none
integer iunit,im,jm,nrec1,nrec2
_RL field1(im,jm)
_RL field2(im,jm)
integer i,j
real*4 f1(im,jm), f2(im,jm)
C--------- Read file -----------------------------------------------
read(iunit,rec=nrec1) f1
read(iunit,rec=nrec2) f2
#ifdef _BYTESWAPIO
call MDS_BYTESWAPR4( im*jm, f1)
call MDS_BYTESWAPR4( im*jm, f2)
#endif
do j=1,jm
do i=1,im
field1(i,j) = f1(i,j)
field2(i,j) = f2(i,j)
enddo
enddo
return
end
subroutine BCHEADER (iunit, ndmax, nrec,
. cname, cdscrip, im, jm, lat0, lon0, ndatebc,
. nymdbc, nhmsbc, undef, error)
C************************************************************************
C
C!ROUTINE: BCHEADER
C!DESCRIPTION: BCHEADER reads the header from a file and returns the info.
C
C!INPUT PARAMETERS:
C! iunit unit number assigned to the data file
C! ndmax maximum number of date/times of the data
C! nrec record number of the header info (or assume 1 ?)
C
C!OUTPUT PARAMETERS:
C! cname name of the data in the file header
C! cdscrip description of the data in the file header
C! im number of x points
C! jm number of y points
C! lat0 starting latitude for the data grid
C! lon0 starting longitude for the data grid
C! ndatebc number of date/times of the data in the file
C! nymdbc(ndmax) array of dates for the data including century
C! nhmsbc(ndmax) array of times for the data
C! undef value for undefined values in the data
C! error logical TRUE if dataset problem
C
C--------------------------------------------------------------------------
implicit none
integer iunit, ndmax, nrec
character*8 cname
character*80 cdscrip
character*112 dummy112
integer im,jm,ndatebc,nymdbc(ndmax),nhmsbc(ndmax)
_RL lat0,lon0,undef
logical error
integer i
integer*4 im_32,jm_32
integer*4 ndatebc_32,nhmsbc_32(ndmax),nymdbc_32(ndmax)
real*4 lat0_32,lon0_32,undef_32
C--------- Read file -----------------------------------------------
read(iunit,rec=nrec,err=500) cname, cdscrip,
. im_32, jm_32, lat0_32, lon0_32,
. ndatebc_32, undef_32
#ifdef _BYTESWAPIO
call MDS_BYTESWAPI4( 1, im_32)
call MDS_BYTESWAPI4( 1, jm_32)
call MDS_BYTESWAPR4( 1, lat0_32)
call MDS_BYTESWAPR4( 1, lon0_32)
call MDS_BYTESWAPI4( 1, ndatebc_32)
call MDS_BYTESWAPR4( 1, undef_32)
#endif
read(iunit,rec=nrec,err=500) dummy112,
. (nymdbc_32(i), nhmsbc_32(i), i=1,ndatebc_32)
im = im_32
jm = jm_32
lat0 = lat0_32
lon0 = lon0_32
undef = undef_32
ndatebc = ndatebc_32
do i=1,ndatebc
#ifdef _BYTESWAPIO
call MDS_BYTESWAPI4( 1, nymdbc_32(i))
call MDS_BYTESWAPI4( 1, nhmsbc_32(i))
#endif
nymdbc(i) = nymdbc_32(i)
nhmsbc(i) = nhmsbc_32(i)
enddo
return
500 continue
print *, 'Error reading boundary condition from unit ',iunit
error = .true.
return
end