pro read_mrr_ave_txt_sub,file,julian_day,height,transfer_function,$
  frequency,diameter,number,path_attenuation,reflectivity_attenuated,$
  reflectivity,rain_rate,liquid_water_content,fall_velocity

;  Undefine the variable
;pres_o=!NULL
;temp_o=!NULL
;dp_o=!NULL
;rh_o=!NULL
;wspd_o=!NULL
;deg_o=!NULL
;alt_o=!NULL
U=''
SVS=''
Frequency=replicate(!values.f_nan,[64,31])
Diameter=replicate(!values.f_nan,[64,31])
Number=replicate(!values.f_nan,[64,31])
nans=replicate(!values.f_nan,31)
dummy=replicate(!values.F_NAN,[64,31])
dummyd=replicate(!values.F_NAN,[64,31])
dummyn=replicate(!values.F_NAN,[64,31])
first=0

;  Get the data from the txt file
;file='/Users/u0029340/Documents/data/Capricorn/MRR_DATA/AveData/201603/0314.ave'
;file='/uufs/chpc.utah.edu/common/home/mace-group2/capricorn/mrr/0316.ave'
openr,lun,file,/get_lun

while not eof(lun) do begin
  line=''
  readf,lun,line
  ;print,line
  indicator=''
  reads,line,indicator,format='(A1)'
  ;*** M *** TIME header line
  if indicator eq 'M' then begin
    u=''
    reads,line,y,m,d,h,mm,s,U,ave,stp,seal,smp,svs,dvs,dsn,cc,mdq,$
      format='(4x,I2,I2,I2,I2,I2,I2,1x,A3,4x,I6,4x,I6,4x,I6,4x,E6,4x,A8,4x,F5,4x,I11,3x,I8,4x,I4)'    
    if y lt 100 then y=y+2000
    if (first eq 0) then begin
      year=y & month=m & day=d & hour=h & minute=mm & seconds=s & UTC=U
      avg=ave & resolution=stp & asl=seal & sample_rate=smp
      ;Noise_0=nf0 & Noise_1=nf1
      Version_service=svs & Version_firmware=dvs & serial=dsn
      calibration = cc & valid_spectra=mdq          
    endif else begin
      year=[year,y] & month=[month,m] & day=[day,d] & hour=[hour,h] & minute=[minute,mm] & seconds=[seconds,s] & UTC=[UTC,U]
      avg=[avg,ave] & resolution=[resolution,stp] & asl=[asl,seal] & sample_rate=[sample_rate,smp]
      ;Noise_0=[Noise_0,nf0] & Noise_1=[Noise_1,nf1]
      Version_service=[Version_service,svs] & Version_firmware=[Version_firmware,dvs] & serial=[serial,dsn]
      calibration=[calibration,cc] & valid_spectra=[valid_spectra,mdq] 
    endelse  
    ;print,y,m,d,h,mm,s,U,ave,stp,seal,smp,svs,dvs,dsn,cc,mdq     
  ;*** H ***  HEIGHT header line
  endif else if indicator eq 'H' then begin
    he=fltarr(31)
    reads,line,he,format='(3x,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7,I7)'
    if (first eq 0) then begin
      ;height = he
      height = transpose(he)
    endif else begin
      ;height = [[height],[he]]
      height = [height,transpose(he)]
    endelse   
  ;*** T ***  TRANSFER FUNCTION header line     
  endif else if indicator eq 'T' then begin
    trans=fltarr(31)
    parts=strsplit(line,' ',/extract)
    if n_elements(parts) gt 4 then begin
      reads,line,trans,format='(3x,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7)'
    endif else begin
      trans[*]=-9999
    endelse
    if (first eq 0) then begin
      ;Transfer_function = trans
      Transfer_function = transpose(trans)
    endif else begin
      ;Transfer_function = [[Transfer_function],[trans]]  
      Transfer_function = [Transfer_function,transpose(trans)]      
    endelse 
  ;*** F ***  Spectral reflectivities (dB) at each height             
  endif else if indicator eq 'F' then begin
    u=0
    reads,line,u,format='(1x,I2)'
    num=strlen(line)
    num=fix((num-3.)/7.)
    ;num=31
    if (num gt 0) then freq=replicate(!values.f_nan,num) else freq=!values.f_nan
    if (num gt 0) then reads,line,freq,format='(3x,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7)'
    nan_replace=where(freq eq 0.0,count)
    if (count gt 0) then freq[nan_replace]=!values.F_NAN
    ;if (num le 30) then freq=[freq,nans[num:30]]
    ;if (num eq 0) then freq=nans
    if (first eq 0) then begin
      Frequency[u,*]=freq
    endif else begin
      dummy[u,*]=freq       
    endelse  
  ;*** D ***  DROP SIZE (mm)
  endif else if indicator eq 'D' then begin
    reads,line,u,format='(1x,I2)'
    num=strlen(line)
    num=fix((num-3.)/7.)
    if (num gt 0) then diam=replicate(!values.f_nan,num) else diam=!values.f_nan
    if (num gt 0) then reads,line,diam,format='(3x,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7)'
    nan_replace=where(diam eq 0.0,count)
    if (count gt 0) then diam[nan_replace]=!values.F_NAN
    ;if (num le 30) then diam=[diam,nans[num:30]]
    ;if (num eq 0) then diam=nans
    if (first eq 0) then begin
      Diameter[u,*] = diam
    endif else begin
      dummyd[u,*] = diam       
    endelse  
  ;*** N *** Number     
  endif else if indicator eq 'N' then begin
    reads,line,u,format='(1x,I2)'
    num=strlen(line)
    num=fix((num-3.)/7.)
    if (num gt 0) then nums=replicate(!values.f_nan,num) else nums=!values.f_nan
    if (num gt 0) then reads, line, nums, format='(3x,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7)'
    nan_replace=where(nums eq 0.0,count)
    if (count gt 0) then nums[nan_replace]=!values.F_NAN
    if (num le 30) then nums=[nums,nans[num:30]]
    if (num eq 0) then nums=nans
    if (first eq 0) then begin
      Number[u,*] = nums
    endif else begin
      dummyn[u,*] = nums       
    endelse   
  ;*** P *** Path Integrated Attenuation (dB)
  endif else if indicator eq 'P' then begin 
    num=strlen(line)
    num=fix((num-3.)/7.)
    if (num gt 0) then path=replicate(!values.f_nan,num) else path=!values.f_nan
    if (num gt 0) then reads, line, path, format='(3x,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7)'
    if (num le 30) then path=[path,nans[num:30]]
    if (num eq 0) then path=nans
    if (first eq 0) then begin
      ;Path_attenuation = path
      Path_attenuation = transpose(path)
    endif else begin
      ;Path_attenuation = [[Path_attenuation],[path]]  
      Path_attenuation = [Path_attenuation,transpose(path)]      
    endelse  
  ;*** z *** Attenuated Radar Reflectivity (dBZ)     
  endif else if indicator eq 'z' then begin
    num=strlen(line)
    num=fix((num-3.)/7.)
    if (num gt 0) then a_reflect=replicate(!values.f_nan,num) else a_reflect=!values.f_nan
    if (num gt 0) then reads, line, a_reflect, format='(3x,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7)'
    if (num le 30) then a_reflect=[a_reflect,nans[num:30]]
    if (num eq 0) then a_reflect=nans
    if (first eq 0) then begin
      ;Reflectivity_attenuated = a_reflect
      Reflectivity_attenuated = transpose(a_reflect)
    endif else begin
      ;Reflectivity_attenuated = [[Reflectivity_attenuated],[a_reflect]]    
      Reflectivity_attenuated = [Reflectivity_attenuated,transpose(a_reflect)]    
    endelse   
  ;*** Z *** Radar reflectivity 
  endif else if indicator eq 'Z' then begin
    num=strlen(line)
    num=fix((num-3.)/7.)
    if (num gt 0) then reflect=replicate(!values.f_nan,num) else reflect=!values.f_nan
    if (num gt 0) then reads, line, reflect, format='(3x,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7)'
    if (num le 30) then reflect=[reflect,nans[num:30]]
    if (num eq 0) then reflect=nans
    if (first eq 0) then begin
      ;Reflectivity = reflect
      Reflectivity = transpose(reflect)
    endif else begin
      ;Reflectivity = [[Reflectivity],[reflect]]  
      Reflectivity = [Reflectivity,transpose(reflect)]        
    endelse    
  ;*** R ***  rain rate             
  endif else if indicator eq 'R' then begin
    num=strlen(line)
    num=fix((num-3.)/7.)
    if (num gt 0) then rate=replicate(!values.f_nan,num) else rate=!values.f_nan
    if (num gt 0) then reads, line, rate, format='(3x,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7)'
    if (num le 30) then rate=[rate,nans[num:30]]
    if (num eq 0) then rate=nans
    if (first eq 0) then begin
      ;Rain_rate = rate
      Rain_rate=transpose(rate)
    endif else begin
      ;Rain_rate = [[Rain_rate],[rate]]  
      Rain_rate=[Rain_rate,transpose(rate)]     
    endelse    
  ;*** L *** liquid water content
  endif else if indicator eq 'L' then begin
    num=strlen(line)
    num=fix((num-3.)/7.)
    if (num gt 0) then lwc=replicate(!values.f_nan,num) else lwc=!values.f_nan
    if (num gt 0) then reads, line, lwc, format='(3x,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7)'
    if (num le 30) then lwc=[lwc,nans[num:30]]
    if (num eq 0) then lwc=nans
    if (first eq 0) then begin
      ;Liquid_water_content = lwc
      Liquid_water_content=transpose(lwc)
    endif else begin
      ;Liquid_water_content = [[Liquid_water_content],[lwc]]      
      Liquid_water_content = [Liquid_water_content,transpose(lwc)]   
    endelse  
  ;*** W *** Fall Velocity        
  endif else if indicator eq 'W' then begin   
    num=strlen(line)
    num=fix((num-3.)/7.)
    if (num gt 0) then velo=replicate(!values.f_nan,num) else velo=!values.f_nan
    if (num gt 0) then reads, line, velo, format='(3x,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7,F7)'
    if (num le 30) then velo=[velo,nans[num:30]]
    if (num eq 0) then velo=nans
    if (first eq 0) then begin
      ;Fall_velocity = velo
      Fall_velocity = transpose(velo)
    endif else begin 
      ;Fall_velocity = [[Fall_velocity],[velo]]   
      Fall_velocity = [Fall_velocity,transpose(velo)]     
    endelse     
    ;readf, lun,scrap   
    if (first ne 0) then Frequency = [[[Frequency]],[[dummy]]]
    if (first ne 0) then Diameter=[[[Diameter]],[[dummyd]]]
    if (first ne 0) then Number=[[[Number]],[[dummyn]]]
    first=first+1 
    dummy=replicate(!values.F_NAN,[64,31])
    dummyd=replicate(!values.F_NAN,[64,31])
    dummyn=replicate(!values.F_NAN,[64,31])         
  endif    
endwhile
free_lun,lun

julian_day=julday(month,day,year,hour,minute,seconds)

;  I don't need to rotate anymore.  The transposes fix this
;rot_num=6
;reflectivity_attenuated=rotate(reflectivity_attenuated,rot_num)
;reflectivity=rotate(reflectivity,rot_num)
;height=rotate(height,rot_num)
;rain_rate=rotate(rain_rate,rot_num)
;liquid_water_content=rotate(liquid_water_content,rot_num)
;fall_velocity=rotate(fall_velocity,rot_num)

end