Newer
Older
//--------------------------------
// Includes
//--------------------------------
#include "MantidDataHandling/SetSampleMaterial.h"
#include "MantidAPI/ExperimentInfo.h"
#include "MantidAPI/Workspace.h"
#include "MantidGeometry/Crystal/OrientedLattice.h"
#include "MantidKernel/MandatoryValidator.h"
#include "MantidKernel/EnabledWhenProperty.h"
#include "MantidKernel/BoundedValidator.h"
#include "MantidKernel/PhysicalConstants.h"
using namespace Mantid::PhysicalConstants;
namespace Mantid {
namespace DataHandling {
// Register the algorithm into the AlgorithmFactory
DECLARE_ALGORITHM(SetSampleMaterial)
SetSampleMaterial::SetSampleMaterial() : Mantid::API::Algorithm() {}
SetSampleMaterial::~SetSampleMaterial() {}
const std::string SetSampleMaterial::name() const {
return "SetSampleMaterial";
}
int SetSampleMaterial::version() const { return (1); }
const std::string SetSampleMaterial::category() const {
return "Sample;DataHandling";
}
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
using namespace Mantid::DataHandling;
using namespace Mantid::API;
using namespace Kernel;
/**
* Initialize the algorithm
*/
void SetSampleMaterial::init() {
using namespace Mantid::Kernel;
declareProperty(
new WorkspaceProperty<Workspace>("InputWorkspace", "", Direction::InOut),
"The workspace with which to associate the sample ");
declareProperty("ChemicalFormula", "",
"ChemicalFormula or AtomicNumber must be given.");
declareProperty("AtomicNumber", 0,
"ChemicalFormula or AtomicNumber must be given");
declareProperty("MassNumber", 0, "Mass number if ion (default is 0)");
auto mustBePositive = boost::make_shared<BoundedValidator<double>>();
mustBePositive->setLower(0.0);
declareProperty("SampleNumberDensity", EMPTY_DBL(), mustBePositive,
"Optional: This number density of the sample in number of "
"formulas per cubic angstrom will be used instead of "
"calculated");
declareProperty("ZParameter", EMPTY_DBL(), mustBePositive,
"Number of atoms in the unit cell");
declareProperty("UnitCellVolume", EMPTY_DBL(), mustBePositive,
"Unit cell volume in Angstoms^3. Will be calculated from the "
"OrientedLattice if not supplied.");
declareProperty("CoherentXSection", EMPTY_DBL(), mustBePositive,
"Optional: This coherent cross-section for the sample "
"material in barns will be used instead of tabulated");
declareProperty("IncoherentXSection", EMPTY_DBL(), mustBePositive,
"Optional: This incoherent cross-section for the sample "
"material in barns will be used instead of tabulated");
declareProperty("AttenuationXSection", EMPTY_DBL(), mustBePositive,
"Optional: This absorption cross-section for the sample "
"material in barns will be used instead of tabulated");
declareProperty("ScatteringXSection", EMPTY_DBL(), mustBePositive,
"Optional: This total scattering cross-section (coherent + "
"incoherent) for the sample material in barns will be used "
"instead of tabulated");
// Perform Group Associations.
std::string formulaGrp("By Formula or Atomic Number");
setPropertyGroup("ChemicalFormula", formulaGrp);
setPropertyGroup("AtomicNumber", formulaGrp);
setPropertyGroup("MassNumber", formulaGrp);
std::string densityGrp("Sample Density");
setPropertyGroup("SampleNumberDensity", densityGrp);
setPropertyGroup("ZParameter", densityGrp);
setPropertyGroup("UnitCellVolume", densityGrp);
std::string specificValuesGrp("Override Cross Section Values");
setPropertyGroup("CoherentXSection", specificValuesGrp);
setPropertyGroup("IncoherentXSection", specificValuesGrp);
setPropertyGroup("AttenuationXSection", specificValuesGrp);
setPropertyGroup("ScatteringXSection", specificValuesGrp);
// Extra property settings
setPropertySettings(
"AtomicNumber",
new Kernel::EnabledWhenProperty("ChemicalFormula", Kernel::IS_DEFAULT));
setPropertySettings("MassNumber", new Kernel::EnabledWhenProperty(
"ChemicalFormula", Kernel::IS_DEFAULT));
setPropertySettings("UnitCellVolume",
new Kernel::EnabledWhenProperty("SampleNumberDensity",
Kernel::IS_DEFAULT));
setPropertySettings("ZParameter",
new Kernel::EnabledWhenProperty("SampleNumberDensity",
Kernel::IS_DEFAULT));
// output properties
declareProperty(
"SampleNumberDensityResult", EMPTY_DBL(),
"The provided or calculated sample number density in atoms/Angstrom^3",
Direction::Output);
declareProperty("ReferenceWavelength", EMPTY_DBL(),
"The reference wavelength in Angstroms", Direction::Output);
declareProperty("TotalXSectionResult", EMPTY_DBL(),
"The provided or calculated total cross-section for the "
"sample material in barns.",
Direction::Output);
declareProperty("IncoherentXSectionResult", EMPTY_DBL(),
"The provided or calculated incoherent cross-section for the "
"sample material in barns.",
Direction::Output);
declareProperty("CoherentXSectionResult", EMPTY_DBL(),
"The provided or calculated coherent cross-section for the "
"sample material in barns.",
Direction::Output);
declareProperty("AbsorptionXSectionResult", EMPTY_DBL(),
"The provided or calculated Absorption cross-section for the "
"sample material in barns.",
Direction::Output);
}
std::map<std::string, std::string> SetSampleMaterial::validateInputs() {
std::map<std::string, std::string> result;
const std::string chemicalSymbol = getProperty("ChemicalFormula");
const int z_number = getProperty("AtomicNumber");
const int a_number = getProperty("MassNumber");
if (chemicalSymbol.empty()) {
if (z_number <= 0) {
result["ChemicalFormula"] = "Need to specify the material";
} else {
if (z_number > 0)
result["AtomicNumber"] =
"Cannot specify both ChemicalFormula and AtomicNumber";
}
if (a_number > 0 && z_number <= 0)
result["AtomicNumber"] = "Specified MassNumber without AtomicNumber";
/**
* Add the cross sections to the neutron atom if they are not-empty
* numbers. All values are in barns.
*
* @param neutron The neutron to update
* @param coh_xs Coherent cross section
* @param inc_xs Incoherent cross section
* @param abs_xs Absorption cross section
* @param tot_xs Total scattering cross section
*/
void SetSampleMaterial::fixNeutron(NeutronAtom &neutron, double coh_xs,
double inc_xs, double abs_xs,
double tot_xs) {
if (!isEmpty(coh_xs))
neutron.coh_scatt_xs = coh_xs;
if (!isEmpty(inc_xs))
neutron.inc_scatt_xs = inc_xs;
if (!isEmpty(abs_xs))
neutron.abs_scatt_xs = abs_xs;
if (!isEmpty(tot_xs))
neutron.tot_scatt_xs = tot_xs;
}
/**
* Execute the algorithm
*/
void SetSampleMaterial::exec() {
// Get the input workspace
Workspace_sptr workspace = getProperty("InputWorkspace");
// an ExperimentInfo object has a sample
ExperimentInfo_sptr expInfo =
boost::dynamic_pointer_cast<ExperimentInfo>(workspace);
if (!bool(expInfo)) {
throw std::runtime_error("InputWorkspace does not have a sample object");
}
// determine the sample number density
double rho = getProperty("SampleNumberDensity"); // in Angstroms-3
double zParameter = getProperty("ZParameter"); // number of atoms
if (isEmpty(rho)) {
double unitCellVolume = getProperty("UnitCellVolume"); // in Angstroms^3
// get the unit cell volume from the workspace if it isn't set
if (isEmpty(unitCellVolume) && expInfo->sample().hasOrientedLattice()) {
unitCellVolume = expInfo->sample().getOrientedLattice().volume();
g_log.notice() << "found unit cell volume " << unitCellVolume
<< " Angstrom^-3\n";
// density is just number of atoms in the unit cell
// ...but only calculate it if you have both numbers
if ((!isEmpty(zParameter)) && (!isEmpty(unitCellVolume)))
rho = zParameter / unitCellVolume;
}
// get the scattering information - this will override table values
double coh_xs = getProperty("CoherentXSection"); // in barns
double inc_xs = getProperty("IncoherentXSection"); // in barns
double sigma_atten = getProperty("AttenuationXSection"); // in barns
double sigma_s = getProperty("ScatteringXSection"); // in barns
// determine the material
const std::string chemicalSymbol = getProperty("ChemicalFormula");
const int z_number = getProperty("AtomicNumber");
const int a_number = getProperty("MassNumber");
boost::scoped_ptr<Material> mat;
if (!chemicalSymbol.empty()) {
// Use chemical formula if given by user
Material::ChemicalFormula CF;
try {
CF = Material::parseChemicalFormula(chemicalSymbol);
} catch(std::runtime_error &ex) {
UNUSED_ARG(ex);
std::stringstream msg;
msg << "Could not parse chemical formula: " << chemicalSymbol
<< std::endl;
throw std::runtime_error(msg.str());
}
g_log.notice() << "Found " << CF.atoms.size() << " types of atoms in \""
<< chemicalSymbol << "\"\n";
NeutronAtom neutron(0, 0., 0., 0., 0., 0.,
0.); // starting thing for neutronic information
if (CF.atoms.size() == 1 && isEmpty(zParameter) && isEmpty(rho)) {
mat.reset(new Material(chemicalSymbol, CF.atoms[0]->neutron,
CF.atoms[0]->number_density));
} else {
double numAtoms = 0.; // number of atoms in formula
for (size_t i = 0; i < CF.atoms.size(); i++) {
neutron = neutron + CF.numberAtoms[i] * CF.atoms[i]->neutron;
g_log.information() << CF.atoms[i] << ": " << CF.atoms[i]->neutron
<< "\n";
numAtoms += static_cast<double>(CF.numberAtoms[i]);
// normalize the accumulated number by the number of atoms
neutron = (1. / numAtoms) *
neutron; // funny syntax b/c of operators in neutron atom
fixNeutron(neutron, coh_xs, inc_xs, sigma_atten, sigma_s);
// create the material
mat.reset(new Material(chemicalSymbol, neutron, rho));
} else {
// try using the atomic number
Atom atom = getAtom(static_cast<uint16_t>(z_number),
static_cast<uint16_t>(a_number));
NeutronAtom neutron = atom.neutron;
fixNeutron(neutron, coh_xs, inc_xs, sigma_atten, sigma_s);
// create the material
mat.reset(new Material(chemicalSymbol, neutron, rho));
}
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
// set the material but leave the geometry unchanged
auto shapeObject = expInfo->sample().getShape(); // copy
shapeObject.setMaterial(*mat);
expInfo->mutableSample().setShape(shapeObject);
g_log.notice() << "Sample number density ";
if (isEmpty(mat->numberDensity())) {
g_log.notice() << "was not specified\n";
} else {
g_log.notice() << "= " << mat->numberDensity() << " atoms/Angstrom^3\n";
setProperty("SampleNumberDensityResult",
mat->numberDensity()); // in atoms/Angstrom^3
}
g_log.notice() << "Cross sections for wavelength = "
<< NeutronAtom::ReferenceLambda << " Angstroms\n"
<< " Coherent " << mat->cohScatterXSection() << " barns\n"
<< " Incoherent " << mat->incohScatterXSection()
<< " barns\n"
<< " Total " << mat->totalScatterXSection() << " barns\n"
<< " Absorption " << mat->absorbXSection() << " barns\n";
setProperty("CoherentXSectionResult", mat->cohScatterXSection()); // in barns
setProperty("IncoherentXSectionResult",
mat->incohScatterXSection()); // in barns
setProperty("TotalXSectionResult", mat->totalScatterXSection()); // in barns
setProperty("AbsorptionXSectionResult", mat->absorbXSection()); // in barns
setProperty("ReferenceWavelength",
NeutronAtom::ReferenceLambda); // in Angstroms
if (isEmpty(rho)) {
g_log.notice("Unknown value for number density");
} else {
double smu = mat->totalScatterXSection(NeutronAtom::ReferenceLambda) * rho;
double amu = mat->absorbXSection(NeutronAtom::ReferenceLambda) * rho;
g_log.notice() << "Anvred LinearScatteringCoef = " << smu << " 1/cm\n"
<< "Anvred LinearAbsorptionCoef = " << amu << " 1/cm\n";
// Done!
progress(1);
}
}