; --------------  LOAD FUNCTIONS AND PROCEDURES ----------------

    load "/usr/share/ncarg/nclscripts/csm/gsn_code.ncl"
    load "/usr/share/ncarg/nclscripts/csm/gsn_csm.ncl"
    load "/usr/share/ncarg/nclscripts/csm/contributed.ncl"
    load "/usr/share/ncarg/nclscripts/wrf/WRFUserARW.ncl"

    ; --------------  BEGINING OF NCL SCRIPT ----------------

    begin
    ;********************************************************
    ; read in netCDF file and make a loop for all time steps
    ;********************************************************
      in     = addfile("/home/isakhar/WRF/WRFV3/test/em_real/wrfout_d01_2014-01-02_18:00:00.nc","r")


    times = wrf_user_list_times(in)     	; get times in the file
    ntimes = dimsizes(times)            	; number of times in the file

    radt0 = new(ntimes,float)            	; creation of a radt vector at each time step

    do it = 0,ntimes-1				;Loop for the time: it= starting time
    time = it


 ;************************************************
    ;  - Select lon & lat of the point of interest -
    ;************************************************
    ; - The function wrf_user_ll_to_ij find the nearest grid point for a specific lat and lon

    Latitude =   -8.550740
    Longitude = 125.658959
    res = True      
    res@returnInt = True        		               ; False : return real values, True: return interger values
    point = wrf_user_ll_to_ij(in,Longitude,Latitude,res)       ; wrf_user_ll_to_ij(nc_file,lon,lat,opt)
     
     x = point(0)
     y = point(1)


    ;*************************************************************************************
    ;  - Extract the radiation-  *
    ;*************************************************************************************

   cz=in->SWDDIR
  	
   rad = 4.0*atan(1)/180
   cza = acos(cz)/rad
   radt0(it)=cza(x,y)
 
  end do
   

;************************************************************
;  - Print the variable -                    *
;************************************************************

print("Time(hour)        Radiation (W/m2)")

do it = 0,ntimes-1
 
        print (sprintf("%5.0f",it)    +" " \
	    +sprintf ("%0.4f", radt0(it)) +"  " )   
    
  end do

asciiwrite (fName , data)


end