;
; This program processes a one dimensional variable and writes it to the output
; netcdf file.  Determines the number of data points, number of good data
; points,mode, minimum, and maximum value for each time interval
;

pro one_dimensional_mode2,dataname

common ncdf_file, oldcdfname,twoD,time_dim_name,height_dim_name,$
  height_var_name,avgfilename,curdate

common old_grid, num_times,num_heights,base_time,time_offset,$
  height,site,files,num_files,num_days

common new_grid, numtimes, newtimes, numheights, newheights,$
 timespacing, heightspacing

common bracket, startbst

common data, avg_data

;
;  Loop through all the matching files in the directory
;
m=-1  ;counts number of cdf files without the variable
p=-1  ;counts number of cdf files with variable
j=-1  ;counts number of cdf files
for i = 0, num_files-1 do begin
;
;  Don't try to process any gif files that may be in the directory
;
  if strpos(files[i],'gif') eq -1 and $
  strpos(files[i],'png') eq -1 then begin 
    j=j+1
    print,i,j,files[i]
;
;  Open the netcdf file
;
    cdf_id=ncdf_open(files[i])
;
;  Get the dimension id's
;
    time_did=ncdf_dimid(cdf_id,time_dim_name)
;
;  Get the dimensions
;
    ncdf_diminq,cdf_id,time_did,charstring,ntimes
;
;  Create array for nodata value
;
    dummy=make_array(ntimes,/float,value=-8888)
;
;  Get the variable id
;
    data_id=ncdf_varid(cdf_id,dataname)
;
;  Get the variable if it exists
;
    if data_id ne -1 then begin
      print,'Processing 1-D ',dataname
      ncdf_varget, cdf_id, data_id, data
      p=p+1
    endif else begin
      print, 'Variable not in this file'
      m=m+1
    endelse
;
;  If it is the first real variable
;
    if p eq 0 and data_id ne -1 then begin
      print, 'first variable, getting the attributes'
      datainfo=ncdf_varinq(cdf_id, data_id)  ;inquire about the variable
      numatts=datainfo.natts  ;number of attributes assigned to the variable
      parts=strsplit(files[i],'/',/extract)
      part=strsplit(parts[5],'.',/extract)
      origfile=part[0]+'.'+part[1]
;
;  Store the attribute names
;
      dataatr=strarr(numatts)  ;string array to hold attribute names
      for k=0,numatts-1 do begin
;  Returns the name of the attribute
         dataatr[k]=ncdf_attname(cdf_id, data_id, k)  
      endfor
    endif
;
;  If it is not the first real variable
;
    if p eq 0 and m lt 0 and data_id ne -1 then begin
      print,'first file and variable exists' 
      combined_data=data
    endif else if p ge 0 and data_id ne -1 then begin
      print,'adding on existing variable'
      combined_data=[combined_data,data]
    endif else if p lt 0 and m eq 0 and data_id eq -1 then begin
      print,'first file and variable doesnt exist'
      combined_data=dummy
    endif else if m ge 0 and data_id eq -1 then begin
      print,'adding on dummy variable'
      combined_data=[combined_data,dummy]
    endif
;
;  Close the netcdf file
;
    ncdf_close,cdf_id
  endif  ;end of not a gif file
endfor  ;end of looping through the matching files

;
;  Continue processing data because some data exists
;

if j gt -1 and p gt -1 then begin

data=combined_data

;
;  Remove bad values
;

result=where(data ge 1e18, count)
if count ne 0 then data[result]=-9999

;
;  Convert cloudtop and cloudbase from msl to agl
;

if dataname eq 'RadarFirstTop' or dataname eq 'CloudBaseBestEstimate' then begin
  if site eq 'sgp' then begin
    result=where(data ne -9999 and data ne -8888, count)
    if count ne 0 then data[result]=data[result]-318.
  endif
  if site eq 'nsa' then begin
    result=where(data ne -9999 and data ne -8888, count)
    if count ne 0 then data[result]=data[result]-7.
  endif
  if site eq 'twpc1' then begin
    result=where(data ne -9999 and data ne -8888, count)
    if count ne 0 then data[result]=data[result]-3.
  endif
  if site eq 'twpc2' then begin
    result=where(data ne -9999 and data ne -8888, count)
    if count ne 0 then data[result]=data[result]-6.
  endif
endif

;
;  Give all bad or missing data vales the same bad flag of -9999
;

if dataname eq 'CloudBaseBestEstimate' or dataname eq 'RadarFirstTop' then begin
  result=where(data eq 0, count)
  if count ne 0 then data[result]=-9999
endif

if dataname eq 'liq' then begin
  result=where(data eq 0, count)
  if count ne 0 then data[result]=-9999
endif

;
;  Find the mode value
;

possible_points=fltarr(numtimes)  ;all possible data points within the time interval, both useable and nonuseable data values
data_points=fltarr(numtimes)  ;all valid data points within the time interval, only useable data values
avg_data=fltarr(numtimes)  ;the mode for the time interval
hi_ext=fltarr(numtimes)  ;the highest data value in the interval
lo_ext=fltarr(numtimes)  ;the lowest data value in the interval

for t=0,numtimes-2 do begin
  result=where(time_offset ge newtimes[t] and $
               time_offset lt newtimes[t+1] and $
               data ne -8888, count_all)  ;possible data within time interval
  possible_points[t]=count_all  ;possible points within the time interval
  result1=where(time_offset ge newtimes[t] and $
                time_offset lt newtimes[t+1] and $
                data ne -9999 and data ne -8888, count1)  ;good data in the time interval
  data_points[t]=count1  ;possible points with a valid data point

  if count1 ne 0 then begin  
    sub_data=data[result1]  ;good data
    hist=histogram(sub_data,omax=hi_extreme,omin=lo_extreme)
    bins=findgen(n_elements(hist))+min(sub_data)
    hist_max=max(hist,ind)   
    if hist_max eq 1 and count1 gt 1 then avg_data[t]=-9999 else $
    avg_data[t]=bins(ind)  ;mode value
    hi_ext[t]=hi_extreme
    lo_ext[t]=lo_extreme
  endif else begin
    avg_data[t]=-9999
    hi_ext[t]=-9999
    lo_ext[t]=-9999
  endelse
endfor

;
;  Calculate the average for the last data point
;

result_all=where(time_offset ge newtimes[numtimes-1] and $
                 time_offset lt (newtimes[numtimes-1]+300) and $
                 data ne -8888, count_all)
possible_points[numtimes-1]=count_all
result=where(time_offset ge newtimes[numtimes-1] and $
             time_offset lt (newtimes[numtimes-1]+300) and $
             data ne -9999 and data ne -8888, count)
data_points[numtimes-1]=count
if count ne 0 then begin
	sub_data=data[result]  ;good data
	hist=histogram(sub_data,omax=hi_extreme,omin=lo_extreme)
	bins=findgen(n_elements(hist))+min(sub_data)
	hist_max=max(hist,ind)   
	if hist_max eq 1 and count gt 1 then avg_data[numtimes-1]=-9999 else avg_data[numtimes-1]=bins(ind)  ;mode value
	hi_ext[numtimes-1]=hi_extreme
	lo_ext[numtimes-1]=lo_extreme
endif else begin
	avg_data[numtimes-1]=-9999
	hi_ext[numtimes-1]=-9999
	lo_ext[numtimes-1]=-9999
endelse

;
;  Plot to check data
;

plot1='no'  ;flag this yes if you want to see plots to the screen, no if you don't
if plot1 eq 'yes' then begin

;
;  Find the minimum data value
;
	result=where(lo_ext le -8888, count)
	if count ne 0 then lo_ext[result]=1e8
	dminv=min(lo_ext)
	result=where(lo_ext eq 1e8, count)
	if count ne 0 then lo_ext[result]=-9999

;
;  Find the data maximum
;

	dmaxv=max(hi_ext)

	loadct,5
	window,2,title='avg  '+dataname
	plot,newtimes,avg_data,$
		psym=6,$;color=255,$
		yrange=[dminv,dmaxv]
	oplot,newtimes,hi_ext,psym=4,color=355
	oplot,newtimes,lo_ext,psym=5,color=255

	window,3,title='raw '+dataname
	plot,time_offset,data,$
		psym=6,$
		yrange=[dminv,dmaxv],$
		xrange=[newtimes[0],newtimes[numtimes-1]]
endif  ;plots to the screen

;
;  Open the netcdf file
;

avg_id=ncdf_open(avgfilename, /write)

;
;  Get the variable dimension
;

time_id=ncdf_dimid(avg_id,'time')

;
;  Put the netcdf file into define mode
;

ncdf_control, avg_id, /redef

;
;  Define variable
;

possible_points_id=ncdf_vardef(avg_id,dataname+'PossiblePoints',time_id,/short)
data_points_id=ncdf_vardef(avg_id,dataname+'UseablePoints',time_id,/short)
avg_data_id=ncdf_vardef(avg_id, dataname+'Mode', time_id)
hi_ext_id=ncdf_vardef(avg_id, dataname+'MaxExtreme', time_id)
lo_ext_id=ncdf_vardef(avg_id, dataname+'MinExtreme', time_id)
;
;  Go back to an old file and get the information
;  for a variable.   Necessary for attribute copying
;

for i=0,num_files-1 do begin
  if strpos(files[i],'gif') eq -1 and $
  strpos(files[i],'png') eq -1 then begin 
    cdf_id=ncdf_open(files[i])
    data_id=ncdf_varid(cdf_id,dataname)
    if data_id ne -1 then goto, found_var
    ncdf_close,cdf_id
  endif
endfor
found_var:  ;skip out of loop because found a file containing the variable
;
;  Define attributes
;

for i=0,numatts-1 do begin
  nowatt=dataatr[i]
  mytry=ncdf_attcopy(cdf_id, data_id, nowatt, avg_id, avg_data_id)
  mytry=ncdf_attcopy(cdf_id, data_id, nowatt, avg_id, hi_ext_id)
  mytry=ncdf_attcopy(cdf_id, data_id, nowatt, avg_id, lo_ext_id)
endfor
ncdf_close,cdf_id
if dataname eq 'RadarFirstTop' or dataname eq 'CloudBaseBestEstimate' then begin
  ncdf_attdel,avg_id,avg_data_id,'units'
  ncdf_attput,avg_id,avg_data_id,'units ','meters AGL'
  ncdf_attdel,avg_id,hi_ext_id,'units'
  ncdf_attput,avg_id,hi_ext_id,'units ','meters AGL'
  ncdf_attdel,avg_id,lo_ext_id,'units'
  ncdf_attput,avg_id,lo_ext_id,'units ','meters AGL'
endif
ncdf_attput, avg_id, avg_data_id, 'original data source file ',origfile
ncdf_attput, avg_id, possible_points_id, 'long_name','number of possible data points within the time interval'
ncdf_attput, avg_id, data_points_id, 'long_name','number of useable points withing the time interval'

;
;  Change from define mode to data mode
;

ncdf_control, avg_id, /endef

;
;  Write variable
;

ncdf_varput, avg_id, avg_data_id, avg_data
ncdf_varput, avg_id, hi_ext_id, hi_ext
ncdf_varput, avg_id, lo_ext_id, lo_ext
ncdf_varput, avg_id, possible_points_id, possible_points
ncdf_varput, avg_id, data_points_id, data_points

;
;  Close the netcdf file
;

ncdf_close, avg_id

endif  ;variable does exist in at least one file

return
end