diff --git a/Framework/Crystal/inc/MantidCrystal/FindUBUsingLatticeParameters.h b/Framework/Crystal/inc/MantidCrystal/FindUBUsingLatticeParameters.h
index ae0a6de8e88df35bf913096fa01f221dcf2f8893..dc8eebab33015d1b5cf2c85e76db846ea74efa4d 100644
--- a/Framework/Crystal/inc/MantidCrystal/FindUBUsingLatticeParameters.h
+++ b/Framework/Crystal/inc/MantidCrystal/FindUBUsingLatticeParameters.h
@@ -1,8 +1,8 @@
#ifndef MANTID_CRYSTAL_FIND_UB_USING_LATTICE_PARAMETERS_H_
#define MANTID_CRYSTAL_FIND_UB_USING_LATTICE_PARAMETERS_H_
-#include "MantidKernel/System.h"
#include "MantidAPI/Algorithm.h"
+#include "MantidKernel/System.h"
namespace Mantid {
namespace Crystal {
diff --git a/Framework/Crystal/src/FindUBUsingLatticeParameters.cpp b/Framework/Crystal/src/FindUBUsingLatticeParameters.cpp
index 12d301bc14f030b9c72c87bf59bc0ca836176eb0..03b93f8e3dd35477005d5a165cae055c30c44d9b 100644
--- a/Framework/Crystal/src/FindUBUsingLatticeParameters.cpp
+++ b/Framework/Crystal/src/FindUBUsingLatticeParameters.cpp
@@ -44,6 +44,8 @@ void FindUBUsingLatticeParameters::init() {
"Lattice parameter gamma");
this->declareProperty("NumInitial", 15, moreThan2Int,
"Number of Peaks to Use on First Pass(15)");
+ this->declareProperty("FixAll", false,
+ "Do not optimise the UB lattice parameters");
this->declareProperty("Tolerance", 0.15, mustBePositive,
"Indexing Tolerance (0.15)");
}
@@ -59,6 +61,7 @@ void FindUBUsingLatticeParameters::exec() {
double gamma = this->getProperty("gamma");
int num_initial = this->getProperty("NumInitial");
double tolerance = this->getProperty("Tolerance");
+ auto fixAll = this->getProperty("FixAll");
int base_index = -1; // these "could" be properties if need be
double degrees_per_step = 1.5;
@@ -79,7 +82,7 @@ void FindUBUsingLatticeParameters::exec() {
Matrix<double> UB(3, 3, false);
double error = IndexingUtils::Find_UB(UB, q_vectors, a, b, c, alpha, beta,
gamma, tolerance, base_index,
- num_initial, degrees_per_step);
+ num_initial, degrees_per_step, fixAll);
std::cout << "Error = " << error << '\n';
std::cout << "UB = " << UB << '\n';
@@ -93,15 +96,17 @@ void FindUBUsingLatticeParameters::exec() {
{ // and save the UB in the sample
std::vector<double> sigabc(7);
- std::vector<V3D> miller_ind;
- std::vector<V3D> indexed_qs;
- double fit_error;
- miller_ind.reserve(q_vectors.size());
- indexed_qs.reserve(q_vectors.size());
- IndexingUtils::GetIndexedPeaks(UB, q_vectors, tolerance, miller_ind,
- indexed_qs, fit_error);
-
- IndexingUtils::Optimize_UB(UB, miller_ind, indexed_qs, sigabc);
+// std::vector<V3D> miller_ind;
+// std::vector<V3D> indexed_qs;
+// double fit_error;
+// miller_ind.reserve(q_vectors.size());
+// indexed_qs.reserve(q_vectors.size());
+// IndexingUtils::GetIndexedPeaks(UB, q_vectors, tolerance, miller_ind,
+// indexed_qs, fit_error);
+
+//// if(!fixAll) {
+//// IndexingUtils::Optimize_UB(UB, miller_ind, indexed_qs, sigabc, fixAngles);
+//// }
char logInfo[200];
int num_indexed = IndexingUtils::NumberIndexed(UB, q_vectors, tolerance);
diff --git a/Framework/Crystal/test/FindUBUsingLatticeParametersTest.h b/Framework/Crystal/test/FindUBUsingLatticeParametersTest.h
index 22082cccdaf8cb1b6c55917b8e002cba82976518..fb547cf7da481e235a7c027b21b85946d9ff2e45 100644
--- a/Framework/Crystal/test/FindUBUsingLatticeParametersTest.h
+++ b/Framework/Crystal/test/FindUBUsingLatticeParametersTest.h
@@ -27,28 +27,12 @@ public:
}
void test_exec() {
- // Name of the output workspace.
- std::string WSName("peaks");
- LoadIsawPeaks loader;
- TS_ASSERT_THROWS_NOTHING(loader.initialize());
- TS_ASSERT(loader.isInitialized());
- loader.setPropertyValue("Filename", "TOPAZ_3007.peaks");
- loader.setPropertyValue("OutputWorkspace", WSName);
-
- TS_ASSERT(loader.execute());
- TS_ASSERT(loader.isExecuted());
+ auto ws = loadPeaksWorkspace();
- PeaksWorkspace_sptr ws;
- TS_ASSERT_THROWS_NOTHING(
- ws = boost::dynamic_pointer_cast<PeaksWorkspace>(
- AnalysisDataService::Instance().retrieve(WSName)));
- TS_ASSERT(ws);
- if (!ws)
- return;
FindUBUsingLatticeParameters alg;
TS_ASSERT_THROWS_NOTHING(alg.initialize())
TS_ASSERT(alg.isInitialized())
- TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("PeaksWorkspace", WSName));
+ TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("PeaksWorkspace", ws->name()));
TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("a", "14.131"));
TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("b", "19.247"));
TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("c", "8.606"));
@@ -75,8 +59,85 @@ public:
TS_ASSERT_DELTA(correct_UB[i], UB_calculated[i], 5e-4);
}
+ TS_ASSERT_DELTA(latt.a(), 14.131, 5e-4);
+ TS_ASSERT_DELTA(latt.b(), 19.247, 5e-4);
+ TS_ASSERT_DELTA(latt.c(), 8.606, 5e-4);
+
+ TS_ASSERT_DELTA(latt.alpha(), 90.0, 5e-1);
+ TS_ASSERT_DELTA(latt.beta(), 105.071, 5e-1);
+ TS_ASSERT_DELTA(latt.gamma(), 90.0, 5e-1);
+
+ // Remove workspace from the data service.
+ AnalysisDataService::Instance().remove(ws->name());
+ }
+
+ void test_fixAll() {
+ auto ws = loadPeaksWorkspace();
+
+ FindUBUsingLatticeParameters alg;
+ TS_ASSERT_THROWS_NOTHING(alg.initialize())
+ TS_ASSERT(alg.isInitialized())
+ TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("PeaksWorkspace", ws->name()));
+ TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("a", "14.131"));
+ TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("b", "19.247"));
+ TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("c", "8.606"));
+ TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("alpha", "90.0"));
+ TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("beta", "105.071"));
+ TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("gamma", "90.0"));
+ TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("NumInitial", "15"));
+ TS_ASSERT_THROWS_NOTHING(alg.setPropertyValue("Tolerance", "0.12"));
+ TS_ASSERT_THROWS_NOTHING(alg.setProperty("FixAll", true));
+ TS_ASSERT_THROWS_NOTHING(alg.execute());
+ TS_ASSERT(alg.isExecuted());
+
+ // Check that we set an oriented lattice
+ TS_ASSERT(ws->mutableSample().hasOrientedLattice());
+ // Check that the UB matrix is the same as in TOPAZ_3007.mat
+ OrientedLattice latt = ws->mutableSample().getOrientedLattice();
+
+ double correct_UB[] = {0.04542050, 0.040619900, 0.0127661,
+ -0.00198382, -0.00264404, -0.1165450,
+ -0.05749760, 0.03223800, -0.0257623};
+
+ std::vector<double> UB_calculated = latt.getUB().getVector();
+
+ for (size_t i = 0; i < 9; i++) {
+ TS_ASSERT_DELTA(correct_UB[i], UB_calculated[i], 5e-4);
+ }
+
+ TS_ASSERT_DELTA(latt.a(), 14.131, 5e-10);
+ TS_ASSERT_DELTA(latt.b(), 19.247, 5e-10);
+ TS_ASSERT_DELTA(latt.c(), 8.606, 5e-10);
+
+ TS_ASSERT_DELTA(latt.alpha(), 90.0, 5e-10);
+ TS_ASSERT_DELTA(latt.beta(), 105.071, 5e-10);
+ TS_ASSERT_DELTA(latt.gamma(), 90.0, 5e-10);
+
// Remove workspace from the data service.
- AnalysisDataService::Instance().remove(WSName);
+ AnalysisDataService::Instance().remove(ws->name());
+ }
+
+private:
+ /*
+ * Load a peaks workspace to use as input data
+ */
+ PeaksWorkspace_sptr loadPeaksWorkspace() const {
+ std::string WSName("peaks");
+ LoadIsawPeaks loader;
+ TS_ASSERT_THROWS_NOTHING(loader.initialize());
+ TS_ASSERT(loader.isInitialized());
+ loader.setPropertyValue("Filename", "TOPAZ_3007.peaks");
+ loader.setPropertyValue("OutputWorkspace", WSName);
+
+ TS_ASSERT(loader.execute());
+ TS_ASSERT(loader.isExecuted());
+
+ PeaksWorkspace_sptr ws;
+ TS_ASSERT_THROWS_NOTHING(
+ ws = boost::dynamic_pointer_cast<PeaksWorkspace>(
+ AnalysisDataService::Instance().retrieve(WSName)));
+ TS_ASSERT(ws);
+ return ws;
}
};
diff --git a/Framework/Geometry/inc/MantidGeometry/Crystal/IndexingUtils.h b/Framework/Geometry/inc/MantidGeometry/Crystal/IndexingUtils.h
index 0756868ab09a659e5a216157c9fd53b7a9707052..67cecfc0fabdf6b00798c97808e2cd3cb7c2bf0e 100644
--- a/Framework/Geometry/inc/MantidGeometry/Crystal/IndexingUtils.h
+++ b/Framework/Geometry/inc/MantidGeometry/Crystal/IndexingUtils.h
@@ -54,7 +54,8 @@ public:
const std::vector<Kernel::V3D> &q_vectors, double a,
double b, double c, double alpha, double beta,
double gamma, double required_tolerance, int base_index,
- size_t num_initial, double degrees_per_step);
+ size_t num_initial, double degrees_per_step,
+ bool fixAll = false);
/// Find the UB matrix that most nearly indexes the specified qxyz values
/// given the range of possible real space unit cell edge lengths.
diff --git a/Framework/Geometry/src/Crystal/IndexingUtils.cpp b/Framework/Geometry/src/Crystal/IndexingUtils.cpp
index 4e423f22dcc8285bd72c4aa4a3d8dc03ecf3cd39..f0e49f6e355502ffa63a200640d0fcbeff4e9dca 100644
--- a/Framework/Geometry/src/Crystal/IndexingUtils.cpp
+++ b/Framework/Geometry/src/Crystal/IndexingUtils.cpp
@@ -90,7 +90,8 @@ double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors,
double a, double b, double c, double alpha,
double beta, double gamma,
double required_tolerance, int base_index,
- size_t num_initial, double degrees_per_step) {
+ size_t num_initial, double degrees_per_step,
+ bool fixAll) {
if (UB.numRows() != 3 || UB.numCols() != 3) {
throw std::invalid_argument("Find_UB(): UB matrix NULL or not 3X3");
}
@@ -159,11 +160,21 @@ double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors,
ScanFor_UB(UB, some_qs, a, b, c, alpha, beta, gamma, degrees_per_step,
required_tolerance);
+
double fit_error = 0;
std::vector<V3D> miller_ind;
std::vector<V3D> indexed_qs;
miller_ind.reserve(q_vectors.size());
indexed_qs.reserve(q_vectors.size());
+
+ // If the user has requested we fix all parameters then return here with the
+ // found UB and do not optimise it.
+ if (fixAll) {
+ GetIndexedPeaks(UB, q_vectors, required_tolerance, miller_ind, indexed_qs,
+ fit_error);
+ return fit_error;
+ }
+
// now gradually bring in the remaining
// peaks and re-optimize the UB to index
// them as well
@@ -180,10 +191,11 @@ double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors,
some_qs.push_back(sorted_qs[i]);
try {
+ std::vector<double> sigabc(7);
GetIndexedPeaks(UB, some_qs, required_tolerance, miller_ind, indexed_qs,
fit_error);
Matrix<double> temp_UB(3, 3, false);
- fit_error = Optimize_UB(temp_UB, miller_ind, indexed_qs);
+ fit_error = Optimize_UB(temp_UB, miller_ind, indexed_qs, sigabc);
UB = temp_UB;
} catch (...) {
// failed to fit using these peaks, so add some more and try again
@@ -193,10 +205,11 @@ double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors,
if (q_vectors.size() >= 3) // try one last refinement using all peaks
{
try {
+ std::vector<double> sigabc(7);
GetIndexedPeaks(UB, q_vectors, required_tolerance, miller_ind, indexed_qs,
fit_error);
Matrix<double> temp_UB(3, 3, false);
- fit_error = Optimize_UB(temp_UB, miller_ind, indexed_qs);
+ fit_error = Optimize_UB(temp_UB, miller_ind, indexed_qs, sigabc);
UB = temp_UB;
} catch (...) {
// failed to improve UB using these peaks, so just return the current UB
@@ -583,7 +596,9 @@ double IndexingUtils::Optimize_UB(DblMatrix &UB,
const std::vector<V3D> &hkl_vectors,
const std::vector<V3D> &q_vectors,
std::vector<double> &sigabc) {
- double result = Optimize_UB(UB, hkl_vectors, q_vectors);
+ double result = 0;
+ result = Optimize_UB(UB, hkl_vectors, q_vectors);
+
if (sigabc.size() < 6) {
sigabc.clear();
return result;
@@ -611,11 +626,10 @@ double IndexingUtils::Optimize_UB(DblMatrix &UB,
for (int c = 0; c < 3; c++) {
UB[r][c] += SMALL;
-
GetLatticeParameters(UB, latNew);
UB[r][c] -= SMALL;
- for (int l = 0; l < 7; l++)
+ for (size_t l = 0; l < 7; l++)
derivs[c][l] = (latNew[l] - latOrig[l]) / SMALL;
}