Mantid
Loading...
Searching...
No Matches
FindUBUsingLatticeParameters.cpp
Go to the documentation of this file.
1// Mantid Repository : https://github.com/mantidproject/mantid
2//
3// Copyright © 2018 ISIS Rutherford Appleton Laboratory UKRI,
4// NScD Oak Ridge National Laboratory, European Spallation Source,
5// Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
6// SPDX - License - Identifier: GPL - 3.0 +
9#include "MantidAPI/Sample.h"
15
16namespace Mantid::Crystal {
17// Register the algorithm into the AlgorithmFactory
18DECLARE_ALGORITHM(FindUBUsingLatticeParameters)
19
20using namespace Mantid::Kernel;
21using namespace Mantid::API;
22using namespace Mantid::DataObjects;
23using namespace Mantid::Geometry;
24
28 this->declareProperty(std::make_unique<WorkspaceProperty<IPeaksWorkspace>>("PeaksWorkspace", "", Direction::InOut),
29 "Input Peaks Workspace");
30
31 auto mustBePositive = std::make_shared<BoundedValidator<double>>();
32 mustBePositive->setLower(0.0);
33
34 auto moreThan2Int = std::make_shared<BoundedValidator<int>>();
35 moreThan2Int->setLower(2);
36
37 auto reasonable_angle = std::make_shared<BoundedValidator<double>>();
38 reasonable_angle->setLower(5.0);
39 reasonable_angle->setUpper(175.0);
40
41 // use negative values, force user to input all parameters
42 this->declareProperty("a", -1.0, mustBePositive, "Lattice parameter a");
43 this->declareProperty("b", -1.0, mustBePositive, "Lattice parameter b");
44 this->declareProperty("c", -1.0, mustBePositive, "Lattice parameter c");
45 this->declareProperty("alpha", -1.0, reasonable_angle, "Lattice parameter alpha");
46 this->declareProperty("beta", -1.0, reasonable_angle, "Lattice parameter beta");
47 this->declareProperty("gamma", -1.0, reasonable_angle, "Lattice parameter gamma");
48 this->declareProperty("NumInitial", 15, moreThan2Int, "Number of Peaks to Use on First Pass(15)");
49 this->declareProperty("Tolerance", 0.15, mustBePositive, "Indexing Tolerance (0.15)");
50 this->declareProperty("FixParameters", false, "Do not optimise the UB after finding the orientation");
51 this->declareProperty("Iterations", 1, "Iterations to refine UB (1)");
52}
53
57 double a = this->getProperty("a");
58 double b = this->getProperty("b");
59 double c = this->getProperty("c");
60 double alpha = this->getProperty("alpha");
61 double beta = this->getProperty("beta");
62 double gamma = this->getProperty("gamma");
63 int num_initial = this->getProperty("NumInitial");
64 double tolerance = this->getProperty("Tolerance");
65 auto fixAll = this->getProperty("FixParameters");
66 auto iterations = this->getProperty("Iterations");
67
68 int base_index = -1; // these "could" be properties if need be
69 double degrees_per_step = 1.5;
70
71 IPeaksWorkspace_sptr ws = this->getProperty("PeaksWorkspace");
72 if (!ws)
73 throw std::runtime_error("Could not read the peaks workspace");
74
75 const int n_peaks = ws->getNumberPeaks();
76
77 std::vector<V3D> q_vectors;
78 q_vectors.reserve(n_peaks);
79
80 for (int i = 0; i < n_peaks; i++)
81 q_vectors.emplace_back(ws->getPeak(i).getQSampleFrame());
82
83 Matrix<double> UB(3, 3, false);
84 OrientedLattice lattice(a, b, c, alpha, beta, gamma);
85 double error = IndexingUtils::Find_UB(UB, q_vectors, lattice, tolerance, base_index, num_initial, degrees_per_step,
86 fixAll, iterations);
87
88 g_log.notice() << "Error = " << error << '\n';
89 g_log.notice() << "UB = " << UB << '\n';
90
91 if (!IndexingUtils::CheckUB(UB)) // UB not found correctly
92 {
93 g_log.notice(std::string("Found Invalid UB...peaks used might not be linearly independent"));
94 g_log.notice(std::string("UB NOT SAVED."));
95 } else // tell user how many would be indexed
96 { // and save the UB in the sample
97 char logInfo[200];
98 const int num_indexed = IndexingUtils::NumberIndexed(UB, q_vectors, tolerance);
99 sprintf(logInfo,
100 "New UB will index %1d Peaks out of %1d with tolerance "
101 "%5.3f",
102 num_indexed, n_peaks, tolerance);
103 g_log.notice(std::string(logInfo));
104
105 double calc_a = lattice.a();
106 double calc_b = lattice.b();
107 double calc_c = lattice.c();
108 double calc_alpha = lattice.alpha();
109 double calc_beta = lattice.beta();
110 double calc_gamma = lattice.gamma();
111 // Show the modified lattice parameters
112 g_log.notice() << lattice << "\n";
113
114 sprintf(logInfo,
115 "Lattice Parameters (Refined - Input): %11.6f "
116 "%11.6f %11.6f %11.6f %11.6f %11.6f",
117 calc_a - a, calc_b - b, calc_c - c, calc_alpha - alpha, calc_beta - beta, calc_gamma - gamma);
118 g_log.notice(std::string(logInfo));
119 ws->mutableSample().setOrientedLattice(std::make_unique<OrientedLattice>(lattice));
120 }
121}
122
123} // namespace Mantid::Crystal
#define DECLARE_ALGORITHM(classname)
Definition: Algorithm.h:576
double error
Definition: IndexPeaks.cpp:133
double tolerance
void declareProperty(std::unique_ptr< Kernel::Property > p, const std::string &doc="") override
Add a property to the list of managed properties.
Definition: Algorithm.cpp:1913
TypedValue getProperty(const std::string &name) const override
Get the value of a property.
Definition: Algorithm.cpp:2076
Kernel::Logger & g_log
Definition: Algorithm.h:451
A property class for workspaces.
void init() override
Initialise the properties.
static double Find_UB(Kernel::DblMatrix &UB, 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, int iterations=1)
Find the UB matrix that most nearly indexes the specified qxyz values given the lattice parameters.
static bool CheckUB(const Kernel::DblMatrix &UB)
Check that the specified UB is reasonable for an orientation matrix.
static int NumberIndexed(const Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &q_vectors, double tolerance)
Calculate the number of Q vectors that are mapped to integer indices by UB.
Class to implement UB matrix.
double alpha() const
Get lattice parameter.
Definition: UnitCell.cpp:133
double a(int nd) const
Get lattice parameter a1-a3 as function of index (0-2)
Definition: UnitCell.cpp:94
double c() const
Get lattice parameter.
Definition: UnitCell.cpp:128
double beta() const
Get lattice parameter.
Definition: UnitCell.cpp:138
double b() const
Get lattice parameter.
Definition: UnitCell.cpp:123
double gamma() const
Get lattice parameter.
Definition: UnitCell.cpp:143
void notice(const std::string &msg)
Logs at notice level.
Definition: Logger.cpp:95
Numerical Matrix class.
Definition: Matrix.h:42
std::shared_ptr< IPeaksWorkspace > IPeaksWorkspace_sptr
shared pointer to Mantid::API::IPeaksWorkspace
@ InOut
Both an input & output workspace.
Definition: Property.h:55