Mantid
|
This class contains static utility methods for indexing peaks and finding the UB matrix. More...
#include <IndexingUtils.h>
Static Public Member Functions | |
static Kernel::V3D | CalculateMillerIndices (const Kernel::DblMatrix &inverseUB, const Kernel::V3D &q_vector) |
Given a UB, calculate the miller indices for given q vector. More... | |
static bool | CalculateMillerIndices (const Kernel::DblMatrix &inverseUB, const Kernel::V3D &q_vector, double tolerance, Kernel::V3D &miller_indices) |
Given a UB, calculate the miller indices for given q vector. More... | |
static int | CalculateMillerIndices (const Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &q_vectors, double tolerance, std::vector< Kernel::V3D > &miller_indices, double &ave_error) |
Given a UB, get list of Miller indices for specifed Qs and tolerance. More... | |
static bool | CheckUB (const Kernel::DblMatrix &UB) |
Check that the specified UB is reasonable for an orientation matrix. More... | |
static void | DiscardDuplicates (std::vector< Kernel::V3D > &new_list, std::vector< Kernel::V3D > &directions, const std::vector< Kernel::V3D > &q_vectors, double required_tolerance, double len_tol, double ang_tol) |
Construct a sublist of the specified list of a,b,c directions, by removing all directions that seem to be duplicates. More... | |
static size_t | FFTScanFor_Directions (std::vector< Kernel::V3D > &directions, const std::vector< Kernel::V3D > &q_vectors, double min_d, double max_d, double required_tolerance, double degrees_per_step) |
Use FFT to get list of possible directions and lengths for real space unit cell. More... | |
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, int iterations=4) |
Find the UB matrix that most nearly indexes the specified qxyz values using FFTs, given the range of possible real space unit cell edge lengths. More... | |
static double | Find_UB (Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &q_vectors, double min_d, double max_d, double required_tolerance, int base_index, size_t num_initial, double degrees_per_step) |
Find the UB matrix that most nearly indexes the specified qxyz values given the range of possible real space unit cell edge lengths. More... | |
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. More... | |
static bool | FormUB_From_abc_Vectors (Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &directions, const std::vector< Kernel::V3D > &q_vectors, double req_tolerance, double min_vol) |
Form a UB matrix by choosing three vectors from list of possible a,b,c to maximize the number of peaks indexed and minimize cell volume. More... | |
static bool | FormUB_From_abc_Vectors (Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &directions, size_t a_index, double min_d, double max_d) |
Try to form a UB matrix from three vectors chosen from list of possible a,b,c directions, starting with the a vector at the given index. More... | |
static double | GetFirstMaxIndex (const double magnitude_fft[], size_t N, double threshold) |
Get the location of the first maximum (beyond the DC term) in the |FFT| that exceeds the specified threshold. More... | |
static int | GetIndexedPeaks (const Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &q_vectors, double required_tolerance, std::vector< Kernel::V3D > &miller_indices, std::vector< Kernel::V3D > &indexed_qs, double &fit_error) |
Get lists of indices and Qs for peaks indexed by the specified UB matrix. More... | |
static int | GetIndexedPeaks_1D (const Kernel::V3D &direction, const std::vector< Kernel::V3D > &q_vectors, double required_tolerance, std::vector< int > &index_vals, std::vector< Kernel::V3D > &indexed_qs, double &fit_error) |
Get lists of indices and Qs for peaks indexed in the specified direction. More... | |
static int | GetIndexedPeaks_3D (const Kernel::V3D &direction_1, const Kernel::V3D &direction_2, const Kernel::V3D &direction_3, const std::vector< Kernel::V3D > &q_vectors, double required_tolerance, std::vector< Kernel::V3D > &miller_indices, std::vector< Kernel::V3D > &indexed_qs, double &fit_error) |
Get lists of indices and Qs for peaks indexed in three given directions. More... | |
static bool | GetLatticeParameters (const Kernel::DblMatrix &UB, std::vector< double > &lattice_par) |
Get the lattice parameters for the specified orientation matrix. More... | |
static std::string | GetLatticeParameterString (const Kernel::DblMatrix &UB) |
Get a formatted string listing the lattice parameters and cell volume. More... | |
static double | GetMagFFT (const std::vector< Kernel::V3D > &q_vectors, const Kernel::V3D ¤t_dir, const size_t N, double projections[], double index_factor, double magnitude_fft[]) |
Get the magnitude of the FFT of the projections of the q_vectors on the current direction vector. More... | |
static bool | GetModulationVector (const Kernel::DblMatrix &UB, const Kernel::DblMatrix &ModUB, Kernel::V3D &ModVec, const int j) |
static int | GetModulationVectors (const Kernel::DblMatrix &UB, const Kernel::DblMatrix &ModUB, Kernel::V3D &ModVec1, Kernel::V3D &ModVec2, Kernel::V3D &ModVec3) |
static double | IndexingError (const Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &hkls, const std::vector< Kernel::V3D > &q_vectors) |
Find the average indexing error for UB with the specified q's and hkls. More... | |
static Kernel::V3D | makeCDir (const Kernel::V3D &a_dir, const Kernel::V3D &b_dir, const double c, const double cosAlpha, const double cosBeta, const double cosGamma, const double sinGamma) |
Get the vector in the direction of "c" given other unit cell information. More... | |
static std::vector< Kernel::V3D > | MakeCircleDirections (int n_steps, const Kernel::V3D &axis, double angle_degrees) |
Make list of the circle of direction vectors that form a fixed angle with the specified axis. More... | |
static std::vector< Kernel::V3D > | MakeHemisphereDirections (int n_steps) |
Make list of direction vectors uniformly distributed over a hemisphere. More... | |
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. More... | |
static int | NumberIndexed_1D (const Kernel::V3D &direction, const std::vector< Kernel::V3D > &q_vectors, double tolerance) |
Calculate the number of Q vectors that map to an integer index in one direction. More... | |
static int | NumberIndexed_3D (const Kernel::V3D &a_dir, const Kernel::V3D &b_dir, const Kernel::V3D &c_dir, const std::vector< Kernel::V3D > &q_vectors, double tolerance) |
Calculate the number of Q vectors that map to integer indices simutlaneously in three directions. More... | |
static int | NumberOfValidIndexes (const std::vector< Kernel::V3D > &hkls, double tolerance, double &average_error) |
Find number of valid HKLs and average error, in list of HKLs. More... | |
static double | Optimize_6dUB (Kernel::DblMatrix &UB, Kernel::DblMatrix &ModUB, const std::vector< Kernel::V3D > &hkl_vectors, const std::vector< Kernel::V3D > &mnp_vectors, const int &ModDim, const std::vector< Kernel::V3D > &q_vectors) |
STATIC method Optimize_6dUB: Calculates the 6-dimensional matrix that most nearly maps the specified hkl_vectors and mnp_vectors to the specified q_vectors. More... | |
static double | Optimize_6dUB (Kernel::DblMatrix &UB, Kernel::DblMatrix &ModUB, const std::vector< Kernel::V3D > &hkl_vectors, const std::vector< Kernel::V3D > &mnp_vectors, const int &ModDim, const std::vector< Kernel::V3D > &q_vectors, std::vector< double > &sigabc, std::vector< double > &sigq) |
STATIC method Optimize_UB: Calculates the matrix that most nearly maps the specified hkl_vectors to the specified q_vectors. More... | |
static double | Optimize_Direction (Kernel::V3D &best_vec, const std::vector< int > &index_values, const std::vector< Kernel::V3D > &q_vectors) |
Find the vector that best corresponds to plane normal, given 1-D indices. More... | |
static double | Optimize_UB (Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &hkl_vectors, const std::vector< Kernel::V3D > &q_vectors) |
Find the UB matrix that most nearly maps hkl to qxyz for 3 or more peaks. More... | |
static double | Optimize_UB (Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &hkl_vectors, const std::vector< Kernel::V3D > &q_vectors, std::vector< double > &sigabc) |
Find the UB matrix that most nearly maps hkl to qxyz for 3 or more peaks. More... | |
static void | RoundHKL (Kernel::V3D &hkl) |
Round all the components of a HKL objects to the nearest integer. More... | |
static void | RoundHKLs (std::vector< Kernel::V3D > &hkl_list) |
Round all the components of a list of V3D objects, to the nearest integer. More... | |
static size_t | ScanFor_Directions (std::vector< Kernel::V3D > &directions, const std::vector< Kernel::V3D > &q_vectors, double min_d, double max_d, double required_tolerance, double degrees_per_step) |
Get list of possible directions and lengths for real space unit cell. More... | |
static double | ScanFor_UB (Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &q_vectors, const UnitCell &cell, double degrees_per_step, double required_tolerance) |
Scan rotations to find UB that indexes peaks given lattice parameters. More... | |
static int | SelectDirection (Kernel::V3D &best_direction, const std::vector< Kernel::V3D > &q_vectors, const std::vector< Kernel::V3D > &direction_list, double plane_spacing, double required_tolerance) |
Choose the direction in a list of directions, that is most nearly perpendicular to planes with the specified spacing in reciprocal space. More... | |
static bool | ValidIndex (const Kernel::V3D &hkl, double tolerance) |
Check is hkl is within tolerance of integer (h,k,l) non-zero values. More... | |
Static Private Attributes | |
static Mantid::Kernel::Logger & | g_Log |
Static reference to the logger class. More... | |
This class contains static utility methods for indexing peaks and finding the UB matrix.
Definition at line 32 of file IndexingUtils.h.
|
static |
Given a UB, calculate the miller indices for given q vector.
Calculate the Miller Indices for the specified Q vector, using the inverse of the specified UB matrix.
inverseUB | A 3x3 matrix of doubles holding the inverse UB matrix. The matrix is not checked for validity |
q_vector | V3D object containing Q vector in sample frame |
Definition at line 2366 of file IndexingUtils.cpp.
|
static |
Given a UB, calculate the miller indices for given q vector.
Calculate the Miller Indices for the specified Q vector, using the inverse of the specified UB matrix.
If the peak could not be indexed it is set to (0,0,0)
inverseUB | A 3x3 matrix of doubles holding the inverse UB matrix. The matrix is not checked for validity |
q_vector | std::vector of V3D objects that contains the list of q_vectors that are to be indexed. |
tolerance | The maximum allowed distance between each component of UB^(-1)*Q and the nearest integer value, required to to count the peak as indexed by UB. |
miller_indices | This vector returns a list of Miller Indices, with one entry for each given Q vector. |
Definition at line 2344 of file IndexingUtils.cpp.
References CalculateMillerIndices(), tolerance, and ValidIndex().
|
static |
Given a UB, get list of Miller indices for specifed Qs and tolerance.
Calculate the Miller Indices for each of the specified Q vectors, using the inverse of the specified UB matrix.
The Miller Indices will be set to 0, 0, 0 for any peak for which h, k or l differs from an intenger by more than the specified tolerance. Also (h,k,l) = (0,0,0) the peak will NOT be counted as indexed, since (0,0,0) is not a valid index of any peak.
UB | A 3x3 matrix of doubles holding the UB matrix |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are to be indexed. |
tolerance | The maximum allowed distance between each component of UB^(-1)*Q and the nearest integer value, required to to count the peak as indexed by UB. |
miller_indices | This vector returns a list of Miller Indices, with one entry for each given Q vector. |
ave_error | The average error from all lattice directions. |
std::invalid_argument | exception if the UB matrix is not a 3X3 matrix, or if UB is singular. |
Definition at line 2294 of file IndexingUtils.cpp.
References CalculateMillerIndices(), CheckUB(), count, Mantid::Kernel::V3D::hklError(), Mantid::Kernel::Matrix< T >::Invert(), and tolerance.
Referenced by Mantid::Crystal::SelectCellWithForm::ApplyTransform(), CalculateMillerIndices(), Mantid::MDAlgorithms::ConvertCWSDMDtoHKL::convertFromQSampleToHKL(), and Mantid::Crystal::FindUBUsingIndexedPeaks::exec().
|
static |
Check that the specified UB is reasonable for an orientation matrix.
Check whether or not the specified matrix is reasonable for an orientation matrix.
In particular, check that it is a 3x3 matrix without any nan or infinite values and that its determinant is within a reasonable range, for an orientation matrix.
UB | A 3x3 matrix of doubles holding the UB matrix |
Definition at line 2136 of file IndexingUtils.cpp.
References Mantid::Kernel::Matrix< T >::determinant(), fabs, Mantid::Kernel::Matrix< T >::numCols(), and Mantid::Kernel::Matrix< T >::numRows().
Referenced by CalculateMillerIndices(), Mantid::Crystal::CalculatePeaksHKL::exec(), Mantid::Crystal::FindUBUsingFFT::exec(), Mantid::Crystal::FindUBUsingIndexedPeaks::exec(), Mantid::Crystal::FindUBUsingLatticeParameters::exec(), Mantid::Crystal::FindUBUsingMinMaxD::exec(), Mantid::Crystal::SelectCellOfType::exec(), Mantid::Crystal::SelectCellWithForm::exec(), Mantid::Crystal::ShowPossibleCells::exec(), Mantid::Crystal::TransformHKL::exec(), GetIndexedPeaks(), IndexingError(), NumberIndexed(), Optimize_6dUB(), and Optimize_UB().
|
static |
Construct a sublist of the specified list of a,b,c directions, by removing all directions that seem to be duplicates.
If several directions all have the same length (within the specified length tolerance) and have the same direction (within the specified angle tolerange) then only one of those directions will be recorded in the sublist. The one that indexes the most peaks, within the specified tolerance will be kept.
new_list | This vector will be cleared and filled with the vectors from the directions list that are not duplicates of other vectors in the list. |
directions | Input list of possible a,b,c directions, listed in order of increasing vector norm. This list will be cleared by this method. |
q_vectors | List of q_vectors that should be indexed |
required_tolerance | The tolerance for indexing |
len_tol | The tolerance on the relative difference in length for two directions to be considered equal. Eg. if relative differences must be less than 5% for two lengths to be considered the same, pass in .05 for the len_tol. |
ang_tol | The tolerance for the difference in directions, specified in degrees. |
Definition at line 1925 of file IndexingUtils.cpp.
References Mantid::Kernel::V3D::angle(), fabs, Mantid::Kernel::V3D::norm(), and NumberIndexed_1D().
Referenced by FFTScanFor_Directions().
|
static |
Use FFT to get list of possible directions and lengths for real space unit cell.
Get list of possible edge vectors for the real space unit cell.
This method uses FFTs to find directions for which projections of the peaks on those directions have repetitive patterns. This list of directions found will consist of vectors, V, for which V dot Q is essentially an integer for the most Q vectors. The difference between V dot Q and an integer must be less than the required tolerance for it to count as an integer.
directions | Vector that will be filled with the directions that may correspond to unit cell edges. |
q_vectors | Vector of new Vector3D objects that contains the list of q_vectors that are to be indexed. |
min_d | Lower bound on shortest unit cell edge length. This does not have to be specified exactly but must be strictly less than the smallest edge length, in Angstroms. |
max_d | Upper bound on longest unit cell edge length. This does not have to be specified exactly but must be strictly more than the longest edge length in angstroms. |
required_tolerance | The maximum allowed deviation of Miller indices from integer values for a peak to be indexed. |
degrees_per_step | The number of degrees between directions that are checked while scanning for an initial indexing of the peaks with lowest |Q|. |
Definition at line 1415 of file IndexingUtils.cpp.
References Mantid::Kernel::V3D::compareMagnitude(), count, DiscardDuplicates(), GetFirstMaxIndex(), GetIndexedPeaks_1D(), GetMagFFT(), index, MakeHemisphereDirections(), Mantid::Kernel::V3D::norm(), NumberIndexed_1D(), Optimize_Direction(), and position.
Referenced by Find_UB().
|
static |
Find the UB matrix that most nearly indexes the specified qxyz values using FFTs, given the range of possible real space unit cell edge lengths.
STATIC method Find_UB: This method will attempt to calculate the matrix that most nearly indexes the specified q_vectors, using FFTs to find patterns in projected Q-vectors, given only a range of possible unit cell edge lengths.
If successful, the resulting matrix should correspond to the Niggli reduced cell. The resolution of the search through possible orientations is specified by the degrees_per_step parameter. One to two degrees per step is usually adequate. NOTE: The execution time is O(n^3) where n is the number of degrees per step, so decreasing the resolution to 0.5 degree per step will take about 8 times longer than using 1 degree per step. It should not be necessary to decrease this value below 1 degree per step, and users will have to be VERY patient, if it is decreased much below 1 degree per step. The specified q_vectors should correspond to a single crystal, for this to work reliably.
UB | 3x3 matrix that will be set to the UB matrix |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are to be indexed NOTE: There must be at least 4 q_vectors and it really should have at least 10 or more peaks for this to work quite consistently. |
min_d | Lower bound on shortest unit cell edge length. This does not have to be specified exactly but must be strictly less than the smallest edge length, in Angstroms. |
max_d | Upper bound on longest unit cell edge length. This does not have to be specified exactly but must be strictly more than the longest edge length in angstroms. |
required_tolerance | The maximum allowed deviation of Miller indices from integer values for a peak to be indexed. |
degrees_per_step | The number of degrees between different orientations used during the initial scan. |
iterations | Number of refinements of UB |
std::invalid_argument | exception if UB is not a 3X3 matrix, if there are not at least 3 q vectors, if min_d >= max_d or min_d <= 0, if the required_tolerance or degrees_per_step is <= 0, if at least three possible a,b,c directions were not found, or if a valid UB matrix could not be formed from the a,b,c directions that were found. |
Definition at line 451 of file IndexingUtils.cpp.
References Mantid::Kernel::V3D::compareMagnitude(), FFTScanFor_Directions(), FormUB_From_abc_Vectors(), GetIndexedPeaks(), Mantid::Geometry::NiggliCell::MakeNiggliUB(), Mantid::Kernel::Matrix< T >::numCols(), Mantid::Kernel::Matrix< T >::numRows(), and Optimize_UB().
|
static |
Find the UB matrix that most nearly indexes the specified qxyz values given the range of possible real space unit cell edge lengths.
STATIC method Find_UB: This method will attempt to calculate the matrix that most nearly indexes the specified q_vectors, given only a range of possible unit cell edge lengths.
If successful, the matrix should correspond to the Niggli reduced cell. The resolution of the search through possible orientations is specified by the degrees_per_step parameter. Approximately 1-3 degrees_per_step is usually adequate. NOTE: This is an expensive calculation which takes approximately 1 second using 1 degree_per_step. However, the execution time is O(n^3) so decreasing the resolution to 0.5 degree per step will take about 8 seconds, etc. It should not be necessary to decrease this value below 1 degree per step, and users will have to be VERY patient, if it is decreased much below 1 degree per step. The number of peaks used to obtain an initial indexing is specified by the "NumInitial" parameter. Good values for this are typically around 15-25. The specified q_vectors must correspond to a single crystal. If several crystallites are present or there are other sources of "noise" leading to invalid peaks, this method will not work well. The method that uses lattice parameters may be better in such cases. Alternatively, adjust the list of specified q_vectors so it does not include noise peaks or peaks from more than one crystal, by increasing the threshold for what counts as a peak, or by other methods.
UB | 3x3 matrix that will be set to the UB matrix |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are to be indexed NOTE: There must be at least 3 q_vectors. |
min_d | Lower bound on shortest unit cell edge length. This does not have to be specified exactly but must be strictly less than the smallest edge length, in Angstroms. |
max_d | Upper bound on longest unit cell edge length. This does not have to be specified exactly but must be strictly more than the longest edge length in angstroms. |
required_tolerance | The maximum allowed deviation of Miller indices from integer values for a peak to be indexed. |
base_index | The sequence number of the peak that should be used as the central peak. On the first scan for a UB matrix that fits the data, the remaining peaks in the list of q_vectors will be shifted by -base_peak, where base_peak is the q_vector with the specified base index. If fewer than 6 peaks are specified in the q_vectors list, this parameter is ignored. If this parameter is -1, and there are at least five peaks in the q_vector list, then a base index will be calculated internally. In most cases, it should suffice to set this to -1. |
num_initial | The number of low |Q| peaks that should be used to scan for an initial orientation matrix. |
degrees_per_step | The number of degrees between different orientations used during the initial scan. |
std::invalid_argument | exception if UB is not a 3X3 matrix, if there are not at least 3 q vectors, if min_d >= max_d or min_d <= 0 if num_initial is < 3, or if the required_tolerance or degrees_per_step is <= 0. |
Definition at line 277 of file IndexingUtils.cpp.
References Mantid::Kernel::V3D::compareMagnitude(), FormUB_From_abc_Vectors(), GetIndexedPeaks(), Mantid::Geometry::NiggliCell::MakeNiggliUB(), Mantid::Kernel::Matrix< T >::numCols(), Mantid::Kernel::Matrix< T >::numRows(), Optimize_UB(), and ScanFor_Directions().
|
static |
Find the UB matrix that most nearly indexes the specified qxyz values given the lattice parameters.
STATIC method Find_UB: Calculates the matrix that most nearly indexes the specified q_vectors, given the lattice parameters.
The sum of the squares of the residual errors is returned. This method first sorts the specified q_vectors in order of increasing magnitude. It then searches through all possible orientations to find an initial UB that indexes the lowest magnitude peaks. The resolution of the search through possible orientations is specified by the degrees_per_step parameter. Approximately 2-4 degrees_per_step is usually adequate. NOTE: This is an expensive calculation which takes approximately 1 second using 2 degrees_per_step. However, the execution time is O(n^3) so decreasing the resolution to 1 degree per step will take about 8 seconds, etc. It should not be necessary to decrease this value below 1 degree per step, and users will have to be VERY patient, if it is decreased much below 1 degree per step. The number of peaks used to obtain an initial indexing is specified by the "NumInitial" parameter. Good values for this are typically around 10-15, though with accurate peak positions, and good values for the lattice paramters, as few as 2 can be used. Using substantially more than 15 peaks initially typically has no benefit and increases execution time.
UB | 3x3 matrix that will be set to the UB matrix |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are to be indexed NOTE: There must be at least 2 q_vectors. |
lattice | The orientated lattice with the lattice parameters a,b,c and alpha, beta, gamma. The found UB and errors will be set on this lattice. |
required_tolerance | The maximum allowed deviation of Miller indices from integer values for a peak to be indexed. |
base_index | The sequence number of the peak that should be used as the central peak. On the first scan for a UB matrix that fits the data, the remaining peaks in the list of q_vectors will be shifted by -base_peak, where base_peak is the q_vector with the specified base index. If fewer than 5 peaks are specified in the q_vectors list, this parameter is ignored. If this parameter is -1, and there are at least four peaks in the q_vector list, then a base index will be calculated internally. In most cases, it should suffice to set this to -1. |
num_initial | The number of low |Q| peaks that should be used to scan for an initial orientation matrix. |
degrees_per_step | The number of degrees between different orientations used during the initial scan. |
fixAll | Fix the lattice parameters and do not optimise the UB matrix. |
iterations | Number of refinements of UB |
std::invalid_argument | exception if UB is not a 3X3 matrix, if there are not at least 2 q vectors, if num_initial is < 2, or if the required_tolerance or degrees_per_step is <= 0. |
Definition at line 94 of file IndexingUtils.cpp.
References Mantid::Kernel::V3D::compareMagnitude(), GetIndexedPeaks(), Mantid::Kernel::Matrix< T >::numCols(), Mantid::Kernel::Matrix< T >::numRows(), Optimize_UB(), ScanFor_UB(), Mantid::Geometry::UnitCell::setError(), and Mantid::Geometry::OrientedLattice::setUB().
Referenced by Mantid::Crystal::FindUBUsingFFT::exec(), Mantid::Crystal::FindUBUsingLatticeParameters::exec(), and Mantid::Crystal::FindUBUsingMinMaxD::exec().
|
static |
Form a UB matrix by choosing three vectors from list of possible a,b,c to maximize the number of peaks indexed and minimize cell volume.
Form a UB matrix from the given list of possible directions, using the three directions that correspond to a unit cell with the smallest volume (greater than or equal to the specified minimum volume) that indexes at least 80% of the maximum number of peaks indexed by any set of three distinct vectors chosen from the list.
UB | The calculated UB matrix will be returned in this parameter |
directions | List of possible vectors for a, b, c. This list MUST be sorted in order of increasing magnitude. |
q_vectors | The list of q_vectors that should be indexed |
req_tolerance | The required tolerance on h,k,l to consider a peak to be indexed. |
min_vol | The smallest possible unit cell volume. |
std::invalid_argument | exception if the UB matrix is not a 3X3 matrix. |
std::runtime_error | exception if the matrix inversion fails and UB can't be formed |
Definition at line 1809 of file IndexingUtils.cpp.
References Mantid::Kernel::V3D::cross_prod(), fabs, Mantid::Geometry::OrientedLattice::GetUB(), NumberIndexed_3D(), Mantid::Kernel::Matrix< T >::numCols(), Mantid::Kernel::Matrix< T >::numRows(), and Mantid::Kernel::V3D::scalar_prod().
|
static |
Try to form a UB matrix from three vectors chosen from list of possible a,b,c directions, starting with the a vector at the given index.
Form a UB matrix from the given list of possible directions, using the direction at the specified index for the "a" direction.
The "b" and "c" directions are chosen so that 1) |a| < |b| < |c|, 2) the angle between the a, b, c, vectors is at least a minimum angle based on min_d/max_d 3) c is not in the same plane as a and b.
UB | The calculated UB matrix will be returned in this parameter |
directions | List of possible vectors for a, b, c. This list MUST be sorted in order of increasing magnitude. |
a_index | The index to use for the a vector. The b and c vectors will be choosen from LATER positions in the directions list. |
min_d | Minimum possible real space unit cell edge length. |
max_d | Maximum possible real space unit cell edge length. |
std::invalid_argument | exception if the UB matrix is not a 3X3 matrix. |
std::runtime_error | exception if the matrix inversion fails and UB can't be formed |
Definition at line 1709 of file IndexingUtils.cpp.
References Mantid::Kernel::V3D::angle(), Mantid::Kernel::V3D::cross_prod(), Mantid::Geometry::OrientedLattice::GetUB(), index, Mantid::Kernel::normalize(), Mantid::Kernel::Matrix< T >::numCols(), and Mantid::Kernel::Matrix< T >::numRows().
Referenced by Find_UB().
|
static |
Get the location of the first maximum (beyond the DC term) in the |FFT| that exceeds the specified threshold.
Scan the FFT array for the first maximum that exceeds the specified threshold and is beyond the initial DC term/interval.
magnitude_fft | The array containing the magnitude of the FFT values. |
N | The size of the FFT array. |
threshold | The required threshold for the first peak. This must be positive. |
Definition at line 1646 of file IndexingUtils.cpp.
Referenced by FFTScanFor_Directions().
|
static |
Get lists of indices and Qs for peaks indexed by the specified UB matrix.
Given a list of peak positions and a UB matrix, get the list of Miller indices and corresponding peak positions for the peaks that are indexed to within a specified tolerance, by the UB matrix.
This method is similar to GetIndexedPeaks_3D, but directly uses the inverse of the UB matrix to map Q -> hkl.
UB | The UB matrix that determines the indexing of the peaks. |
q_vectors | List of V3D peaks in reciprocal space |
required_tolerance | The maximum allowed error (as a faction of the corresponding Miller index) for a peak q_vector to be counted as indexed. |
miller_indices | List of the Miller indices (h,k,l) of peaks that were indexed in all specified directions. |
indexed_qs | List of Qxyz value for the peaks that were indexed indexed in all specified directions. |
fit_error | The sum of the squares of the distances from integer values for the UB*Q for the specified UB matrix and the specified q_vectors. This is a measure of the error in HKL space. |
Definition at line 2534 of file IndexingUtils.cpp.
References CheckUB(), error, Mantid::Kernel::Matrix< T >::Invert(), and ValidIndex().
Referenced by Mantid::Crystal::SelectCellWithForm::DetermineErrors(), Mantid::Crystal::FindUBUsingFFT::exec(), Mantid::Crystal::FindUBUsingMinMaxD::exec(), and Find_UB().
|
static |
Get lists of indices and Qs for peaks indexed in the specified direction.
Given one plane normal direction for a family of parallel planes in reciprocal space, find the peaks that lie on these planes to within the specified tolerance.
The direction is specified as a vector with length "a" if the plane spacing in reciprocal space is 1/a. In that way, the dot product of a peak Qxyz with the direction vector will be an integer if the peak lies on one of the planes.
direction | Direction vector in the direction of the normal vector for a family of parallel planes in reciprocal space. The length of this vector must be the reciprocal of the plane spacing. |
q_vectors | List of V3D peaks in reciprocal space |
required_tolerance | The maximum allowed error (as a faction of the corresponding Miller index) for a peak q_vector to be counted as indexed. |
index_vals | List of the one-dimensional Miller indices peaks that were indexed in the specified direction. |
indexed_qs | List of Qxyz value for the peaks that were indexed indexed in the specified direction. |
fit_error | The sum of the squares of the distances from integer values for the projections of the indexed q_vectors on the specified direction. This is a measure of the error in HKL space. |
Definition at line 2398 of file IndexingUtils.cpp.
References error, fabs, Mantid::Kernel::V3D::norm(), and Mantid::Kernel::V3D::scalar_prod().
Referenced by FFTScanFor_Directions(), and ScanFor_Directions().
|
static |
Get lists of indices and Qs for peaks indexed in three given directions.
Given three plane normal directions for three families of parallel planes in reciprocal space, find the peaks that lie on these planes to within the specified tolerance.
The three directions are specified as vectors with lengths that are the reciprocals of the corresponding plane spacings. In that way, the dot product of a peak Qxyz with one of the direction vectors will be an integer if the peak lies on one of the planes corresponding to that direction. If the three directions are properly chosen to correspond to the unit cell edges, then the resulting indices will be proper Miller indices for the peaks. This method is similar to GetIndexedPeaks_1D, but checks three directions simultaneously and requires that the peak lies on all three families of planes simultaneously and does NOT index as (0,0,0).
direction_1 | Direction vector in the direction of the normal vector for the first family of parallel planes. |
direction_2 | Direction vector in the direction of the normal vector for the second family of parallel planes. |
direction_3 | Direction vector in the direction of the normal vector for the third family of parallel planes. |
q_vectors | List of V3D peaks in reciprocal space |
required_tolerance | The maximum allowed error (as a faction of the corresponding Miller index) for a peak q_vector to be counted as indexed. |
miller_indices | List of the Miller indices (h,k,l) of peaks that were indexed in all specified directions. |
indexed_qs | List of Qxyz value for the peaks that were indexed indexed in all specified directions. |
fit_error | The sum of the squares of the distances from integer values for the projections of the indexed q_vectors on the specified directions. This is a measure of the error in HKL space. |
Definition at line 2464 of file IndexingUtils.cpp.
References fabs, Mantid::Kernel::V3D::scalar_prod(), and ValidIndex().
|
static |
Get the lattice parameters for the specified orientation matrix.
Get the lattice parameters, a, b, c, alpha, beta, gamma and cell volume for the specified orientation matrix.
UB | A non-singular matrix representing an orientation matrix. |
lattice_par | std::vector of doubles that will contain the lattice parameters and cell volume as it's first seven entries. |
Definition at line 2780 of file IndexingUtils.cpp.
References Mantid::Geometry::UnitCell::a(), Mantid::Geometry::UnitCell::alpha(), Mantid::Geometry::UnitCell::b(), Mantid::Geometry::UnitCell::beta(), Mantid::Geometry::UnitCell::c(), Mantid::Geometry::UnitCell::gamma(), Mantid::Geometry::OrientedLattice::setUB(), and Mantid::Geometry::UnitCell::volume().
Referenced by Mantid::Geometry::ConventionalCell::ConventionalCell(), Mantid::Crystal::OptimizeLatticeForCellType::exec(), Mantid::Geometry::ScalarUtils::GetCellForForm(), Mantid::Geometry::ScalarUtils::GetCellsUBOnly(), GetLatticeParameterString(), Mantid::Geometry::ConventionalCell::GetSumOfSides(), Optimize_6dUB(), and Optimize_UB().
|
static |
Get a formatted string listing the lattice parameters and cell volume.
Get a nicely formatted string listing the lattice parameters and cell volume.
Definition at line 2838 of file IndexingUtils.cpp.
References GetLatticeParameters().
Referenced by Mantid::Crystal::SelectCellOfType::exec(), Mantid::Crystal::SelectCellWithForm::exec(), and Mantid::Crystal::ShowPossibleCells::exec().
|
static |
Get the magnitude of the FFT of the projections of the q_vectors on the current direction vector.
Fill an array with the magnitude of the FFT of the projections of the specified q_vectors on the specified direction.
The largest value in the magnitude FFT that occurs at index 5 or more is returned as the value of the function.
q_vectors | The list of Q vectors to project on the specified direction. |
current_dir | The direction the Q vectors will be projected on. |
N | The size of the projections[] array. This MUST BE a power of 2. The magnitude_fft[] array must be half the size. |
projections | Array to hold the projections of the Q vectors. This must be long enough so that all projected values map map to a valid index, after they are multiplied by the index_factor. |
index_factor | Factor that when multiplied by a projected Q vector will give a valid index into the projections array. |
magnitude_fft | Array that will be filled out with the magnitude of the FFT of the projections. |
Definition at line 1600 of file IndexingUtils.cpp.
References fabs, index, and Mantid::Kernel::V3D::scalar_prod().
Referenced by FFTScanFor_Directions().
|
static |
Definition at line 2819 of file IndexingUtils.cpp.
References Mantid::Geometry::UnitCell::getModVec(), Mantid::Geometry::OrientedLattice::setModUB(), and Mantid::Geometry::OrientedLattice::setUB().
|
static |
Definition at line 2798 of file IndexingUtils.cpp.
References Mantid::Geometry::UnitCell::getdh(), Mantid::Geometry::UnitCell::getdk(), Mantid::Geometry::UnitCell::getdl(), Mantid::Geometry::UnitCell::getModVec(), Mantid::Geometry::OrientedLattice::setModUB(), and Mantid::Geometry::OrientedLattice::setUB().
|
static |
Find the average indexing error for UB with the specified q's and hkls.
Definition at line 2096 of file IndexingUtils.cpp.
References CheckUB(), fabs, and Mantid::Kernel::Matrix< T >::Invert().
Referenced by Mantid::Crystal::TransformHKL::exec().
|
static |
Get the vector in the direction of "c" given other unit cell information.
For a rotated unit cell, calculate the vector in the direction of edge "c" given two vectors a_dir and b_dir in the directions of edges "a" and "b", with lengths a and b, and the cell angles.
The calculation is explained in the Mantid document UBMatriximplementationnotes.pdf, pg 3, Andre Savici.
a_dir | V3D object with length "a" in the direction of the rotated cell edge "a" |
b_dir | V3D object with length "b" in the direction of the rotated cell edge "b" |
c | The length of the third cell edge, c. |
cosAlpha | cos angle between edges b and c in radians. |
cosBeta | cos angle between edges c and a in radians. |
cosGamma | cos angle between edges a and b in radians. |
sinGamma | sin angle between edges a and b in radians. |
Definition at line 1884 of file IndexingUtils.cpp.
References Mantid::Kernel::toV3D(), and Mantid::Kernel::toVector3d().
Referenced by ScanFor_UB().
|
static |
Make list of the circle of direction vectors that form a fixed angle with the specified axis.
Make a list of directions, uniformly distributed around a circle, all of which form the specified angle with the specified axis.
n_steps | The number of vectors to generate around the circle. |
axis | The specified axis |
angle_degrees | The specified angle |
std::invalid_argument | exception if the number of steps is <= 0, or if the axix length is 0. |
Definition at line 2645 of file IndexingUtils.cpp.
References Mantid::Kernel::V3D::cross_prod(), fabs, Mantid::Kernel::V3D::normalize(), Mantid::Kernel::normalize(), Mantid::Kernel::Quat::rotate(), and Mantid::Kernel::Quat::setAngleAxis().
Referenced by ScanFor_UB().
|
static |
Make list of direction vectors uniformly distributed over a hemisphere.
Make a list of directions, approximately uniformly distributed over a hemisphere, with the angular separation between direction vectors approximately 90 degrees/n_steps.
NOTE: This method provides a list of possible directions for plane normals for reciprocal lattice planes. This facilitates a brute force search for lattice planes with a specific spacing between planes. This will be used for finding the UB matrix, given the lattice parameters.
n_steps | The number of subdivisions in latitude in the upper hemisphere. |
std::invalid_argument | exception if the number of steps is <= 0. |
Definition at line 2594 of file IndexingUtils.cpp.
References fabs.
Referenced by Mantid::Crystal::GoniometerAnglesFromPhiRotation::exec(), FFTScanFor_Directions(), ScanFor_Directions(), and ScanFor_UB().
|
static |
Calculate the number of Q vectors that are mapped to integer indices by UB.
Calculate the number of Q vectors that are mapped to integer h,k,l values by UB.
Each of the Miller indexes, h, k and l must be within the specified tolerance of an integer, in order to count the peak as indexed. Also, if (h,k,l) = (0,0,0) the peak will NOT be counted as indexed, since (0,0,0) is not a valid index of any peak.
UB | A 3x3 matrix of doubles holding the UB matrix. The UB matrix must not be singular. |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors. |
tolerance | The maximum allowed distance between each component of UB^(-1)*Q and the nearest integer value, required to to count the peak as indexed by UB. |
std::invalid_argument | exception if the UB matrix is not a 3X3 matrix, or if UB is singular. |
Definition at line 2175 of file IndexingUtils.cpp.
References CheckUB(), count, Mantid::Kernel::Matrix< T >::Invert(), tolerance, and ValidIndex().
Referenced by Mantid::Crystal::SelectCellWithForm::DetermineErrors(), Mantid::Crystal::FindUBUsingFFT::exec(), Mantid::Crystal::FindUBUsingIndexedPeaks::exec(), Mantid::Crystal::FindUBUsingLatticeParameters::exec(), and Mantid::Crystal::FindUBUsingMinMaxD::exec().
|
static |
Calculate the number of Q vectors that map to an integer index in one direction.
Calculate the number of Q vectors that are mapped to a integer index value by taking the dot product with the specified direction vector.
The direction vector represents a possible unit cell edge vector in real space. The dot product must be within the specified tolerance of an integer, in order to count as indexed.
direction | A V3D specifying a possible edge vector in real space. |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors. |
tolerance | The maximum allowed distance to an integer from the dot products of peaks with the specified direction. |
Definition at line 2213 of file IndexingUtils.cpp.
References count, error, fabs, Mantid::Kernel::V3D::norm(), Mantid::Kernel::V3D::scalar_prod(), and tolerance.
Referenced by DiscardDuplicates(), and FFTScanFor_Directions().
|
static |
Calculate the number of Q vectors that map to integer indices simutlaneously in three directions.
Calculate the number of Q vectors for which the dot product with three specified direction vectors is an integer triple, NOT equal to (0,0,0) to within the specified tolerance.
This give the number of peaks that would be indexed by the UB matrix formed from the specified those three real space unit cell edge vectors. NOTE: This method assumes that the three edge vectors are linearly independent and could be used to form a valid UB matrix.
a_dir | A Vector3D representing unit cell edge vector a |
b_dir | A Vector3D representing unit cell edge vector b |
c_dir | A Vector3D representing unit cell edge vector c |
q_vectors | Vector of Vector3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors. |
tolerance | The maximum allowed distance to an integer from the dot products of peaks with the specified direction. |
Definition at line 2251 of file IndexingUtils.cpp.
References count, Mantid::Kernel::V3D::norm(), Mantid::Kernel::V3D::scalar_prod(), tolerance, and ValidIndex().
Referenced by FormUB_From_abc_Vectors().
|
static |
Find number of valid HKLs and average error, in list of HKLs.
Find number of valid HKLs and average error for the valid Miller indices, in a list of HKLs.
hkls | List of V3D objects containing hkl values |
tolerance | The maximum acceptable deviation from integer values for the Miller indices. |
average_error | This is set to the average error in the hkl values for the hkl values that are valid Miller indices. |
Definition at line 2070 of file IndexingUtils.cpp.
References count, fabs, tolerance, and ValidIndex().
|
static |
STATIC method Optimize_6dUB: Calculates the 6-dimensional matrix that most nearly maps the specified hkl_vectors and mnp_vectors to the specified q_vectors.
The calculated UB minimizes the sum squared differences between UB|ModUB*(h,k,l,m,n,p) and the corresponding (qx,qy,qz) for all of the specified hklmnp and Q vectors. The sum of the squares of the residual errors is returned. This method is used to optimize the UB matrix and ModUB matrix once an initial indexing has been found.
ModUB matrix is a 3x3 defines the modulation vectors in Q-sample. each colomn corresponds to a modulation vector.
UB | 3x3 matrix that will be set to the UB matrix |
ModUB | 3x3 matrix that will be set to the ModUB matrix |
hkl_vectors | std::vector of V3D objects that contains the list of hkl values |
mnp_vectors | std::vector of V3D objects that contains the list of mnp values |
ModDim | int value from 1 to 3, defines the number of dimensions of modulation. |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors. NOTE: The number of hkl_vectors and mnp_vectors and q_vectors must be the same, and must be at least 4. |
std::invalid_argument | exception if there are not at least 3 hkl and q vectors, or if the numbers of hkl and q vectors are not the same, or if the UB matrix is not a 3x3 matrix. or if the ModDim is not 1 or 2 or 3. |
std::runtime_error | exception if the QR factorization fails or the UB matrix can't be calculated or if UB is a singular matrix. |
Created by Shiyun Jin on 7/16/18.
Definition at line 886 of file IndexingUtils.cpp.
References CheckUB(), Mantid::Kernel::Matrix< T >::numCols(), Mantid::Kernel::Matrix< T >::numRows(), Mantid::Kernel::Matrix< T >::setRow(), and value.
|
static |
STATIC method Optimize_UB: Calculates the matrix that most nearly maps the specified hkl_vectors to the specified q_vectors.
The calculated UB minimizes the sum squared differences between UB*(h,k,l) and the corresponding (qx,qy,qz) for all of the specified hkl and Q vectors. The sum of the squares of the residual errors is returned. This method is used to optimize the UB matrix once an initial indexing has been found.
UB | 3x3 matrix that will be set to the UB matrix |
ModUB | 3x3 matrix that will be set to the ModUB matrix |
hkl_vectors | std::vector of V3D objects that contains the list of hkl values |
mnp_vectors | std::vector of V3D objects that contains the |
ModDim | int value from 1 to 3, defines the number of dimensions of modulation. |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors. |
sigabc | error in the crystal lattice parameter values if length is at least 6. NOTE: Calculation of these errors is based on SCD FORTRAN code base at IPNS. Contributors to the least squares application(1979) are J.Marc Overhage, G.Anderson, P. C. W. Leung, R. G. Teller, and A. J. Schultz NOTE: The number of hkl_vectors and q_vectors must be the same, and must be at least 3. |
sigq | error in the modulation vectors. |
std::invalid_argument | exception if there are not at least 3 hkl and q vectors, or if the numbers of hkl and q vectors are not the same, or if the UB matrix is not a 3x3 matrix. |
std::runtime_error | exception if the QR factorization fails or the UB matrix can't be calculated or if UB is a singular matrix. |
Definition at line 770 of file IndexingUtils.cpp.
References Mantid::Kernel::delta, GetLatticeParameters(), Mantid::Kernel::Matrix< T >::Invert(), Mantid::Geometry::m, n, and Optimize_6dUB().
Referenced by Mantid::Crystal::FindUBUsingIndexedPeaks::exec(), Mantid::MDAlgorithms::IntegrateEllipsoidsV1::exec(), and Optimize_6dUB().
|
static |
Find the vector that best corresponds to plane normal, given 1-D indices.
STATIC method Optimize_Direction: Calculates the vector for which the dot product of the the vector with each of the specified Qxyz vectors is most nearly the corresponding integer index.
The calculated best_vec minimizes the sum squared differences between best_vec dot (qx,qy,z) and the corresponding index for all of the specified Q vectors and indices. The sum of the squares of the residual errors is returned. NOTE: This method is similar the Optimize_UB method, but this method only optimizes the plane normal in one direction. Also, this optimizes the mapping from (qx,qy,qz) to one index (Q to index), while the Optimize_UB method optimizes the mapping from three (h,k,l) to (qx,qy,qz) (3 indices to Q).
best_vec | V3D vector that will be set to a vector whose direction most nearly corresponds to the plane normal direction and whose magnitude is d. The corresponding plane spacing in reciprocal space is 1/d. |
index_values | std::vector of ints that contains the list of indices |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are indexed in one direction by the corresponding index values. NOTE: The number of index_values and q_vectors must be the same, and must be at least 3. |
std::invalid_argument | exception if there are not at least 3 indices and q vectors, or if the numbers of indices and q vectors are not the same. |
std::runtime_error | exception if the QR factorization fails or the best direction can't be calculated. |
Definition at line 1027 of file IndexingUtils.cpp.
References value, and Mantid::Geometry::x.
Referenced by FFTScanFor_Directions(), and ScanFor_Directions().
|
static |
Find the UB matrix that most nearly maps hkl to qxyz for 3 or more peaks.
STATIC method Optimize_UB: Calculates the matrix that most nearly maps the specified hkl_vectors to the specified q_vectors.
The calculated UB minimizes the sum squared differences between UB*(h,k,l) and the corresponding (qx,qy,qz) for all of the specified hkl and Q vectors. The sum of the squares of the residual errors is returned. This method is used to optimize the UB matrix once an initial indexing has been found.
UB | 3x3 matrix that will be set to the UB matrix |
hkl_vectors | std::vector of V3D objects that contains the list of hkl values |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors. NOTE: The number of hkl_vectors and q_vectors must be the same, and must be at least 3. |
std::invalid_argument | exception if there are not at least 3 hkl and q vectors, or if the numbers of hkl and q vectors are not the same, or if the UB matrix is not a 3x3 matrix. |
std::runtime_error | exception if the QR factorization fails or the UB matrix can't be calculated or if UB is a singular matrix. |
Definition at line 644 of file IndexingUtils.cpp.
References CheckUB(), Mantid::Kernel::Matrix< T >::numCols(), Mantid::Kernel::Matrix< T >::numRows(), Mantid::Kernel::Matrix< T >::setRow(), and value.
|
static |
Find the UB matrix that most nearly maps hkl to qxyz for 3 or more peaks.
STATIC method Optimize_UB: Calculates the matrix that most nearly maps the specified hkl_vectors to the specified q_vectors.
The calculated UB minimizes the sum squared differences between UB*(h,k,l) and the corresponding (qx,qy,qz) for all of the specified hkl and Q vectors. The sum of the squares of the residual errors is returned. This method is used to optimize the UB matrix once an initial indexing has been found.
UB | 3x3 matrix that will be set to the UB matrix |
hkl_vectors | std::vector of V3D objects that contains the list of hkl values |
q_vectors | std::vector of V3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors. |
sigabc | error in the crystal lattice parameter values if length is at least 6. NOTE: Calculation of these errors is based on SCD FORTRAN code base at IPNS. Contributors to the least squares application(1979) are J.Marc Overhage, G.Anderson, P. C. W. Leung, R. G. Teller, and A. J. Schultz NOTE: The number of hkl_vectors and q_vectors must be the same, and must be at least 3. |
std::invalid_argument | exception if there are not at least 3 hkl and q vectors, or if the numbers of hkl and q vectors are not the same, or if the UB matrix is not a 3x3 matrix. |
std::runtime_error | exception if the QR factorization fails or the UB matrix can't be calculated or if UB is a singular matrix. |
Definition at line 562 of file IndexingUtils.cpp.
References Mantid::Kernel::delta, GetLatticeParameters(), Mantid::Kernel::Matrix< T >::Invert(), Mantid::Geometry::m, n, and Optimize_UB().
Referenced by Mantid::Crystal::SelectCellWithForm::DetermineErrors(), Mantid::Crystal::FindUBUsingFFT::exec(), Mantid::Crystal::FindUBUsingMinMaxD::exec(), Mantid::MDAlgorithms::IntegrateEllipsoidsTwoStep::exec(), Find_UB(), and Optimize_UB().
|
static |
Round all the components of a HKL objects to the nearest integer.
Round all of the components of the V3D to the nearest integer.
hkl | V3D object whose components will be rounded. |
Definition at line 2005 of file IndexingUtils.cpp.
Referenced by RoundHKLs().
|
static |
Round all the components of a list of V3D objects, to the nearest integer.
Round all of the components of all vectors to the nearest integer.
This is useful when the vectors in the list represent Miller indices. Since the PeaksWorkspace records the Miller indices as sets of three doubles, there is no guarantee that the Miller indices will be integers.
hkl_list | Vector of V3D objects whose components will be rounded. |
Definition at line 2020 of file IndexingUtils.cpp.
References RoundHKL().
|
static |
Get list of possible directions and lengths for real space unit cell.
Get list of possible edge vectors for the real space unit cell.
This list will consist of vectors, V, for which V dot Q is essentially an integer for the most Q vectors. The difference between V dot Q and an integer must be less than the required tolerance for it to count as an integer.
directions | Vector that will be filled with the directions that may correspond to unit cell edges. |
q_vectors | Vector of new Vector3D objects that contains the list of q_vectors that are to be indexed. |
min_d | Lower bound on shortest unit cell edge length. This does not have to be specified exactly but must be strictly less than the smallest edge length, in Angstroms. |
max_d | Upper bound on longest unit cell edge length. This does not have to be specified exactly but must be strictly more than the longest edge length in angstroms. |
required_tolerance | The maximum allowed deviation of Miller indices from integer values for a peak to be indexed. |
degrees_per_step | The number of degrees between directions that are checked while scanning for an initial indexing of the peaks with lowest |Q|. |
Definition at line 1297 of file IndexingUtils.cpp.
References error, fabs, GetIndexedPeaks_1D(), MakeHemisphereDirections(), Mantid::Kernel::V3D::norm(), Optimize_Direction(), and Mantid::Kernel::V3D::scalar_prod().
Referenced by Find_UB().
|
static |
Scan rotations to find UB that indexes peaks given lattice parameters.
The method uses two passes to scan across all possible directions and orientations to find the direction and orientation for the unit cell that best fits the specified list of peaks.
On the first pass, only those sets of directions that index the most peaks are kept. On the second pass, the directions that minimize the sum-squared deviations from integer indices are selected from that smaller set of directions. This method should be most useful if number of peaks is on the order of 10-20, and most of the peaks belong to the same crystallite.
UB | This will be set to the UB matrix that best indexes the supplied list of q_vectors. |
q_vectors | List of locations of peaks in "Q". |
cell | Unit cell defining the parameters a,b,c and alpha, beta, gamma. |
degrees_per_step | The number of degrees per step used when scanning through all possible directions and orientations for the unit cell. NOTE: The work required rises very rapidly as the number of degrees per step decreases. A value of 1 degree leads to about 10 seconds of compute time. while a value of 2 only requires a bit more than 1 sec. The required time is O(n^3) where n = 1/degrees_per_step. |
required_tolerance | The maximum distance from an integer that the calculated h,k,l values can have if a peak is to be considered indexed. |
std::invalid_argument | exception if the UB matrix is not a 3X3 matrix. |
std::runtime_error | exception if the matrix inversion fails and UB can't be formed |
Definition at line 1134 of file IndexingUtils.cpp.
References Mantid::Geometry::UnitCell::a(), Mantid::Geometry::UnitCell::alpha(), Mantid::Geometry::UnitCell::b(), Mantid::Geometry::UnitCell::beta(), Mantid::Geometry::UnitCell::c(), error, fabs, Mantid::Geometry::UnitCell::gamma(), Mantid::Geometry::OrientedLattice::GetUB(), makeCDir(), MakeCircleDirections(), MakeHemisphereDirections(), Mantid::Kernel::Matrix< T >::numCols(), Mantid::Kernel::Matrix< T >::numRows(), and Mantid::Kernel::V3D::scalar_prod().
Referenced by Find_UB().
|
static |
Choose the direction in a list of directions, that is most nearly perpendicular to planes with the specified spacing in reciprocal space.
Choose the direction vector that most nearly corresponds to a family of planes in the list of q_vectors, with spacing equal to the specified plane_spacing.
The direction is chosen from the specified direction_list.
best_direction | This will be set to the direction that minimizes the sum squared distances of projections of peaks from integer multiples of the specified plane spacing. |
q_vectors | List of peak positions, specified according to the convention that |q| = 1/d. (i.e. Q/2PI) |
direction_list | List of possible directions for plane normals. Initially, this will be a long list of possible directions from MakeHemisphereDirections(). |
plane_spacing | The required spacing between planes in reciprocal space. |
required_tolerance | The maximum deviation of the component of a peak Qxyz in the direction of the best_direction vector for that peak to count as being indexed. NOTE: The tolerance is specified in terms of Miller Index. That is, the distance between adjacent planes is effectively normalized to one for measuring the error in the computed index. |
invalid_argument | exception of no Q vectors or directions are specified. |
Definition at line 2726 of file IndexingUtils.cpp.
References error, Mantid::Kernel::V3D::normalize(), and Mantid::Kernel::V3D::scalar_prod().
|
static |
Check is hkl is within tolerance of integer (h,k,l) non-zero values.
Check whether or not the components of the specified vector are within the specified tolerance of integer values, other than (0,0,0).
hkl | A V3D object containing what may be valid Miller indices for a peak. |
tolerance | The maximum acceptable deviation from integer values for the Miller indices. |
Definition at line 2053 of file IndexingUtils.cpp.
References tolerance, and withinTol().
Referenced by CalculateMillerIndices(), Mantid::Crystal::FindUBUsingIndexedPeaks::exec(), Mantid::Crystal::OptimizeCrystalPlacement::exec(), Mantid::MDAlgorithms::IntegrateEllipsoidsTwoStep::exec(), Mantid::MDAlgorithms::IntegrateEllipsoidsV1::exec(), Mantid::MDAlgorithms::Integrate3DEvents::getHklMnpKey(), Mantid::MDAlgorithms::Integrate3DEvents::getHklMnpKey2(), GetIndexedPeaks(), GetIndexedPeaks_3D(), Mantid::Crystal::FindUBUsingIndexedPeaks::isPeakIndexed(), NumberIndexed(), NumberIndexed_3D(), and NumberOfValidIndexes().
|
staticprivate |
Static reference to the logger class.
Definition at line 199 of file IndexingUtils.h.