Program filter_TC
! for: filter_TC based on the output of tcdetect scheme
! by: hapt
use netcdf
implicit none
include 'netcdf.inc'

Integer:: total_num
Integer, parameter:: max_tcstep=100
Integer:: inum, istep, itc, ntc
Integer:: nstep
Integer,dimension(300):: tc_step
Integer:: life_length
Real,parameter:: flag=-99.99
Real, dimension (:,:), allocatable:: tc_lon, tc_lat
Integer, dimension (:,:), allocatable:: tc_year, tc_mon, tc_day, tc_hour
Character:: ctmp
Character*100:: iDir, oDir
Character*100:: cfile, ofile
Character(len=4):: cyear
Character(len=2):: cmon, cday
Integer,parameter:: nmon=12
Integer, dimension(nmon):: ntc_mon
Integer:: imon
! Namelist
namelist/nml/iDir,cfile,oDir,life_length

!iDir="/work/users/huyhd/Detect/for_CoKim/2013/test_CCAM.105.0/exp15/ccam_cfsnl+atmos/"

!iDir=""
!cfile="TCs.txt"

Open (111,file='info.txt',status='unknown')
      Read(111,nml)
!      Write(*,*) iDir
Close(111)
   
!Write(*,*) iDir
itc=0
ntc_mon=0
Open (11, file=trim(iDir)//trim(cfile)//".txt", status="unknown")
      Read (11,*) ctmp, total_num
      Allocate (tc_lon(total_num,max_tcstep))
      Allocate (tc_lat(total_num,max_tcstep))
      Allocate (tc_year(total_num,max_tcstep))
      Allocate (tc_mon(total_num,max_tcstep))
      Allocate (tc_day(total_num,max_tcstep))
      Allocate (tc_hour(total_num,max_tcstep))
      !Write (*,*) total_num
      Do inum=1,total_num
         Read(11,*)
         Read(11,*)
         Read(11,*) ctmp, nstep
         !Write(*,*)'inum', inum,nstep
         If (nstep.ge.life_length) then
             itc=itc+1
             tc_step(itc)=nstep
             !Write(*,*) itc, nstep
             Do istep=1,nstep
                Read(11,*) tc_year(itc,istep), tc_mon(itc,istep), tc_day(itc,istep), tc_hour(itc,istep), tc_lon(itc,istep), tc_lat(itc,istep)
                If (istep.eq.1) then
                   ntc_mon(tc_mon(itc,istep))=ntc_mon(tc_mon(itc,istep))+1
                End if
             End do
         Else
             Do istep=1,nstep
                Read(11,*)
                !Write(*,*) inum
             End do
         End if
      End do
Close(11)
ntc=itc
!Write(*,*) itc

! Write output data 
oDir=iDir
Open(4,file=trim(oDir)//"output_"//trim(cfile)//".txt",status="unknown")
     Write(4,*) "TC# ", ntc
     Do itc=1,ntc
        Write(4,*) " ---------------------------------- "
        Write(4,*) "TC#", itc
        Write(4,*) "Num_Obs:", tc_step(itc)
        Do istep=1,tc_step(itc)
           Write(4,"(2(I2,3x),I4,2x,I4,2x,3(F6.1,2x))") tc_day(itc,istep), tc_mon(itc,istep), tc_year(itc,istep), tc_hour(itc,istep), tc_lon(itc,istep), tc_lat(itc,istep)
        End do
     End do
Close(4)

Open(5,file=trim(oDir)//"firstget_"//trim(cfile)//".txt",status="unknown")
     Do itc=1,ntc
        istep=1
        Write(5,"(2(I2,3x),I4,2x,I4,2x,3(F6.1,2x))") tc_day(itc,istep),tc_mon(itc,istep),tc_year(itc,istep), tc_hour(itc,istep), tc_lon(itc,istep),tc_lat(itc,istep)
     End do
Close(5)

Do itc=1,ntc
   Write (cyear,'(I4)') tc_year(itc,1)
   Write (cmon,'(I2)') tc_mon(itc,1)
   Write (cday,'(I2)') tc_day(itc,1)
   !write(*,*) cmon
   If (tc_mon(itc,1).lt.10) then
       Write (cmon,'(I1)') tc_mon(itc,1)
       cmon=trim("0"//trim(cmon))
   End if

   If (tc_day(itc,1).lt.10) then
       Write (cday,'(I1)') tc_day(itc,1)
       cday=trim("0"//trim(cday))
   End if
   Open (6,file=trim(oDir)//"TC_"//trim(cyear)//trim(cmon)//trim(cday)//".txt")
          Do istep=1,tc_step(itc)
           Write(6,"(2(I2,3x),I4,2x,I4,2x,3(F6.1,2x))") tc_day(itc,istep), tc_mon(itc,istep), tc_year(itc,istep), tc_hour(itc,istep), tc_lon(itc,istep),tc_lat(itc,istep)
          End do
   Close (6)   
End do

Open(3,file=trim(Odir)//"TCmon_"//trim(cfile)//".txt",status="unknown")
     Write(3,'(A3,2x,I4,a3)') "Mon", tc_year(1,1), "_dt"
     Do imon=1,nmon
        Write(3,'(I3,2x,I4)') imon, ntc_mon(imon)
     End do
Close(3)

Open(22,file=trim(Odir)//"TC_track_"//trim(cfile)//".txt",status="unknown")
     Write(22,*) "CLatitude     CLongitude      MaxSpeed        Pressure     Category        Name"
     Do itc=1,ntc
         Do istep=1,tc_step(itc)
           Write(22,"(2(F6.1,2x),3(a7,2x),I4)") tc_lat(itc,istep), tc_lon(itc,istep),'55','-999','TS',itc
        End do
     End do
Close(22)
End