; Copyright (c) 1998-2012 A.P. Hitchcock  All rights reserved
;+
;NAME:
;	READ_XRF_STACK
;
;LAST CHANGED: ----------------------------------- 01-Sep-12 (aph)
;
;PURPOSE:
;	This function reads in a single file containing X-ray fluorescence spectra
; at each pixel, written by STXM_control or by Oxford Inca (*.raw, *.rpl) files
;
;CATEGORY:
;	STAND ALONE: read in utility
;
;CALLING SEQUENCE:
;	Result = READ_XRF_STACK(file=file, NOF=nof, inca=inca, verbose=verbose, group=group, $
;            no_write=no_write, silent=silent, energy = energy, one_image=one_image, extra = extra)
;
;CALLED FROM AXIS:
;	Read->Stacks
;
;INPUTS: none. All input parameters are passed as keywords.
;
;KEYWORDS:
;	FILE 	filename to read
;	NOF		no filter
;	GROUP	group ID (aXis_ID)
;	INCA	if set, read in inca format stacks (default is STXM_control)
;	ENERGY	energy of incident excitation (eleectron or X-ray)
;	one_image = index of the image from a stack that is requested to read to define E
;	VERBOSE print log information about processing
;	NO_WRITE  if set, do not write *.ncb (used to load image_stack when reading XAS-XRF 4d stack)
;	SILENT	if set, do not display parameters
;	EXTRA 	passed on parameters
;
;OUTPUTS:
;	RESULT = AXIS (*,dat, *,ncb) stack
;
;COMMON BLOCKS: none
;
;MODIFICATION HISTORY:
; (24-Jan-10 aph) first version
; (16-Feb-10 aph) add inca option
; (07-mar-10 aph) adapt to XRF map files written by STXM_Control from XGLabs data
; (17-apr-10 aph) desensitize to lower/upper case of structure label from read_ascii
;                 make one_image read work for dimensions
; (20-apr-10 aph) add energy call parameter to pass photon energy
; (25-may-10 aph) fix no-rerad when size of xrf map does not match EI image
; (01-Aug-10 aph) add upper limit selection
; (02-Oct-10 aph) remove 'replace -ve' code;  added correction for missing last point
;                 change read_Ascii method to give correct size of xrf_stack
; (10-dec-11 aph) add ax_wait hourglass for Oxford (Inca) stack read-in
; (01-Sep-12 aph) correct error on 'get_num'
;-

Function read_xrf_stack, file=file, NOF=nof, inca=inca, verbose=verbose, group=group, $
            no_write=no_write, silent=silent, energy = energy, one_image=one_image, extra = extra
@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='*.xrf
	if keyword_set(inca) then fltr = '*.raw'
	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)

widget_control, /hourglass

; -------------- READ in STXM_Control xrf stack (XRF at each pixel) --------------------
if NOT keyword_set(inca) then begin
	axis_log, 'XRF STXM stacks are large ....'
	axis_log, '  . . . be patient'
	ax_wait
;	a = read_ascii(file, data_start=0,num_records=npts, delimiter='	')
;	a = read_ascii(file, data_start=0,TEMPLATE=ASCII_TEMPLATE(file))

	a = read_ascii(file, data_start=0) ; ,num_records=npts, delimiter='	')
	xrf_stack = fix(a.(0))
	IF NOT keyword_set(silent) then print, ' size after reading with read_ascii'
	IF NOT keyword_set(silent) then help, xrf_stack


; --------- kluge for when STXM_control throws -ves into file ---------------
; ----- remove -ve values
;	index = where(xrf_stack LT 0, count)
;	if count NE 0 then begin
;		t=size(xrf_stack)		; first index is size of each spectrum; second is number of images
;		print, 'before changing  -ve ', t
;		xrf_stack(index) = 0
;		print, count, ' negative values replaced with zero'
;		t=size(xrf_stack)		; first index is size of each spectrum; second is number of images
;		print,  'after removing -ve ',t
;	endif

	; ---------- check that the second dimension agrees with the image size
	t = ax_name(file)
	test = t(1)
	if strmid(test,strlen(test)-4,1) EQ '_' then $
		t(1)=strmid(t(1),0, strlen(t(1))-4)			; remove last 3 characters from the name
	file_image=t(0)+t(1) + '.hdr'
	IF NOT keyword_set(silent) then print, 'Reading STXM_control XRF stack ', file_image
	IF NOT keyword_set(silent) then print, ' reading image dimensions from ', file_image
	if not keyword_set(one_image) then one_image=1
	tmp = read_sdf(file_image, channel=1, one_image=one_image, /no_display, /silent)

	energy = tmp.e
	nx = n_elements(tmp.x)
	ny = n_elements(tmp.y)
	x_start = tmp.x(0)  & x_stop = tmp.x(nx-1)  & x_step = (x_stop - x_start)/nx
	y_start = tmp.y(0)  & y_stop = tmp.y(ny-1)  & y_step = (y_stop - y_start)/ny

	t=size(xrf_stack)		; first index is size of each spectrum; second is number of images
	npixels=t(2)
	npts = t(1) & n_data = npts

	IF NOT keyword_set(silent) then  $
	     print, 'nx, ny, nx*ny, #pixels ', nx, ny, nx*ny, npixels

	IF nx*ny NE npixels then begin
		n_lines_diff = fix(npixels - nx*ny)
; ---------- sometimes last point is missed (e.g. 101001032)
		IF  n_lines_diff EQ '-1' then begin
			last = XRF_stack(*,npixels-1)
			txrf = fltarr(t(1), t(2)+1)
			txrf(*,0:t(2)-1)=XRF_stack
			txrf(*,t(2))=last(*,0)
			XRF_stack = txrf
			IF NOT keyword_set(silent) then help, xrf_stack
			axis_log, 'XRF map data set extended by 1 pixel'
		ENDIF ELSE BEGIN
; ---------- sometimes there are more XRF spectral than pixels
; kluge to force correct number of spectra
; error was fixed but this code left for compatibility with early data
		     axis_log, 'WARNING: data size from ' + fileshort + 'disagrees with image size'
		     n_lines_diff = strtrim(string(fix(npixels - nx*ny)))
		     axis_log, 'XRF map has ' + n_lines_diff + ' more lines'
			if keyword_set(group) then begin
				nx = get_num(prompt='Number of columns (x)', val = nx, group = axis_ID)
			endif else	nx = get_num(prompt='Number of columns (x)', val = nx)
			ny = fix(npixels/nx)-1
			if keyword_set(group) then begin
				ny = get_num(prompt='Number of rows (y)', val = ny, group = axis_ID)
			endif else	ny = get_num(prompt='Number of rows (y)', val = ny)

			xrf_stack = xrf_stack(*,0:nx*ny-1)
			IF NOT keyword_set(silent)then print, size(xrf_stack)
		ENDELSE
	endif

	n_cols=nx & n_rows=ny
;	help, xrf_stack
	xrf_stack = transpose(xrf_stack)
	image_stack=reform(xrf_stack, n_cols, n_rows, npts)

ENDIF ELSE BEGIN
	ax_wait
; ---------------- READ in INCA format stacks --------------------

; --------  get array dimensions from the *.rpl file
	info_file = file_path + fileshort + '.rpl'
	openr, lun, info_file, /get_lun
	text = '' & test = 'q'
	while test NE 'DEPTH' DO BEGIN
		if eof(lun) then goto, jump_out
		readf, lun, text
;		print, text
		len = strpos(text,'::')
		if len gt 0 then begin
			l1 = strpos(text, ')')
			text2 = strmid(text, len+2, l1-len-2)
;			print, text2
			l2 = strpos(text2,' ')
			par_label = strmid(text2, 0, l2)
;			print, par_label
			CASE par_label OF
			   'WIDTH' : nx = fix(strcompress(strmid(text2,l2+1), /remove_all))
			   'HEIGHT' : ny = fix(strcompress(strmid(text2,l2+1), /remove_all))
			   'DEPTH' : npts = fix(strcompress(strmid(text2,l2+1), /remove_all))
			ELSE : goto, keep_on
			ENDCASE
			keep_on:
			test = par_label
		endif
	endwhile
	jump_out:
	close, lun & free_lun, lun
	IF NOT keyword_set(silent) then $
	       print, 'Reading INCA stack (*.raw) of dimensions (x,y,#E) ', nx, ny, npts

; ------- allow user to modify npts -------
	if keyword_set(group) then begin
 		npts = get_num(prompt='# of XRF energies', val=npts, group=group)
 		nx = get_num(prompt='# columns (nx)', val=nx, group=group)
 		ny = get_num(prompt='# rows (ny)', val=ny, group=group)
 		bw = dialog_message ('Byte or word (2-byte) ? (Yes = Byte)', /question)
	endif else begin
		npts = get_num(prompt='# of XRF energies', val=npts)
 		nx = get_num(prompt='# columns (nx)', val=nx)
 		ny = get_num(prompt='# rows (ny)', val=ny)
 		bw = dialog_message ('Byte or word (2-byte) ? (Yes = Byte)', /question)
	endelse
; ------------- use IDL READ_BINARY to read-in stack
	if bw EQ ' Yes' then image_stack=bytarr(nx,ny,npts) $
	           else image_stack=intarr(nx,ny,npts)
	openr, lun, file, /get_lun
	on_ioerror, no_read
	readu,lun, image_stack
	close, lun & free_lun, lun
	goto, continue
	no_read:
		axis_log, !ERROR_STATE.MSG
		ax_wait, /off
		return, s
	continue:

; ---------- set up x,y scales - in absence of other info = pixels
	x_start = 0  & x_stop = nx  & x_step = (x_stop - x_start)/nx
	y_start = 0  & y_stop = ny  & y_step = (y_stop - y_start)/ny
	ax_wait, /off
; ----- get incident electron energy (eV) -----
	energy = 20000.
	if keyword_set(group) then begin
		energy = get_num(prompt='Electron beam energy ', val = energy, group = axis_ID)
	endif else	energy = get_num(prompt='Electron beam energy ', val = energy)

ENDELSE

; ------------- prepare & write aXsi2000 stack file
n_cols=nx & n_rows=ny			; needed to pass elsewhere via bsif_com
n_data = npts

; -------- Create aXis2000 stack
ev_max= npts & ev_min = 0
slope = (ev_max - ev_min)/npts
ev = ev_min + slope*findgen(npts)

filename_ev_msec_list = strarr(npts)
for i = 0, npts-1 do filename_ev_msec_list(i) = string(ev(i))

s=fileshort

; --------- ask user what file name they want - default name is same as files without counter
if not keyword_set(no_write) then begin
	fileout = file_path + fileshort + '.ncb'
	sname=pickfile2(/write, file= fileout, filter='*.ncb', title = 'name of binary stack file')
	stack_wb, sname, /silent
endif

IF NOT keyword_set(silent) then begin
	E_string = ' at E = ' + string(energy, format="(F8.2)") + ' eV'
	axis_log, 'XRF stack read from ' + fileshort  + E_string
ENDIF

return, s

end

