; $Id: dcs_sumidenticals.pro,v 1.11 2005/02/14 17:50:14 azuah Exp $
; Written by J.R.D. Copley.
;************************************************************************************************
;###############################################################################
;
; LICENSE:
;  The software in this file is written by an employee of 
;  National Institute of Standards and Technology 
;  as part of the DAVE software project.
;
;  The DAVE software package is not subject to copyright protection
;  and is in the public domain. It should be considered as an
;  experimental neutron scattering data reduction, visualization, and
;  analysis system. As such, the authors assume no responsibility
;  whatsoever for its use, and make no guarantees, expressed or
;  implied, about its quality, reliability, or any other
;  characteristic. The use of certain trade names or commercial
;  products does not imply any endorsement of a particular product,
;  nor does it imply that the named product is necessarily the best
;  product for the stated purpose. We would appreciate acknowledgment
;  if the DAVE software is used or if the code in this file is
;  included in another product.
;
;###############################################################################
pro dcs_sumidenticals_dothesum,dcs_base,sumfile,filelist,nfiles,files_selected,ntchan,nhedet,$
	mask,normalization_values,fbnormalization_values,bitflips,maxbitflips,davePtr=davePtr,$
	logentry=logentry
;************************************************************************************************
;
compile_opt strictarr
;
omithisto=0
printinfo=0
filetype="binary"
;
unitvec_ntchan=bytarr(ntchan)+1
unitvec_nhedet=bytarr(nhedet)+1
;
lines=strarr(nfiles+4)
lines[0]="Summing identical files, applying masks, correcting bitflips,"
lines[1]="dividing each channel by its width in time, and saving to disk."
bases=intarr(nfiles)
;
; Quantities from raw data files (other than histodata) that are to be summed.
squantity=["badbrd","baddet","duration","nframes","ncycles",$
	"detsum","grandsum","histohigh","repeats","resets",$
	"timsum","totals"]
nsq=n_elements(squantity)
;
for j=0,nfiles-1 do begin
	index=filelist[j]
	filepath=files_selected[index]
	lines[j+3]="Reading file "+strcompress(index+1)+":"+filepath
;	widget_control,/hourglass
	dcs_message_display,lines,base,group_leader=dcs_base,alignment=0
	bases[j]=base
	dcs_readdatafile,filepath,omithisto,printinfo,filetype,davePtr=davePtr
;
; If necessary perform bitflip correction	to data.
	for k=0,maxbitflips-1 do begin
		if (bitflips[1,k,index] lt 1000) then begin
			bftim=bitflips[0,k,index]
			bfdet=bitflips[1,k,index]
			bfamt=bitflips[2,k,index]
			if (!debug) then print,"A bitflip has been corrected in dcs_sumidenticals."
			if (!debug) then print,j,bftim,bfdet,bfamt
			if (bfdet lt nhedet) then begin
				(*(*(*davePtr).dataStrPtr).commonStr.histPtr).qty[bftim,bfdet] = $
					(*(*(*davePtr).dataStrPtr).commonStr.histPtr).qty[bftim,bfdet]-bfamt
			endif else begin
				(*(*(*davePtr).dataStrPtr).specificPtr).histohigh[bftim,bfdet-nhedet] = $
					(*(*(*davePtr).dataStrPtr).specificPtr).histohigh[bftim,bfdet-nhedet]-bfamt
			endelse
		endif
	endfor
;
; Next line creates matrix with dimensions [ntchan,nhedet].
	maskj=unitvec_ntchan#mask[*,index]
;
	if (j eq 0) then begin
;
; These quantities are to be summed.
		pdata_qty=ptr_new(dcs_operate((*(*(*davePtr).dataStrPtr).commonStr.histPtr).qty,maskj,"*"))
		normalization_sum=dcs_operate(normalization_values[index],mask[*,index],"*")
		fbnormalization_sum=dcs_operate(fbnormalization_values[index],mask[*,index],"*")
	endif
;
	if (j eq 0 and nfiles gt 1) then begin
;
; These quantities are to be summed.
		for k=0,nsq-1 do begin
			com=squantity[k]+"sum=(*(*(*davePtr).dataStrPtr).specificPtr)."+squantity[k]
			res=execute(com)
		endfor
;
; These quantities, from the first file, will be placed in the output sum file.
		temp_setpoint=(*(*(*davePtr).dataStrPtr).specificPtr).temp_setpoint
		command=(*(*(*davePtr).dataStrPtr).specificPtr).command
		comments=(*(*(*davePtr).dataStrPtr).specificPtr).comments
		sample_desc=(*(*(*davePtr).dataStrPtr).specificPtr).sample_desc
		start_date=(*(*(*davePtr).dataStrPtr).specificPtr).start_date
		stop_date=(*(*(*davePtr).dataStrPtr).specificPtr).stop_date
		user=(*(*(*davePtr).dataStrPtr).specificPtr).user
;
; These quantities are to be concatenated.
		duration_csum=(*(*(*davePtr).dataStrPtr).specificPtr).duration_csum
		he3_csum=(*(*(*davePtr).dataStrPtr).specificPtr).he3_csum
		high_csum=(*(*(*davePtr).dataStrPtr).specificPtr).high_csum
		temp_control=(*(*(*davePtr).dataStrPtr).specificPtr).temp_control
		temp_sample=(*(*(*davePtr).dataStrPtr).specificPtr).temp_sample
;
	endif
;
	if (j gt 0) then begin
;
; Quantities to be summed.
		*pdata_qty=*pdata_qty+(*(*(*davePtr).dataStrPtr).commonStr.histPtr).qty*maskj
		normalization_sum=normalization_sum+normalization_values[index]*mask[*,index]
		fbnormalization_sum=fbnormalization_sum+fbnormalization_values[index]*mask[*,index]
		for k=0,nsq-1 do begin
			com=squantity[k]+"sum="+squantity[k]+$
				"sum+(*(*(*davePtr).dataStrPtr).specificPtr)."+squantity[k]
			res=execute(com)
		endfor
;
; Quantities to be concatenated.
;
		void=where(duration_csum ge 0,nt0)
		void=where((*(*(*davePtr).dataStrPtr).specificPtr).duration_csum ge 0,nt)
		if (nt gt 0) then begin
			duration_csum[nt0:nt0+nt-1]=(*(*(*davePtr).dataStrPtr).specificPtr).duration_csum[0:nt-1]
			he3_csum[nt0:nt0+nt-1]=(*(*(*davePtr).dataStrPtr).specificPtr).he3_csum[0:nt-1]
			high_csum[*,nt0:nt0+nt-1]=(*(*(*davePtr).dataStrPtr).specificPtr).high_csum[*,0:nt-1]
		endif
;
		void=where(temp_control gt 0,ntc0)
		void=where((*(*(*davePtr).dataStrPtr).specificPtr).temp_control gt 0,ntc)
		if (ntc gt 0) then $
			temp_control[ntc0:ntc0+ntc-1]=(*(*(*davePtr).dataStrPtr).specificPtr).temp_control[0:ntc-1]
;
		void=where(temp_sample gt 0,nts0)
		void=where((*(*(*davePtr).dataStrPtr).specificPtr).temp_sample gt 0,nts)
		if (nts gt 0) then $
			temp_sample[nts0:nts0+nts-1]=(*(*(*davePtr).dataStrPtr).specificPtr).temp_sample[0:nts-1]
;
	endif
endfor
;
; Summed quantities.
(*(*(*davePtr).dataStrPtr).commonStr.histPtr).qty=*pdata_qty
(*(*(*davePtr).dataStrPtr).commonStr.histPtr).err=sqrt(*pdata_qty)
(*(*(*davePtr).dataStrPtr).specificPtr).normalization_sum=normalization_sum
(*(*(*davePtr).dataStrPtr).specificPtr).fbnormalization_sum=fbnormalization_sum

if (nfiles gt 1) then begin
	for k=0,nsq-1 do begin
		com="(*(*(*davePtr).dataStrPtr).specificPtr)."+squantity[k]+"="+squantity[k]+"sum"
		res=execute(com)
	endfor
;
; Maximum quantities.
	(*(*(*davePtr).dataStrPtr).specificPtr).datamax=max(*pdata_qty)
	(*(*(*davePtr).dataStrPtr).specificPtr).highmax=max(histohighsum)
;
; Quantities from first file.
	(*(*(*davePtr).dataStrPtr).specificPtr).temp_setpoint=temp_setpoint
	(*(*(*davePtr).dataStrPtr).specificPtr).command=command
	(*(*(*davePtr).dataStrPtr).specificPtr).comments=comments
	(*(*(*davePtr).dataStrPtr).specificPtr).sample_desc=sample_desc
	(*(*(*davePtr).dataStrPtr).specificPtr).start_date=start_date
	(*(*(*davePtr).dataStrPtr).specificPtr).stop_date=stop_date
	(*(*(*davePtr).dataStrPtr).specificPtr).user=user
;
; Concatenated quantities.
	(*(*(*davePtr).dataStrPtr).specificPtr).duration_csum=duration_csum
	(*(*(*davePtr).dataStrPtr).specificPtr).he3_csum=he3_csum
	(*(*(*davePtr).dataStrPtr).specificPtr).high_csum=high_csum
	(*(*(*davePtr).dataStrPtr).specificPtr).temp_control=temp_control
	(*(*(*davePtr).dataStrPtr).specificPtr).temp_sample=temp_sample
;
endif
;
lines[nfiles+3]="Writing output file "+sumfile
dcs_message_display,lines,base,group_leader=dcs_base,alignment=0
;;
;widget_control,/hourglass
;
; Dereference histPtr, placing (most of) qty, (most of) err and x and y
;	into new pointers; qty and err are first divided by the appropriate time
; channel width.
; The x array is redefined as times of flight from sample to detector.
ngtchan0=(*(*(*davePtr).dataStrPtr).specificPtr).ngtchan0
times=(*(*(*(*davePtr).dataStrPtr).specificPtr).time_propsPtr)[0:ngtchan0]
twidPtr=ptr_new((times[1:ngtchan0]-times[0:ngtchan0-1])#unitvec_nhedet)
qtyPtr=ptr_new((*(*(*davePtr).dataStrPtr).commonStr.histPtr).qty[0:ngtchan0-1,*]/*twidPtr)
errPtr=ptr_new((*(*(*davePtr).dataStrPtr).commonStr.histPtr).err[0:ngtchan0-1,*]/*twidPtr)
xPtr=ptr_new(times)
yPtr=ptr_new((*(*(*davePtr).dataStrPtr).commonStr.histPtr).y)
;
; Free and then redefine histPtr. Free the temporary pointers.
ptr_free,(*(*davePtr).dataStrPtr).commonStr.histPtr
(*(*davePtr).dataStrPtr).commonStr.histPtr=ptr_new({qty:*qtyPtr,err:*errPtr,x:*xPtr,y:*yPtr})
ptr_free,qtyPtr
ptr_free,errPtr
ptr_free,xPtr
ptr_free,yPtr
ptr_free,twidPtr
;
(*(*davePtr).dataStrPtr).commonStr.xlabel = 'Time-of-flight from sample to detector'
(*(*davePtr).dataStrPtr).commonStr.xunits = 'time:us'
(*(*davePtr).dataStrPtr).commonStr.xtype = 'histogram'
;
(*(*davePtr).dataStrPtr).commonStr.histlabel = 'Summed intensity'
;
(*(*(*davePtr).dataStrPtr).commonStr.treatmentPtr)=strsplit(logentry,string(10B),/extract)
;
save,davePtr,filename=sumfile
;
ptr_free,pdata_qty
;
widget_control,base,/destroy
for j=nfiles-1,0,-1 do widget_control,bases[j],/destroy
;
end


;************************************************************************************************
pro dcs_sumidenticals,dcs_base,files_selected,nfiles_selected,nruns,normalization_values,$
	fbnormalization_values,mask,assignments,ntchan,nhedet,sumfiles,n_sumfiles,$
	bitflips,maxbitflips,davePtr=davePtr,logentry=logentry
;************************************************************************************************
; This routine takes raw data input files and creates summed data output files.
;
compile_opt strictarr
;
n_sample_conditions=max(assignments); this is the number of different sample conditions
number=where(nfiles_selected gt 0, count); this tells us which file types are present
n_summed_total=n_sample_conditions+count-1; this is the number of output files
;
n_sumfiles=intarr(5); this will contain the number of output files for each file type (S,C,V,A,B).
sumfiles=strarr(n_summed_total); this will contain the names of the output files
;
; Get path and base filename.
pathbase=!dcs_outp_directory+!dcs_rootname+"_"
name=["sam","can","van","abs","bgd"]
ext="sum.dave"
;
unitvec_ntchan=bytarr(ntchan)+1
;
logentry=logentry+string(10B)+"Nominally identical files have been summed."
if (not lmgr(/demo)) then begin
	openw,unit,!logfile,/get_lun
	printf,unit,logentry
	free_lun,unit
endif
;
; Loop over run types treating sample runs first.
n_sumfiles[0]=n_sample_conditions
for k=0,n_sample_conditions-1 do begin; loop over sample conditions
	sumfiles[k]=pathbase+name[0]+strcompress(k+1,/remove_all)+ext; output file name
	filelist=where(assignments eq k+1,nfiles); identify input files
	dcs_sumidenticals_dothesum,dcs_base,sumfiles[k],filelist,nfiles,files_selected,ntchan,nhedet,$
		mask,normalization_values,fbnormalization_values,bitflips,maxbitflips,davePtr=davePtr,$
		logentry=logentry
endfor
nfiles_cum=nfiles_selected[0]
k=n_sample_conditions-1
	;
; Now treat runs other than sample runs.
for itype=1,4 do begin
	n_sumfiles[itype]=0
	if (nfiles_selected[itype] gt 0) then begin
		n_sumfiles[itype]=1
		k=k+1
		sumfiles[k]=pathbase+name[itype]+ext; output file name
		nfiles=nfiles_selected[itype]; number of files to be summed
		filelist=indgen(nfiles)+nfiles_cum
		nfiles_cum=nfiles_cum+nfiles
		dcs_sumidenticals_dothesum,dcs_base,sumfiles[k],filelist,nfiles,files_selected,ntchan,nhedet,$
			mask,normalization_values,fbnormalization_values,bitflips,maxbitflips,davePtr=davePtr,$
			logentry=logentry
		endif
endfor
;
free_all_ptrs,davePtr
;
end
