load "$NCARG_ROOT/nclscripts/csm/contributed.ncl"


begin

print("Beginning the program...")

  dir_in  = (/"/archive/data/cfsr/1deg/wind/"/)

  months  = sprinti("%0.2i",ispan(1,12,1))
  years   = sprinti("%0.4i",ispan(1995,2009,1))

  nyears  = dimsizes(years)
  nmonths = dimsizes(months)

  nfiles = nmonths*nyears

  all_fil_u = new((/nfiles/),"string")   ; Specify input file strings

  do yr = 0,nyears-1,1         
     do m = 0,nmonths-1,1
       all_fil_u(nmonths*yr+m) = dir_in+"wnd10mx0.5.gdas."+years(yr)+months(m)+".grb2"
     end do
  end do

  print("File strings created.")


print("Reading in all times now...")
inf = addfiles(all_fil_u,"r")
ListSetType(inf, "cat")

print("Times have been read in.")
time = inf[:]->initial_time0_hours
ut_time = ut_calendar(time, 0)
year = floattointeger(ut_time(:,0))
mn = floattointeger(ut_time(:,1))
dy = floattointeger(ut_time(:,2))
hour =  floattointeger(ut_time(:,3))
doy = day_of_year(year, mn, dy)


do hr=0,18,6 ; this will step through every initial hour. 
   do fhr=0,5;	this will step through the forecasts, which are the equiv times of (hr+fhr) UTC

      climo = new( (/ 366, dimsizes(inf[0]->lat_0), dimsizes(inf[0]->lon_0) /), "float" )
      climo!0 = "time"
      climo!1 = "lat"
      climo!2 = "lon"

         do d=1,365
            print("hr: "+hr+ ", fhr: "+fhr+ ", d: "+d)
            doyind = ind(doy.eq.d .and. hour.eq. hr )   ; doyind now contains indices for e.g. Jan 1st’s at 0z so 31 indices. 
            data = inf[:]->UGRD_P0_L103_GLL0(doyind, fhr, ::-1, :)  
            climo(d-1,:,:) = dim_avg_n_Wrap(data,0)
         end do

      printVarSummary(climo)
      climo(365,:,:) = dim_avg_n_Wrap( climo( (/0,365/) ,:,:) ,0)
      climoHarm =  smthClmDayTLL(climo, 4)  ;first 4 harmonics
      copy_VarAtts(climo, climoHarm )

     print( "Writing the file..." )
     fileName = "/archive/data/cfsr/1deg/climos/u_10m_"+(hr+fhr)+"UTC.nc"  
     system( "rm " + fileName )
     setfileoption( "nc", "Format", "LargeFile" )
     outFile = addfile( fileName, "c" )
     outFile->climotology = climoHarm 

    delete(data)
    delete(climoHarm)
    delete(climo)
   delete(doyind)

   end do
end do

end