! periodic.f90
! f90 program written by Tim Mitchell on 20.12.00
! last modified on 12.11.01
! main program with which to manipulate data from .per files
! pgf90 -Mstandard -Minfo -fast -Mscalarsse -Mvect=sse -Mflushz 
! 	-o ./../obs/periodic time.f90 filenames.f90 loadperfiles.f90 
!	saveperfiles.f90 dayfiles.f90 annfiles.f90 linfiles.f90 sortmod.f90 
!	usesort.f90 basicfun.f90 plainperm.f90 regress.f90 cetgeneral.f90 
!	./../obs/periodic.f90 2> /tyn1/tim/scratch/stderr.txt

program Periodic

use FileNames
use Time
use LoadPERFiles
use SavePERFiles
use DAYFiles
use ANNFiles
use LinFiles
use SortMod
use UseSort
use BasicFun
use PlainPerm
use Regress
use CETGeneral

implicit none

real, pointer, dimension (:,:,:)		:: Data, FileDaily
real, pointer, dimension (:,:)			:: FileMonthly, FileSeasonal, FileAll
real, pointer, dimension (:)			:: FileAnnual, OpVec, YearVec, ResultVecA, ResultVecB
real, allocatable, dimension (:,:), target	:: TargFileAll

integer, pointer, dimension (:,:)		:: MonthLengths
integer, pointer, dimension (:)			:: YearAD, FileYearAD
integer, pointer, dimension (:)			:: VariCode
integer, allocatable, dimension (:), target	:: YearTarg

character (len=80), pointer, dimension (:)		:: FileCols,FileSet
character (len=80), allocatable, dimension (:), target	:: TargFileCols
character (len= 9), pointer, dimension (:)		:: AnnLabels

character (len=80), dimension (17)		:: CalendarNames

real, parameter :: MissVal = -999.0

logical :: Verbose, Talkative

real :: MissAccept,MissThresh
real :: OpTot,OpMiss,OpEn
real :: MonMiss,MonEn,MonFrac,SeaMiss,SeaEn,SeaFrac,AnnMiss,AnnEn,AnnFrac
real :: Base,Constant,Difference,Precision,Aye,Bee, DataMin,DataMax

integer :: ReadStatus, AllocStat
integer :: MenuChoice, NextFree, NonBlank 
integer :: ChosenCol, ChosenColA, ChosenColB, ChosenPeriod, ChosenOp, ChosenSumm, ChosenPer
integer :: StartYear, EndYear, YearN, ColumnN, FileColN, SegN, SavePerN, FileN, SaveColN
integer :: XYear, XMonth, XFileYear, XMasterYear, XSeason, XPeriod, XStartYear, XEndYear, XDestineYear
integer :: XCol, XSeg, XOrigYear, XDestYear, XPer,XFile,XSaveCol
integer :: XAnnual, XFileCol
integer :: FileYear0,FileYear1,MasterYear0,MasterYear1,FileYearN,Month0,Month1
integer :: StartYearAD,EndYearAD,DestineYearAD
integer :: FileVariCode
integer :: THalf,BaseLen,DisagreeN,Iter,SegLen, TotA,TotB
integer :: QPerMon,QGradCept,QTimeCol

character (len=80) :: LoadFile,SaveFile,Blank
character (len=20) :: String20

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

print*
print*, "  > ##### periodic.f90 manipulates .per files #####"
print*

call Intro

do
  print*, "  > Main menu. Make your choice. (0=list)"
  do
	read (*,*,iostat=ReadStatus), MenuChoice
	if (ReadStatus.LE.0) exit
  end do
	
  if      (MenuChoice.EQ. 1)  then
    call Conclude
    call Intro
  else if (MenuChoice.EQ. 3)  then
    call ReplicateSeg
  else if (MenuChoice.EQ. 4)  then
    call CopyPeriodCol
  else if (MenuChoice.EQ. 5)  then
    call OperateColA
  else if (MenuChoice.EQ. 6)  then
    call OperateColAB
  else if (MenuChoice.EQ. 7)  then
    call CombineColAB
  else if (MenuChoice.EQ. 8)  then
    call ReplicateYear
  else if (MenuChoice.EQ.10) then
    call LoadPerSingle
  else if (MenuChoice.EQ.12) then
    call LoadDayFile
  else if (MenuChoice.EQ.13) then
    call LoadLinFile
  else if (MenuChoice.EQ.15) then
    call LoadPerBulk
  else if (MenuChoice.EQ.20) then
    QPerMon = 1 ; call SavePerMon
  else if (MenuChoice.EQ.23) then
    call SaveLinFile
  else if (MenuChoice.EQ.24) then
    call SaveAnnCol
  else if (MenuChoice.EQ.30) then
    call CheckMaster
  else if (MenuChoice.EQ.40)  then
    call DeTrend
  else if (MenuChoice.EQ.41)  then
    call TruncRange
  else if (MenuChoice.EQ.42)  then
    call AnomalisePeriod
  else if (MenuChoice.EQ.43)  then
    call SmoothCol
  else if (MenuChoice.EQ.50)  then
    call SeasonaliseA
  else if (MenuChoice.NE.99) then
    print*, "  >  1. Reinitialise"
    print*, "  >  3. Replicate segment"
    print*, "  >  4. Copy one period into entire column"
    print*, "  >  5. Operate column  A"
    print*, "  >  6. Operate columns A.B"
    print*, "  >  7. Combine columns A.B (check disputes)" 
    print*, "  >  8. Copy one year into entire column"
    print*, "  > 10. Load .per file (single)"
    print*, "  > 12. Load .day file"
    print*, "  > 13. Load .lin file" 
    print*, "  > 15. Load .per file (batch)"
    print*, "  > 20. Save .per file"
    print*, "  > 23. Save .lin file"
    print*, "  > 24. Save .ann file"
    print*, "  > 30. Check data array"
    print*, "  > 40. Detrend (v. time or col)"
    print*, "  > 41. Truncate data range"
    print*, "  > 42. Anomalise v. period"
    print*, "  > 43. Smooth with Gaussian filter"
    print*, "  > 50. Seasonalise A --> .ann"
    print*, "  > 99. Exit"
  end if
  
  if (MenuChoice.EQ.99) exit
end do

call Conclude
print*

contains

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

subroutine Intro

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

print*, "  > Enter the number of series to allow in the master array:"
do
	read (*,*,iostat=ReadStatus), ColumnN
	if (ReadStatus.LE.0 .AND. ColumnN.GE.1) exit
end do

print*, "  > Enter the start and end years AD of the master array:"
do
	read (*,*,iostat=ReadStatus), StartYear, EndYear
	if (ReadStatus.LE.0 .AND. StartYear.LE.EndYear) exit
end do

Blank = "                                                                                 "
YearN      = EndYear - StartYear + 1
NextFree   = 1
MissAccept = 10				! % missing values that is acceptable
Talkative = .TRUE.

allocate ( Data    (ColumnN, YearN, 17), &
	   YearTarg(YearN),     	 &
	   VariCode(ColumnN),		 stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: Intro: Allocation failure #####"
Data     = MissVal
VariCode = MissVal
YearAD   => YearTarg

do XYear = 1, YearN
  YearAD (XYear) = StartYear + XYear - 1
end do

CalendarNames (1)  = "Jan" ; CalendarNames (2)  = "Feb" ; CalendarNames (3)  = "Mar"
CalendarNames (4)  = "Apr" ; CalendarNames (5)  = "May" ; CalendarNames (6)  = "Jun"
CalendarNames (7)  = "Jul" ; CalendarNames (8)  = "Aug" ; CalendarNames (9)  = "Sep"
CalendarNames (10) = "Oct" ; CalendarNames (11) = "Nov" ; CalendarNames (12) = "Dec"
CalendarNames (13) = "MAM" ; CalendarNames (14) = "JJA" ; CalendarNames (15) = "SON"
CalendarNames (16) = "DJF" ; CalendarNames (17) = "Ann"

print*

end subroutine Intro

!*******************************************************************************
! select column

subroutine SelectColumn

print*, "  > Select the column: "
do
	read (*,*,iostat=ReadStatus), ChosenCol
	if (ReadStatus.LE.0 .AND. ChosenCol.GE.1 .AND. ChosenCol.LE.ColumnN) exit
end do

end subroutine SelectColumn

!*******************************************************************************
! anomalise against period from within column

subroutine AnomalisePeriod

print*, "  > Select the column from which to calculate the base."
call SelectColumn

allocate (OpVec (17), stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: AnomalisePeriod: Allocation failure #####"
OpVec = MissVal

print*, "  > Select the start,end years AD of the base period: "
do
	read (*,*,iostat=ReadStatus), StartYear, EndYear
	if (StartYear.LT.YearAD(1)) print*, "  > Too early a start. Try again."
	if (EndYear.GT.YearAD(YearN)) print*, "  > Too late a finish. Try again."
	if (StartYear.GT.EndYear) print*, "  > Try getting it the right way round. Try again."
	if (ReadStatus.LE.0.AND.StartYear.GE.YearAD(1).AND.EndYear.LE.YearAD(YearN).AND. &
		EndYear.GE.StartYear) exit
end do

XStartYear = StartYear - YearAD(1) + 1
XEndYear   = XStartYear + EndYear - StartYear
BaseLen    = EndYear - StartYear + 1
MissThresh = (MissAccept/100.0) * BaseLen

do XPeriod = 1, 17
  OpMiss = 0.0
  OpTot  = 0.0
  
  do XYear = XStartYear, XEndYear
    if (Data(ChosenCol,XYear,XPeriod).NE.MissVal) then
      OpTot  = OpTot  + Data(ChosenCol,XYear,XPeriod)
    else
      OpMiss = OpMiss + 1
    end if
  end do
  
  if (OpMiss.LE.MissThresh) then
  	OpVec(XPeriod) = OpTot / (BaseLen-OpMiss)
  end if
end do

do XPeriod = 1, 17
    do XYear = 1, YearN
    	Data(NextFree,XYear,XPeriod) = OpVec(XPeriod)
    end do
end do

VariCode(NextFree) = VariCode(ChosenCol)
call StoreMoveOn

print*, "  > Select the column to anomalise."
call SelectColumn

do XPeriod = 1, 17
  if (OpVec(XPeriod).NE.MissVal) then
    do XYear = 1, YearN
      if (Data(ChosenCol,XYear,XPeriod).NE.MissVal) &
    			Data(NextFree,XYear,XPeriod) = Data(ChosenCol,XYear,XPeriod) - OpVec(XPeriod)
    end do
  end if
end do

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

VariCode(NextFree) = VariCode(ChosenCol)
call StoreMoveOn

end subroutine AnomalisePeriod

!*******************************************************************************
! truncate data range

subroutine TruncRange

call SelectColumn

print*, "  > Select the min,max (-999=no truncation): "
do
	read (*,*,iostat=ReadStatus), DataMin,DataMax
	if (DataMin.NE.MissVal.AND.DataMax.NE.MissVal.AND.DataMin.GT.DataMax) then
	  print*, "  > Range reversed. Try again."
	  ReadStatus=1
	end if
	if (ReadStatus.LE.0) exit
end do

TotA=0 ; TotB=0
do XPeriod = 1, 17
  do XYear = 1, YearN
    if (Data(ChosenCol,XYear,XPeriod).NE.MissVal) then
      if      (DataMin.NE.MissVal.AND.Data(ChosenCol,XYear,XPeriod).LT.DataMin) then
      		Data(NextFree,XYear,XPeriod) = DataMin
      		TotA=TotA+1
      else if (DataMax.NE.MissVal.AND.Data(ChosenCol,XYear,XPeriod).GT.DataMax) then
      		Data(NextFree,XYear,XPeriod) = DataMax
      		TotB=TotB+1
      else
      		Data(NextFree,XYear,XPeriod) = Data(ChosenCol,XYear,XPeriod)
      end if
    end if
  end do
end do
print*, "  > Data truncated to min, max: ", TotA,TotB

VariCode(NextFree) = VariCode(ChosenCol)

call StoreMoveOn

end subroutine TruncRange

!*******************************************************************************
! smooth column with Gaussian filter

subroutine SmoothCol

call SelectColumn

print*, "  > Select the periodicity at which to smooth: "
do
	read (*,*,iostat=ReadStatus), THalf
	if (ReadStatus.LE.0 .AND. THalf.GE.1 .AND. THalf.LE.YearN) exit
end do

allocate (OpVec      (YearN), &
	  ResultVecA (YearN), &
	  ResultVecB (YearN), stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SmoothCol: Allocation failure #####"
OpVec = MissVal

do XPeriod = 1, 17
  do XYear = 1, YearN
    OpVec (XYear) = Data (ChosenCol,XYear,XPeriod)
  end do
  
  call GaussSmooth (YearN,THalf,1,OpVec,ResultVecA,ResultVecB)
  
  do XYear = 1, YearN
    Data (NextFree,XYear,XPeriod) = ResultVecA (XYear)
  end do
end do

deallocate (OpVec,ResultVecA,ResultVecB, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SmoothCol: Deallocation failure #####"

VariCode(NextFree) = VariCode(ChosenCol)

call StoreMoveOn

end subroutine SmoothCol

!*******************************************************************************
! average or summate over months into new season --> .ann

subroutine SeasonaliseA

allocate (TargFileAll  (YearN,1), & 
	  AnnLabels (1)      , &
	  MonthLengths (YearN,12), stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SeasonaliseA: Allocation failure #####"
FileAll => TargFileAll
FileAll=MissVal 

call SelectColumn
call GetMonthLengths (YearAD,MonthLengths)

print*, "  > Select the first, last months in season (overlaps=OK): "
do
	read (*,*,iostat=ReadStatus), Month0,Month1
	if (ReadStatus.LE.0.AND.Month0.GE.1.AND.Month1.GE.1.AND.Month0.LE.12.AND.Month1.LE.12) exit
end do

AnnLabels (1) = trim(CalendarNames(Month0)) // "-" // trim(CalendarNames(Month1))

if (Month1.GE.Month0) then
  do XYear = 1, YearN
    OpTot  = 0.0
    OpMiss = 0.0
    OpEn   = 0.0
    
    do XMonth = Month0, Month1
      if (Data(ChosenCol,XYear,XMonth).NE.MissVal) then
        OpTot  = OpTot + Data(ChosenCol,XYear,XMonth) * MonthLengths (XYear,XMonth)
        OpEn   = OpEn  + MonthLengths (XYear,XMonth)
      else
        OpMiss = OpMiss + 1
      end if
    end do
    
    if (OpMiss.EQ.0) FileAll(XYear,1) = OpTot / OpEn
  end do
else
  do XYear = 1, (YearN-1)
    OpTot  = 0.0
    OpMiss = 0.0
    OpEn   = 0.0
    
    do XMonth = Month0, 12
      if (Data(ChosenCol,XYear,XMonth).NE.MissVal) then
        OpTot  = OpTot + Data(ChosenCol,XYear,XMonth) * MonthLengths (XYear,XMonth)
        OpEn   = OpEn  + MonthLengths (XYear,XMonth)
      else
        OpMiss = OpMiss + 1
      end if
    end do
    
    do XMonth = 1, Month1
      if (Data(ChosenCol,(XYear+1),XMonth).NE.MissVal) then
        OpTot  = OpTot + Data(ChosenCol,(XYear+1),XMonth) * MonthLengths ((XYear+1),XMonth)
        OpEn   = OpEn  + MonthLengths ((XYear+1),XMonth)
      else
        OpMiss = OpMiss + 1
      end if
    end do
    
    if (OpMiss.EQ.0) FileAll(XYear,1) = OpTot / OpEn
  end do
end if

call SaveANN (Blank, YearAD, AnnLabels, FileAll)

deallocate (FileAll,AnnLabels,MonthLengths, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SeasonaliseA: Deallocation failure #####"

print*

end subroutine SeasonaliseA

!*******************************************************************************
! replicate segments from other columns to make new column

subroutine ReplicateSeg

print*, "  > Enter the number of segments with which to make: ", NextFree
do
	read (*,*,iostat=ReadStatus), SegN
	if (ReadStatus.LE.0 .AND. SegN.GE.1) exit
end do

do XSeg = 1, SegN
  print "(a,i3)", "   > Segment ", XSeg
  
  call SelectColumn
  
  print "(a,i3)", "   > Select the start,end years AD from origin: ", ChosenCol
  do
	read (*,*,iostat=ReadStatus), StartYearAD,EndYearAD
	if (ReadStatus.LE.0 .AND. StartYearAD.GE.YearAD(1) .AND. EndYearAD.LE.YearAD(YearN) &
		.AND. StartYearAD.LE.EndYearAD) exit
  end do
  SegLen = EndYearAD - StartYearAD + 1
  
  print "(a,i3)", "   > Select the start year AD in destination:   ", NextFree
  do
	read (*,*,iostat=ReadStatus), DestineYearAD
	if (ReadStatus.LE.0 .AND. DestineYearAD.GE.YearAD(1) .AND. &
		(DestineYearAD+SegLen-1).LE.YearAD(YearN)) exit
  end do
  print "(a,i4)", "   > End year AD in destination: ", (DestineYearAD+SegLen-1)
  
  XStartYear   = StartYearAD   - YearAD(1) + 1 
  XDestineYear = DestineYearAD - YearAD(1) + 1
  
  do XYear = 1, SegLen
    XOrigYear = XStartYear   + XYear - 1
    XDestYear = XDestineYear + XYear - 1
    
    do XPeriod = 1, 17
      Data (NextFree,XDestYear,XPeriod) = Data (ChosenCol,XOrigYear,XPeriod)
    end do
  end do 
end do

call StoreMoveOn

end subroutine ReplicateSeg

!*******************************************************************************
! replicate single year from old column to make new column

subroutine ReplicateYear

call SelectColumn
  
print*, "  > Enter the year to replicate: "
do
	read (*,*,iostat=ReadStatus), StartYearAD
	XStartYear = StartYearAD - YearAD(1) + 1
	if (ReadStatus.LE.0 .AND. XStartYear.GE.1 .AND. XStartYear.LE.YearN) exit
end do

do XYear = 1, YearN
    do XPeriod = 1, 17
      Data (NextFree,XYear,XPeriod) = Data (ChosenCol,XStartYear,XPeriod)
    end do
end do

call StoreMoveOn

end subroutine ReplicateYear

!*******************************************************************************
! copy a single period (e.g. Jan) into an entire column

subroutine CopyPeriodCol

call SelectColumn

print*, "  > Select the period (1-12=months,13-16=MAM-DJF,17=annual): "
do
	read (*,*,iostat=ReadStatus), ChosenPeriod
	if (ReadStatus.LE.0 .AND. ChosenPeriod.GE.1 .AND. ChosenPeriod.LE.17) exit
end do

do XYear = 1, YearN
  do XPeriod = 1, 17
    Data (NextFree,XYear,XPeriod) = Data (ChosenCol,XYear,ChosenPeriod)
  end do
end do

VariCode(NextFree) = VariCode(ChosenCol)

call StoreMoveOn

end subroutine CopyPeriodCol

!*******************************************************************************
! operate on column A

subroutine OperateColA

print*, "  > Select the operation (A/k=1, A*k=2, A-k=3, A+k=4, sqrt=5, **k=6, abs=7): "
do
	read (*,*,iostat=ReadStatus), ChosenOp
	if (ReadStatus.LE.0 .AND. ChosenOp.GE.1 .AND. ChosenOp.LE.7) exit
end do

call SelectColumn

if (ChosenOp.NE.5.AND.ChosenOp.NE.7) then
  print*, "  > Select the constant: "
  do
	read (*,*,iostat=ReadStatus), Constant
	if (ReadStatus.LE.0) exit
  end do
end if

do XPeriod = 1, 17
  do XYear = 1, YearN
    Data(NextFree,XYear,XPeriod) = OpAwithB (Data(ChosenCol,XYear,XPeriod),Constant,ChosenOp)
  end do
end do

call StoreMoveOn

end subroutine OperateColA

!*******************************************************************************
! operate on the combination of columns A and B

subroutine OperateColAB

print*, "  > Select the operation (A/B=1, A*B=2, A-B=3, A+B=4): "
do
	read (*,*,iostat=ReadStatus), ChosenOp
	if (ReadStatus.LE.0 .AND. ChosenOp.GE.1 .AND. ChosenOp.LE.4) exit
end do

print*, "  > Select columns A and B in turn."
call SelectColumn
ChosenColA = ChosenCol
call SelectColumn
ChosenColB = ChosenCol

do XPeriod = 1, 17
  do XYear = 1, YearN
    Data(NextFree,XYear,XPeriod) = OpAwithB (Data(ChosenColA,XYear,XPeriod), &
    		Data(ChosenColB,XYear,XPeriod),ChosenOp)
  end do
end do

call StoreMoveOn

end subroutine OperateColAB

!*******************************************************************************
! combines columns A and B, checking that they don't disagree

subroutine CombineColAB

DisagreeN = 0

print*, "  > Select columns A and B in turn."
call SelectColumn
ChosenColA = ChosenCol
call SelectColumn
ChosenColB = ChosenCol

print*, "  > Enter allowable difference between columns (e.g. 0.1): "
do
	read (*,*,iostat=ReadStatus), Precision
	if (ReadStatus.LE.0) exit
end do  

do XYear = 1, YearN
  do XPeriod = 1, 17
    if (Data(ChosenColA,XYear,XPeriod).NE.MissVal) then
      if (Data(ChosenColB,XYear,XPeriod).NE.MissVal) then
        Difference = Data(ChosenColA,XYear,XPeriod) - Data(ChosenColB,XYear,XPeriod)
        
        if (abs(Difference).GT.Precision) then
          DisagreeN = DisagreeN + 1
          if (DisagreeN.EQ.1) write (99,"(a36)"), "  Year   Per    Column A    Column B"          
          write (99,"(2i6,3f12.4)"), & 
            YearAD(XYear), XPeriod, Data(ChosenColA,XYear,XPeriod), Data(ChosenColB,XYear,XPeriod), Difference
        else
          Data(NextFree,XYear,XPeriod) = Data(ChosenColA,XYear,XPeriod)
        end if
      else
        Data(NextFree,XYear,XPeriod) = Data(ChosenColA,XYear,XPeriod)
      end if
    else
      Data(NextFree,XYear,XPeriod) = Data(ChosenColB,XYear,XPeriod)
    end if
  end do
end do

if (DisagreeN.GE.1) then
  print*, "  > Columns are not in complete agreement. See log file for details."
  print*
  
  do XYear = 1, YearN
    do XPeriod = 1, 17
      Data(NextFree,XYear,XPeriod) = MissVal
    end do
  end do
else
  call StoreMoveOn
end if

end subroutine CombineColAB

!*******************************************************************************
! detrend column

subroutine DeTrend

print*, "  > Enter the column to detrend."
call SelectColumn
ChosenColA = ChosenCol

print*, "  > Detrend against time (=1) or another column (=2) ? "
do
	read (*,*,iostat=ReadStatus), QTimeCol
	if (ReadStatus.LE.0 .AND. QTimeCol.GE.1 .AND. QTimeCol.LE.2) exit
end do

if (QTimeCol.EQ.2) then
  call SelectColumn
  ChosenColB = ChosenCol
end if

print*, "  > Remove gradient only (=1) or intercept also (=2): "
do
	read (*,*,iostat=ReadStatus), QGradCept
	if (ReadStatus.LE.0 .AND. QGradCept.GE.1 .AND. QGradCept.LE.2) exit
end do

allocate (OpVec   (YearN), &
	  YearVec (YearN), stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: DeTrend: Allocation failure #####"
OpVec = MissVal ; YearVec = MissVal

do XPeriod = 1, 17
  do XYear = 1, YearN
    OpVec (XYear) = Data(ChosenColA,XYear,XPeriod)
    if (QTimeCol.EQ.1) YearVec (XYear) = real (YearAD(XYear))
    if (QTimeCol.EQ.2) YearVec (XYear) = Data(ChosenColB,XYear,XPeriod)
  end do
  
  call RobustLSR (YearN,YearVec,OpVec,ResultVecA,Aye,Bee,Iter)
  
  write (99,"(i4,2f12.4,i12)"), XPeriod,Aye,Bee,Iter		! records LSR info to file
    
  deallocate (ResultVecA, stat=AllocStat)
  if (AllocStat.NE.0) print*, "  > ##### ERROR: DeTrend: Deallocation failure #####"
  
  if (Aye.NE.MissVal.AND.Bee.NE.MissVal) then
   if (QGradCept.EQ.1) then
    do XYear = 1, YearN
      if (OpVec(XYear).NE.MissVal) Data(NextFree,XYear,XPeriod) = OpVec(XYear) - (Bee * YearVec(XYear)) 
    end do
   else
    do XYear = 1, YearN
      if (OpVec(XYear).NE.MissVal) Data(NextFree,XYear,XPeriod) = OpVec(XYear) - (Aye+(Bee * YearVec(XYear))) 
    end do
   end if
  end if
end do

deallocate (OpVec,YearVec, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: DeTrend: Deallocation failure #####"

call StoreMoveOn

end subroutine DeTrend

!*******************************************************************************
! load .per file (single)

subroutine LoadPerSingle

LoadFile=Blank ; Verbose=.TRUE. ; call LoadPerFile

end subroutine LoadPerSingle

!*******************************************************************************
! load .per file (bulk)

subroutine LoadPerBulk

print*, "  > Select the .per files using a filter."
call GetBatch (Blank,FileSet)
FileN = size (FileSet, 1) ; Verbose=.FALSE. ; Talkative=.FALSE.

do XFile=1,FileN
  LoadFile=FileSet(XFile)
  call LoadPerFile
end do

deallocate (FileSet, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: LoadPerBulk: Deallocation failure #####"
Talkative=.TRUE.

end subroutine LoadPerBulk

!*******************************************************************************
! load .per file

subroutine LoadPerFile

call PreLoadPer (LoadFile,FileYearN)

allocate (FileYearAD(FileYearN),FileMonthly(FileYearN,12), &
	  FileSeasonal(FileYearN,4),FileAnnual(FileYearN), stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: LoadPerFile: alloc #####"

call LoadPER (LoadFile,FileVariCode,FileYearAD, &
		FileMonthly,FileSeasonal,FileAnnual,Silent=1)
call CommonVecPer (FileYearAD,YearAD,FileYear0,FileYear1,MasterYear0,MasterYear1)

if (FileYear0.EQ.MissVal) then
  if (Verbose) print*, "  > The loaded file has no period in common with the master array."
else
  if (Verbose) print*, "  > Common period loaded: ", FileYearAD(FileYear0), FileYearAD(FileYear1)

  do XFileYear = FileYear0, FileYear1
    XMasterYear = MasterYear0 + XFileYear - FileYear0
  
    do XMonth = 1, 12
      Data(NextFree,XMasterYear,XMonth) = FileMonthly(XFileYear,XMonth)
    end do
  
    do XSeason = 1, 4
      Data(NextFree,XMasterYear,(XSeason+12)) = FileSeasonal(XFileYear,XSeason)
    end do
  
    Data(NextFree,XMasterYear,17) = FileAnnual(XFileYear)
  end do
end if

deallocate (FileYearAD,FileMonthly,FileSeasonal,FileAnnual, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: LoadPerFile: Deallocation failure #####"

VariCode(NextFree) = FileVariCode

call StoreMoveOn

end subroutine LoadPerFile

!*******************************************************************************
! load .day file
! this loads a .day file and summaries it by period
! for more complicated manipulations use daily.f90

subroutine LoadDayFile

call LoadDAY (Blank,FileVariCode,FileYearAD,FileDaily)
FileYearN = size(FileDaily,1)

print*, "  > Select the summarisation (-1=min,0=mean,1=max,2=sum): "
do
	read (*,*,iostat=ReadStatus), ChosenSumm
	if (ReadStatus.LE.0 .AND. ChosenSumm.GE.-1 .AND. ChosenSumm.LE.2) exit
end do

call CommonVecPer (FileYearAD,YearAD,FileYear0,FileYear1,MasterYear0,MasterYear1)

allocate (FileMonthly  (FileYearN,12), &
	  FileSeasonal (FileYearN,4),  &
	  FileAnnual   (FileYearN),    stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: LoadDayFile: Allocation failure #####"
FileMonthly  = MissVal
FileSeasonal = MissVal
FileAnnual   = MissVal

if (FileYear0.EQ.MissVal) then
  print*, "  > The loaded file has no period in common with the master array."
else
  print*, "  > Common period loaded: ", FileYearAD(FileYear0), FileYearAD(FileYear1)
  
  if      (ChosenSumm.EQ.-1) then
    call FillDailyMin   (MissAccept,FileYearAD,FileDaily,FileMonthly)
    call FillSeaAnnMin  (FileYearAD,FileMonthly,FileSeasonal,FileAnnual)  
  else if (ChosenSumm.EQ. 0) then
    call FillDailyMean  (MissAccept,FileYearAD,FileDaily,FileMonthly)
    call FillSeaAnnMean (FileYearAD,FileMonthly,FileSeasonal,FileAnnual)  
  else if (ChosenSumm.EQ. 1) then
    call FillDailyMax   (MissAccept,FileYearAD,FileDaily,FileMonthly)
    call FillSeaAnnMax  (FileYearAD,FileMonthly,FileSeasonal,FileAnnual)  
  else if (ChosenSumm.EQ. 2) then
    call FillDailySum   (MissAccept,FileYearAD,FileDaily,FileMonthly)
    call FillSeaAnnSum  (FileYearAD,FileMonthly,FileSeasonal,FileAnnual)  
  end if

  do XFileYear = FileYear0, FileYear1
    XMasterYear = MasterYear0 + XFileYear - FileYear0
  
    do XMonth = 1, 12
      Data(NextFree,XMasterYear,XMonth) = FileMonthly(XFileYear,XMonth)
    end do
  
    do XSeason = 1, 4
      Data(NextFree,XMasterYear,(XSeason+12)) = FileSeasonal(XFileYear,XSeason)
    end do
  
    Data(NextFree,XMasterYear,17) = FileAnnual(XFileYear)
  end do
end if

deallocate (FileYearAD,FileDaily,FileMonthly,FileSeasonal,FileAnnual, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: LoadPerFile: Deallocation failure #####"

VariCode(NextFree) = FileVariCode

call StoreMoveOn

end subroutine LoadDayFile

!*******************************************************************************
! load .lin file

subroutine LoadLinFile

call LoadLinAuto  (FileYearAD, FileCols, FileAll)

call CommonVecPer (FileYearAD,YearAD,FileYear0,FileYear1,MasterYear0,MasterYear1)

if (FileYear0.EQ.MissVal) then
  print*, "  > The loaded file has no period in common with the master array."
else
  print*, "  > Common period: ", FileYearAD(FileYear0), FileYearAD(FileYear1)

  FileColN = size (FileCols)

  print*, "  > COL NAME"
  do XCol = 1, FileColN
    print "(i8,a1,a)", XCol, " ", adjustl(trim(FileCols(XCol)))
  end do
  
  do
    if (FileColN.NE.1) then
      print*, "  > Enter the file column to fill next master column (-1=end): "
      do
	read (*,*,iostat=ReadStatus), XFileCol
	if (ReadStatus.LE.0 .AND. XFileCol.GE.-1 .AND. XFileCol.LE.FileColN) exit
      end do
    else
      XFileCol = 1
    end if
    
    if (XFileCol.GE.1 .AND. XFileCol.LE.FileColN) then
      do XFileYear = FileYear0, FileYear1
        XMasterYear = MasterYear0 + XFileYear - FileYear0
  
        do XCol = 1, 17
          Data(NextFree,XMasterYear,XCol) = FileAll(XFileCol,XFileYear)
        end do
      end do
      
      call StoreMoveOn
    end if
    
    if (XFileCol.EQ.-1.OR.FileColN.EQ.1) exit
  end do    
end if

deallocate (FileYearAD,FileCols,FileAll, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: LoadLinFile: Deallocation failure #####"

end subroutine LoadLinFile

!*******************************************************************************
! save .per file

subroutine SavePerMon

call SelectColumn

allocate (FileMonthly  (YearN,12), &
	  FileSeasonal (YearN,4),  &
	  FileAnnual   (YearN),    stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SavePerMon: Allocation failure #####"
FileMonthly  = MissVal
FileSeasonal = MissVal
FileAnnual   = MissVal

do XYear = 1, YearN
  do XMonth = 1, 12
    FileMonthly(XYear,XMonth) = Data(ChosenCol,XYear,XMonth)
  end do
  
  do XSeason = 1, 4
    FileSeasonal(XYear,XSeason) = Data(ChosenCol,XYear,(XSeason+12))
  end do
  
  FileAnnual(XYear) = Data(ChosenCol,XYear,17)
end do

if (QPerMon.EQ.1) call SavePER    (Blank,YearAD,VariCode(ChosenCol),Monthly=FileMonthly,&
		Seasonal=FileSeasonal,Annual=FileAnnual)

deallocate (FileMonthly,FileSeasonal,FileAnnual, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SavePerMon: Deallocation failure #####"

print*

end subroutine SavePerMon

!*******************************************************************************
! save .lin file

subroutine SaveLinFile

call SelectColumn

allocate (FileAll   (17,YearN), & 
	  FileCols  (17),       stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SaveLinFile: Allocation failure #####"
FileAll  = MissVal
FileCols = ""

do XYear = 1, YearN
  do XMonth = 1, 12
    FileAll(XMonth,XYear) = Data(ChosenCol,XYear,XMonth)
  end do
  
  do XSeason = 1, 4
    FileAll((XSeason+12),XYear) = Data(ChosenCol,XYear,(XSeason+12))
  end do
  
  FileAll(17,XYear) = Data(ChosenCol,XYear,17)
end do

call SaveLin (Blank, Blank, FileCols, YearAD, FileAll)

deallocate (FileAll,FileCols, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SaveLinFile: Deallocation failure #####"

print*

end subroutine SaveLinFile

!*******************************************************************************
! save .ann file from various columns

subroutine SaveAnnCol

print*, "  > Enter the number of cols in the .ann file: "
do
	read (*,*,iostat=ReadStatus), SaveColN
	if (ReadStatus.LE.0 .AND. SaveColN.GE.1) exit
end do

allocate (FileAll   (YearN,SaveColN), & 
	  AnnLabels  (SaveColN),       stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SaveAnnCol: Alloc failure #####"
FileAll=MissVal ; AnnLabels=""

do XSaveCol=1,SaveColN
  print "(a,i4)", "   > Enter the col, per (1-17) to save to .ann col:", XSaveCol
  do
	read (*,*,iostat=ReadStatus), ChosenCol,ChosenPer
	if (ReadStatus.LE.0 .AND. ChosenPer.GE.1 .AND. ChosenPer.LE.17 .AND. &
		ChosenCol.GE.1.AND.ChosenCol.LE.ColumnN) exit
  end do

  do XYear = 1, YearN
        FileAll(XYear,XSaveCol) = Data(ChosenCol,XYear,ChosenPer)
  end do
  
  String20 = GetTextFromInt (ChosenCol)
  AnnLabels(XSaveCol) = trim(String20) // ":" // trim(CalendarNames(ChosenPer))
  AnnLabels(XSaveCol) = adjustr(AnnLabels(XSaveCol))
end do

SaveFile=Blank
call SaveANN (SaveFile, YearAD, AnnLabels, FileAll)

deallocate (FileAll,AnnLabels, stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: SaveAnnCol: Deallocation failure #####"
print*

end subroutine SaveAnnCol

!*******************************************************************************
! store data in column and move on

subroutine StoreMoveOn

call CheckBlank

if (NonBlank.EQ.0) then
  if (Talkative) print*, "  > Operation has no valid results. No results stored."
  VariCode(NextFree) = MissVal
else
  if (Talkative) print*, "  > Codes for column, variable: ", NextFree, VariCode(NextFree)
  
  NextFree = NextFree + 1
  if (NextFree.GT.ColumnN) then
    NextFree = NextFree - 1
    if (Talkative) print*, "  > All columns are being used. Next free set to: ", NextFree
  end if
end if

if (Talkative) print*

end subroutine StoreMoveOn

!*******************************************************************************
! check whether column is (=0) or is not (=1) blank

subroutine CheckBlank

NonBlank = 0

XYear    = 0
do
  XYear = XYear + 1
  
  XPeriod = 0
  do
    XPeriod = XPeriod + 1
    
    if (Data(NextFree,XYear,XPeriod).NE.MissVal) NonBlank = 1
    if (NonBlank.EQ.1) exit
    if (XPeriod.EQ.17) exit
  end do
  
  if (NonBlank.EQ.1) exit
  if (XYear.EQ.YearN) exit
end do

end subroutine CheckBlank

!*******************************************************************************
! check master array

subroutine CheckMaster

print*, "  > Column  Months Seasons Annuals (missing%)"

do XCol = 1, ColumnN
 MonMiss = 0.0
 MonEn   = 0.0
 SeaMiss = 0.0
 SeaEn   = 0.0
 AnnMiss = 0.0
 AnnEn   = 0.0

 do XYear = 1, YearN
  do XMonth  = 1, 12
      if (Data (XCol,XYear,XMonth)  .NE. MissVal) then
        MonEn   = MonEn   + 1.0
      else
        MonMiss = MonMiss + 1.0
      end if
  end do

  do XSeason = 13,16
      if (Data (XCol,XYear,XSeason) .NE. MissVal) then
        SeaEn   = SeaEn   + 1.0
      else
        SeaMiss = SeaMiss + 1.0
      end if
  end do

  do XAnnual = 17,17
      if (Data (XCol,XYear,XAnnual) .NE. MissVal) then
        AnnEn   = AnnEn   + 1.0
      else
        AnnMiss = AnnMiss + 1.0
      end if
  end do
 end do

 MonFrac = 100.0 * MonMiss / (MonMiss + MonEn)
 SeaFrac = 100.0 * SeaMiss / (SeaMiss + SeaEn)
 AnnFrac = 100.0 * AnnMiss / (AnnMiss + AnnEn)

 print "(i11,3f8.2)", XCol, MonFrac, SeaFrac, AnnFrac
end do

print*

end subroutine CheckMaster

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

subroutine Conclude

deallocate (Data,YearAD,stat=AllocStat)
if (AllocStat.NE.0) print*, "  > ##### ERROR: Conclude: Deallocation failure #####"

close (99)

end subroutine Conclude

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

end program Periodic