pro get_vceil25k_list,varname,data,site,date,missing,num_files

;***************************************************************
;  Inputs - set here or passed in through the command line
;***************************************************************

;  Site to process
;site='sgp'

;  Date to process
;date='20000713'

;  Variables to get
;varname=['base_time','time_offset','first_cbh']

;Missing data flag
;missing=-8888

;**************************************************************
;  Establish some constants
;**************************************************************

;  Save a copy of the original variables
varname_orig=varname

;  Print flag
dp='yes'

;  Main arm site directory
arm_dir='/uufs/chpc.utah.edu/common/home/mace_grp/data/arm/'+site+'/'

;  Base_time and time_offset names are the same
base_timename='base_time'
time_offsetname='time_offset'

;  Arm data directory, datastream, and varnames depends on site
if site eq 'sgp' then begin
  limit_date_12=julday(5,24,2000,0,0,0) 
  limit_date_23=julday(10,13,2003,0,0,0)   
  datastream_1='sgpblcC1.a1' 
  datastream_2='sgpvceil25kC1.a1' 
  datastream_3='sgpvceil25kC1.b1'
 
  ;  These will have to be picked up separately for each file
  ;  datastream 1
  ;first_cbh_name='cloud1'      ;m above ground level 
  ;second_cbh_name='cloud2'
  ;third_cbh_name='cloud3'
  ;detection_status_name='none'
  ;status_flag_name='status_flags'
  ;vertical_visibility_name='visibility'
  ;alt_highest_signal_name='none'
  ;  datastream 2&3
  ;first_cbh_name='first_cbh'   ;m above ground level
  ;second_cbh_name='second_cbh'
  ;third_cbh_name='third_cbh'
  ;detection_status_name='detection_status'
  ;status_flag_name='status_flag'
  ;vertical_visibility_name='vertical_visibility'
  ;alt_highest_signal_name='alt_highest_signal'
endif else if site eq 'nsa' then begin
  datastream='nsavceil25kC1.b1'
  first_cbh_name='first_cbh'   ;m above ground level
  second_cbh_name='second_cbh'
  third_cbh_name='third_cbh'
  detection_status_name='detection_status'
  status_flag_name='status_flag'
  vertical_visibility_name='vertical_visibility'
  alt_highest_signal_name='alt_highest_signal'
endif else if site eq 'twpc1' then begin
  datastream='twpvceil25kC1.b1'
  first_cbh_name='first_cbh'   ;m above ground level
  second_cbh_name='second_cbh'
  third_cbh_name='third_cbh'
  detection_status_name='detection_status'
  status_flag_name='status_flag'
  vertical_visibility_name='vertical_visibility'
  alt_highest_signal_name='alt_highest_signal'
endif else if site eq 'twpc2' then begin
  datastream='twpvceil25kC2.b1'
  first_cbh_name='first_cbh'   ;m above ground level
  second_cbh_name='second_cbh'
  third_cbh_name='third_cbh'
  detection_status_name='detection_status'
  status_flag_name='status_flag'
  vertical_visibility_name='vertical_visibility'
  alt_highest_signal_name='alt_highest_signal'
endif else if site eq 'twpc3' then begin
  datastream='twpvceil25kC3.b1'
  first_cbh_name='first_cbh'   ;m above ground level
  second_cbh_name='second_cbh'
  third_cbh_name='third_cbh'
  detection_status_name='detection_status'
  status_flag_name='status_flag'
  vertical_visibility_name='vertical_visibility'
  alt_highest_signal_name='alt_highest_signal'
endif 

;***************************************************************
;  Make sure the dates are in an array
;***************************************************************

;  Number of different dates to pick up
numdates=n_elements(date)

;  Put the search dates in an array of dates called curdate
curdate=make_array(/string,numdates)
if numdates gt 1 then begin
  for i=0,numdates-1 do begin
    curdate[i]=strcompress(string(date[i]),/remove_all)
  endfor
endif else begin
  curdate[0]=strcompress(string(date),/remove_all)
endelse

;***************************************************************
;  Loop through the dates and make an array of file names
;***************************************************************

;  Loop through the dates
for i=0,numdates-1 do begin

  ;  Separate the date strings
  yy=strmid(curdate[i],0,4)
  mm=strmid(curdate[i],4,2)
  dd=strmid(curdate[i],6,2)

  ;  sgp had datastream name changes
  check_date=julday(mm,dd,yy,0,0,0)
  if site eq 'sgp' then begin
    if check_date ge limit_date_23 then begin
      datastream=datastream_3
    endif else if check_date ge limit_date_12 then begin
      datastream=datastream_2
    endif else begin
      datastream=datastream_1
    endelse
  endif 
  
  ;  Create the filename
  arm_file=arm_dir+datastream+'/'+yy+'/'+datastream+'.'+curdate[i]+'*cdf'

  ;  Unzip them if necessary
  unix_command='gunzip -fq '+arm_file+'*gz'
  spawn,unix_command

  ;  Get a list of the files for the day
  filesd=file_search(arm_file,count=num_files)  

  ;  Put these files in one big list if files exist for this date
  if num_files gt 0 then begin
    if n_elements(files) eq 0 then begin
      files=filesd
    endif else begin
      files=[files,filesd]
    endelse
  endif
endfor  ;end of loop through the dates


;  Number of files to read and cat together
num_files=n_elements(files)

;  Continue on if files were found
if num_files gt 0 then begin

;***************************************************************
;  Loop through all the matching files and read in the data
;***************************************************************

for i=0,num_files-1 do begin  

;  Print the filename
  if dp eq 'yes' then print,files[i]

;  Separate out the date parts
  parts=strsplit(files[i],datastream,/regex)
  year=fix(strmid(files[i],parts[2]+1,4))
  month=fix(strmid(files[i],parts[2]+5,2))
  day=fix(strmid(files[i],parts[2]+7,2))

;  Calculate check date to determine var names
  check_date=julday(month,day,year,0,0,0)

;  Choose the varname depending on date

  if site eq 'sgp' then begin
    if check_date ge limit_date_12 then begin
      first_cbh_name='first_cbh'   ;m above ground level
      second_cbh_name='second_cbh'
      third_cbh_name='third_cbh'
      detection_status_name='detection_status'
      status_flag_name='status_flag'
      vertical_visibility_name='vertical_visibility'
      alt_highest_signal_name='alt_highest_signal'
    endif else begin
      first_cbh_name='cloud1'      ;m above ground level 
      second_cbh_name='cloud2'
      third_cbh_name='cloud3'
      detection_status_name='none'
      status_flag_name='status_flags'
      vertical_visibility_name='visibility'
      alt_highest_signal_name='none'
    endelse
  endif

;  List of 1d variables in the file and their dimensions
;  These variables are (time)
  var1=[first_cbh_name,$
        second_cbh_name,$
        third_cbh_name,$
        detection_status_name,$
        status_flag_name,$
        vertical_visibility_name,$
        alt_highest_signal_name]
  vardim1=make_array(n_elements(var1),/int,value=1)

;  Variables that aren't to be combined across files
  var0=[base_timename,time_offsetname]
  vardim0=make_array(n_elements(var0),/int,value=0)

;  Combine these all into one big array
  var=[var1,var0]
  vardim=[vardim1,vardim0]

;  Put the new varname in the varname array
  result=where(varname_orig eq 'cloud1' or $
               varname_orig eq 'first_cbh',count)
  if count gt 0 then varname[result]=first_cbh_name
  result=where(varname_orig eq 'cloud2' or $
               varname_orig eq 'second_cbh',count)
  if count gt 0 then varname[result]=second_cbh_name
  result=where(varname_orig eq 'cloud3' or $
               varname_orig eq 'third_cbh',count)
  if count gt 0 then varname[result]=third_cbh_name
  result=where(varname_orig eq 'status_flag' or $
               varname_orig eq 'status_flags',count)
  if count gt 0 then varname[result]=status_flag_name
  result=where(varname_orig eq 'vertical_visibility' or $
               varname_orig eq 'visibility',count)
  if count gt 0 then varname[result]=vertical_visibility_name

;  Open the netcdf file
  cdf_id=ncdf_open(files[i])

;  Get the dimension id's
  time_did=ncdf_dimid(cdf_id,'time')

;  Get the dimensions
  ncdf_diminq,cdf_id,time_did,charstring,numtimes

;  Create the dummy variables
  dummy1d=make_array(numtimes,/float,value=missing)

;  Get the variable ids
  for v=0,n_elements(varname)-1 do begin
    xstr="v"+varname[v]+"_id=ncdf_varid(cdf_id,'"+varname[v]+"')"
    result=execute(xstr)
    if dp eq 'yes' then print,varname[v]
  endfor
;
;  Loop through the variable names.  Get the variable if it exists.
;  If it doesn't exist, set the variable equal to the dummy array.
;
  for v=0,n_elements(varname)-1 do begin

;  Determine the dimension of this variable
    xstr="idx=where(var eq '"+varname[v]+"',count)"
    result=execute(xstr)
    if count eq 0 then begin 
      print,varname[v],'variable not in list of acceptable variables'
      stop
    endif else begin
      dim=vardim[idx]
      if dim eq 1 then dumname='dummy1d' $
      else if dim eq 2 then dumname='dummy2d' $
      else if dim eq 3 then dumname='dummy3d'
    endelse
    if dp eq 'yes' then print,'processing ',varname[v],'dimension',dim

;  Get the variable
    if dim ne 0 then begin
      xstr='if v'+varname[v]+'_id gt -1 then ncdf_varget,cdf_id,v'+$
        varname[v]+'_id,v'+varname_orig[v]+' else v'+varname_orig[v]+'='+dumname
    endif else begin
      xstr='if v'+varname[v]+'_id gt -1 then ncdf_varget,cdf_id,v'+$
        varname[v]+'_id,'+varname[v]   
    endelse
    result=execute(xstr)

;  Transpose 2d variables
    if dim gt 1 then begin
      if dp eq 'yes' then print,'transposing'
      xstr='v'+varname_orig[v]+'=transpose(v'+varname_orig[v]+')'
      result=execute(xstr)
    endif
  endfor

;
;  Combine the data 
;

    for v=0,n_elements(varname)-1 do begin

      if varname[v] eq 'time_offset' then begin
        if n_elements(com_btto) eq 0 then begin
          com_btto=long(base_time)+long(time_offset)
        endif else begin
          com_btto=[long(com_btto),long(base_time)+long(time_offset)]
        endelse
      endif
    
      if varname[v] ne 'base_time' and varname[v] ne 'time_offset' then begin
        xstr='if n_elements(com_'+varname_orig[v]+') eq 0 then com_'+$
              varname_orig[v]+'=v'+varname_orig[v]+' else com_'+varname_orig[v]+$
              '=[com_'+varname_orig[v]+',v'+varname_orig[v]+']'
        result=execute(xstr)
      endif

    endfor


;  Close the netcdf file
    ncdf_close,cdf_id

  endfor  ;loop through number of matching files

;***************************************************************
;  Remove the 'com_' prefix from the data array names
;***************************************************************

  for v=0,n_elements(varname)-1 do begin

    if varname[v] eq 'base_time' then begin
      base_time=long(com_btto[0])  
    endif else if varname[v] eq 'time_offset' then begin
      time_offset=long(com_btto)-long(base_time)        
    endif else begin
      xstr=varname_orig[v]+'=com_'+varname_orig[v]
      result=execute(xstr)
    endelse

  endfor  ;end of loop through variables


;*****************************************
;  Put the data in a structure to pass back out
;*****************************************

  data={data,filename:'string'}
  for v=0,n_elements(varname)-1 do begin
    xstr='data=create_struct(data,"'+varname_orig[v]+'",'+varname_orig[v]+')'
    result=execute(xstr)
  endfor  

;  Add in a variable to hold the files
  xstr='data=create_struct(data,"files",files)'
  result=execute(xstr)

;***************************************************************
;  Zip up the data files
;***************************************************************

;for i=0,num_files-1 do begin
;  unix_command='gzip -fq  '+files[i]  ;zip the data files
;  spawn,unix_command
;endfor

;***************************************************************
;  Plot the data
;***************************************************************

do_plot='yes'
if do_plot eq 'yes' then begin

;
;  Find the size of the arrays
;
  numtimes=n_elements(time_offset)  ;number of times
  numheights=n_elements(height)     ;number of heights
;
;  Calculate time array for plotting
;
;  Number of seconds in a day
  sec_per_day=24D*60D*60D

;  Get IDL julian day of day 1
  day1=julday(1,1,1970,0,0,0)
 
;  Get IDL julian day of the entire time array
  julian_day=double(day1 + (base_time+time_offset)/sec_per_day)

;  Check the dates
  caldat,julian_day,mm,dd,yy,hh,mi,ss
  print,yy[0],mm[0],dd[0],hh[0],mi[0],ss[0],'start time'
  print,yy[numtimes-1],mm[numtimes-1],dd[numtimes-1],hh[numtimes-1],$
    mi[numtimes-1],ss[numtimes-1],'end time'

;  Set up the plot
  set_plot,'Z'
  device,set_resolution=[900,600]
  loadct,39

  !p.background=!d.n_colors-1
  dummy=label_date(date_format=['%H:%I !C %N/%D'])

;  Postion plots
  xl=0.11 & xr=0.85
  yb=0.07 & yt=0.9
  sx=0.05 & sy=0.08
  numplots_x=1 & numplots_y=1
  position_plots,xl,xr,yb,yt,sx,sy,numplots_x,numplots_y,pos

  result=where(first_cbh gt missing,count)
  dmin=min(first_cbh[result])
  dmax=max(first_cbh[result])

  plot,julian_day,first_cbh,position=pos[0,*],xtickformat='label_date',$
    color=0,psym=4,yrange=[dmin,dmax],xstyle=1,ystyle=1

  image_name=datastream+'.'+date+'.gif'
  write_gif,image_name,tvrd()
  

;
;  Scale the data
;
  ;image=bytscl(CloudMaskMplCloth,max=dmax,min=dmin)

  ;contour,CloudMaskMplCloth,julian_day,heights/1000.0,/nodata,xstyle=4,ystyle=4,$
  ;  position=pos[0,*]
  ;px=!x.window*!d.x_vsize
  ;py=!y.window*!d.y_vsize
  ;sx=px[1]-px[0]+1
  ;sy=py[1]-py[0]+1
  ;tv,congrid(image,sx,sy),px[0],py[0]
  ;contour,CloudMaskMplCloth,julian_day,heights/1000.0,/noerase,/nodata,$
  ;  xtickformat='label_date',xstyle=1,ystyle=1,$
  ;  xtickunits='hour',xticklen=-0.02,yticklen=-0.02,$
  ;  color=0,xtitle='signal',$
  ;  ytitle='Height (km)',position=pos[0,*]

  ;colorbar_fanning,maxrange=dmax,minrange=dmin,$
  ;title='signal',$
  ;ncolors=256,format='(F6.2)',$
  ;vertical=1,color=0,position=cbpos[0,*]

  ;image_name=datastream+'.gif'
  ;write_gif,image_name,tvrd()
  ;stop
endif  ;end of do_plot

endif  ;end of found files

return
end