; Copyright (c) 1998-2010 A.P. Hitchcock  All rights reserved
;+
;NAME:
;	READ_DM_DATACUBE
;
;LAST CHANGED: ----------------------------------- 22-Jul-10 (aph)
;
;PURPOSE:
;	This function reads in a binary datacube I(x,y,E) file
; exported from Digital Micrograph (Gatan)
; and writes out an axis format ncb file
; The image and spectral dimensions are either 
; in an associated header file 
; or are supplied by the user
; data is 4-bytes

;
;CATEGORY:
;	STAND ALONE: read in utility
;CALLED FROM AXIS: Read~stacks~DM datacube
;
;CALLING SEQUENCE:
;	Result = READ_DM_DATACUBE( [file=file, NOF=nof, $
;	    silent=silent, verbose=verbose, group=group])
;
;INPUTS: none. All input parameters are passed as keywords.
;
;KEYWORDS:
;	FILE 	filename to read
;	GROUP	group ID (aXisID)
; NOF   no filter (*.*)
;	SILENT  if set, do not print variables etc
;	VERBOSE print log information about processing
;	EXTRA 	passed on parameters
;
;OUTPUTS:
;	RESULT = AXIS stack (*.ncb,*.dat) written to disk
;
;COMMON BLOCKS: 
; axis_com
; stack_process_com
; bsif_com
; volume_data
;
;MODIFICATION HISTORY:
; (14-Feb-10 aph) first version
; (20-Feb-10 aph) change name to avoid overwriting *.dat files
; (23-feb-10 aph) remove real for now as it gave problems with other stacks
; (22-jul-10 aph) adapt to read *.txt if it exists

Function get_DM_variable, lun, string
; ---- reads value of a variable from an OPEN file (lun)
 line = ' '
  test = -1
  while test LT 0 do begin 
    readf,lun, line
    test = strpos(line,string) 
  endwhile  
  test = strpos(line,':')
  value = float(strmid(line, test+1))
;  print, string, ' : ', value
return, value
End  
  

Function read_dm_datacube, file=file, NOF=nof, $
     silent=silent, verbose=verbose, group=group
@axis_com
@stack_process_com
COMMON volume_data, image_stack
@bsif_com
on_error,2

if NOT keyword_set(group) then group = 0 else begin
	if n_elements(axis_ID) EQ 1 then group = axis_ID
 endelse

IF NOT keyword_set(file) then begin
	if keyword_set(nof) then fltr='*' else  fltr='*.dat'
	file = pickfile2(/READ, FILTER=fltr, /LPATH, DEFPATH=defpath)
ENDIF
s = ''
if file EQ '' THEN RETURN, s  ; bail-out if no filename
check = findfile(file)		 ; bail-out if non-existent file
if strlen(check(0)) EQ 0 then begin
	axis_log, 'Cannot find ' + file
	return, s
endif
t = ax_name(file)
file_path = t(0)
fileshort = t(1)

; ---------- define default parameters of the datacube
nx =24 & ny = 24 & nen = 1340
dwell = 1. &  x_step = 1. & y_step = 1. & e_step = 0.05
; ---- try to read from an associated *.txt header file
header_file = file_path+fileshort+'.txt'
; check if file exists - only proceed if it does
file_exist = findfile(header_file, count=count)
if count GT 0 then begin
  get_lun,lun
  openr, lun, header_file
  dwell = get_DM_variable(lun, 'Dwell')
  nx = get_DM_variable(lun,'DimensionX')
  ny = get_DM_variable(lun,'DimensionY')
  nen = get_DM_variable(lun,'DimensionE')
  x_step = get_DM_variable(lun,'ScaleX')
  y_step = get_DM_variable(lun,'ScaleY')
  E_step = get_DM_variable(lun,'Dispersion')
 close, lun & free_lun, lun
  
endif else begin
; ---- if the header file does not exist, ask the user ---  
  if keyword_set(group) then begin
  	nx = get_num(prompt='Number of columns (x)', val = nx, group = axis_ID)
  	ny = get_num(prompt='Number of rows (y)', val = ny, group = axis_ID)
    nen = get_num(prompt='Number of energy channels', val = nen, group = axis_ID)
    dwell = get_num(prompt='Dwell time (s)', val = dwell, group = axis_ID)
    e_step = get_num(prompt='Energy (eV) per channel', val = evc, group = axis_ID)
  endif else	begin
    nx = get_num(prompt='Number of columns (x)', val = nx)
    ny = get_num(prompt='Number of rows (y)', val = ny)
  	nen = get_num(prompt='Number of energy channels', val = nen)
  	dwell = get_num(prompt='Dwell time (s)', val = dwell)
    e_step = get_num(prompt='Energy (eV) per channel', val = evc)
  endelse
endelse

; ------ check that parameter read-in worked
ncols = nx
nrows = ny
npts = nen
if not keyword_set(silent) then begin
  print, 'Parameters for ', fileshort
  axis_log, ' #x pts ' + strtrim(string(fix(nx)),2)
  axis_log, ' X-step (um) ' + string(X_step, format='(F8.4)')
  axis_log, ' #y pts ' + strtrim(string(fix(ny)),2)
  axis_log, ' Y-step (um) ' + string(Y_step, format='(F8.4)')
  axis_log, ' #E pts ' + strtrim(string(fix(npts)),2)
  axis_log, ' E-step (eV) ' + string(E_step, format='(F6.3)')
  axis_log, ' dwell (s) ' + string(dwell, format='(F4.2)')
endif

dc = fltarr(nx,ny,nen)

on_ioerror, read_error
get_lun, iunit
openu, iunit, file
readu, iunit, dc
close, iunit & free_lun, iunit

; -------- Create aXis2000 stack
image_stack=dc              ;reform(dc, ncols, nrows, npts)
ev_max= npts*e_step & ev_min = 0
slope = (ev_max - ev_min)/npts
ev = ev_min + slope*findgen(npts)
x_start = 0.
x_stop = float(ncols-1)*X_step
y_start = 0.
y_stop = float(nrows-1)*Y_step
filename_ev_msec_list = strarr(npts)
for i = 0, npts-1 do filename_ev_msec_list(i) = $
         fileshort + ' ' + string(ev(i), format='(F8.2)') + ' ' + string(dwell, format='(F5.2)')

; --------- ask user what file name they want - default name is same as files without counter
fileout = file_path + fileshort + '_axis2000.ncb'
sname=pickfile2(/write, file= fileout, filter='*.ncb', title = 'name of binary stack file')
stack_wb, sname , /real		; use REAL keyword to avoid loss of z-resolution

axis_log, 'DM data cube read from ' + fileshort
return, s

read_error :
axis_log, 'read_dm_datacube: read error - perhaps data cube dimensions are wrong'
end

