/[MITgcm]/MITgcm/pkg/fizhi/update_ocean_exports.F
ViewVC logotype

Diff of /MITgcm/pkg/fizhi/update_ocean_exports.F

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

revision 1.21 by molod, Wed May 4 22:13:46 2005 UTC revision 1.30 by jmc, Tue May 12 19:56:35 2009 UTC
# Line 18  c--------------------------------------- Line 18  c---------------------------------------
18  #include "EEPARAMS.h"  #include "EEPARAMS.h"
19  #include "chronos.h"  #include "chronos.h"
20  #ifdef ALLOW_EXCH2  #ifdef ALLOW_EXCH2
21    #include "W2_EXCH2_SIZE.h"
22  #include "W2_EXCH2_TOPOLOGY.h"  #include "W2_EXCH2_TOPOLOGY.h"
 #include "W2_EXCH2_PARAMS.h"  
23  #endif /* ALLOW_EXCH2 */  #endif /* ALLOW_EXCH2 */
24    
25         integer myIter, myThid         integer myIter, myThid
26         _RL myTime         _RL myTime
27    
28           INTEGER xySize
29    #if defined(ALLOW_EXCH2)
30           PARAMETER ( xySize = W2_ioBufferSize )
31    #else
32           PARAMETER ( xySize = Nx*Ny )
33    #endif
34         integer i, j, bi, bj, bislot, bjslot         integer i, j, bi, bj, bislot, bjslot
35         integer im1, im2, jm1, jm2, idim1, idim2, jdim1, jdim2         integer im1, im2, jm1, jm2, idim1, idim2, jdim1, jdim2
36         integer xsize, ysize         integer xsize, ysize
37         integer ksst,kice         _RL        sstmin
        _RL        sstmin  
38         parameter ( sstmin = 273.16 )         parameter ( sstmin = 273.16 )
39    
40  #if defined(ALLOW_EXCH2)         _RL sst1 (xySize), sst2 (xySize)
41         PARAMETER ( xsize = exch2_domain_nxt * sNx )         _RL sice1(xySize), sice2(xySize)
42         PARAMETER ( ysize = exch2_domain_nyt * sNy )  c      _RL sst1(xsize,ysize),sst2(xsize,ysize)
43  #else  c      _RL sice1(xsize,ysize),sice2(xsize,ysize)
44         PARAMETER ( xsize = Nx )         integer nymd1sst(nSx,nSy),nymd2sst(nSx,nSy)
45         PARAMETER ( ysize = Ny )         integer nymd1sice(nSx,nSy),nymd2sice(nSx,nSy)
46  #endif         integer nhms1sst(nSx,nSy),nhms2sst(nSx,nSy)
47         _RL sst1(xsize,ysize),sst2(xsize,ysize)         integer nhms1sice(nSx,nSy),nhms2sice(nSx,nSy)
48         _RL sice1(xsize,ysize),sice2(xsize,ysize)         integer sstdates(370,nSx,nSy),sicedates(370,nSx,nSy)
49           integer ssttimes(370,nSx,nSy),sicetimes(370,nSx,nSy)
50           logical first(nSx,nSy)
51           integer nSxnSy
52           parameter(nSxnSy = nSx*nSy)
53           data first/nSxnSy*.true./
54    
55           save nymd1sst,nymd2sst,nymd1sice,nymd2sice
56           save nhms1sst,nhms2sst,nhms1sice,nhms2sice
57         save sst1, sst2, sice1, sice2         save sst1, sst2, sice1, sice2
58           save sstdates, sicedates
59           save ssttimes, sicetimes
60    
61    #if defined(ALLOW_EXCH2)
62           xsize = exch2_global_Nx
63           ysize = exch2_global_Ny
64    #else
65           xsize = Nx
66           ysize = Ny
67    #endif
68         idim1 = 1-OLx         idim1 = 1-OLx
69         idim2 = sNx+OLx         idim2 = sNx+OLx
70         jdim1 = 1-OLy         jdim1 = 1-OLy
# Line 53  c--------------------------------------- Line 74  c---------------------------------------
74         jm1 = 1         jm1 = 1
75         jm2 = sNy         jm2 = sNy
76    
        call mdsfindunit( ksst, myThid )  
        call mdsfindunit( kice, myThid )  
   
77  C***********************************************************************  C***********************************************************************
78    
79         DO BJ = myByLo(myThid),myByHi(myThid)         DO BJ = myByLo(myThid),myByHi(myThid)
# Line 68  C*************************************** Line 86  C***************************************
86         bjslot = myYGlobalLo-1+(bj-1)*sNy         bjslot = myYGlobalLo-1+(bj-1)*sNy
87  #endif  #endif
88    
89         call getsst(ksst,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSx,         call getsst(ksst,sstclim,idim1,idim2,jdim1,jdim2,im1,im2,
90       .  nSy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,sst1,sst2,sst)       .  jm1,jm2,nSx,nSy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,
91         call getsice(kice,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSx,       .  sst1,sst2,first(bi,bj),nymd1sst(bi,bj),nymd2sst(bi,bj),
92       .  nSy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,sice1,sice2,sice)       .  nhms1sst(bi,bj),nhms2sst(bi,bj),sstdates(1,bi,bj),
93         .  ssttimes(1,bi,bj),sst,myThid)
94           call getsice(kice,siceclim,idim1,idim2,jdim1,jdim2,im1,im2,
95         .  jm1,jm2,nSx,nSy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,
96         .  sice1,sice2,first(bi,bj),nymd1sice(bi,bj),nymd2sice(bi,bj),
97         .  nhms1sice(bi,bj),nhms2sice(bi,bj),sicedates(1,bi,bj),
98         .  sicetimes(1,bi,bj),sice,myThid)
99    
100  c Check for Minimum Open-Water SST  c Check for Minimum Open-Water SST
101  c --------------------------------  c --------------------------------
# Line 84  c -------------------------------- Line 108  c --------------------------------
108    
109         ENDDO         ENDDO
110         ENDDO         ENDDO
111         _EXCH_XYZ_R8(sst,myThid)         _EXCH_XY_RL(sst,myThid)
112         _EXCH_XYZ_R8(sice,myThid)         _EXCH_XY_RL(sice,myThid)
113    
114         return         return
115         end         end
116    
117         subroutine getsice(iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,         subroutine getsice(iunit,clim,idim1,idim2,jdim1,jdim2,im1,im2,
118       .            nSumx,nSumy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,       .   jm1,jm2,nSumx,nSumy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,
119       .                                             sicebc1,sicebc2,sice)       .   sicebc1,sicebc2,first,nymdbc1,nymdbc2,nhmsbc1,nhmsbc2,
120         .   nymdbc,nhmsbc,sice,mythid)
121  C***********************************************************************  C***********************************************************************
122  C  C
123  C!ROUTINE:      GETSICE  C!ROUTINE: GETSICE
124  C!DESCRIPTION:  GETSICE returns the sea ice depth.  C!DESCRIPTION:  GETSICE returns the sea ice depth.
125  C!              This routine is adaptable for any frequency  C!              This routine is adaptable for any frequency
126  C!              data upto a daily frequency.    C!              data upto a daily frequency.
127  C!              note: for diurnal data ndmax should be increased.  C!              note: for diurnal data ndmax should be increased.
128  C  C
129  C!INPUT PARAMETERS:  C!INPUT PARAMETERS:
130  C!      iunit     Unit number assigned to the sice data file  C!      iunit     Unit number assigned to the sice data file
131  C!      idim1     Start dimension in x-direction  C!      idim1     Start dimension in x-direction
132  C!      idim2     End dimension in x-direction  C!      idim2     End dimension in x-direction
133  C!      jdim1     Start dimension in y-direction  C!      jdim1     Start dimension in y-direction
# Line 119  C!      bi        Processor number in x- Line 144  C!      bi        Processor number in x-
144  C!      bj        Processor number in y-direction (local to processor)  C!      bj        Processor number in y-direction (local to processor)
145  C!      bislot  Processor number in x-direction (global)  C!      bislot  Processor number in x-direction (global)
146  C!      bjslot  Processor number in y-direction (global)  C!      bjslot  Processor number in y-direction (global)
147  C!      nymd      YYMMDD of the current model timestep  C!      nymd    YYMMDD of the current model timestep
148  C!      nhms      HHMMSS of the model time  C!      nhms    HHMMSS of the model time
149  C  C
150  C!OUTPUT PARAMETERS:  C!OUTPUT PARAMETERS:
151  C!      sice(idim1:idim2,jdim1:jdim2,nSumx,nSumy) Sea ice depth in meters  C!      sice(idim1:idim2,jdim1:jdim2,nSumx,nSumy) Sea ice depth in meters
152  C  C
153  C!ROUTINES CALLED:  C!ROUTINES CALLED:
154  C  C
155  C!      bcdata       Reads the data for a given unit number  C!      bcdata       Reads the data for a given unit number
156  C!      bcheader     Reads the header info for a given unit number  C!      bcheader     Reads the header info for a given unit number
157  C!      interp_time  Returns weights for linear interpolation  C!      interp_time  Returns weights for linear interpolation
158  C  C
159  C--------------------------------------------------------------------------  C--------------------------------------------------------------------------
# Line 137  C--------------------------------------- Line 162  C---------------------------------------
162  #include "SIZE.h"  #include "SIZE.h"
163    
164        integer iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSumx,nSumy        integer iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSumx,nSumy
165        integer xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms        integer xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,mythid
166          logical clim
167    
168        _RL sicebc1(xsize,ysize)        _RL sicebc1(xsize,ysize)
169        _RL sicebc2(xsize,ysize)        _RL sicebc2(xsize,ysize)
170        _RL sice(idim1:idim2,jdim1:jdim2,nSumx,nSumy)        _RL sice(idim1:idim2,jdim1:jdim2,nSumx,nSumy)
171          integer nhmsbc1,nhmsbc2,nymdbc1,nymdbc2
172          logical first
173    
174  C Maximum number of dates in one year for the data  C Maximum number of dates in one year for the data
175        integer   ndmax        integer ndmax
176        parameter (ndmax = 370)        parameter (ndmax = 370)
177          integer nymdbc(ndmax),nhmsbc(ndmax)
178    
179        character*8  cname        character*8  cname
180        character*80 cdscrip        character*80 cdscrip
181        character*40 sicedata        character*22 sicedata
182        _RL fac1, fac2, lat0, lon0, timebc1, timebc2, timemod, undef        _RL fac1, fac2, lat0, lon0, timebc1, timebc2, timemod, undef
183        logical first, found, error        logical found, error
184        integer i,j,n,nn,iyear,iyearbc,nd,ndby3,imbc,jmbc        integer i,j,n,nn,iyear,iyearbc,nd,ndby3,imbc,jmbc
185        integer ndatebc,nhmsbc(ndmax), nhmsbc1, nhmsbc2,nrec        integer ndatebc,nrec
186        integer nymdbc(ndmax),nymdbc1,nymdbc2,nymdmod        integer nymdmod
187    
188    
189  C--------- Variable Initialization ---------------------------------  C--------- Variable Initialization ---------------------------------
190    
       data first /.true./  
191        data error /.false./        data error /.false./
192    
193  c  save header info  c  save header info
194        save imbc,jmbc,lat0,lon0,ndatebc,undef,nymdbc,nhmsbc        save imbc,jmbc,lat0,lon0,ndatebc,undef
       save nymdbc1, nymdbc2, nhmsbc1, nhmsbc2  
       save first  
195    
196  c  this only works for between 1950-2050  c  this only works for between 1950-2050
197        if (nymd .lt. 500101) then        if (nymd .lt. 500101) then
# Line 176  c  this only works for between 1950-2050 Line 202  c  this only works for between 1950-2050
202          nymdmod = nymd          nymdmod = nymd
203        endif        endif
204    
205        sicedata='sice19232.weekly.clim'        iyear   = nymdmod/10000
206          if(clim) then
207           if(xsize.eq.192)sicedata='sice19232.weekly.clim'
208           if(xsize.eq.612)sicedata='sice612102.weekly.clim'
209          else
210           if(xsize.eq.192)
211         .           WRITE(sicedata,'(A,I4)')'sice19232.weekly.y',iyear
212           if(xsize.eq.612)
213         .           WRITE(sicedata,'(A,I4)')'sice612102.weekly.y',iyear
214          endif
215    
216  c  initialize so that first time through they have values for the check  c  initialize so that first time through they have values for the check
217  c  these vaules make the iyear .ne. iyearbc true anyways for  c  these values make the iyear .ne. iyearbc true anyways for
218  c  for the first time so first isnt checked below.  c  for the first time so first isnt checked below.
219    
220        if (first) then        if (first) then
# Line 188  c  for the first time so first isnt chec Line 223  c  for the first time so first isnt chec
223          nymdbc2   = 0          nymdbc2   = 0
224          nhmsbc1   = 0          nhmsbc1   = 0
225          nhmsbc2   = 0          nhmsbc2   = 0
226  c       first = .false.          first = .false.
227        endif        endif
228    
229  C---------- Read in Header file ----------------------------------  C---------- Read in Header file ----------------------------------
230    
       iyear   = nymdmod/10000  
231        iyearbc = nymdbc(2)/10000        iyearbc = nymdbc(2)/10000
232    
233        if( iyear.ne.iyearbc ) then        if( iyear.ne.iyearbc ) then
234    
235           close(iunit)
236         open (iunit,file=sicedata,form='unformatted',access='direct',         open (iunit,file=sicedata,form='unformatted',access='direct',
237       .                                         recl=xsize*ysize*4)       .                                         recl=xsize*ysize*4)
238         nrec = 1         nrec = 1
239         call bcheader (iunit, ndmax, nrec,         call bcheader (iunit, ndmax, nrec,
240       .          cname, cdscrip, imbc, jmbc, lat0, lon0,       .          cname, cdscrip, imbc, jmbc, lat0, lon0,
241       .          ndatebc, nymdbc, nhmsbc, undef, error)       .          ndatebc, nymdbc, nhmsbc, undef, error)
242    
243  C--------- Check data for Compatibility ------------------------------  C--------- Check data for Compatibility ------------------------------
# Line 221  C Check Horizontal Resolution Line 256  C Check Horizontal Resolution
256          error = .true.          error = .true.
257         endif         endif
258    
259  C Check Year  C Check Year
260         iyearbc = nymdbc(2)/10000         iyearbc = nymdbc(2)/10000
261         if (.not.error .and. iyear.ne.iyearbc .and. iyearbc.ne.0) then         if (.not.error .and. iyear.ne.iyearbc .and. iyearbc.ne.0) then
262          write(6,*)'     B.C. Year DOES NOT match REQUESTED Year!'          write(6,*)'     B.C. Year DOES NOT match REQUESTED Year!'
# Line 232  C Check Year Line 267  C Check Year
267    
268         if (.not.error)   then         if (.not.error)   then
269  C if climatology, fill dates for data with current model year  C if climatology, fill dates for data with current model year
270          if (iyearbc.eq.0) then                    if (iyearbc.eq.0) then
271           write(6,*)           write(6,*)
272           write(6,*) 'Climatological Dataset is being used.'             write(6,*) 'Climatological Dataset is being used.'
273           write(6,*) 'Current model year to be used to fill Header Dates'           write(6,*) 'Current model year to be used to fill Header Dates'
274           do n = 2, ndatebc-1           do n = 2, ndatebc-1
275            nymdbc(n) = nymdbc(n) +(nymdmod/10000)*10000            nymdbc(n) = nymdbc(n) +(nymdmod/10000)*10000
276           enddo           enddo
277  C  For the first date subtract 1 year from the current model NYMD  C  For the first date subtract 1 year from the current model NYMD
278           n = 1           n = 1
279           nymdbc(n) = nymdbc(n) +(nymdmod/10000-1)*10000           nymdbc(n) = nymdbc(n) +(nymdmod/10000-1)*10000
280  C  For the last date add 1 year to the current model NYMD  C  For the last date add 1 year to the current model NYMD
281           n = ndatebc           n = ndatebc
282           nymdbc(n) = nymdbc(n) +(nymdmod/10000+1)*10000           nymdbc(n) = nymdbc(n) +(nymdmod/10000+1)*10000
283          endif          endif
284    
285  C  Write out header info  C  Write out header info
286            _BEGIN_MASTER( myThid )
287          write(6,*) ' Updated boundary condition data'          write(6,*) ' Updated boundary condition data'
288          write(6,*) ' ---------------------------------'          write(6,*) ' ---------------------------------'
289          write(6,*) ' Variable: ',cname          write(6,*) ' Variable: ',cname
# Line 262  C  Write out header info Line 298  C  Write out header info
298   1000    format(3(2x,i3,':',i8,2x,i8))   1000    format(3(2x,i3,':',i8,2x,i8))
299          enddo          enddo
300          write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=ndby3*3+1,ndatebc)          write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=ndby3*3+1,ndatebc)
301         endif            _END_MASTER( myThid )
302           endif
303    
304        endif        endif
305    
# Line 271  C---------- Read sice data if necessary Line 308  C---------- Read sice data if necessary
308        found = .false.        found = .false.
309        nd = 2        nd = 2
310    
311  c  If model time is not within the times of saved sice data  c  If model time is not within the times of saved sice data
312  c  from previous call to getsice then read new data  c  from previous call to getsice then read new data
313    
314        timemod = float(nymdmod) + float(nhms)   /1000000        timemod = dfloat(nymdmod) + dfloat(nhms)   /1000000
315        timebc1 = float(nymdbc1) + float(nhmsbc1)/1000000        timebc1 = dfloat(nymdbc1) + dfloat(nhmsbc1)/1000000
316        timebc2 = float(nymdbc2) + float(nhmsbc2)/1000000        timebc2 = dfloat(nymdbc2) + dfloat(nhmsbc2)/1000000
317    
318        if (timemod .lt. timebc1 .or. timemod .ge. timebc2) then        if (timemod .lt. timebc1 .or. timemod .ge. timebc2) then
319    
320         do while (.not.found .and. nd .le. ndatebc)         do while (.not.found .and. nd .le. ndatebc)
321          timebc2 = float(nymdbc(nd)) + float(nhmsbc(nd))/1000000          timebc2 = dfloat(nymdbc(nd)) + dfloat(nhmsbc(nd))/1000000
322          if (timebc2 .gt. timemod) then            if (timebc2 .gt. timemod) then
323           nymdbc1 = nymdbc(nd-1)           nymdbc1 = nymdbc(nd-1)
324           nymdbc2 = nymdbc(nd)           nymdbc2 = nymdbc(nd)
325           nhmsbc1 = nhmsbc(nd-1)           nhmsbc1 = nhmsbc(nd-1)
# Line 294  c  from previous call to getsice then re Line 331  c  from previous call to getsice then re
331          endif          endif
332         enddo         enddo
333    
334  c  Otherwise the data from the last time in getsice surrounds the  c  Otherwise the data from the last time in getsice surrounds the
335  c  current model time.  c  current model time.
336    
337        else        else
# Line 302  c  current model time. Line 339  c  current model time.
339        endif        endif
340    
341        if (.not.found) then        if (.not.found) then
342         print *, 'STOP: Could not find SICE dates for model time.'         print *, 'STOP: Could not find SICE dates for model time.'
343         call my_finalize         call my_finalize
344         call my_exit (101)         call my_exit (101)
345        endif        endif
# Line 314  C---------- Interpolate sice data ------ Line 351  C---------- Interpolate sice data ------
351    
352        do j = jm1,jm2        do j = jm1,jm2
353        do i = im1,im2        do i = im1,im2
354         sice(i,j,bi,bj) = sicebc1(i+bislot,j+bjslot)*fac1         sice(i,j,bi,bj) = sicebc1(i+bislot,j+bjslot)*fac1
355       .                 + sicebc2(i+bislot,j+bjslot)*fac2       .                 + sicebc2(i+bislot,j+bjslot)*fac2
356  c average to 0 or 1  c average to 0 or 1
357  c -----------------  c -----------------
358         if (sice(i,j,bi,bj) .ge. 0.5) then         if (sice(i,j,bi,bj) .ge. 0.5) then
# Line 329  c ----------------- Line 366  c -----------------
366  C---------- Fill sice with depth of ice ------------------------------------  C---------- Fill sice with depth of ice ------------------------------------
367        do j = jm1,jm2        do j = jm1,jm2
368        do i = im1,im2        do i = im1,im2
369         if (sice(i,j,bi,bj) .eq. 1.) then                         if (sice(i,j,bi,bj) .eq. 1.) then
370          sice(i,j,bi,bj) = 3.          sice(i,j,bi,bj) = 3.
371         endif         endif
372        enddo        enddo
373        enddo        enddo
374  C---------------------------------------------------------------------------  C---------------------------------------------------------------------------
375        close(iunit)  
   
376        return        return
377        end        end
378        subroutine getsst(iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,        subroutine getsst(iunit,clim,idim1,idim2,jdim1,jdim2,im1,im2,
379       .            nSumx,nSumy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,       .   jm1,jm2,nSumx,nSumy,xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,
380       .                                                sstbc1,sstbc2,sst)       .   sstbc1,sstbc2,first,nymdbc1,nymdbc2,nhmsbc1,nhmsbc2,
381         .   nymdbc,nhmsbc,sst,mythid)
382  C***********************************************************************  C***********************************************************************
383  C  C
384  C!ROUTINE:      GETSST  C!ROUTINE: GETSST
385  C!DESCRIPTION:  GETSST gets the SST data.  C!DESCRIPTION:  GETSST gets the SST data.
386  C!              This routine is adaptable for any frequency  C!              This routine is adaptable for any frequency
387  C!              data upto a daily frequency.    C!              data upto a daily frequency.
388  C!              note: for diurnal data ndmax should be increased.  C!              note: for diurnal data ndmax should be increased.
389  C  C
390  C!INPUT PARAMETERS:  C!INPUT PARAMETERS:
391  C!      iunit     Unit number assigned to the sice data file  C!      iunit     Unit number assigned to the sice data file
392  C!      idim1     Start dimension in x-direction  C!      idim1     Start dimension in x-direction
393  C!      idim2     End dimension in x-direction  C!      idim2     End dimension in x-direction
394  C!      jdim1     Start dimension in y-direction  C!      jdim1     Start dimension in y-direction
# Line 368  C!      bi        Processor number in x- Line 405  C!      bi        Processor number in x-
405  C!      bj        Processor number in y-direction (local to processor)  C!      bj        Processor number in y-direction (local to processor)
406  C!      bislot    Slot number into global array in x-direction (global)  C!      bislot    Slot number into global array in x-direction (global)
407  C!      bjslot    Slot number into global array in y-direction (global)  C!      bjslot    Slot number into global array in y-direction (global)
408  C!      nymd      YYMMDD of the current model timestep  C!      nymd      YYMMDD of the current model timestep
409  C!      nhms      HHMMSS of the model time  C!      nhms      HHMMSS of the model time
410  C  C
411  C!OUTPUT PARAMETERS:  C!OUTPUT PARAMETERS:
412  C!      sst(idim1:idim2,jdim1:jdim2,nSumx,nSumy) Sea surface temperature (K)  C!     sst(idim1:idim2,jdim1:jdim2,nSumx,nSumy) Sea surface temperature (K)
413  C  C
414  C!ROUTINES CALLED:  C!ROUTINES CALLED:
415  C  C
416  C!      bcdata          Reads the data for a given unit number  C!     bcdata     Reads the data for a given unit number
417  C!      bcheader        Reads the header info for a given unit number  C!     bcheader   Reads the header info for a given unit number
418  C!     interp_time      Returns weights for linear interpolation  C!     interp_time  Returns weights for linear interpolation
419  C  C
420  C--------------------------------------------------------------------------  C--------------------------------------------------------------------------
421    
# Line 386  C--------------------------------------- Line 423  C---------------------------------------
423  #include "SIZE.h"  #include "SIZE.h"
424    
425        integer iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSumx,nSumy        integer iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSumx,nSumy
426        integer xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms        integer xsize,ysize,bi,bj,bislot,bjslot,nymd,nhms,mythid
427          logical clim
428    
429        _RL sstbc1(xsize,ysize)        _RL sstbc1(xsize,ysize)
430        _RL sstbc2(xsize,ysize)        _RL sstbc2(xsize,ysize)
431        _RL sst(idim1:idim2,jdim1:jdim2,nSumx,nSumy)        _RL sst(idim1:idim2,jdim1:jdim2,nSumx,nSumy)
432          integer nhmsbc1,nhmsbc2,nymdbc1,nymdbc2
433          logical first
434    
435  C Maximum number of dates in one year for the data  C Maximum number of dates in one year for the data
436        integer   ndmax        integer ndmax
437        parameter (ndmax = 370)        parameter (ndmax = 370)
438          integer nymdbc(ndmax),nhmsbc(ndmax)
439    
440        character*8  cname        character*8  cname
441        character*80 cdscrip        character*80 cdscrip
442        character*20 sstdata        character*21 sstdata
443        _RL fac1, fac2, lat0, lon0, timebc1, timebc2, timemod, undef        _RL fac1, fac2, lat0, lon0, timebc1, timebc2, timemod, undef
444        logical first, found, error        logical found, error
445        integer i,j,n,nn,iyear,iyearbc,nd,ndby3,imbc,jmbc        integer i,j,n,nn,iyear,iyearbc,nd,ndby3,imbc,jmbc
446        integer ndatebc,nhmsbc(ndmax), nhmsbc1, nhmsbc2,nrec        integer ndatebc,nrec
447        integer nymdbc(ndmax),nymdbc1,nymdbc2,nymdmod        integer nymdmod
448    
449    
450  C--------- Variable Initialization ---------------------------------  C--------- Variable Initialization ---------------------------------
451    
       data first /.true./  
452        data error /.false./        data error /.false./
453    
454  c  save header info  c  save header info
455        save imbc,jmbc,lat0,lon0,ndatebc,undef,nymdbc,nhmsbc        save imbc,jmbc,lat0,lon0,ndatebc,undef
       save nymdbc1, nymdbc2, nhmsbc1, nhmsbc2  
       save first  
456    
457  c  this only works for between 1950-2050  c  this only works for between 1950-2050
458        if (nymd .lt. 500101) then        if (nymd .lt. 500101) then
# Line 425  c  this only works for between 1950-2050 Line 463  c  this only works for between 1950-2050
463          nymdmod = nymd          nymdmod = nymd
464        endif        endif
465    
466        sstdata='sst19232.weekly.clim'        iyear   = nymdmod/10000
467          if(clim) then
468           if(xsize.eq.192)sstdata='sst19232.weekly.clim'
469           if(xsize.eq.612)sstdata='sst612102.weekly.clim'
470          else
471           if(xsize.eq.192)
472         .           WRITE(sstdata,'(A,I4)')'sst19232.weekly.y',iyear
473           if(xsize.eq.612)
474         .           WRITE(sstdata,'(A,I4)')'sst612102.weekly.y',iyear
475          endif
476    
477  c  initialize so that first time through they have values for the check  c  initialize so that first time through they have values for the check
478  c  these vaules make the iyear .ne. iyearbc true anyways for  c  these vaules make the iyear .ne. iyearbc true anyways for
# Line 436  c  for the first time so first isnt chec Line 483  c  for the first time so first isnt chec
483          nymdbc2   = 0          nymdbc2   = 0
484          nhmsbc1   = 0          nhmsbc1   = 0
485          nhmsbc2   = 0          nhmsbc2   = 0
486  c       first = .false.          first = .false.
487        endif        endif
488    
489  C---------- Read in Header file ----------------------------------  C---------- Read in Header file ----------------------------------
490    
       iyear   = nymdmod/10000  
491        iyearbc = nymdbc(2)/10000        iyearbc = nymdbc(2)/10000
492    
493        if( iyear.ne.iyearbc ) then        if( iyear.ne.iyearbc ) then
494    
495           close(iunit)
496         open (iunit,file=sstdata,form='unformatted',access='direct',         open (iunit,file=sstdata,form='unformatted',access='direct',
497       .                                        recl=xsize*ysize*4)       .                                        recl=xsize*ysize*4)
498         nrec = 1         nrec = 1
499         call bcheader (iunit, ndmax, nrec,         call bcheader (iunit, ndmax, nrec,
500       .          cname, cdscrip, imbc, jmbc, lat0, lon0,       .          cname, cdscrip, imbc, jmbc, lat0, lon0,
501       .          ndatebc, nymdbc, nhmsbc, undef, error)       .          ndatebc, nymdbc, nhmsbc, undef, error)
502    
503  C--------- Check data for Compatibility  C--------- Check data for Compatibility
# Line 469  C Check Horizontal Resolution Line 516  C Check Horizontal Resolution
516          error = .true.          error = .true.
517         endif         endif
518    
519  C Check Year  C Check Year
520         iyearbc = nymdbc(2)/10000         iyearbc = nymdbc(2)/10000
521         if (.not.error .and. iyear.ne.iyearbc .and. iyearbc.ne.0) then         if (.not.error .and. iyear.ne.iyearbc .and. iyearbc.ne.0) then
522          write(6,*)'     B.C. Year DOES NOT match REQUESTED Year!'          write(6,*)'     B.C. Year DOES NOT match REQUESTED Year!'
# Line 480  C Check Year Line 527  C Check Year
527    
528         if (.not.error)   then         if (.not.error)   then
529  C if climatology, fill dates for data with current model year  C if climatology, fill dates for data with current model year
530          if (iyearbc.eq.0) then                    if (iyearbc.eq.0) then
531           write(6,*)           write(6,*)
532           write(6,*)'Climatological Dataset is being used.'             write(6,*)'Climatological Dataset is being used.'
533           write(6,*)'Current model year is used to fill Header Dates'           write(6,*)'Current model year is used to fill Header Dates'
534           do n = 2, ndatebc-1           do n = 2, ndatebc-1
535            nymdbc(n) = nymdbc(n) +(nymdmod/10000)*10000            nymdbc(n) = nymdbc(n) +(nymdmod/10000)*10000
536           enddo           enddo
537  C  For the first date subtract 1 year from the current model NYMD  C  For the first date subtract 1 year from the current model NYMD
538           n = 1           n = 1
539           nymdbc(n) = nymdbc(n) +(nymdmod/10000-1)*10000           nymdbc(n) = nymdbc(n) +(nymdmod/10000-1)*10000
540  C  For the last date add 1 year to the current model NYMD  C  For the last date add 1 year to the current model NYMD
541           n = ndatebc           n = ndatebc
542           nymdbc(n) = nymdbc(n) +(nymdmod/10000+1)*10000           nymdbc(n) = nymdbc(n) +(nymdmod/10000+1)*10000
543          endif          endif
544    
545  C  Write out header info  C  Write out header info
546            _BEGIN_MASTER( myThid )
547          write(6,*) ' Updated boundary condition data'          write(6,*) ' Updated boundary condition data'
548          write(6,*) ' ---------------------------------'          write(6,*) ' ---------------------------------'
549          write(6,*) ' Variable: ',cname          write(6,*) ' Variable: ',cname
# Line 510  C  Write out header info Line 558  C  Write out header info
558   1000    format(3(2x,i3,':',i8,2x,i8))   1000    format(3(2x,i3,':',i8,2x,i8))
559          enddo          enddo
560          write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=ndby3*3+1,ndatebc)          write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=ndby3*3+1,ndatebc)
561            _END_MASTER( myThid )
562         endif         endif
563    
564         if( error ) call my_exit (101)         if( error ) call my_exit (101)
# Line 521  C---------- Read SST data if necessary - Line 570  C---------- Read SST data if necessary -
570        found = .false.        found = .false.
571        nd = 2        nd = 2
572    
573  c  If model time is not within the times of saved sst data  c  If model time is not within the times of saved sst data
574  c  from previous call to getsst then read new data  c  from previous call to getsst then read new data
575    
576          timemod = dfloat(nymdmod) + dfloat(nhms)   /1000000
577          timebc1 = dfloat(nymdbc1) + dfloat(nhmsbc1)/1000000
578          timebc2 = dfloat(nymdbc2) + dfloat(nhmsbc2)/1000000
579    
       timemod = float(nymdmod) + float(nhms)   /1000000  
       timebc1 = float(nymdbc1) + float(nhmsbc1)/1000000  
       timebc2 = float(nymdbc2) + float(nhmsbc2)/1000000  
580        if (timemod .lt. timebc1 .or. timemod .ge. timebc2) then        if (timemod .lt. timebc1 .or. timemod .ge. timebc2) then
581    
582         do while (.not.found .and. nd .le. ndatebc)         do while (.not.found .and. nd .le. ndatebc)
583          timebc2 = float(nymdbc(nd)) + float(nhmsbc(nd))/1000000          timebc2 = dfloat(nymdbc(nd)) + dfloat(nhmsbc(nd))/1000000
584          if (timebc2 .gt. timemod) then            if (timebc2 .gt. timemod) then
585           nymdbc1 = nymdbc(nd-1)           nymdbc1 = nymdbc(nd-1)
586           nymdbc2 = nymdbc(nd)           nymdbc2 = nymdbc(nd)
587           nhmsbc1 = nhmsbc(nd-1)           nhmsbc1 = nhmsbc(nd-1)
# Line 543  c  from previous call to getsst then rea Line 593  c  from previous call to getsst then rea
593          endif          endif
594         enddo         enddo
595    
596  c  Otherwise the data from the last time in getsst surrounds the  c  Otherwise the data from the last time in getsst surrounds the
597  c  current model time.  c  current model time.
598    
599        else        else
600         found = .true.         found = .true.
601        endif        endif
602    
603        if (.not.found) then        if (.not.found) then
604         print *, 'STOP: Could not find SST dates for model time.'         print *, 'STOP: Could not find SST dates for model time.'
605         call my_finalize         call my_finalize
# Line 563  C---------- Interpolate SST data ------- Line 613  C---------- Interpolate SST data -------
613    
614        do j = jm1,jm2        do j = jm1,jm2
615        do i = im1,im2        do i = im1,im2
616         sst(i,j,bi,bj) = sstbc1(i+bislot,j+bjslot)*fac1         sst(i,j,bi,bj) = sstbc1(i+bislot,j+bjslot)*fac1
617       .                + sstbc2(i+bislot,j+bjslot)*fac2       .                + sstbc2(i+bislot,j+bjslot)*fac2
618        enddo        enddo
619        enddo        enddo
620    
       close(iunit)  
621    
622        return        return
623        end        end
# Line 579  C Line 628  C
628  C!ROUTINE:      BCDATA  C!ROUTINE:      BCDATA
629  C!DESCRIPTION:  BCDATA reads the data from the file assigned to the  C!DESCRIPTION:  BCDATA reads the data from the file assigned to the
630  C!              passed unit number and returns data from the two times  C!              passed unit number and returns data from the two times
631  C!              surrounding the current model time.  The two record  C!              surrounding the current model time.  The two record
632  C!              postitions are not assumed to be next to each other.  C!              postitions are not assumed to be next to each other.
633  C  C
634  C!INPUT PARAMETERS:  C!INPUT PARAMETERS:
# Line 621  C--------- Read file ------------------- Line 670  C--------- Read file -------------------
670        return        return
671        end        end
672        subroutine bcheader (iunit, ndmax, nrec,        subroutine bcheader (iunit, ndmax, nrec,
673       .           cname, cdscrip, im, jm, lat0, lon0, ndatebc,       .           cname, cdscrip, im, jm, lat0, lon0, ndatebc,
674       .           nymdbc, nhmsbc, undef, error)       .           nymdbc, nhmsbc, undef, error)
675  C************************************************************************  C************************************************************************
676  C  C
677  C!ROUTINE:     BCHEADER  C!ROUTINE:     BCHEADER
678  C!DESCRIPTION: BCHEADER reads the header from a file and returns the info.    C!DESCRIPTION: BCHEADER reads the header from a file and returns the info.
679  C  C
680  C!INPUT PARAMETERS:  C!INPUT PARAMETERS:
681  C!      iunit    unit number assigned to the data file  C!      iunit    unit number assigned to the data file
# Line 635  C!      nrec     record number of the he Line 684  C!      nrec     record number of the he
684  C  C
685  C!OUTPUT PARAMETERS:  C!OUTPUT PARAMETERS:
686  C!      cname         name of the data in the file header  C!      cname         name of the data in the file header
687  C!      cdscrip       description of the data in the file header  C!      cdscrip       description of the data in the file header
688  C!      im            number of x points  C!      im            number of x points
689  C!      jm            number of y points  C!      jm            number of y points
690  C!      lat0          starting latitude for the data grid  C!      lat0          starting latitude for the data grid
691  C!      lon0          starting longitude for the data grid  C!      lon0          starting longitude for the data grid
692  C!      ndatebc       number of date/times of the data in the file  C!      ndatebc       number of date/times of the data in the file
693  C!      nymdbc(ndmax) array of dates for the data including century  C!      nymdbc(ndmax) array of dates for the data including century
694  C!      nhmsbc(ndmax) array of times for the data  C!      nhmsbc(ndmax) array of times for the data
695  C!      undef         value for undefined values in the data  C!      undef         value for undefined values in the data
# Line 665  C--------------------------------------- Line 714  C---------------------------------------
714    
715  C--------- Read file -----------------------------------------------  C--------- Read file -----------------------------------------------
716    
717        read(iunit,rec=nrec,err=500) cname, cdscrip,        read(iunit,rec=nrec,err=500) cname, cdscrip,
718       .     im_32, jm_32, lat0_32, lon0_32,       .     im_32, jm_32, lat0_32, lon0_32,
719       .     ndatebc_32, undef_32       .     ndatebc_32, undef_32
720    
721  #ifdef _BYTESWAPIO  #ifdef _BYTESWAPIO
Legend:
Removed from v.1.21  
changed lines
  Added in v.1.30
  ViewVC Help
Powered by ViewVC 1.1.22