! gsodgrab.f90
! f90 program written by Tim Mitchell on 30.11.01
! last modified on 22.04.02
! tool to obtain subset of all stations in Global Summary Of Stn Stn sets from NCEP
! pgf90 -Mstandard -Minfo -fast -Mscalarsse -Mvect=sse -Mflushz 
! 	-o ./../obs/gsodgrab filenames.f90 time.f90 ./../obs/gsodgrab.f90 
!	2> /tyn1/tim/scratch/stderr.txt

program GSoDGrab

use FileNames
use Time

implicit none

integer, pointer, dimension (:,:) 		:: MonthLengths
integer, pointer, dimension (:) 		:: YearAD, GetStnWMO

character (len=150), pointer, dimension (:,:) 	:: GetStnText
character (len=80), pointer, dimension (:,:) 	:: Files
character (len=80), pointer, dimension (:) 	:: Batch, BatchName, DumpFile

character (len=3), dimension (12), parameter :: MonthNames = &
	(/'jan','feb','mar','apr','may','jun','jul','aug','sep','oct','nov','dec'/)

real, parameter :: MissVal = -999.0

character (len=80), parameter :: Blank = ""

integer :: ReadStatus, AllocStat
integer :: BegYearAD,EndYearAD,BegMonth,EndMonth
integer :: XGetStn,StnN,BatchN, LineN, YearN
integer :: GetStnN,XStn,XBatch, XLine, XYear, XMonth, XYearAD, XDay
integer :: QZip
integer :: Month0,Month1, GotYear,GotMonth, LastSlash,LastDot, LoadFileLen, ThisYear,ThisMonth,ThisDay
integer :: Station,DateCode,LatDeg,LatMin,LonDeg,LonMin,Elev,StnWMOTextLen

character (len=150) :: RestOfLine
character (len=80) :: GivenFile, SaveFile, LoadFile, InfoFile, GenericFile, Trash, DataDir
character (len=22) :: FullName
character (len=20) :: YearStr20,StnWMOText
character (len=12) :: DateText,TimeText
character (len= 6) :: CallName
character (len= 4) :: Year,YearStr04
character (len= 2) :: Month,Stn,Hour,Minute
character (len= 1) :: Ch,LatDir,LonDir

!*******************************************************************************
! main

call Intro
call Specifics
call DumpHeadersFirst
call LoadStations
call Conclude 					

contains

!*******************************************************************************
! intro

subroutine Intro

open (99,file="/tyn1/tim/scratch/log-gsodgrab.dat",status="replace",action="write")

print*
print*, "  > ***** GSoDGrab.f90 : obtains subset from GSoD set *****"
print*

end subroutine Intro

!*******************************************************************************
! make specificiations

subroutine Specifics

DataDir="/cru/tyn1/f709762/daily/gsod/data/????.??.txt*"
call GetBatch (DataDir,Batch,Silent=1)
BatchN = size (Batch,1)
BegYearAD=1994 ; BegMonth=1
EndMonth=nint(mod(real(BatchN),12.0)) 
if (EndMonth.EQ.0) EndMonth=12
EndYearAD=nint(real(BatchN-EndMonth)/12.0)+1994
YearN = EndYearAD - BegYearAD + 1
print "(2(a,i4,x,i2))", "   > Range: ", BegYearAD, BegMonth, &
			" to ", EndYearAD, EndMonth

allocate (Files     (YearN,12), &
	  BatchName (BatchN),   stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SortBatch: Allocation failure #####"
Files = Blank ; BatchName = Blank

do XBatch = 1, BatchN
  GivenFile = Batch(XBatch)
  LastSlash = index(GivenFile,"/",.TRUE.)
  BatchName(XBatch) = adjustl(trim(GivenFile((LastSlash+1):80)))
end do

XBatch=0
do XYear=1,YearN
  do XMonth=1,12
    XBatch=XBatch+1
    if (XBatch.LE.BatchN) Files(XYear,XMonth) = Batch(XBatch)
  end do
end do


print*, "  > Enter the number of stations to obtain:"
do
	read (*,*,iostat=ReadStatus), GetStnN	! get the number of stations
	if (ReadStatus.LE.0.AND.GetStnN.GE.1) exit
end do

allocate (GetStnWMO (GetStnN),    &
	  GetStnText(GetStnN,31), &
	  YearAD    (YearN),      &
	  DumpFile  (GetStnN),    &
	  MonthLengths(YearN,12), stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: Specifics: Allocation failure #####"

StnN = 999999					! get the station WMO codes
print*, "  > Enter each station's WMO code in turn: "
do XGetStn = 1, GetStnN
  do
	read (*,*,iostat=ReadStatus), GetStnWMO(XGetStn)
	if (GetStnWMO(XGetStn).LT.1.OR.GetStnWMO(XGetStn).GT.StnN) then
	  ReadStatus = 1 ; print*, "  > WMO code out of range. Try again."	
	end if
	if (ReadStatus.LE.0) exit
  end do
end do

do XYear = 1, YearN				! get MonthLengths
  YearAD(XYear) = BegYearAD + XYear - 1
end do
call GetMonthLengths (YearAD,MonthLengths)

print*, "  > Enter the generic .txt file to save (station code added auto):"
do
	read (*,*,iostat=ReadStatus), GivenFile
	if (ReadStatus.LE.0.AND.GivenFile.NE."") exit
end do
GenericFile = SavePath (GivenFile,".txt")
LastDot = index(GenericFile,".",.TRUE.)

do XGetStn = 1, GetStnN
  StnWMOText = GetTextFromInt (GetStnWMO(XGetStn))
  StnWMOTextLen = len_trim(StnWMOText)
  DumpFile (XGetStn) = GenericFile(1:LastDot) // StnWMOText(1:StnWMOTextLen) // ".txt"
end do

end subroutine Specifics

!*******************************************************************************
! load station data

subroutine LoadStations

print*, "  > Loading station data from..."

do XYearAD = BegYearAD, EndYearAD
  XYear = XYearAD-BegYearAD+1
  
  Month0 = 1  ; if (XYearAD.EQ.BegYearAD) Month0 = BegMonth
  Month1 = 12 ; if (XYearAD.EQ.EndYearAD) Month1 = EndMonth
  
  do XMonth = Month0, Month1
    print "(2a)", "   > ", trim(Files(XYear,XMonth))
    
    GetStnText = ""
    
    LoadFile = Files(XYear,XMonth)
    LoadFileLen = len_trim(LoadFile)
    if (LoadFileLen.GT.1.AND.LoadFile((LoadFileLen-1):LoadFileLen).EQ.".Z") then
      QZip = 2							! file is zipped
      call system ('uncompress ' // LoadFile)
      LoadFile ((LoadFileLen-1):LoadFileLen) = "  "		! change filename to the unzipped file
    else
      QZip = 1							! file already unzipped
    end if
    
    call system ('wc -l ' // LoadFile // ' > /tyn1/tim/scratch/trashme.txt')
    open  (3,file='/tyn1/tim/scratch/trashme.txt',status="old",access="sequential",form="formatted",action="read")
    read  (3,"(i8)"), LineN
    close (3)
    call system ('rm /tyn1/tim/scratch/trashme.txt')
    
    open  (2,file=LoadFile,status="old",access="sequential",form="formatted",action="read")    
    read  (2,*), Trash
    do XLine = 2, LineN
      read  (2,"(i6,2x,i4,i2,i2,a150)"), Station, ThisYear,ThisMonth,ThisDay, RestOfLine
      
      if (ThisMonth.EQ.XMonth) then
        do XGetStn = 1, GetStnN
          if (Station.EQ.GetStnWMO(XGetStn)) GetStnText(XGetStn,ThisDay) = RestOfLine
        end do
      end if
    end do    
    close (2)
    
    if (QZip.EQ.2) call system ('compress ' // LoadFile // ' &')
    
    call DumpMonth
  end do
end do

end subroutine LoadStations

!*******************************************************************************
! dump headers first

subroutine DumpHeadersFirst

call date_and_time (DateText, TimeText)
Year  = DateText (1:4)
Month = DateText (5:6)
Stn   = DateText (7:8)
Hour  = TimeText (1:2)
Minute= TimeText (3:4)

do XGetStn = 1, GetStnN
  open  (2,file=DumpFile(XGetStn),status="replace",access="sequential",form="formatted",action="write")    
  
  write (2,"(a,a2,a1,a2,a1,a4,a4,a2,a1,a2,a20)"), "Tyndall Centre file created on ", &
		Stn, ".", Month, ".", Year, " at ", Hour, ":", Minute, " by Dr. Tim Mitchell"
  write (2,"(2(3a,i4))"), "Global Summary of the Day: ",  MonthNames(BegMonth), " ", BegYearAD, " to ", &
							MonthNames(EndMonth), " ", EndYearAD
  write (2,"(a,i6)"), "WMO station code: ", GetStnWMO(XGetStn)
  write (2,"(a)"), "WMO     YYYYMMDD    TEMP       DEWP      SLP        STP       VISIB      WDSP     MXSPD   GUST    MAX     MIN   PRCP   SNDP   FRSHTT"  

  close (2)
end do

end subroutine DumpHeadersFirst

!*******************************************************************************
! dump month

subroutine DumpMonth

do XGetStn = 1, GetStnN
  open  (3,file=DumpFile(XGetStn),status="old",access="sequential",form="formatted",&
  							position="append",action="write")    
  
  do XDay = 1, MonthLengths(XYear,XMonth)
    write (3,"(i6,a2,i4,i2,i2,a)"), GetStnWMO(XGetStn),"  ",XYearAD,XMonth,XDay, &
    					trim(GetStnText(XGetStn,XDay))
  end do
  
  close (3)
end do

end subroutine DumpMonth

!*******************************************************************************
! conclude

subroutine Conclude

deallocate (Files,Batch,BatchName,DumpFile,YearAD,MonthLengths,GetStnWMO,GetStnText,stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: Conclude: Deallocation failure #####"

close (99)

print*

end subroutine Conclude

!*******************************************************************************

end program GSoDGrab