! timfiles.f90
! module procedure written by Tim Mitchell in Dec 99
! last modification on 12.10.01
! contains all the .tim file routines
!	LoadTim, SaveTim, TimToGlo, TimAvLin

module TimFiles

implicit none

contains

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

subroutine LoadTim (JobYearN, JobADYear, LinLineN, LinNames, LinAnaSeries)

integer, pointer, dimension (:)			:: JobADYear

! note change to len=20, 9.8.00, to harmonise with region names
character (len=20), pointer, dimension (:) 	:: LinNames
					! names of individual lines
real, pointer, dimension (:,:)			:: LinAnaSeries
					! datum: line, year

integer, intent (in)	:: JobYearN	! no. years in job
integer, intent (out)	:: LinLineN	! no. lines

integer, allocatable, dimension (:) :: FileADYears

real, parameter :: MissVal = -999.0

real :: TimMin, TimMax
real :: TotalMiss, PerCentMiss

character (len=80) :: GivenFile, LinFilePath, LinTitle
character (len=10) :: LinFormat

integer :: LinTimeN			! no. time steps in file
integer :: XLine, XTime, XTen
integer :: Time0, Time1
integer :: TimeTenN
integer :: AllocStat			! status of allocation statement
integer :: ReadStatus			! status of user input
integer :: TimeMisMatch
integer :: QComp			! 0=uncompressed +=pos of '.' in '.Z'

print*, "  > Enter the filepath of the .tim file to load (include any .Z suffix): "
do
	do
		read (*,*,iostat=ReadStatus), GivenFile
		if (ReadStatus.LE.0) exit
	end do
	
	inquire (file=GivenFile, name=LinFilePath)
	open (1, file=LinFilePath, status="old", iostat=ReadStatus)
	if (ReadStatus .EQ. 0) close (1)
	if (ReadStatus .EQ. 0) exit
end do

QComp = index(LinFilePath,'.Z')
if (QComp.GT.0) then
  print*, "  > Temporarily uncompressing .z file..."
  call system ('uncompress ' // LinFilePath)
  LinFilePath (QComp:(QComp+1)) = "  " 	
end if
	
open (2, file=LinFilePath, status="old", access="sequential", form="formatted", &
		action="read")

read (2, fmt="(A80)"), LinTitle
print*, adjustl(trim(LinTitle))

read (2, fmt="(A10)"), LinFormat
read (2, fmt="(2I6)"), LinLineN, LinTimeN

TimeTenN = LinTimeN / 10

allocate (LinNames    (LinLineN), &
	  LinAnaSeries(LinLineN,LinTimeN), &
	  FileADYears (LinTimeN), stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: LoadTim: Allocation failure #####"

do XLine = 1, LinLineN
  read (2, fmt="(a20)"), LinNames (XLine)
end do

do XTime = 1, LinTimeN
  read (2, fmt="(i4)"), FileADYears (XTime)
end do

TimeMisMatch = 0
do XTime = 1, LinTimeN
  if (FileADYears(XTime).NE.MissVal) then
    if (JobADYear(XTime).NE.MissVal) then
      if (FileADYears(XTime).NE.JobADYear(XTime)) then
      	TimeMisMatch = TimeMisMatch + 1
      end if
    end if
  end if
end do
if (TimeMisMatch.NE.0) print*, &
	"  > ##### ERROR: LoadTim: Job and file time mismatches, totalling: ", TimeMisMatch

do XLine = 1, LinLineN
 do XTen = 1, TimeTenN
  Time0 = (XTen-1)*10 + 1
  Time1 = Time0 + 9
  read (2, fmt=LinFormat), (LinAnaSeries (XLine, XTime), XTime=Time0,Time1)
 end do
end do

close (2)

TimMin    =  100000.0
TimMax    = -100000.0
TotalMiss = 0.0
do XLine = 1, LinLineN
    do XTime = 1, LinTimeN
      if      (LinAnaSeries (XLine, XTime).EQ.MissVal) then
      	TotalMiss = TotalMiss + 1.0
      else if (LinAnaSeries (XLine, XTime).GT.TimMax)  then
        TimMax    = LinAnaSeries (XLine, XTime)
      else if (LinAnaSeries (XLine, XTime).LT.TimMin)  then
        TimMin    = LinAnaSeries (XLine, XTime)
      end if
    end do
end do
PerCentMiss = 100.0 * TotalMiss / (LinLineN*LinTimeN)
print "(a33,3f10.2)", "   > % missing, min val, max val: ", PerCentMiss, TimMin, TimMax

if (QComp.GT.0) then
  print*, "  > Recompressing .z file..."
  call system ('compress ' // LinFilePath)
end if

end subroutine LoadTim

!*******************************************************************************
! extract any number of .glo from a typical tim structure

subroutine TimToGlo (LongN, LatN, RegN, YearN, ModelFilePath, ADTimes, MapIDLReg, GivenTim)

real, pointer, dimension (:)			:: GloToSave
real, pointer, dimension (:,:)			:: GivenTim
integer, pointer, dimension (:) 		:: ADTimes
integer, pointer, dimension (:,:) 		:: MapIDLReg

integer, intent (in) 		:: LongN, LatN, RegN, YearN
character (len=80), intent (in)	:: ModelFilePath

real, parameter :: MissVal = -999.0

integer :: ChosenADYear, SelectADYear
integer :: ReadStatus, AllocStat
integer :: XYear, XReg

character (len=80) :: GloTitle, Blank

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

Blank = ""

allocate (GloToSave  (RegN), stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: TimToGlo: Allocation failure #####"
GloToSave  = MissVal

ChosenADYear = 0

do
  do
    print*, "  > Enter the year AD to extract (-99=end): "
    do
        read (*,*,iostat=ReadStatus), SelectADYear
	if (ReadStatus.LE.0) exit
    end do
    
    if (SelectADYear.NE.-99) then
      do XYear = 1, YearN
        if (ADTimes(XYear).EQ.SelectADYear) ChosenADYear = XYear
      end do
    else
      ChosenADYear = 1
    end if
    
    if (ChosenADYear.EQ.0) print*, "  > Year out of range. Try again."
    if (ChosenADYear.GE.1) exit
  end do
  
  if (SelectADYear.NE.-99) then
    do XReg = 1, RegN
      GloToSave(XReg) = GivenTim(XReg,ChosenADYear)
    end do
    
    call SaveGlo (LongN,LatN,RegN,ModelFilePath,Blank,Blank,GloToSave,MapIDLReg)
  end if
  
  if (SelectADYear.EQ.-99) exit
end do

end subroutine TimToGlo

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

subroutine TimAvLin (RegN, YearN, MissAccept, ADYear, GivenTim)

real, pointer, dimension (:,:) 			:: GivenTim, GloAvLin
integer, pointer, dimension (:) 		:: ADYear
character (len=80), pointer, dimension (:)	:: LinNames

real, intent(in)	:: MissAccept
integer, intent(inout)	:: RegN, YearN

real, parameter 	:: MissVal = -999.0

real			:: RegThresh, OpTotal, OpEn
integer 		:: XYear, XReg
integer 		:: AllocStat, ReadStatus
character (len=80) 	:: LinTitle

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

allocate (GloAvLin (1,YearN), &
	  LinNames (1)      , stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: TimAvLin: Allocation failure #####"

GloAvLin     = MissVal
LinNames (1) = "globally averaged array in .tim form"

RegThresh = RegN * (100.0 - MissAccept) / 100.0

do XYear = 1, YearN
  OpTotal = 0.0
  OpEn    = 0.0
  
  do XReg = 1, RegN
    if (GivenTim(XReg,XYear).NE.MissVal) then
      OpTotal = OpTotal + GivenTim(XReg,XYear)
      OpEn    = OpEn    + 1.0
    end if
  end do
  
  if (OpEn.GT.RegThresh) GloAvLin (1,XYear) = OpTotal / OpEn
end do

call SaveLin (1, YearN, LinNames, ADYear, GloAvLin)

deallocate (GloAvLin, LinNames, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: TimAvLin: Deallocation failure #####"

end subroutine TimAvLin

!*******************************************************************************
! altered on 2.10.00 to prompt for title and filepath, which can be blank

subroutine SaveTim (LinLineN, LinTimeN, CallTitle, CallFile, LinNames, LinADTimes, LinAnaSeries)

! note change to len=20, on 9.8.00, to harmonise with calls to save region names in .tim files
character (len=20), pointer, dimension (:) 	:: LinNames
					! names of individual lines
integer, pointer, dimension (:) 		:: LinADTimes
					! year AD of each time step
real, pointer, dimension (:,:)			:: LinAnaSeries
					! datum: line, year

integer, intent (in)				:: LinLineN, LinTimeN
					! no. lines, time steps
character (len=80), intent (in)			:: CallTitle,CallFile		! can be left blank

real, parameter :: MissVal = -999.0

character (len=80) :: GivenFile, LinFilePath, TimTitle
character (len=10) :: LinFormat

integer :: XLine, XTime, XTen
integer :: Time0, Time1
integer :: TimeTenN
integer :: AllocStat			! status of allocation statement
integer :: ReadStatus			! status of user input

if (CallTitle.EQ."") then
  print*, "  > Enter the .tim file title: "
  do
	read (*,*,iostat=ReadStatus), TimTitle
	if (ReadStatus.LE.0.AND.TimTitle.NE."") exit
  end do
else
  TimTitle = CallTitle
end if

if (CallFile.EQ."") then
  print*, "  > Enter the filepath of the .tim file to save: "
  do
	do
		read (*,*,iostat=ReadStatus), GivenFile
		if (ReadStatus.LE.0) exit
	end do
	
	inquire (file=GivenFile, name=LinFilePath)
	open (1, file=LinFilePath, status="new", iostat=ReadStatus)
	if (ReadStatus .EQ. 0) close (1)
	if (ReadStatus .EQ. 0) exit
  end do
else
  LinFilePath = CallFile
end if

if (LinLineN.GE.100) print*, "  > Saving to .tim ..."

LinFormat = "(10E12.4)"

TimeTenN = LinTimeN / 10

open (2, file=LinFilePath, status="replace", access="sequential", form="formatted", &
		action="write")

write (2, fmt="(A80)"), TimTitle
write (2, fmt="(A10)"), LinFormat
write (2, fmt="(2I6)"), LinLineN, LinTimeN

do XLine = 1, LinLineN
  write (2, fmt="(a20)"), LinNames (XLine)
end do

do XTime = 1, LinTimeN
  write (2, fmt="(i4)"), LinADTimes (XTime)
end do

do XLine = 1, LinLineN
 do XTen = 1, TimeTenN
  Time0 = (XTen-1)*10 + 1
  Time1 = Time0 + 9
  write (2, fmt=LinFormat), (LinAnaSeries (XLine, XTime), XTime=Time0,Time1)
 end do
end do

close (2)

! print*, "  > Compressing into a .z file..."
! call system ('compress ' // LinFilePath)

end subroutine SaveTim

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

end module TimFiles
