;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; Description: Filter data for each equatorial wave type
; These files still have to be loaded manually;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
load "$NCARG_ROOT/lib/ncarg/nclscripts/contrib/kf_filter.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/contrib/time_axis_labels.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/contrib/cd_string.ncl"
;<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<;;;;;

;begin
   
  print_clock( "Here we go! " )
  f = addfile("u_v_2011_2020.nc", "r")

; These are some parameters that could be useful to have up top
  mis = -999
  mis@_FillValue = -999

  ymin    =  -30
  ymax    =   30
  
  timeUnits = f->time@units
	print(timeUnits)
	startDate = ut_inv_calendar( 2020, 05, 01, 12, 0, 0, timeUnits, 0 )
	;print (startDate)
	endDate   = ut_inv_calendar( 2020, 06, 30, 12, 0, 0, timeUnits, 0 )
	;print (endDate)
  ;if( .not.isvar("filtName") ) then
    ;filtName = "k09"
  ;end if
  obsPerDay = 1
  if( .not.isvar("varName") ) then
    varName = "u"
	;varName = "olr"
  end if
  print( varName )

  ;basePath = "~/data/trmm/"
   basePath = "/mnt/e/Super_cyclone/Datasets/wind_daily/"
; basePath = "~/data/nasa_ymc/fcst_verif/"
; basePath = "~/data/olr/current/"
; pathIn = basePath + ".std.nc"
; pathOut = basePath + ".waves.std.nc"
  ;pathIn  = basePath + "olr.day.1.nc"
  pathOut = basePath + varName + ".gonu.anom.waves.nc"
  ;pathOut = basePath + varName + ".kyarr.waves.nc"
; pathOut = basePath + "trmm3b42." + filtName + ".nc"
; pathIn  = basePath + varName + ".anom.nc"
; pathOut = basePath + varName + "." + filtName + ".nc"
; pathIn   = basePath + "anom/" + varName + ".anom.nc"
; pathOut  = basePath + "waves/" + varName + ".anom.waves.nc"

  
  ;calcMrg = True
  calcMjo = True
  calcKelvin = True
  calcEr = True
  calcMtd = True
  ;calcTd = True
  ;calcWig = True
  makeNewFile = True

; Open the input files
  ;fin  = addfile( pathIn, "r" )
  time = f->time({startDate:endDate})
  ;print(time)
  lat  = f->latitude({ymin:ymax})
  lon  = f->longitude

  lat@actual_range = (/ min(lat), max(lat) /)

; Open the output files
  setfileoption("nc","Format","LargeFile")
; setfileoption("nc","Format","NetCDF4")
; setfileoption("nc","CompressionLevel",1)
  if( makeNewFile ) then
    system( "rm " + pathOut )
    fout  = addfile( pathOut, "c" )
    fout->time = time
    fout->lat  = lat
    fout->lon  = lon
  else
    fout  = addfile( pathOut, "w" )
  end if

  print_clock( "Reading the input data... " )

  data = f->$varName$({startDate:endDate},0,0,{ymin:ymax},:)
  data&time@beginning_date = cd_string( data&time(0), "" )
  data&time@ending_date    = cd_string( data&time(dimsizes(data&time)-1), "" )
   data!1 = "lat"
   data!2 = "lon"
  ;if( any( ismissing(data) ) ) then
    ;print_clock( "WARNING: Setting missing data to zero" )
    ;data = where( ismissing(data), 0, data )
  ;end if

  printVarSummary(data)
  ;fout->unfilt = data
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 if( calcMjo ) then
    print_clock( "Filtering MJO..." )

    mjo = data
    mjo@long_name     = "Madden-Julian Oscillation in " + str_upper(varName)
    mjo@filter        = "Landu et al (2020)"
    mjo@wavenumber    = (/   1,   6 /)
    mjo@period        = (/  30,  70 /)
    mjo@depth         = (/ mis, mis /)
   

    do y = 0, ( dimsizes(lat) - 1 )
  ;    print_clock( y + " " + ( dimsizes(lat) - 1 ) + " " )

    mjo(time|:,{lat|lat(y)},lon|:) = (/ kf_filter(  \\
       data(time|:,{lat|lat(y)},lon|:), obsPerDay,  \\
       mjo@period(0),     mjo@period(1),     \\
       mjo@wavenumber(0), mjo@wavenumber(1), \\
       mjo@depth(0),      mjo@depth(1),      \\
       "none" ) /)

    end do

    ;print_clock( "Writing MJO... " )
    ;fout->mjo = mjo
 end if	

;end ; filter_waves