;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

PRO read_sm, filename, sm_par, khz, col_dist_um, row_dist_um, $
             pixel_channel_names,i_main,i_clock,$
             header_only=header_only, help=help, print=print,$
             display=display,xinch=xinch
  
  ;; First version, September 2, 1998 CJJ
  ;; Deal correctly with the pixel_channels, and fix the
  ;;   clock.  CJJ, 30 Nov. 1999
  ;; Corrected bug in display calculation, which lead to always 
  ;;   only looking at channel ADC0 for ADC data, MF, 12 Oct. 2000
  ;; Adapted code from sm_image.cpp to properly calculate the khz
  ;;   array for proportional counter data, MF, 09-Feb-2001
  ;; Made sm_fixclock get called in all cases.  CJ, 22-May-2001
  
  
@bsif_common
  
  if (keyword_set(help) or (n_elements(filename) eq 0)) then begin
      print,'read_sm,filename,sm_par,khz,col_dist_um,row_dist_um'
      print,'  Reads SM data from a netCDF file to khz, "image_data" in BSIF_COMMON'
      print,'  Options: /header_only'
      return
  endif
  
  IF (NOT keyword_set(xinch)) THEN xinch = 3.
  
  ;; For dealing with older files
  got_dispinfo = 1
  
  ;; Open input file
  cdfid = ncdf_open(filename, /NOWRITE)
  
  ;; Go into inquire mode
  ncdf_inq = ncdf_inquire( cdfid )
  
  ;; These dimensions were written in this order in write_sm.pro
  ncdf_diminq, cdfid, 0, dummy, n_cols
  ncdf_diminq, cdfid, 1, dummy, n_rows
  ncdf_diminq, cdfid, 2, dummy, n_data
  
  sm_ipar,sm_par
  sm_par.n_cols = n_cols
  sm_par.n_rows = n_rows
  sm_par.n_data = n_data
  
  ;; Now we can get the variables
  image_data_id = ncdf_varid(cdfid, 'image_data')
  col_dist_um_id = ncdf_varid(cdfid, 'col_dist_um')
  row_dist_um_id = ncdf_varid(cdfid, 'row_dist_um')
  
  ;; Now we can read the global header info
  temp_long = long(0)
  ncdf_attget, cdfid, /GLOBAL, "filever", temp_long
  sm_par.filever = temp_long
  ncdf_attget, cdfid, /GLOBAL, "progver", temp_long
  sm_par.progver = temp_long
  ncdf_attget, cdfid, /GLOBAL, "col_device", temp_long
  sm_par.col_device = temp_long
  ncdf_attget, cdfid, /GLOBAL, "row_device", temp_long
  sm_par.row_device = temp_long
  ncdf_attget, cdfid, /GLOBAL, "strlen", temp_long
  sm_par.strlen = temp_long
  ncdf_attget, cdfid, /GLOBAL, "adc_ports", temp_long
  sm_par.adc_ports = temp_long
  ncdf_attget, cdfid, /GLOBAL, "scaler_ports", temp_long
  sm_par.scaler_ports = temp_long
  ncdf_attget, cdfid, /GLOBAL, "dac_ports", temp_long
  sm_par.dac_ports = temp_long
  ncdf_attget, cdfid, /GLOBAL, "maxpointlist", temp_long
  sm_par.maxpointlist = temp_long
  ;; Note that we skip n_cols, n_rows, n_data
  ncdf_attget, cdfid, /GLOBAL, "prescan_pixels", temp_long
  sm_par.prescan_pixels = temp_long
  ncdf_attget, cdfid, /GLOBAL, "shutopen", temp_long
  sm_par.shutopen = temp_long
  
  ;; For dealing with older files
  attstruct = ncdf_attinq(cdfid,/global,'bidirectional')
  IF (attstruct.datatype NE 'UNKNOWN') THEN BEGIN
      ncdf_attget, cdfid, /GLOBAL, "bidirectional", temp_long
      sm_par.bidirectional = temp_long
  ENDIF
  
  ncdf_attget, cdfid, /GLOBAL, "clock_data_channel", temp_long
  sm_par.clock_data_channel = temp_long
  ncdf_attget, cdfid, /GLOBAL, "pixel_channels", temp_long
  sm_par.pixel_channels = temp_long
  ncdf_attget, cdfid, /GLOBAL, "n_oneread_channels", temp_long
  sm_par.n_oneread_channels = temp_long
  ncdf_attget, cdfid, /GLOBAL, "oneread_channels", temp_long
  sm_par.oneread_channels = temp_long
  ncdf_attget, cdfid, /GLOBAL, "n_pointlist", temp_long
  sm_par.n_pointlist = temp_long
  
  ;; For dealing with older files
  attstruct = ncdf_attinq(cdfid,/global,'disp_num_channels')
  IF (attstruct.datatype NE 'UNKNOWN') THEN BEGIN
      ncdf_attget, cdfid, /GLOBAL, "disp_num_channels", temp_long
      sm_par.disp_num_channels = temp_long
  ENDIF ELSE BEGIN
      got_dispinfo = 0
  ENDELSE
  attstruct = ncdf_attinq(cdfid,/global,'disp_denom_channels')
  IF (attstruct.datatype NE 'UNKNOWN') THEN BEGIN
      ncdf_attget, cdfid, /GLOBAL, "disp_denom_channels", temp_long
      sm_par.disp_denom_channels = temp_long
  ENDIF ELSE BEGIN
      got_dispinfo = 0
  ENDELSE
  attstruct = ncdf_attinq(cdfid,/global,'disp_raw')
  IF (attstruct.datatype NE 'UNKNOWN') THEN BEGIN
      ncdf_attget, cdfid, /GLOBAL, "disp_raw", temp_long
      sm_par.disp_raw = temp_long
  ENDIF ELSE BEGIN
      got_dispinfo = 0
  ENDELSE
  
  temp_float = float(0)
  ncdf_attget, cdfid, /GLOBAL, "col_pixel_um", temp_float
  sm_par.col_pixel_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "row_pixel_um", temp_float
  sm_par.row_pixel_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "x_center_um", temp_float
  sm_par.x_center_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "y_center_um", temp_float
  sm_par.y_center_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "z_center_um", temp_float
  sm_par.z_center_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "x_pzt_um", temp_float
  sm_par.x_pzt_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "y_pzt_um", temp_float
  sm_par.y_pzt_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "z_pzt_um", temp_float
  sm_par.z_pzt_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "x_stg_um", temp_float
  sm_par.x_stg_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "y_stg_um", temp_float
  sm_par.y_stg_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "z_stg_um", temp_float
  sm_par.z_stg_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "x_det_um", temp_float
  sm_par.x_det_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "y_det_um", temp_float
  sm_par.y_det_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "z_det_um", temp_float
  sm_par.z_det_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "z_osa_um", temp_float
  sm_par.z_osa_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "col_start_um", temp_float
  sm_par.col_start_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "col_stop_um", temp_float
  sm_par.col_stop_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "row_start_um", temp_float
  sm_par.row_start_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "row_stop_um", temp_float
  sm_par.row_stop_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "start_ev", temp_float
  sm_par.start_ev = temp_float
  ncdf_attget, cdfid, /GLOBAL, "stop_ev", temp_float
  sm_par.stop_ev = temp_float
  ncdf_attget, cdfid, /GLOBAL, "dwell_msec", temp_float
  sm_par.dwell_msec = temp_float
  ncdf_attget, cdfid, /GLOBAL, "ring_gev", temp_float
  sm_par.ring_gev = temp_float
  ncdf_attget, cdfid, /GLOBAL, "ring_ma", temp_float
  sm_par.ring_ma = temp_float
  ncdf_attget, cdfid, /GLOBAL, "undulator_gap_mm", temp_float
  sm_par.undulator_gap_mm = temp_float
  ncdf_attget, cdfid, /GLOBAL, "undulator_period_mm", temp_float
  sm_par.undulator_period_mm = temp_float
  ncdf_attget, cdfid, /GLOBAL, "ens_slit_um", temp_float
  sm_par.ens_slit_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "exs_slit_um", temp_float
  sm_par.exs_slit_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "exsy_slit_um", temp_float
  sm_par.exsy_slit_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "sgm_zero_order_steps", temp_float
  sm_par.sgm_zero_order_steps = temp_float
  ncdf_attget, cdfid, /GLOBAL, "sgm_steps_per_angstrom", temp_float
  sm_par.sgm_steps_per_angstrom = temp_float
  ncdf_attget, cdfid, /GLOBAL, "zp_d_um", temp_float
  sm_par.zp_d_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "zp_delta_nm", temp_float
  sm_par.zp_delta_nm = temp_float
  ncdf_attget, cdfid, /GLOBAL, "zp_stop_um", temp_float
  sm_par.zp_stop_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "zp_osa_um", temp_float
  sm_par.zp_osa_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "zp_osa_z_um", temp_float
  sm_par.zp_osa_z_um = temp_float
  ncdf_attget, cdfid, /GLOBAL, "tilt_degrees", temp_float
  sm_par.tilt_degrees = temp_float
  ncdf_attget, cdfid, /GLOBAL, "xpzt_nm_per_volt", temp_float
  sm_par.xpzt_nm_per_volt = temp_float
  ncdf_attget, cdfid, /GLOBAL, "ypzt_nm_per_volt", temp_float
  sm_par.ypzt_nm_per_volt = temp_float
  ncdf_attget, cdfid, /GLOBAL, "zpzt_nm_per_volt", temp_float
  sm_par.zpzt_nm_per_volt = temp_float
  ncdf_attget, cdfid, /GLOBAL, "xum_per_zum", temp_float
  sm_par.xum_per_zum = temp_float
  ncdf_attget, cdfid, /GLOBAL, "yum_per_zum", temp_float
  sm_par.yum_per_zum = temp_float
  ncdf_attget, cdfid, /GLOBAL, "clock_hertz", temp_float
  sm_par.clock_hertz = temp_float
  ncdf_attget, cdfid, /GLOBAL, "det1_volt1", temp_float
  sm_par.det1_volt1 = temp_float
  ncdf_attget, cdfid, /GLOBAL, "det1_volt2", temp_float
  sm_par.det1_volt2 = temp_float
  ncdf_attget, cdfid, /GLOBAL, "det2_volt1", temp_float
  sm_par.det2_volt1 = temp_float
  ncdf_attget, cdfid, /GLOBAL, "det2_volt2", temp_float
  sm_par.det2_volt2 = temp_float
  ncdf_attget, cdfid, /GLOBAL, "temp1_celsius", temp_float
  sm_par.temp1_celsius = temp_float
  ncdf_attget, cdfid, /GLOBAL, "temp2_celsius", temp_float
  sm_par.temp2_celsius = temp_float
  ncdf_attget, cdfid, /GLOBAL, "disp_min", temp_float
  sm_par.disp_min = temp_float
  ncdf_attget, cdfid, /GLOBAL, "disp_max", temp_float
  sm_par.disp_max = temp_float
  ncdf_attget, cdfid, /GLOBAL, "disp_gamma", temp_float
  sm_par.disp_gamma = temp_float
  temp_floats = fltarr(sm_par.adc_ports)
  ncdf_attget, cdfid, /GLOBAL, "adc_voltsperbit", temp_floats
  sm_par.adc_voltsperbit = temp_floats
  ncdf_attget, cdfid, /GLOBAL, "adc_voltsoffset", temp_floats
  sm_par.adc_voltsoffset = temp_floats
  temp_floats = fltarr(sm_par.dac_ports)
  ncdf_attget, cdfid, /GLOBAL, "dac_voltsperbit", temp_floats
  sm_par.dac_voltsperbit = temp_floats
  ncdf_attget, cdfid, /GLOBAL, "dac_voltsoffset", temp_floats
  sm_par.dac_voltsoffset = temp_floats
  
  temp_floats = fltarr(sm_par.adc_ports+sm_par.scaler_ports+1)
  ;; For dealing with older files
  attstruct = ncdf_attinq(cdfid,/global,'disp_num_factors')
  IF (attstruct.datatype NE 'UNKNOWN') THEN BEGIN
      ncdf_attget, cdfid, /GLOBAL, "disp_num_factors", temp_floats
      sm_par.disp_num_factors = temp_floats
  ENDIF ELSE BEGIN
      got_dispinfo = 0
  ENDELSE
  attstruct = ncdf_attinq(cdfid,/global,'disp_denom_factors')
  IF (attstruct.datatype NE 'UNKNOWN') THEN BEGIN
      ncdf_attget, cdfid, /GLOBAL, "disp_denom_factors", temp_floats
      sm_par.disp_denom_factors = temp_floats
  ENDIF ELSE BEGIN
      got_dispinfo = 0
  ENDELSE
  
  temp_string = bytarr(sm_par.strlen)
  ncdf_attget, cdfid, /GLOBAL, "col_title", temp_string
  sm_par.col_title = temp_string
  ncdf_attget, cdfid, /GLOBAL, "row_title", temp_string
  sm_par.row_title = temp_string
  ncdf_attget, cdfid, /GLOBAL, "beamline", temp_string
  sm_par.beamline = temp_string
  ncdf_attget, cdfid, /GLOBAL, "systime", temp_string
  sm_par.systime = temp_string
  ncdf_attget, cdfid, /GLOBAL, "scientist", temp_string
  sm_par.scientist = temp_string
  ncdf_attget, cdfid, /GLOBAL, "sample", temp_string
  sm_par.sample = temp_string
  ncdf_attget, cdfid, /GLOBAL, "comments1", temp_string
  sm_par.comments1 = temp_string
  ncdf_attget, cdfid, /GLOBAL, "comments2", temp_string
  sm_par.comments2 = temp_string
  ncdf_attget, cdfid, /GLOBAL, "original_filename", temp_string
  sm_par.original_filename = temp_string
  ncdf_attget, cdfid, /GLOBAL, "zp_name", temp_string
  sm_par.zp_name = temp_string
  ncdf_attget, cdfid, /GLOBAL, "temp1_name", temp_string
  sm_par.temp1_name = temp_string
  ncdf_attget, cdfid, /GLOBAL, "temp2_name", temp_string
  sm_par.temp2_name = temp_string
  temp_bytes = bytarr(sm_par.scaler_ports)
  ncdf_attget, cdfid, /GLOBAL, "scaler_divideby", temp_bytes
  sm_par.scaler_divideby = temp_bytes
  temp_string = bytarr(sm_par.strlen*sm_par.adc_ports)
  ncdf_attget, cdfid, /GLOBAL, "adc_channel_names", temp_string
  sm_par.adc_channel_names = temp_string
  temp_string = bytarr(sm_par.strlen*sm_par.scaler_ports)
  ncdf_attget, cdfid, /GLOBAL, "scaler_channel_names", temp_string
  sm_par.scaler_channel_names = temp_string
  temp_string = bytarr(sm_par.maxpointlist*sm_par.strlen)
  ncdf_attget, cdfid, /GLOBAL, "point_namestr", temp_string
  sm_par.point_namestr = temp_string
  
  ;; For dealing with older files
  attstruct = ncdf_attinq(cdfid,/global,'disp_name')
  IF (attstruct.datatype NE 'UNKNOWN') THEN BEGIN
      ncdf_attget, cdfid, /GLOBAL, "disp_name", temp_string
      sm_par.disp_name = temp_string
  ENDIF ELSE BEGIN
      got_dispinfo = 0
  ENDELSE
  
  temp_ints = intarr(sm_par.adc_ports+sm_par.scaler_ports)
  ncdf_attget, cdfid, /GLOBAL, "oneread_data", temp_ints
  sm_par.oneread_data = temp_ints
  
  ;; We always get the following in case they are needed for
  ;; x_dist, y_dist
  col_dist_um = fltarr(sm_par.n_cols)
  ncdf_varget, cdfid, col_dist_um_id, col_dist_um
  row_dist_um = fltarr(sm_par.n_rows)
  ncdf_varget, cdfid, row_dist_um_id, row_dist_um
  
  x_dist = col_dist_um
  y_dist = row_dist_um
  
  ;; Now set BSIF variables
  version = 1                   ; no longer used
  first_data_record = 0         ; no longer used
  ;;n_rows = sm_par.n_rows
  ;;n_cols = sm_par.n_cols
  ;;n_data = sm_par.n_data
  x_normal = 1
  y_normal = 1
  rotated = 0
  x_start = sm_par.col_start_um
  x_stop = sm_par.col_stop_um
  y_start = sm_par.row_start_um
  y_stop = sm_par.row_stop_um
  data_type = 1 ;; not sure of this anymore, but not used by anything
  compression_type = 0
  data_min = sm_par.disp_min
  data_max = sm_par.disp_max
  image_title = string(sm_par.original_filename)
  x_title = string(sm_par.col_title)
  y_title = string(sm_par.row_title)
  data_title = ' '
  
  pixel_channel_names = strarr(n_data)
  
  IF (NOT keyword_set(header_only)) THEN BEGIN
      image_data = intarr(n_cols,n_rows,n_data)
      ncdf_varget, cdfid, image_data_id, image_data
      image_data = reform(image_data, n_cols, n_rows, n_data)
      IF (got_dispinfo EQ 0) THEN BEGIN
          ;; If got_dispinfo is zero, then this is an older file and we
          ;; have to make a guess as to which channel is the clock channel.
          IF (n_data EQ 1) THEN BEGIN
              i_clock = -1
              i_data = 0
              i_main = 0
          ENDIF ELSE IF (n_data EQ 3) THEN BEGIN
              ;; Take a guess
              i_data = 1
              i_main = 1
              i_clock = 2
          ENDIF ELSE BEGIN
              ncdf_close, cdfid
              message,'Not sure what to do with this older file'
              return
          ENDELSE
          IF (i_clock GE 0) THEN BEGIN
              sm_fixclock,image_data,i_clock
              khz = float(sm_par.clock_hertz)*$
                float(image_data(*,*,i_main))/$
                float(image_data(*,*,i_clock))
          ENDIF ELSE BEGIN
              khz = float(sm_par.clock_hertz)*$
                float(image_data(*,*,i_main))
          ENDELSE
      ENDIF ELSE BEGIN
          ;; This is a newer file, which has all the proper display
          ;; parameters.
          ;; Add the constant offset to the denominator image
          denom_image = fltarr(n_cols,n_rows)+sm_par.disp_denom_factors(0)
          ;; offset for numerator image
          num_image = fltarr(n_cols,n_rows)+sm_par.disp_num_factors(0)
          
          i_clock = -1
          i_data = 0
          clock_present = 0
          ;; get i_clock: index in image_data that corresponds to the
          ;; clock signal
          FOR i_bit=0,(sm_par.adc_ports+sm_par.scaler_ports-1) DO BEGIN
              bitmask = long(2)^long(i_bit)
              IF ((bitmask AND long(sm_par.pixel_channels)) NE 0) THEN BEGIN
                  IF ((bitmask AND sm_par.clock_data_channel) NE 0) $
                    THEN BEGIN
                      i_clock = i_data
                      clock_present =  1
                  ENDIF
                  i_data =  i_data+1
              ENDIF
          ENDFOR
          IF ( clock_present NE 1 ) THEN BEGIN
              message, 'read_sm: No clock has been found. Setting clock to 0'
              i_clock = 0
          ENDIF ELSE BEGIN
              sm_fixclock,image_data,i_clock
          ENDELSE
          
          i_data = 0
          FOR i_bit=0,(sm_par.adc_ports+sm_par.scaler_ports-1) DO BEGIN
              bitmask = long(2)^long(i_bit)
              channel_was_recorded = 0
              IF ((bitmask AND long(sm_par.pixel_channels)) NE 0) THEN BEGIN
                  channel_was_recorded = 1
                  IF ((bitmask AND long(sm_par.disp_num_channels)) NE 0) $
                    THEN BEGIN
                      IF (i_bit LT 8) THEN BEGIN 
                          ;; this is ADC data
                          IF (sm_par.disp_raw EQ 0) THEN BEGIN 
                              ;; ADC data to be scaled
                              num_image(0:(n_cols-1),0:(n_rows-1)) = $
                                num_image(0:(n_cols-1),0:(n_rows-1))+$
                                sm_par.disp_num_factors(i_bit+1)*$
                                (sm_par.adc_voltsoffset(i_bit)+$
                                 sm_par.adc_voltsperbit(i_bit)*$
                                 float(image_data(0:(n_cols-1),$
                                                  0:(n_rows-1),i_data)))
                          ENDIF ELSE BEGIN ;; ADC data raw!
                              num_image(0:(n_cols-1),0:(n_rows-1)) = $
                                num_image(0:(n_cols-1),0:(n_rows-1))+$
                                sm_par.disp_num_factors(i_bit+1)*$
                                float(image_data(0:(n_cols-1),$
                                                 0:(n_rows-1),i_data))
                          ENDELSE
                      ENDIF ELSE BEGIN ;; scaler data
                          IF (sm_par.scaler_divideby(i_bit-8) EQ 0) $
                            THEN BEGIN
                              undo_divideby = 1.
                          ENDIF ELSE BEGIN
                              undo_divideby = 8.
                          ENDELSE
                          ;; we do not want to divide the clock by itself!
                          IF ( (sm_par.disp_raw EQ 0) AND $
                               (i_bit LT (sm_par.adc_ports+sm_par.scaler_ports-2)) ) THEN BEGIN 
                              ;; scaler data to be scaled
                              num_image(0:(n_cols-1),0:(n_rows-1)) = $
                                num_image(0:(n_cols-1),0:(n_rows-1))+$
                                sm_par.disp_num_factors(i_bit+1)*$
                                (undo_divideby*$
                                 float(uint(image_data(0:(n_cols-1),$
                                                       0:(n_rows-1),$
                                                       i_data))))/$
                                float(uint(image_data(0:(n_cols-1), $
                                                      0:(n_rows-1), $
                                                      i_clock )))*$
                                sm_par.clock_hertz/1000.
                          ENDIF ELSE BEGIN ;; scaler data raw!
                              num_image(0:(n_cols-1),0:(n_rows-1)) = $
                                num_image(0:(n_cols-1),0:(n_rows-1))+$
                                sm_par.disp_num_factors(i_bit+1)*$
                                (undo_divideby*$
                                 float(uint(image_data(0:(n_cols-1),$
                                                       0:(n_rows-1),$
                                                       i_data))))
                          ENDELSE
                      ENDELSE
                  ENDIF
                  IF ((bitmask AND long(sm_par.disp_denom_channels)) NE 0) $
                    THEN BEGIN
                      IF (i_bit LT 8) THEN BEGIN 
                          ;; this is ADC data
                          IF (sm_par.disp_raw EQ 0) THEN BEGIN 
                              ;; ADC data to be scaled
                              denom_image(0:(n_cols-1),0:(n_rows-1)) = $
                                denom_image(0:(n_cols-1),0:(n_rows-1))+$
                                sm_par.disp_denom_factors(i_bit+1)*$
                                (sm_par.adc_voltsoffset(i_bit)+$
                                 sm_par.adc_voltsperbit(i_bit)*$
                                 image_data(0:(n_cols-1),0:(n_rows-1),$
                                            i_data))
                          ENDIF ELSE BEGIN ;; ADC data raw!
                              denom_image(0:(n_cols-1),0:(n_rows-1)) = $
                                denom_image(0:(n_cols-1),0:(n_rows-1))+$
                                sm_par.disp_denom_factors(i_bit+1)*$
                                float(image_data(0:(n_cols-1),0:(n_rows-1),i_data))
                          ENDELSE
                      ENDIF ELSE BEGIN ;; scaler data
                          IF (sm_par.scaler_divideby(i_bit-8) EQ 0) $
                            THEN BEGIN
                              undo_divideby = 1.
                          ENDIF ELSE BEGIN
                              undo_divideby = 8.
                          ENDELSE
                          ;; we do not want to divide the clock by itself!
                          IF ( (sm_par.disp_raw EQ 0) AND $
                               (i_bit LT (sm_par.adc_ports+sm_par.scaler_ports-2)) ) THEN BEGIN 
                              ;; scaler data to be scaled
                              denom_image(0:(n_cols-1),0:(n_rows-1)) = $
                                denom_image(0:(n_cols-1),0:(n_rows-1))+$
                                sm_par.disp_denom_factors(i_bit+1)*$
                                (undo_divideby*$
                                 float(UINT(image_data(0:(n_cols-1),$
                                                  0:(n_rows-1),i_data))))/$
                                float(UINT(image_data(0:(n_cols-1), $
                                                 0:(n_rows-1), i_clock )))*$
                                sm_par.clock_hertz/1000.
                          ENDIF ELSE BEGIN ;; scaler data raw!
                              denom_image(0:(n_cols-1),0:(n_rows-1)) = $
                                denom_image(0:(n_cols-1),0:(n_rows-1))+$
                                sm_par.disp_denom_factors(i_bit+1)*$
                                (undo_divideby*$
                                 (float(UINT(image_data(0:(n_cols-1),0:(n_rows-1),$
                                             i_data)))))
                          ENDELSE
                      ENDELSE
                  ENDIF
              ENDIF
              IF (channel_was_recorded NE 0) THEN BEGIN
                  i_data = i_data+1
              ENDIF
          ENDFOR
          denom_zeroes = where(denom_image EQ 0.,n_denom_zeroes)
          IF (n_denom_zeroes NE 0) THEN BEGIN
              khz = fltarr(n_cols,n_rows)
              denom_nonzeroes = where(denom_image NE 0.,n_denom_nonzeroes)
              IF (n_denom_nonzeroes NE 0) THEN BEGIN
                  khz(denom_nonzeroes) = $
                    num_image(denom_nonzeroes)/denom_image(denom_nonzeroes)
              ENDIF
          ENDIF ELSE BEGIN
              khz = num_image/denom_image
          ENDELSE
      ENDELSE
  ENDIF
  
  IF (NOT keyword_set(header_only)) THEN BEGIN
      IF keyword_set(display) THEN BEGIN
          IF ((sm_devicename(sm_par.col_device) EQ 'EV') AND $
              (n_rows EQ 1)) THEN BEGIN
              plot,(1.239852/col_dist_um),khz,xtitle='eV',ytitle='kHz',$
                title=filename
          ENDIF ELSE BEGIN
              img_scl,khz,/scalebar
          ENDELSE
      ENDIF
  ENDIF
  
  ncdf_close, cdfid
  
  IF keyword_set(print) THEN BEGIN
      IF ((sm_devicename(sm_par.col_device) EQ 'EV') AND $
          (n_rows EQ 1)) THEN BEGIN
          old_plot = !d.name
          old_font = !p.font
          set_plot,'ps'
          psfilename = randomfile()
          device,file=psfilename,/inch,xsize=5.,ysize=3.5
          !p.font = 0
          plot,(1.239852/col_dist_um),khz,xtitle='eV',ytitle='kHz',$
            title=filename
          device,/close
          spawn,'lpr -h -r '+psfilename
          !p.font = old_font
          set_plot,old_plot
      ENDIF ELSE BEGIN
;+$
;          '!C!CXSTG='+strtrim(string(sm_par.x_stg_um,format='(f10.1)'),2)+$
;          ', YSTG='+strtrim(string(sm_par.y_stg_um,format='(f10.1)'),2)+$
;          ', ZSTG='+strtrim(string(sm_par.z_stg_um,format='(f10.1)'),2)+$
;          ',!CXDET='+strtrim(string(sm_par.x_det_um,format='(f10.1)'),2)+$
;          ', YDET='+strtrim(string(sm_par.y_det_um,format='(f10.1)'),2)+$
;          ', ZDET='+strtrim(string(sm_par.z_det_um,format='(f10.1)'),2)+$
;          ',!CColumn: '+strtrim(string(sm_par.col_title),2)+$
;          ', '+strtrim(string(sm_par.n_cols),2)+$
;          ' pixels, '+strtrim(string(sm_par.col_pixel_um,format='(f10.3)'),2)+$
;          ' um steps,!C  from '+strtrim(string(x_dist(0),format='(f10.1)'),2)+$
;          ' to '+strtrim(string(x_dist(sm_par.n_cols-1),format='(f10.1)'),2)+$
;          ',!CRow:  '+strtrim(string(sm_par.row_title),2)+$
;          ', '+strtrim(string(sm_par.n_rows),2)+$
;          ' pixels, '+strtrim(string(sm_par.row_pixel_um,format='(f10.3)'),2)+$
;          ' um steps,!C  from '+strtrim(string(y_dist(0),format='(f10.1)'),2)+$
;          ' to '+strtrim(string(y_dist(sm_par.n_rows-1),format='(f10.1)'),2)+$
;          ',!CDwell='+strtrim(string(sm_par.msec_dwell,format='(f10.1)'),2)+$
;          ' msec'
          psfilename = randomfile()
          img_ps,khz,psfilename,title=filename,xinch=xinch
          spawn,'lpr -h -r '+psfilename
      ENDELSE
  ENDIF
  
  return
END




