diff --git a/Framework/Algorithms/src/CalculatePlaczekSelfScattering.cpp b/Framework/Algorithms/src/CalculatePlaczekSelfScattering.cpp
index 68d1342a98f7ac56c1f931442ff970487d01ae9c..d834c27ded871c0ad68cf1b4187dd7b2c4b07f82 100644
--- a/Framework/Algorithms/src/CalculatePlaczekSelfScattering.cpp
+++ b/Framework/Algorithms/src/CalculatePlaczekSelfScattering.cpp
@@ -34,13 +34,16 @@ getSampleSpeciesInfo(const API::MatrixWorkspace_const_sptr ws) {
for (const Kernel::Material::Material::FormulaUnit element : formula) {
const std::map<std::string, double> atomMap{
{"mass", element.atom->mass},
- {"concentration", element.multiplicity},
+ {"stoich", element.multiplicity},
{"bSqrdBar", bSqrdBar}};
atomSpecies[element.atom->symbol] = atomMap;
totalStoich += element.multiplicity;
}
- for (auto atom : atomSpecies) {
- atom.second["concentration"] = atom.second["stoich"] / totalStoich;
+ std::map<std::string, std::map<std::string, double>>::iterator atom =
+ atomSpecies.begin();
+ while (atom != atomSpecies.end()) {
+ atom->second["concentration"] = atom->second["stoich"] / totalStoich;
+ atom++;
}
return atomSpecies;
}
@@ -113,8 +116,8 @@ void CalculatePlaczekSelfScattering::exec() {
// set detector law term (eps1)
const double lambdaD = 1.44;
std::vector<double> eps1;
- for (size_t i = 0; i < xLambda.size(); i++) {
- double xTerm = -xLambda[i] / lambdaD;
+ for (size_t i = 0; i < xLambda.size() - 1; i++) {
+ double xTerm = -(xLambda[i] + dx) / lambdaD;
eps1.push_back(xTerm * exp(xTerm) / (1.0 - exp(xTerm)));
}
/* Placzek
@@ -148,7 +151,6 @@ void CalculatePlaczekSelfScattering::exec() {
outputWS->setDistribution(true);
outputWS->setYUnit("");
outputWS->setYUnitLabel("Counts");
- int nSpec = 0;
for (size_t specIndex = 0; specIndex < specInfo.size(); specIndex++) {
auto &y = outputWS->mutableY(specIndex);
auto &x = outputWS->mutableX(specIndex);
@@ -156,28 +158,29 @@ void CalculatePlaczekSelfScattering::exec() {
const double pathLength = specInfo.l1() + specInfo.l2(specIndex);
const double f = specInfo.l1() / pathLength;
const double sinThetaBy2 = sin(specInfo.twoTheta(specIndex) / 2.0);
+ Kernel::Units::Wavelength wavelength;
+ wavelength.initialize(specInfo.l1(), specInfo.l2(specIndex),
+ specInfo.twoTheta(specIndex), 0, 1.0, 1.0);
for (size_t xIndex = 0; xIndex < xLambda.size() - 1; xIndex++) {
const double term1 = (f - 1.0) * phi1[xIndex];
- const double term2 = f * eps1[xIndex];
- const double term3 = f - 3;
+ const double term2 = -f * eps1[xIndex];
+ const double term3 = f - 3.0;
const double inelasticPlaczekSelfCorrection =
- 2.0 * (term1 - term2 + term3) * sinThetaBy2 * sinThetaBy2 *
+ 2.0 * (term1 + term2 + term3) * sinThetaBy2 * sinThetaBy2 *
summationTerm;
- x[xIndex] = xLambda[xIndex];
+ x[xIndex] = wavelength.singleToTOF(xLambda[xIndex]);
y[xIndex] = inelasticPlaczekSelfCorrection;
}
- xLambdas.push_back(xLambda.back());
- nSpec += 1;
+ x.back() = wavelength.singleToTOF(xLambda.back());
} else {
for (size_t xIndex = 0; xIndex < xLambda.size() - 1; xIndex++) {
x[xIndex] = xLambda[xIndex];
y[xIndex] = 0;
}
x.back() = xLambda.back();
- nSpec += 1;
}
}
- auto incidentUnit = incidentWS->getAxis(0)->unit();
+ auto incidentUnit = inWS->getAxis(0)->unit();
outputWS->getAxis(0)->unit() = incidentUnit;
setProperty("OutputWorkspace", outputWS);
}
diff --git a/scripts/Diffraction/isis_powder/polaris_routines/polaris_algs.py b/scripts/Diffraction/isis_powder/polaris_routines/polaris_algs.py
index 2368c051045d3be087be3972f396cf5a22238ce3..bdefcd7ba63ce080f3d644b77d6c09d2f9f96589 100644
--- a/scripts/Diffraction/isis_powder/polaris_routines/polaris_algs.py
+++ b/scripts/Diffraction/isis_powder/polaris_routines/polaris_algs.py
@@ -133,33 +133,51 @@ def _generate_grouped_ts_pdf(run_number, focused_ws, q_lims, cal_file_name, samp
min_x = min(np.min(x_data), min_x)
max_x = max(np.max(x_data), max_x)
num_x = max(x_data.size, num_x)
- binning = [min_x, (max_x-min_x)/num_x, max_x]
+ width_x = (max_x-min_x)/num_x
+ binning = [min_x, width_x, max_x]
focused_ws = mantid.Rebin(InputWorkspace=focused_ws, Params=binning)
focused_data_combined = mantid.ConjoinSpectra(InputWorkspaces=focused_ws)
mantid.ConvertFromDistribution(Workspace=focused_data_combined)
+
raw_ws = mantid.LoadRaw(Filename='POL'+str(run_number))
mantid.SetSample(InputWorkspace=raw_ws,
Geometry=common.generate_sample_geometry(sample_details),
Material=common.generate_sample_material(sample_details))
-
monitor = mantid.ExtractSpectra(InputWorkspace=raw_ws, WorkspaceIndexList=[11])
-
- fit_spectra = mantid.FitIncidentSpectrum(InputWorkspace=monitor, FitSpectrumWith="GaussConvCubicSpline")
+ monitor = mantid.ConvertUnits(InputWorkspace=monitor, Target="Wavelength")
+ x_data = monitor.dataX(0)
+ min_x = np.min(x_data)
+ max_x = np.max(x_data)
+ width_x = (max_x-min_x)/x_data.size
+ fit_spectra = mantid.FitIncidentSpectrum(InputWorkspace=monitor,
+ BinningForCalc=[min_x, 1*width_x, max_x],
+ BinningForFit=[min_x, 50*width_x, max_x],
+ FitSpectrumWith="GaussConvCubicSpline")
placzek = mantid.CalculatePlaczekSelfScattering(InputWorkspace=raw_ws, IncidentSpecta=fit_spectra)
mantid.ConvertFromDistribution(Workspace=placzek)
+ cal_workspace = mantid.LoadCalFile(InputWorkspace=placzek,
+ CalFileName=cal_file_name,
+ Workspacename='cal_workspace',
+ MakeOffsetsWorkspace=False,
+ MakeMaskWorkspace=False)
placzek = mantid.AlignDetectors(InputWorkspace=placzek, CalibrationFile=cal_file_name)
- placzek = mantid.DiffractionFocussing(InputWorkspace=placzek, GroupingFileName=cal_file_name)
+ placzek = mantid.DiffractionFocussing(InputWorkspace=placzek, GroupingFilename=cal_file_name)
+ n_pixel = np.zeros(placzek.getNumberHistograms())
+ for i in range(cal_workspace.getNumberHistograms()):
+ grouping = cal_workspace.dataY(i)
+ if grouping[0] > 0:
+ n_pixel[int(grouping[0]-1)] += 1
+ correction_ws = mantid.CreateWorkspace(DataY=n_pixel, DataX=[0, 1], NSpec=placzek.getNumberHistograms())
+ placzek = mantid.Divide(LHSWorkspace=placzek, RHSWorkspace=correction_ws)
placzek = mantid.ConvertUnits(InputWorkspace=placzek, Target="MomentumTransfer", EMode='Elastic')
- mantid.RebinToWorkspace(WorkspaceToRebin=placzek,
- WorkspaceToMatch=focused_data_combined,
- OutputWorkspace=placzek)
+ placzek = mantid.RebinToWorkspace(WorkspaceToRebin=placzek, WorkspaceToMatch=focused_data_combined)
mantid.Subtract(LHSWorkspace=focused_data_combined,
RHSWorkspace=placzek,
OutputWorkspace=focused_data_combined)
mantid.MatchSpectra(InputWorkspace=focused_data_combined,
OutputWorkspace=focused_data_combined,
- ReferenceSpectrum=5)
+ ReferenceSpectrum=1)
if type(q_lims) == str:
q_min = []
q_max = []