Table of Contents
Removes background (constant for now) calculated in TOF units from a matrix workspace, expressed in units, different from TOF
Name | Direction | Type | Default | Description |
---|---|---|---|---|
InputWorkspace | Input | MatrixWorkspace | Mandatory | Workspace containing the input data |
OutputWorkspace | Output | MatrixWorkspace | Mandatory | The name to give the output workspace |
BkgWorkspace | Input | MatrixWorkspace | Mandatory | An optional histogram workspace in the units of TOF defining background for removal during rebinning.The workspace has to have single value or contain the same number of spectra as the “InputWorkspace” and single Y value per each spectra,representing flat background in the background time region. If such workspace is present, the value of the flat background provided by this workspace is removed from each spectra of the rebinned workspace. This works for histogram and event workspace when events are not retained but actually useful mainly for removing background while rebinning an event workspace in the units different from TOF. |
EMode | Input | string | Direct | The energy conversion mode used to define the conversion from the units of the InputWorkspace to TOF. Allowed values: [‘Elastic’, ‘Direct’, ‘Indirect’] |
NullifyNegativeValues | Input | boolean | False | When background is subtracted, signals in some time channels may become negative. If this option is true, signal in such bins is nullified and the module of the removed signalis added to the error. If false, the negative signal and correspondent errors remain unchanged |
Algorithm removes flat background, defined by the a single bin histogram workspace with X-axis in the units of TOF from a histogram workspace in any units with known conversion to TOF.
These options are especially useful during reduction of event workspaces in multirep mode, where different event regions are associated with different incident energies and rebinned into appropriate energy range.
The algorithm used during background removal is equivalent to the proof-of concept one, presented below, except intermediate workspaces are not created.
Errors of the background workspace are currently ignored and their value is calculated as the square root of correspondent background signal.
NOTE: If background is a single value workspace with zero signal and error, and property NullifyNegativeValues is set to True, the algorithm can be used to remove negative signal values and estimate errors during removal.
Proof of concept background removal algorithm:
from mantid.simpleapi import *
from mantid import config
import numpy as np
import sys
import os
maps_dir = '/home/user/InstrumentFiles/let/'
data_dir ='/home/user/results'
ref_data_dir = '/home/user/SystemTests/AnalysisTests/ReferenceResults'
config.setDataSearchDirs('{0};{1};{2}'.format(data_dir,maps_dir,ref_data_dir))
config['defaultsave.directory'] = data_dir # folder to save resulting spe/nxspe files. Defaults are in
# the name of a workspace containing background
filename = 'LET0007438'
groupedFilename = filename+'rings';
#
Ei= 25
e_min = -20
e_max = 20
dE = 0.1
bgRange = [15000,18000]
if not("Tgrid" in mtd):
if not(groupedFilename in mtd):
Load(Filename=filename+'.nxs', OutputWorkspace=filename, LoadMonitors=True)
GroupDetectors(InputWorkspace=filename, OutputWorkspace=groupedFilename , MapFile='LET_one2one_123.map', Behaviour='Average')
wsParent = mtd[groupedFilename];
nHist = wsParent.getNumberHistograms();
print("Parent workspace contains {0:10} histograms".format(nHist))
# Get the energy binning correspondent to the binning produced by rebin function (not to re-implement the same function)
ws1s = ExtractSingleSpectrum(wsParent,0);
ws1s = ConvertUnits(ws1s,'DeltaE','Direct',Ei);
ws1s = Rebin(ws1s,Params=[e_min,dE,e_max]);
e_bins = ws1s.dataX(0);
nBins = e_bins.size;
x=[e_bins[i] for i in xrange(0,nBins)]
y=[0 for xx in xrange(0,len(x)-1)]*nHist
x = x*nHist
DeleteWorkspace(ws1s);
eGrid = CreateWorkspace(DataX=x,DataY=y,UnitX='DeltaE',NSpec=nHist,VerticalAxisUnit='SpectraNumber',ParentWorkspace=wsParent)
Tgrid=ConvertUnits(eGrid,'TOF',Emode='Direct',EFixed=Ei)
else:
Tgrid = mtd['Tgrid'];
eGrid = mtd['eGrid'];
nHist = Tgrid.getNumberHistograms();
nBins = Tgrid.dataX(0).size;
if not('Bg' in mtd):
Bg=Rebin(InputWorkspace=groupedFilename, Params=[bgRange[0],bgRange[1]-bgRange[0],bgRange[1]],PreserveEvents=False)
#Bg=CalculateFlatBackground(InputWorkspace=groupedFilename, StartX=bgRange[0], EndX=bgRange[1], Mode='Mean', OutputMode='Return Background', SkipMonitors=True)
else:
Bg = mtd['Bg']
# Assign constant background to the Time grid workspace, minding different time bin width
for nspec in xrange(0,nHist):
bg = Bg.dataY(nspec)
if bg[0]>0:
bgT = Bg.dataX(nspec)
TimeScale = Tgrid.dataX(nspec);
# Jacobian for the unit conversion
Jac = (TimeScale[1:nBins]-TimeScale[0:nBins-1])*(bg[0]/(bgT[1]-bgT[0]));
error = np.sqrt(Jac);
eGrid.setY(nspec, Jac)
eGrid.setE(nspec, error)
else: # signal and error for background is 0 anyway.
pass
#print(" bg at spectra {0} equal to : {1}".format(nspec,bg[0]))
background = eGrid;
resultEt = ConvertUnits(groupedFilename,'DeltaE',Emode='Direct',EFixed=Ei)
result = Rebin(InputWorkspace=resultEt, Params=[e_min,dE,e_max],PreserveEvents=False)
fr = result-background;
#
sourceSum = SumSpectra(result,0,nHist);
bckgrdSum = SumSpectra(background ,0,nHist);
removedBkgSum = SumSpectra(fr ,0,nHist);
The results of executing this script on workspace contained measured background and the results of the background removal are presented on the following picture:
Blue line on this image represents the results, obtained using Rebin and Background removal after that. The results produced using the script below and shifted by one to show that there is another result plotted on the image, as both results are identical:
from mantid.simpleapi import *
from mantid import config
import numpy as np
import sys
import os
maps_dir = '/home/user/InstrumentFiles/let/'
data_dir ='/home/user/results'
ref_data_dir = '/home/user/SystemTests/AnalysisTests/ReferenceResults'
config.setDataSearchDirs('{0};{1};{2}'.format(data_dir,maps_dir,ref_data_dir))
config['defaultsave.directory'] = data_dir # folder to save resulting spe/nxspe files. Defaults are in
# the name of a workspace containing background
filename = 'LET0007438'
groupedFilename = filename+'rings';
#
Ei= 25
e_min = -20
e_max = 20
dE = 0.1
bgRange = [15000,18000]
if not(groupedFilename in mtd):
Load(Filename=filename+'.nxs', OutputWorkspace=filename, LoadMonitors=True)
GroupDetectors(InputWorkspace=filename, OutputWorkspace=groupedFilename , MapFile='LET_one2one_123.map', Behaviour='Average')
if not('Bg' in mtd):
Bg=Rebin(InputWorkspace=groupedFilename, Params=[bgRange[0],bgRange[1]-bgRange[0],bgRange[1]],PreserveEvents=False)
else:
Bg = mtd['Bg']
if 'resultEtransf' in mtd:
resultEtransf = mtd['resultEtransf']
else:
resultEtransf = ConvertUnits(groupedFilename,'DeltaE',Emode='Direct',EFixed=Ei)
noBgWorkspace = Rebin(InputWorkspace=resultEtransf, Params=[e_min,dE,e_max],PreserveEvents=False)
noBgWorkspace= Rebin(InputWorkspace=noBgWorkspace,BkgWorkspace='Bg',EMode='Direct')
nHist = Bg.getNumberHistograms()
removedBkgSum = SumSpectra(noBgWorkspace ,0,nHist-1);
Example - Background removal from a workspace in energy transfer units
# Create sample workspace with events
Test=CreateSampleWorkspace(WorkspaceType='Event', Function='Flat background')
# Add sample log necessary for unit conversion
AddSampleLog(Test,'Ei',LogText='25.',LogType='Number');
# Calculate background
Bg = Rebin(Test,Params='15000,5000,20000',PreserveEvents=False);
# Convert event's units
Test_BgDE=ConvertUnits(Test,Target='DeltaE',EMode='Direct');
# Calculate histograms for event workspace in energy binning
Sample = Rebin(Test_BgDE,Params='-20,2,20',PreserveEvents=False);
# Calculate histograms for event workspace in energy binning
Result = Rebin(Test_BgDE,Params='-20,2,20',PreserveEvents=False);
# Remove flat background in-place
Result = RemoveBackground(Result,BkgWorkspace='Bg',EMode='Direct');
# Get access to the results
XS = Sample.dataX(0);
XR = Result .dataX(0);
YS = Sample.dataY(0);
YR = Result .dataY(0);
ES = Sample.dataE(0);
ER = Result .dataE(0);
# print first spectra, Note invalid error calculations
print("| x sampl | x result | S sample | S no bg | Err samp | Err no_bg|")
for i in range(0,20):
print("|{0:10}|{1:10}|{2:10.4f}|{3:10.3f}|{4:10.3f}|{5:10.3f}|".format(XS[i],XR[i],YS[i],YR[i],ES[i],ER[i]))
| x sampl | x result | S sample | S no bg | Err samp | Err no_bg|
| -20.0| -20.0| 1.0000| -0.959| 1.000| 1.008|
| -18.0| -18.0| 2.0000| -0.101| 1.414| 1.420|
| -16.0| -16.0| 3.0000| 0.740| 1.732| 1.738|
| -14.0| -14.0| 1.0000| -1.441| 1.000| 1.012|
| -12.0| -12.0| 5.0000| 2.353| 2.236| 2.242|
| -10.0| -10.0| 2.0000| -0.885| 1.414| 1.426|
| -8.0| -8.0| 5.0000| 1.841| 2.236| 2.245|
| -6.0| -6.0| 2.0000| -1.481| 1.414| 1.431|
| -4.0| -4.0| 4.0000| 0.139| 2.000| 2.015|
| -2.0| -2.0| 3.0000| -1.315| 1.732| 1.753|
| 0.0| 0.0| 6.0000| 1.133| 2.449| 2.469|
| 2.0| 2.0| 7.0000| 1.454| 2.646| 2.669|
| 4.0| 4.0| 5.0000| -1.400| 2.236| 2.272|
| 6.0| 6.0| 7.0000| -0.499| 2.646| 2.688|
| 8.0| 8.0| 9.0000| 0.047| 3.000| 3.053|
| 10.0| 10.0| 11.0000| 0.054| 3.317| 3.388|
| 12.0| 12.0| 16.0000| 2.190| 4.000| 4.094|
| 14.0| 14.0| 16.0000| -2.188| 4.000| 4.162|
| 16.0| 16.0| 26.0000| 0.490| 5.099| 5.348|
| 18.0| 18.0| 39.0000| -0.581| 6.245| 6.728|
Categories: Algorithms | CorrectionFunctions\BackgroundCorrections
C++ source: RemoveBackground.cpp (last modified: 2016-12-22)
C++ header: RemoveBackground.h (last modified: 2016-10-10)