begin
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; Part1: Reading .dat file(2_stations.dat)
;  csv_2.ncl

	print("Part1: Reading .dat files (2_stations.dat)")
  
	drst	= "./"
	flst	= "2_stations.dat"
	ptst	= drst + flst

	lnst	= asciiread(ptst,-1,"string")

	delim	= ","

;---Read fields 1-5 from fname_Report
	stid	= tointeger(str_get_field(lnst(1:),1,delim))
	stcd	= str_get_field(lnst(1:),2,delim)
	stnm	= str_get_field(lnst(1:),3,delim)
	stlon	= tofloat(str_get_field(lnst(1:),4,delim))
	stlat	= tofloat(str_get_field(lnst(1:),5,delim))
  
;---Print the information
	print("stid is '" + stid + "', stnm is " + stnm)
;	print("End Part1 Reading 2_stations.dat files")

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; Part2: Creating Coordinates for SMOPS data

	nlat				= 720
	latitude 			= fspan(-89.875, 89.875, nlat)
;	latitude			= latGlobeFo(nlat, "lat", "latitude", "degrees_north")
	latitude@units		= "degrees_north"
	latitude@long_name	= "latitude"
	latitude!0			= "latitude"
	latitude&latitude	=  latitude
	printVarSummary(latitude)
	printMinMax(latitude,0)
;	print(latitude)
	
	mlon				= 1440
	longitude			= lonGlobeFo(mlon, "longitude", "longitude", "degrees_east")  ; 0->360
	longitude			= (/ longitude - 180. /)  ; subtract 180 from all values 
	longitude&longitude	= longitude            ; update coordinates       printVarSummary(lon)
	printVarSummary(longitude)
	printMinMax(longitude,0)
;	print("End Creating Coordinates for SMOPS data")

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; Part3: Reading SMOPS data

	SMOPS_fnames		= systemfunc("ls ./NPR_SMOPS_CMAP_D*.nc")
;	print(SMOPS_fnames)
	SMOPS_fnamesc		= str_get_cols(SMOPS_fnames,18,25)
	print(SMOPS_fnamesc)

	print("dimsizes(SMOPS_fnames)	= "+dimsizes(SMOPS_fnames))
	print("dimsizes(stnm)			= "+dimsizes(stnm))
	SMOPS_fod						= new((/dimsizes(SMOPS_fnames),dimsizes(stnm)/),float)

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;Creating file to write SM_SMOPS over grid points.
	header_SMOPS_G		=(/"   date      lat      lon    sm_SMOPS_Grid"/)
	hlist_SMOPS_G		=[/header_SMOPS_G/]

	grid_output_file	= "20_SMOPS_Grid.txt"
	interp_output_file	= "20_SMOPS_Interp.txt"
;;;	write_table(grid_output_file,"w",hlist_SMOPS_G,"%s ")
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
	do_loop_start_time	= get_cpu_time()
	do i=0, dimsizes(SMOPS_fnames)-1
 
	  setfileoption("nc", "FileStructure", "Advanced")
      SMOPS_fi			= addfile(SMOPS_fnames(i), "r")
      printVarSummary(SMOPS_fi)
      ;print(SMOPS_fi)
    
	  sm				= (/short2flt( SMOPS_fi->Blended_SM )/)	; short2flt( SMOPS_fi->Blended_SM )	; 
	  sm!0				= "latitude"
	  sm!1				= "longitude"
	  sm&latitude		= latitude
	  sm&longitude		= longitude 
	  printVarSummary(sm)

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;-----'fill' regions containing missing values [_FillValue] with interpolated values
;	  https://www.ncl.ucar.edu/Applications/grid_fill.shtml		---> grid_fill_1.ncl
	  nscan = 2000             ; usually *much* fewer
	  eps   = 0.001            ; variable depended
	  gtype = False            ; "gendat" does not generate cyclic fields
	  guess = 0                ; use zonal means
	  relc  = 0.6              ; standard relaxation coef
	  opt   = 0

	  poisson_grid_fill(sm, gtype, guess, nscan, eps, relc, opt)
;;;	  sm	= linmsg(sm,0)

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
	  
;-----Writing sm_SMOPS over grid points in the file
      start_time		= get_cpu_time()
      nlat				= dimsizes(sm&latitude)
      nlon				= dimsizes(sm&longitude)
      dims				= (/nlat,nlon/)
      smops_fnamesc_tmp	:= ndtooned(conform_dims(dims,SMOPS_fnamesc(i),-1))
      lat1d				:= ndtooned(conform_dims(dims,sm&latitude,0))	; or ndtooned(conform_dims(dims,latitude,0))
      lon1d				:= ndtooned(conform_dims(dims,sm&longitude,1))	; or ndtooned(conform_dims(dims,longitude,0))
;;;	  write_table(grid_output_file,"a",[/smops_fnamesc_tmp,lat1d,lon1d, sm/], " %s %7.4f %7.4f    %5.1f")
      end_time = get_cpu_time()
      print("Elapsed time for writing SMOPS_SM_Grid file (FAST) = " + (end_time-start_time) + " CPU seconds")
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
	  SMOPS_fo = linint2_points_Wrap(sm&longitude,sm&latitude,sm, True , stlon, stlat, 0)
	  SMOPS_fod(i,:)=SMOPS_fo

	end do
	do_loop_end_time = get_cpu_time()
	print("Elapsed time for do loop = " + (do_loop_end_time-do_loop_start_time) + " CPU seconds")
	print("Finish Writing sm_SMOPS over grid points in 21_SMOPS_Grid.txt")

	;Generating dates 20171001-20171031
	TIME = yyyymmdd_time(2017, 2017, "integer")
	time = TIME({20171001:20171031})    ; coordinate subscripting
	
	header_SMOPS=(/"   stid  stcd     lat     lon     date    sm_SMOPS"/)
	hlist_SMOPS=[/header_SMOPS/]
	write_table(interp_output_file,"w",hlist_SMOPS,"%s ")

	start_time = get_cpu_time()
	ntim = dimsizes(time)
	nlat = dimsizes(stlat)
	dims = (/nlat,ntim/)
	
	time1d   := ndtooned(conform_dims(dims,time,1))
	stid1d   := ndtooned(conform_dims(dims,stid,0))
	stcd1d   := ndtooned(conform_dims(dims,stcd,0))
	stlat1d  := ndtooned(conform_dims(dims,stlat,0))
	stlon1d  := ndtooned(conform_dims(dims,stlon,0))
	write_table(interp_output_file,"a",[/stid1d,stcd1d,stlat1d,stlon1d,time1d,SMOPS_fod/], \
				"  %5i  %s   %7.4f %7.4f %8i   %9.6f")
	end_time = get_cpu_time()
	print("Elapsed time for writing SMOPS_SM_Interpolated file (FAST) = " + (end_time-start_time) + " CPU seconds")

	print("End Part2 Interpolating SMOPS files")

 end