program rdeof

c  This routine reads in the two daily spatial EOF patterns of 2.5 degree resolution
c  OLR from 20S-20N (144 longitudes X 17 latitudes = 2448 values) and assigns a day
c  then writes out an unformatted (binary) file of the dates and the two patterns 
c  in a pair of 144X17 arrays. These arrays go from 20S-20N so that (1,1) corresponds
c  to 0E, 20S, (2,2) corresponds to 2.5E, 17.5S, etc.
c  
c  This routine is set up to scale the EOFs so that multiplying them by the normalized
c  PCs produces the properly scaled MJO OLR field. Details are given in Kiladis et al.
c  (2014; Monthly Weather Review, in review) and available on the web page:
c  http://www.esrl.noaa.gov/psd/mjo/mjoindex/

        parameter (nlon=144,nlat=17)

        real data(nlon,nlat)
        real in1(nlon,nlat),in2(nlon,nlat)
        real out1(nlon,nlat),out2(nlon,nlat)
        
        character cfile1*70,cfile2*70
        integer nyy,nmm,ndd,nii,ivar,igrid,level
        integer mm,id,ih
    
      open(16,file = '/data8/test/eof.out.b',
     .        form='unformatted')

        nyy=0000

        do 222 id = 1,366
        
c---code the input file:

        if(id .ge. 1 .and. id  .le. 9) then
        write(cfile1,10) '/data1/data/eof/eof1/eof00',id,'.txt'
        write(cfile2,10) '/data1/data/eof/eof2/eof00',id,'.txt'
 10     format(a26,i1, a4)
        end if

        if(id .ge. 10 .and. id  .le. 99) then
        write(cfile1,11) '/data1/data/eof/eof1/eof0',id,'.txt'
        write(cfile2,11) '/data1/data/eof/eof2/eof0',id,'.txt'
 11     format(a25,i2, a4)
        end if

        if(id .ge. 100) then
        write(cfile1,12) '/data1/data/eof/eof1/eof',id,'.txt'
        write(cfile2,12) '/data1/data/eof/eof2/eof',id,'.txt'
 12     format(a24,i3, a4)
        end if

        print*, cfile1
        print*, cfile2

      open(10,file = cfile1,
     .        form='formatted')
      open(11,file = cfile2,
     .        form='formatted')
      open(12,file = 'eof.out',
     .        form='formatted')

      do 91 lat=1,nlat
       do 90 lon=1,nlon
      read(10,*)in1(lon,lat)
      read(11,*)in2(lon,lat)
  90  continue
  91  continue

      call julian(id,nmm,ndd,nii)

      do 191 lat=1,nlat
       do 190 lon=1,nlon
c      out1(lon,lat)=in1(lon,lat)
c      out2(lon,lat)=in2(lon,lat)
       out1(lon,lat)=in1(lon,lat)*212.710159 !for omi
       out2(lon,lat)=in2(lon,lat)*212.710159 !for omi
 190  continue
 191  continue

c write out a binary file of two EOFs per day

c     do 81 lat=1,nlat
c      do 80 lon=1,nlon
c     write(12,*) nyy,nmm,ndd,nii
c     write(12,*)out1(lon,lat)
c     write(12,*)out2(lon,lat)
c 80  continue
c 81  continue
 
      write(16) nyy,nmm,ndd,nii
      write(16) out1
      write(16) out2

        close(10)	! one file for a day
        close(11)	! one file for a day

 222    continue      	!end loop
        
      stop
      end

      subroutine julian(id,nmm,ndd,nii)

      integer id

c   Computes a date given the julian day. 

      if(id .ge. 1 .and. id .le. 31)then 
      nmm=1
      ndd=id
      end if

      if(id .ge. 32 .and. id .le. 60)then
      nmm=2
      ndd=id-31
      end if

      if(id .ge. 61 .and. id .le. 91)then
      nmm=3
      ndd=id-60
      end if

      if(id .ge. 92 .and. id .le. 121)then
      nmm=4
      ndd=id-91
      end if

      if(id .ge. 122 .and. id .le. 152)then
      nmm=5
      ndd=id-121
      end if

      if(id .ge. 153 .and. id .le. 182)then
      nmm=6
      ndd=id-152
      end if

      if(id .ge. 183 .and. id .le. 213)then
      nmm=7
      ndd=id-182
      end if

      if(id .ge. 214 .and. id .le. 244)then
      nmm=8
      ndd=id-213
      end if

      if(id .ge. 245 .and. id .le. 274)then
      nmm=9
      ndd=id-244
      end if

      if(id .ge. 275 .and. id .le. 305)then
      nmm=10
      ndd=id-274
      end if

      if(id .ge. 306 .and. id .le. 335)then
      nmm=11
      ndd=id-305
      end if

      if(id .ge. 336 .and. id .le. 366)then
      nmm=12
      ndd=id-335
      end if

      return
      end