/[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.2 by molod, Mon Jun 7 18:11:38 2004 UTC revision 1.11 by molod, Mon Jul 19 22:06:59 2004 UTC
# Line 5  C $Name$ Line 5  C $Name$
5  c----------------------------------------------------------------------  c----------------------------------------------------------------------
6  c  Subroutine update_ocean_exports - 'Wrapper' routine to update  c  Subroutine update_ocean_exports - 'Wrapper' routine to update
7  c        the fields related to the ocean's surface that are needed  c        the fields related to the ocean's surface that are needed
8  c        by fizhi.  c        by fizhi (sst and sea ice extent).
 c        Also: Set up "bi, bj loop" and some timers and clocks here.  
9  c  c
10  c Call:  getsst  (Return the current sst field-read dataset if needed)  c Call:  getsst  (Return the current sst field-read dataset if needed)
11  c        getsice (Return the current sea ice field-read data if needed)  c        getsice (Return the current sea ice field-read data if needed)
# Line 17  c--------------------------------------- Line 16  c---------------------------------------
16  #include "GRID.h"  #include "GRID.h"
17  #include "fizhi_ocean_coms.h"  #include "fizhi_ocean_coms.h"
18  #include "EEPARAMS.h"  #include "EEPARAMS.h"
19    #include "chronos.h"
20    
21         integer myTime, myIter, myThid         integer myTime, myIter, myThid
22    
23         integer i, j, L, bi, bj         integer i, j, bi, bj, biglobal, bjglobal
24         integer im1, im2, jm1, jm2, idim1, idim2, jdim1, jdim2         integer im1, im2, jm1, jm2, idim1, idim2, jdim1, jdim2
25           integer nSxglobal, nSyglobal
26           integer ksst,kice
27           _RL        sstmin
28           parameter ( sstmin = 273.16 )
29    
30         im1 = 1-OLx         idim1 = 1-OLx
31         im2 = sNx+OLx         idim2 = sNx+OLx
32         jm1 = 1-OLy         jdim1 = 1-OLy
33         jm2 = sNy+OLy         jdim2 = sNy+OLy
34         idim1 = 1         im1 = 1
35         idim2 = sNx         im2 = sNx
36         jdim1 = 1         jm1 = 1
37         jdim2 = sNy         jm2 = sNy
38           nSxglobal = nSx*nPx
39           nSyglobal = nSy*nPy
40    
41         do bj = myByLo(myThid), myByHi(myThid)         call mdsfindunit( ksst, myThid )
42         do bi = myBxLo(myThid), myBxHi(myThid)         call mdsfindunit( kice, myThid )
43    
44    C***********************************************************************
45    
46  c      call getsst         DO BJ = myByLo(myThid),myByHi(myThid)
47  c      call getsice         DO BI = myBxLo(myThid),myBxHi(myThid)
48    
49           biglobal=bi+(myXGlobalLo-1)/im2
50           bjglobal=bj+(myYGlobalLo-1)/jm2
51    
52           call getsst(ksst,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSx,
53         .  nSy,nSxglobal,nSyglobal,bi,bj,biglobal,bjglobal,nymd,nhms,sst)
54           call getsice(kice,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSx,
55         .  nSy,nSxglobal,nSyglobal,bi,bj,biglobal,bjglobal,nymd,nhms,sice)
56    
57    c Check for Minimum Open-Water SST
58    c --------------------------------
59           do j=jm1,jm2
60           do i=im1,im2
61           if(sice(i,j,bi,bj).eq.0.0 .and. sst(i,j,bi,bj).lt.sstmin)
62         .                                          sst(i,j,bi,bj) = sstmin
63         enddo         enddo
64         enddo         enddo
65    
66           ENDDO
67           ENDDO
68    
69         return         return
70         end         end
71    
72           subroutine getsice(iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,
73         .     nSumx,nSumy,nPgx,nPgy,bi,bj,biglobal,bjglobal,nymd,nhms,sice)
74    C************************************************************************
75    C
76    C!ROUTINE:      GETSICE
77    C!DESCRIPTION:  GETSICE returns the sea ice depth.
78    C!              This routine is adaptable for any frequency
79    C!              data upto a daily frequency.  
80    C!              note: for diurnal data ndmax should be increased.
81    C
82    C!INPUT PARAMETERS:
83    C!      iunit     Unit number assigned to the sice data file
84    C!      idim1     Start dimension in x-direction
85    C!      idim2     End dimension in x-direction
86    C!      jdim1     Start dimension in y-direction
87    C!      jdim2     End dimension in y-direction
88    C!      im1       Begin of x-direction span for filling sice
89    C!      im2       End of x-direction span for filling sice
90    C!      jm1       Begin of y-direction span for filling sice
91    C!      jm2       End of y-direction span for filling sice
92    C!      nSumx     Number of processors in x-direction (local processor)
93    C!      nSumy     Number of processors in y-direction (local processor)
94    C!      nPgx      Number of processors in x-direction (global)
95    C!      nPgx      Number of processors in y-direction (global)
96    C!      bi        Processor number in x-direction (local to processor)
97    C!      bj        Processor number in y-direction (local to processor)
98    C!      biglobal  Processor number in x-direction (global)
99    C!      bjglobal  Processor number in y-direction (global)
100    C!      nymd      YYMMDD of the current model timestep
101    C!      nhms      HHMMSS of the model time
102    C
103    C!OUTPUT PARAMETERS:
104    C!      sice(idim1:idim2,jdim1:jdim2,nSumx,nSumy) Sea ice depth in meters
105    C
106    C!ROUTINES CALLED:
107    C
108    C!      bcdata       Reads the data for a given unit number
109    C!      bcheader     Reads the header info for a given unit number
110    C!      interp_time  Returns weights for linear interpolation
111    C
112    C--------------------------------------------------------------------------
113    
114          implicit none
115    #include "CPP_EEOPTIONS.h"
116    #include "SIZE.h"
117    
118          integer iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSumx,nSumy
119          integer nPgx,nPgy,bi,bj,biglobal,bjglobal,nymd,nhms
120    
121          _RL sice(idim1:idim2,jdim1:jdim2,nSumx,nSumy)
122    
123    C Maximum number of dates in one year for the data
124          integer   ndmax
125          parameter (ndmax = 370)
126    
127          character*8  cname
128          character*80 cdscrip
129          real fac1, fac2, lat0, lon0, timebc1, timebc2, timemod, undef
130          logical first, found, error
131          integer i,j,n,nn,iyear,iyearbc,nd,ndby3,imbc,jmbc,npxbc,npybc
132          integer ndatebc,nhmsbc(ndmax), nhmsbc1, nhmsbc2,nrec
133          integer nymdbc(ndmax),nymdbc1,nymdbc2,nymdmod
134    
135          real sicebc1(sNx,sNy,nSx*nPx,nSy*nPy)
136          real sicebc2(sNx,sNy,nSx*nPx,nSy*nPy)
137    
138    C--------- Variable Initialization ---------------------------------
139    
140          data first /.true./
141          data error /.false./
142    
143    c  save header info
144          save imbc,jmbc,npxbc,npybc,lat0,lon0,ndatebc,undef,nymdbc,nhmsbc
145          save nymdbc1, nymdbc2, nhmsbc1, nhmsbc2
146          save first
147          save sicebc1, sicebc2
148    
149    c  this only works for between 1950-2050
150          if (nymd .lt. 500101) then
151            nymdmod = 20000000 + nymd
152          else if (nymd .le. 991231) then
153            nymdmod = 19000000 + nymd
154          else
155            nymdmod = nymd
156          endif
157    
158    c  initialize so that first time through they have values for the check
159    c  these vaules make the iyear .ne. iyearbc true anyways for
160    c  for the first time so first isnt checked below.
161    
162          if (first) then
163            nymdbc(2) = 0
164            nymdbc1   = 0
165            nymdbc2   = 0
166            nhmsbc1   = 0
167            nhmsbc2   = 0
168            first = .false.
169          endif
170    
171    C---------- Read in Header file ----------------------------------
172    
173          iyear   = nymdmod/10000
174          iyearbc = nymdbc(2)/10000
175    
176          if( iyear.ne.iyearbc ) then
177    
178           close(iunit)
179           open (iunit,form='unformatted',access='direct',
180         .                                         recl=im2*jm2*nPgx*nPgy*4)
181           nrec = 1
182           call bcheader (iunit, ndmax, nrec,
183         .          cname, cdscrip, imbc, jmbc, npxbc, npybc, lat0, lon0,
184         .          ndatebc, nymdbc, nhmsbc, undef, error)
185    
186    C--------- Check data for Compatibility ------------------------------
187    
188    C Check for correct data in boundary condition file
189           if (.not.error .and. cname.ne.'SICE') then
190            write(6,*)'Wrong data in SICE boundary condition file => ',cname
191            error = .true.
192           endif
193    
194    C Check Horizontal Resolution
195           if(.not.error.and.imbc*jmbc*npxbc*npybc.ne.im2*jm2*npgx*npgy)then
196            write(6,*) ' B.C. Resolution DOES NOT match Model Resolution!'
197            write(6,*) ' B.C. Resolution:  ',imbc*jmbc*npxbc*npybc
198            write(6,*) 'Model Resolution:  ',im2*jm2*npgx*npgy
199            error = .true.
200           endif
201    
202    C Check Year
203           iyearbc = nymdbc(2)/10000
204           if (.not.error .and. iyear.ne.iyearbc .and. iyearbc.ne.0) then
205            write(6,*)'     B.C. Year DOES NOT match REQUESTED Year!'
206            write(6,*)'     B.C. Year:  ', iyearbc
207            write(6,*)'Requested Year:  ', iyear
208            error = .true.
209           endif
210    
211           if (.not.error)   then
212    C if climatology, fill dates for data with current model year
213            if (iyearbc.eq.0) then          
214             write(6,*)
215             write(6,*) 'Climatological Dataset is being used.'  
216             write(6,*) 'Current model year to be used to fill Header Dates'
217             do n = 2, ndatebc-1
218              nymdbc(n) = nymdbc(n) +(nymdmod/10000)*10000
219             enddo
220    C  For the first date subtract 1 year from the current model NYMD
221             n = 1
222             nymdbc(n) = nymdbc(n) +(nymdmod/10000-1)*10000
223    C  For the last date add 1 year to the current model NYMD
224             n = ndatebc
225             nymdbc(n) = nymdbc(n) +(nymdmod/10000+1)*10000
226            endif
227    
228    C  Write out header info
229            write(6,*) ' Updated boundary condition data'
230            write(6,*) ' ---------------------------------'
231            write(6,*) ' Variable: ',cname
232            write(6,*) ' Description: ',cdscrip
233            write(6,*) ' Resolution: x= ',imbc,' y= ',jmbc,
234         .                                       ' Undefined value = ',undef
235            write(6,*) ' Starting latitude = ',lat0,' longitude =',lon0
236            write(6,*) ' Data valid at these times: '
237            ndby3 = ndatebc/3
238            do n = 1, ndby3*3,3
239             write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=n,n+2)
240     1000    format(3(2x,i3,':',i8,2x,i8))
241            enddo
242            write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=ndby3*3+1,ndatebc)
243           endif  
244    
245          endif
246    
247    C---------- Read sice data if necessary -------------------------------
248    
249          found = .false.
250          nd = 2
251    
252    c  If model time is not within the times of saved sice data
253    c  from previous call to getsice then read new data
254    
255          timemod = float(nymdmod) + float(nhms)   /1000000
256          timebc1 = float(nymdbc1) + float(nhmsbc1)/1000000
257          timebc2 = float(nymdbc2) + float(nhmsbc2)/1000000
258    
259          if (timemod .lt. timebc1 .or. timemod .ge. timebc2) then
260    
261           do while (.not.found .and. nd .le. ndatebc)
262            timebc2 = float(nymdbc(nd)) + float(nhmsbc(nd))/1000000
263            if (timebc2 .gt. timemod) then  
264             nymdbc1 = nymdbc(nd-1)
265             nymdbc2 = nymdbc(nd)
266             nhmsbc1 = nhmsbc(nd-1)
267             nhmsbc2 = nhmsbc(nd)
268             call bcdata (iunit,imbc,jmbc,nPgx,nPgy,nd,nd+1,sicebc1,sicebc2)
269             found = .true.
270            else
271             nd = nd + 1
272            endif
273           enddo
274    
275    c  Otherwise the data from the last time in getsice surrounds the
276    c  current model time.
277    
278          else
279           found = .true.
280          endif
281    
282          if (.not.found) then
283           print *, 'STOP: Could not find SICE dates for model time.'
284           call my_finalize
285           call my_exit (101)
286          endif
287    
288    C---------- Interpolate sice data ------------------------------------
289    
290          call interp_time(nymdmod,nhms,nymdbc1,nhmsbc1,nymdbc2,nhmsbc2,
291         .                                                       fac1,fac2)
292    
293          do j = jm1,jm2
294          do i = im1,im2
295           sice(i,j,bi,bj) = sicebc1(i,j,biglobal,bjglobal)*fac1
296         .                 + sicebc2(i,j,biglobal,bjglobal)*fac2
297    c average to 0 or 1
298    c -----------------
299           if (sice(i,j,bi,bj) .ge. 0.5) then
300            sice(i,j,bi,bj) = 1.
301           else
302            sice(i,j,bi,bj) = 0.
303           endif
304          enddo
305          enddo
306    
307    C---------- Fill sice with depth of ice ------------------------------------
308          do j = jm1,jm2
309          do i = im1,im2
310           if (sice(i,j,bi,bj) .eq. 1.) then                
311            sice(i,j,bi,bj) = 3.
312           endif
313          enddo
314          enddo
315    C---------------------------------------------------------------------------
316    
317          return
318          end
319          subroutine getsst(iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,
320         .      nSumx,nSumy,nPgx,nPgy,bi,bj,biglobal,bjglobal,nymd,nhms,sst)
321    C************************************************************************
322    C
323    C!ROUTINE:      GETSST
324    C!DESCRIPTION:  GETSST gets the SST data.
325    C!              This routine is adaptable for any frequency
326    C!              data upto a daily frequency.  
327    C!              note: for diurnal data ndmax should be increased.
328    C
329    C!INPUT PARAMETERS:
330    C!      iunit     Unit number assigned to the sice data file
331    C!      idim1     Start dimension in x-direction
332    C!      idim2     End dimension in x-direction
333    C!      jdim1     Start dimension in y-direction
334    C!      jdim2     End dimension in y-direction
335    C!      im1       Begin of x-direction span for filling sice
336    C!      im2       End of x-direction span for filling sice
337    C!      jm1       Begin of y-direction span for filling sice
338    C!      jm2       End of y-direction span for filling sice
339    C!      nSumx     Number of processors in x-direction (local processor)
340    C!      nSumy     Number of processors in y-direction (local processor)
341    C!      nPgx      Number of processors in x-direction (global)
342    C!      nPgy      Number of processors in y-direction (global)
343    C!      bi        Processor number in x-direction (local to processor)
344    C!      bj        Processor number in y-direction (local to processor)
345    C!      biglobal  Processor number in x-direction (global)
346    C!      bjglobal  Processor number in y-direction (global)
347    C!      nymd      YYMMDD of the current model timestep
348    C!      nhms      HHMMSS of the model time
349    C
350    C!OUTPUT PARAMETERS:
351    C!      sst(idim1:idim2,jdim1:jdim2,nSumx,nSumy) Sea surface temperature (K)
352    C
353    C!ROUTINES CALLED:
354    C
355    C!      bcdata          Reads the data for a given unit number
356    C!      bcheader        Reads the header info for a given unit number
357    C!     interp_time   Returns weights for linear interpolation
358    C
359    C--------------------------------------------------------------------------
360    
361          implicit none
362    #include "CPP_EEOPTIONS.h"
363    #include "SIZE.h"
364    
365          integer iunit,idim1,idim2,jdim1,jdim2,im1,im2,jm1,jm2,nSumx,nSumy
366          integer nPgx,nPgy,bi,bj,biglobal,bjglobal,nymd,nhms
367    
368          _RL sst(idim1:idim2,jdim1:jdim2,nSumx,nSumy)
369    
370    C Maximum number of dates in one year for the data
371          integer   ndmax
372          parameter (ndmax = 370)
373    
374          character*8  cname
375          character*80 cdscrip
376          real fac1, fac2, lat0, lon0, timebc1, timebc2, timemod, undef
377          logical first, found, error
378          integer i,j,n,nn,iyear,iyearbc,nd,ndby3,imbc,jmbc,npxbc,npybc
379          integer ndatebc,nhmsbc(ndmax), nhmsbc1, nhmsbc2,nrec
380          integer nymdbc(ndmax),nymdbc1,nymdbc2,nymdmod
381    
382          real sstbc1(sNx,sNy,nSx*nPx,nSy*nPy)
383          real sstbc2(sNx,sNy,nSx*nPx,nSy*nPy)
384    
385    C--------- Variable Initialization ---------------------------------
386    
387          data first /.true./
388          data error /.false./
389    
390    c  save header info
391          save imbc,jmbc,npxbc,npybc,lat0,lon0,ndatebc,undef,nymdbc,nhmsbc
392          save nymdbc1, nymdbc2, nhmsbc1, nhmsbc2
393          save first
394          save sstbc1, sstbc2
395    
396    c  this only works for between 1950-2050
397          if (nymd .lt. 500101) then
398            nymdmod = 20000000 + nymd
399          else if (nymd .le. 991231) then
400            nymdmod = 19000000 + nymd
401          else
402            nymdmod = nymd
403          endif
404    
405    c  initialize so that first time through they have values for the check
406    c  these vaules make the iyear .ne. iyearbc true anyways for
407    c  for the first time so first isnt checked below.
408          if (first) then
409            nymdbc(2) = 0
410            nymdbc1   = 0
411            nymdbc2   = 0
412            nhmsbc1   = 0
413            nhmsbc2   = 0
414            first = .false.
415          endif
416    
417    C---------- Read in Header file ----------------------------------
418    
419          iyear   = nymdmod/10000
420          iyearbc = nymdbc(2)/10000
421    
422          if( iyear.ne.iyearbc ) then
423    
424           close(iunit)
425           open (iunit,form='unformatted',access='direct',
426         .                                        recl=im2*jm2*nPgx*nPgy*4)
427           nrec = 1
428           call bcheader (iunit, ndmax, nrec,
429         .          cname, cdscrip, imbc, jmbc, npxbc, npybc, lat0, lon0,
430         .          ndatebc, nymdbc, nhmsbc, undef, error)
431    
432    C--------- Check data for Compatibility
433    
434    C Check for correct data in boundary condition file
435           if (.not.error .and. cname.ne.'SST') then
436            write(6,*)'Wrong data in SST boundary condition file => ',cname
437            error = .true.
438           endif
439    
440    C Check Horizontal Resolution
441           if(.not.error.and.imbc*jmbc*npxbc*npybc.ne.im2*jm2*npgx*npgy)then
442            write(6,*) ' B.C. Resolution DOES NOT match Model Resolution!'
443            write(6,*) ' B.C. Resolution:  ',imbc*jmbc*npxbc*npybc
444            write(6,*) 'Model Resolution:  ',im2*jm2*npgx*npgy
445            error = .true.
446           endif
447    
448    C Check Year
449           iyearbc = nymdbc(2)/10000
450           if (.not.error .and. iyear.ne.iyearbc .and. iyearbc.ne.0) then
451            write(6,*)'     B.C. Year DOES NOT match REQUESTED Year!'
452            write(6,*)'     B.C. Year:  ', iyearbc
453            write(6,*)'Requested Year:  ', iyear
454            error = .true.
455           endif
456    
457           if (.not.error)   then
458    C if climatology, fill dates for data with current model year
459            if (iyearbc.eq.0) then          
460             write(6,*)
461             write(6,*)'Climatological Dataset is being used.'  
462             write(6,*)'Current model year is used to fill Header Dates'
463             do n = 2, ndatebc-1
464              nymdbc(n) = nymdbc(n) +(nymdmod/10000)*10000
465             enddo
466    C  For the first date subtract 1 year from the current model NYMD
467             n = 1
468             nymdbc(n) = nymdbc(n) +(nymdmod/10000-1)*10000
469    C  For the last date add 1 year to the current model NYMD
470             n = ndatebc
471             nymdbc(n) = nymdbc(n) +(nymdmod/10000+1)*10000
472            endif
473    
474    C  Write out header info
475            write(6,*) ' Updated boundary condition data'
476            write(6,*) ' ---------------------------------'
477            write(6,*) ' Variable: ',cname
478            write(6,*) ' Description: ',cdscrip
479            write(6,*) ' Resolution: x= ',imbc,' y= ',jmbc,
480         .                                       ' Undefined value = ',undef
481            write(6,*) ' Starting latitude = ',lat0,' longitude =',lon0
482            write(6,*) ' Data valid at these times: '
483            ndby3 = ndatebc/3
484            do n = 1, ndby3*3,3
485             write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=n,n+2)
486     1000    format(3(2x,i3,':',i8,2x,i8))
487            enddo
488            write(6,1000) (nn,nymdbc(nn),nhmsbc(nn),nn=ndby3*3+1,ndatebc)
489           endif
490    
491           if( error ) call my_exit (101)
492    
493          endif
494    
495    C---------- Read SST data if necessary -------------------------------
496    
497          found = .false.
498          nd = 2
499    
500    c  If model time is not within the times of saved sst data
501    c  from previous call to getsst then read new data
502    
503          timemod = float(nymdmod) + float(nhms)   /1000000
504          timebc1 = float(nymdbc1) + float(nhmsbc1)/1000000
505          timebc2 = float(nymdbc2) + float(nhmsbc2)/1000000
506          if (timemod .lt. timebc1 .or. timemod .ge. timebc2) then
507    
508           do while (.not.found .and. nd .le. ndatebc)
509            timebc2 = float(nymdbc(nd)) + float(nhmsbc(nd))/1000000
510            if (timebc2 .gt. timemod) then  
511             nymdbc1 = nymdbc(nd-1)
512             nymdbc2 = nymdbc(nd)
513             nhmsbc1 = nhmsbc(nd-1)
514             nhmsbc2 = nhmsbc(nd)
515             call bcdata (iunit,imbc,jmbc,nPgx,nPgy,nd,nd+1,sstbc1,sstbc2)
516             found = .true.
517            else
518             nd = nd + 1
519            endif
520           enddo
521    
522    c  Otherwise the data from the last time in getsst surrounds the
523    c  current model time.
524    
525          else
526           found = .true.
527          endif
528    
529          if (.not.found) then
530           print *, 'STOP: Could not find SST dates for model time.'
531           call my_finalize
532           call my_exit (101)
533          endif
534    
535    C---------- Interpolate SST data ------------------------------------
536    
537          call interp_time(nymdmod,nhms,nymdbc1,nhmsbc1,nymdbc2,nhmsbc2,
538         .                                                        fac1,fac2)
539    
540          do j = jm1,jm2
541          do i = im1,im2
542           sst(i,j,bi,bj) = sstbc1(i,j,biglobal,bjglobal)*fac1
543         .                + sstbc2(i,j,biglobal,bjglobal)*fac2
544          enddo
545          enddo
546    
547          return
548          end
549    
550          subroutine bcdata (iunit,im,jm,nPx,nPy,nrec1,nrec2,field1,field2)
551    C************************************************************************
552    C
553    C!ROUTINE:      BCDATA
554    C!DESCRIPTION:  BCDATA reads the data from the file assigned to the
555    C!              passed unit number and returns data from the two times
556    C!              surrounding the current model time.  The two record
557    C!              postitions are not assumed to be next to each other.
558    C
559    C!INPUT PARAMETERS:
560    C!      im      number of x points
561    C!      im      number of x points
562    C!      nPx     number of faces in x-direction
563    C!      nPy     number of faces in y-direction
564    C!      nrec1   record number of the time before the model time
565    C!      nrec2   record number of the time after the model time
566    C
567    C!OUTPUT PARAMETERS:
568    C!      field1(im,jm,nPx,nPy)  data field before the model time
569    C!      field2(im,jm,nPx,nPy)  data field after the model time
570    C
571    C--------------------------------------------------------------------------
572          implicit none
573    #include "CPP_EEOPTIONS.h"
574    
575          integer iunit,im,jm,nPx,nPy,nrec1,nrec2
576    
577          real  field1(im,jm,nPx,nPy)
578          real  field2(im,jm,nPx,nPy)
579    
580          integer i,j,n1,n2
581          real*4 f1(im,jm,nPx,nPy), f2(im,jm,nPx,nPy)
582    
583    C--------- Read file -----------------------------------------------
584          read(iunit,rec=nrec1) f1
585          read(iunit,rec=nrec2) f2
586    
587          do n2=1,nPy
588          do n1=1,nPx
589          do j=1,jm
590          do i=1,im
591           field1(i,j,n1,n2) = f1(i,j,n1,n2)
592           field2(i,j,n1,n2) = f2(i,j,n1,n2)
593          enddo
594          enddo
595          enddo
596          enddo
597    
598          return
599          end
600          subroutine bcheader (iunit, ndmax, nrec,
601         .           cname, cdscrip, im, jm, npx, npy, lat0, lon0, ndatebc,
602         .           nymdbc, nhmsbc, undef, error)
603    C************************************************************************
604    C
605    C!ROUTINE:     BCHEADER
606    C!DESCRIPTION: BCHEADER reads the header from a file and returns the info.  
607    C
608    C!INPUT PARAMETERS:
609    C!      iunit    unit number assigned to the data file
610    C!      ndmax    maximum number of date/times of the data
611    C!      nrec     record number of the header info (or assume 1??)
612    C
613    C!OUTPUT PARAMETERS:
614    C!      cname         name of the data in the file header
615    C!      cdscrip       description of the data in the file header
616    C!      im            number of x points
617    C!      jm            number of y points
618    C!      npx           number of faces (processors) in x-direction
619    C!      npy           number of faces (processors) in x-direction
620    C!      lat0          starting latitude for the data grid
621    C!      lon0          starting longitude for the data grid
622    C!      ndatebc       number of date/times of the data in the file
623    C!      nymdbc(ndmax) array of dates for the data including century
624    C!      nhmsbc(ndmax) array of times for the data
625    C!      undef         value for undefined values in the data
626    C!      error         logical TRUE if dataset problem
627    C
628    C--------------------------------------------------------------------------
629          implicit none
630    #include "CPP_EEOPTIONS.h"
631    
632          integer iunit, ndmax, nrec
633    
634          character*8  cname
635          character*80 cdscrip
636          integer im,jm,npx,npy,ndatebc,nymdbc(ndmax),nhmsbc(ndmax)
637          real lat0,lon0,undef
638          logical error
639    
640          integer i
641          integer*4 im_32,jm_32,npx_32,npy_32
642          integer*4 ndatebc_32,nhmsbc_32(ndmax),nymdbc_32(ndmax)
643          real*4 lat0_32,lon0_32,undef_32
644    
645    C--------- Read file -----------------------------------------------
646    
647          read(iunit,rec=nrec,err=500) cname, cdscrip,
648         .     im_32, jm_32, npx_32, npy_32, lat0_32, lon0_32,
649         .     ndatebc_32, undef_32,
650         .     (nymdbc_32(i), nhmsbc_32(i), i=1,ndatebc_32)
651    
652          im = im_32
653          jm = jm_32
654          npx = npx_32
655          npy = npy_32
656          lat0 = lat0_32
657          lon0 = lon0_32
658          undef = undef_32
659    
660          ndatebc = ndatebc_32
661          do i=1,ndatebc
662          nymdbc(i) = nymdbc_32(i)
663          nhmsbc(i) = nhmsbc_32(i)
664          enddo
665    
666          return
667      500 continue
668          print *, 'Error reading boundary condition from unit ',iunit
669          error = .true.
670          return
671          end

Legend:
Removed from v.1.2  
changed lines
  Added in v.1.11

  ViewVC Help
Powered by ViewVC 1.1.22