diff --git a/Framework/Crystal/src/FindUBUsingFFT.cpp b/Framework/Crystal/src/FindUBUsingFFT.cpp index 1c0d3bd25feb7fb90cb7f9469552203a19664bf2..29789fa0b20e0997df05bfb5a74742f1b492d449 100644 --- a/Framework/Crystal/src/FindUBUsingFFT.cpp +++ b/Framework/Crystal/src/FindUBUsingFFT.cpp @@ -40,6 +40,11 @@ void FindUBUsingFFT::init() { "Upper Bound on Lattice Parameters a, b, c"); this->declareProperty("Tolerance", 0.15, mustBePositive, "Indexing Tolerance (0.15)"); + this->declareProperty("Iterations", 4, "Iterations to refine UB (4)"); + this->declareProperty("DegreesPerStep", 1.5, + "The resolution of the search through possible " + "orientations is specified by this parameter. One to " + "two degrees per step is usually adequate."); } /** Execute the algorithm. @@ -48,8 +53,9 @@ void FindUBUsingFFT::exec() { double min_d = this->getProperty("MinD"); double max_d = this->getProperty("MaxD"); double tolerance = this->getProperty("Tolerance"); + int iterations = this->getProperty("Iterations"); - double degrees_per_step = 1.5; + double degrees_per_step = this->getProperty("DegreesPerStep"); PeaksWorkspace_sptr ws = this->getProperty("PeaksWorkspace"); @@ -63,7 +69,7 @@ void FindUBUsingFFT::exec() { Matrix<double> UB(3, 3, false); double error = IndexingUtils::Find_UB(UB, q_vectors, min_d, max_d, tolerance, - degrees_per_step); + degrees_per_step, iterations); g_log.notice() << "Error = " << error << '\n'; g_log.notice() << "UB = " << UB << '\n'; diff --git a/Framework/Crystal/src/FindUBUsingLatticeParameters.cpp b/Framework/Crystal/src/FindUBUsingLatticeParameters.cpp index eaf220109897dbfb46c5260b8855abd96ab5840d..d7708ab5e4904ff855ebdfde72104fc37caa3c5b 100644 --- a/Framework/Crystal/src/FindUBUsingLatticeParameters.cpp +++ b/Framework/Crystal/src/FindUBUsingLatticeParameters.cpp @@ -48,6 +48,7 @@ void FindUBUsingLatticeParameters::init() { "Indexing Tolerance (0.15)"); this->declareProperty("FixParameters", false, "Do not optimise the UB after finding the orientation"); + this->declareProperty("Iterations", 1, "Iterations to refine UB (1)"); } /** Execute the algorithm. @@ -62,6 +63,7 @@ void FindUBUsingLatticeParameters::exec() { int num_initial = this->getProperty("NumInitial"); double tolerance = this->getProperty("Tolerance"); auto fixAll = this->getProperty("FixParameters"); + auto iterations = this->getProperty("Iterations"); int base_index = -1; // these "could" be properties if need be double degrees_per_step = 1.5; @@ -83,7 +85,7 @@ void FindUBUsingLatticeParameters::exec() { OrientedLattice lattice(a, b, c, alpha, beta, gamma); double error = IndexingUtils::Find_UB(UB, q_vectors, lattice, tolerance, base_index, - num_initial, degrees_per_step, fixAll); + num_initial, degrees_per_step, fixAll, iterations); g_log.notice() << "Error = " << error << '\n'; g_log.notice() << "UB = " << UB << '\n'; diff --git a/Framework/Geometry/inc/MantidGeometry/Crystal/IndexingUtils.h b/Framework/Geometry/inc/MantidGeometry/Crystal/IndexingUtils.h index 59c29fefbe57334c428fab37c1077e56a7fce80f..76fc2b2622fc89d094e19daae6bdf9594b1967c1 100644 --- a/Framework/Geometry/inc/MantidGeometry/Crystal/IndexingUtils.h +++ b/Framework/Geometry/inc/MantidGeometry/Crystal/IndexingUtils.h @@ -55,7 +55,8 @@ public: const std::vector<Kernel::V3D> &q_vectors, OrientedLattice &lattice, double required_tolerance, int base_index, size_t num_initial, - double degrees_per_step, bool fixAll = false); + double degrees_per_step, bool fixAll = false, + int iterations = 1); /// Find the UB matrix that most nearly indexes the specified qxyz values /// given the range of possible real space unit cell edge lengths. @@ -69,7 +70,7 @@ public: static double Find_UB(Kernel::DblMatrix &UB, const std::vector<Kernel::V3D> &q_vectors, double min_d, double max_d, double required_tolerance, - double degrees_per_step); + double degrees_per_step, int iterations = 4); /// Find the UB matrix that most nearly maps hkl to qxyz for 3 or more peaks static double Optimize_UB(Kernel::DblMatrix &UB, diff --git a/Framework/Geometry/src/Crystal/IndexingUtils.cpp b/Framework/Geometry/src/Crystal/IndexingUtils.cpp index 5dda7228290d211940a803b0a3ac4af45d36f2ad..5e4fb0fbe8c6362c2f52fdf1dc3f5d707140ad8d 100644 --- a/Framework/Geometry/src/Crystal/IndexingUtils.cpp +++ b/Framework/Geometry/src/Crystal/IndexingUtils.cpp @@ -77,6 +77,7 @@ const constexpr double RAD_TO_DEG = 180. / M_PI; orientations used during the initial scan. @param fixAll Fix the lattice parameters and do not optimise the UB matrix. + @param iterations Number of refinements of UB @return This will return the sum of the squares of the residual errors. @@ -90,7 +91,7 @@ double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors, OrientedLattice &lattice, double required_tolerance, int base_index, size_t num_initial, double degrees_per_step, - bool fixAll) { + bool fixAll, int iterations) { if (UB.numRows() != 3 || UB.numCols() != 3) { throw std::invalid_argument("Find_UB(): UB matrix NULL or not 3X3"); } @@ -176,15 +177,16 @@ double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors, for (size_t i = some_qs.size(); i < num_initial; i++) some_qs.push_back(sorted_qs[i]); - - try { - 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); - UB = temp_UB; - } catch (...) { - // failed to fit using these peaks, so add some more and try again + for (int counter = 0; counter < iterations; counter++) { + try { + 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); + UB = temp_UB; + } catch (...) { + // failed to fit using these peaks, so add some more and try again + } } } @@ -192,14 +194,16 @@ double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors, if (!fixAll && q_vectors.size() >= 3) // try one last refinement using all peaks { - try { - GetIndexedPeaks(UB, q_vectors, required_tolerance, miller_ind, indexed_qs, - fit_error); - Matrix<double> temp_UB = UB; - 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 + for (int counter = 0; counter < iterations; counter++) { + try { + GetIndexedPeaks(UB, q_vectors, required_tolerance, miller_ind, + indexed_qs, fit_error); + Matrix<double> temp_UB = UB; + 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 + } } } // Regardless of how we got the UB, find the @@ -448,6 +452,7 @@ double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors, from integer values for a peak to be indexed. @param degrees_per_step The number of degrees between different orientations used during the initial scan. + @param iterations Number of refinements of UB @return This will return the sum of the squares of the residual errors. @@ -464,7 +469,7 @@ double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors, double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors, double min_d, double max_d, double required_tolerance, - double degrees_per_step) { + double degrees_per_step, int iterations) { if (UB.numRows() != 3 || UB.numCols() != 3) { throw std::invalid_argument("Find_UB(): UB matrix NULL or not 3X3"); } @@ -529,7 +534,7 @@ double IndexingUtils::Find_UB(DblMatrix &UB, const std::vector<V3D> &q_vectors, GetIndexedPeaks(UB, q_vectors, required_tolerance, miller_ind, indexed_qs, fit_error); - for (int counter = 0; counter < 4; counter++) { + for (int counter = 0; counter < iterations; counter++) { try { fit_error = Optimize_UB(temp_UB, miller_ind, indexed_qs); UB = temp_UB; @@ -1343,7 +1348,7 @@ size_t IndexingUtils::FFTScanFor_Directions(std::vector<V3D> &directions, for (auto &temp_dir : temp_dirs) { current_dir = temp_dir; double length = current_dir.norm(); - if (length >= min_d && length <= max_d) + if (length >= 0.8 * min_d && length <= 1.2 * max_d) temp_dirs_2.push_back(current_dir); } // only keep directions that index at @@ -1508,9 +1513,7 @@ double IndexingUtils::GetFirstMaxIndex(const double magnitude_fft[], size_t N, bool IndexingUtils::FormUB_From_abc_Vectors(DblMatrix &UB, const std::vector<V3D> &directions, size_t a_index, double min_d, - double max_d) - -{ + double max_d) { if (UB.numRows() != 3 || UB.numCols() != 3) { throw std::invalid_argument("Find_UB(): UB matrix NULL or not 3X3"); } @@ -1521,7 +1524,7 @@ bool IndexingUtils::FormUB_From_abc_Vectors(DblMatrix &UB, // the possible range of d-values // implies a bound on the minimum // angle between a,b, c vectors. - double min_deg = (RAD_TO_DEG)*atan(2.0 * min_d / max_d); + double min_deg = (RAD_TO_DEG)*atan(2.0 * std::min(0.2, min_d / max_d)); double epsilon = 5; // tolerance on right angle (degrees) V3D b_dir; diff --git a/Framework/Geometry/test/IndexingUtilsTest.h b/Framework/Geometry/test/IndexingUtilsTest.h index 585d0fbd180f97eebfac603d84739002e6634741..200bd61e8ec9d09bda9f2c7bc6f25e2f5adc4104 100644 --- a/Framework/Geometry/test/IndexingUtilsTest.h +++ b/Framework/Geometry/test/IndexingUtilsTest.h @@ -314,7 +314,7 @@ public: IndexingUtils::FFTScanFor_Directions(directions, q_vectors, d_min, d_max, required_tolerance, degrees_per_step); - TS_ASSERT_EQUALS(5, directions.size()); + TS_ASSERT_EQUALS(8, directions.size()); for (size_t i = 0; i < 3; i++) { V3D vec = directions[i]; diff --git a/docs/source/release/v3.12.0/diffraction.rst b/docs/source/release/v3.12.0/diffraction.rst index 5affcb02f42fd99b4c42b2dbf91baa6c567be595..5f3832109ec2c2e95e88eaf2612b36df93c39c32 100644 --- a/docs/source/release/v3.12.0/diffraction.rst +++ b/docs/source/release/v3.12.0/diffraction.rst @@ -45,6 +45,10 @@ Single Crystal Diffraction - :ref:`IntegratePeaksMDHKL <algm-IntegratePeaksMDHKL>` now has option to specify background shell instead of using default background determination. +- :ref:`FindUBUsingFFT <algm-FindUBUsingFFT>` now has options to specify number of iterations to refine UB and also resolution of the search through possible orientations. Minimum angle between a,b,c vectors reduced for large unit cells. + +- :ref:`FindUBUsingLatticeParameters <algm-FindUBUsingLatticeParameters>` now has option to specify number of iterations to refine UB. + Imaging -------