Mantid
Loading...
Searching...
No Matches
Static Public Member Functions | Static Private Attributes | List of all members
Mantid::Geometry::IndexingUtils Class Reference

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 &current_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::V3DMakeCircleDirections (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::V3DMakeHemisphereDirections (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::Loggerg_Log
 Static reference to the logger class. More...
 

Detailed Description

This class contains static utility methods for indexing peaks and finding the UB matrix.

Author
Dennis Mikkelson
Date
2011-06-14

Definition at line 32 of file IndexingUtils.h.

Member Function Documentation

◆ CalculateMillerIndices() [1/3]

V3D IndexingUtils::CalculateMillerIndices ( const Kernel::DblMatrix inverseUB,
const Kernel::V3D q_vector 
)
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.

Parameters
inverseUBA 3x3 matrix of doubles holding the inverse UB matrix. The matrix is not checked for validity
q_vectorV3D object containing Q vector in sample frame
Returns
The indexes of the given peak. They have not been tested for validity

Definition at line 2366 of file IndexingUtils.cpp.

◆ CalculateMillerIndices() [2/3]

bool IndexingUtils::CalculateMillerIndices ( const Kernel::DblMatrix inverseUB,
const Kernel::V3D q_vector,
double  tolerance,
Kernel::V3D miller_indices 
)
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)

Parameters
inverseUBA 3x3 matrix of doubles holding the inverse UB matrix. The matrix is not checked for validity
q_vectorstd::vector of V3D objects that contains the list of q_vectors that are to be indexed.
toleranceThe 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_indicesThis vector returns a list of Miller Indices, with one entry for each given Q vector.
Returns
True if the peak was index, false otherwise

Definition at line 2344 of file IndexingUtils.cpp.

References CalculateMillerIndices(), tolerance, and ValidIndex().

◆ CalculateMillerIndices() [3/3]

int IndexingUtils::CalculateMillerIndices ( const Kernel::DblMatrix UB,
const std::vector< Kernel::V3D > &  q_vectors,
double  tolerance,
std::vector< Kernel::V3D > &  miller_indices,
double &  ave_error 
)
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.

Parameters
UBA 3x3 matrix of doubles holding the UB matrix
q_vectorsstd::vector of V3D objects that contains the list of q_vectors that are to be indexed.
toleranceThe 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_indicesThis vector returns a list of Miller Indices, with one entry for each given Q vector.
ave_errorThe average error from all lattice directions.
Returns
A non-negative integer giving the number of peaks indexed by UB, within the specified tolerance on h,k,l.
Exceptions
std::invalid_argumentexception 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().

◆ CheckUB()

bool IndexingUtils::CheckUB ( const Kernel::DblMatrix UB)
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.

Parameters
UBA 3x3 matrix of doubles holding the UB matrix
Returns
true if this could be a valid 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().

◆ DiscardDuplicates()

void IndexingUtils::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 
)
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.

Parameters
new_listThis vector will be cleared and filled with the vectors from the directions list that are not duplicates of other vectors in the list.
directionsInput list of possible a,b,c directions, listed in order of increasing vector norm. This list will be cleared by this method.
q_vectorsList of q_vectors that should be indexed
required_toleranceThe tolerance for indexing
len_tolThe 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_tolThe 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().

◆ FFTScanFor_Directions()

size_t IndexingUtils::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 
)
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.

Parameters
directionsVector that will be filled with the directions that may correspond to unit cell edges.
q_vectorsVector of new Vector3D objects that contains the list of q_vectors that are to be indexed.
min_dLower 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_dUpper 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_toleranceThe maximum allowed deviation of Miller indices from integer values for a peak to be indexed.
degrees_per_stepThe 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().

◆ Find_UB() [1/3]

double IndexingUtils::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 
)
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.

Parameters
UB3x3 matrix that will be set to the UB matrix
q_vectorsstd::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_dLower 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_dUpper 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_toleranceThe maximum allowed deviation of Miller indices from integer values for a peak to be indexed.
degrees_per_stepThe number of degrees between different orientations used during the initial scan.
iterationsNumber of refinements of UB
Returns
This will return the sum of the squares of the residual errors.
Exceptions
std::invalid_argumentexception 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().

◆ Find_UB() [2/3]

double IndexingUtils::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 
)
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.

Parameters
UB3x3 matrix that will be set to the UB matrix
q_vectorsstd::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_dLower 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_dUpper 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_toleranceThe maximum allowed deviation of Miller indices from integer values for a peak to be indexed.
base_indexThe 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_initialThe number of low |Q| peaks that should be used to scan for an initial orientation matrix.
degrees_per_stepThe number of degrees between different orientations used during the initial scan.
Returns
This will return the sum of the squares of the residual errors.
Exceptions
std::invalid_argumentexception 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().

◆ Find_UB() [3/3]

double IndexingUtils::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 
)
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.

Parameters
UB3x3 matrix that will be set to the UB matrix
q_vectorsstd::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.
latticeThe 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_toleranceThe maximum allowed deviation of Miller indices from integer values for a peak to be indexed.
base_indexThe 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_initialThe number of low |Q| peaks that should be used to scan for an initial orientation matrix.
degrees_per_stepThe number of degrees between different orientations used during the initial scan.
fixAllFix the lattice parameters and do not optimise the UB matrix.
iterationsNumber of refinements of UB
Returns
This will return the sum of the squares of the residual errors.
Exceptions
std::invalid_argumentexception 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().

◆ FormUB_From_abc_Vectors() [1/2]

bool IndexingUtils::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 
)
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.

Parameters
UBThe calculated UB matrix will be returned in this parameter
directionsList of possible vectors for a, b, c. This list MUST be sorted in order of increasing magnitude.
q_vectorsThe list of q_vectors that should be indexed
req_toleranceThe required tolerance on h,k,l to consider a peak to be indexed.
min_volThe smallest possible unit cell volume.
Returns
true if a UB matrix was set, and false if it not possible to choose a,b,c (i.e. UB) from the list of directions, starting with the specified a_index.
Exceptions
std::invalid_argumentexception if the UB matrix is not a 3X3 matrix.
std::runtime_errorexception 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().

◆ FormUB_From_abc_Vectors() [2/2]

bool IndexingUtils::FormUB_From_abc_Vectors ( Kernel::DblMatrix UB,
const std::vector< Kernel::V3D > &  directions,
size_t  a_index,
double  min_d,
double  max_d 
)
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.

Parameters
UBThe calculated UB matrix will be returned in this parameter
directionsList of possible vectors for a, b, c. This list MUST be sorted in order of increasing magnitude.
a_indexThe index to use for the a vector. The b and c vectors will be choosen from LATER positions in the directions list.
min_dMinimum possible real space unit cell edge length.
max_dMaximum possible real space unit cell edge length.
Returns
true if a UB matrix was set, and false if it not possible to choose a,b,c (i.e. UB) from the list of directions, starting with the specified a_index.
Exceptions
std::invalid_argumentexception if the UB matrix is not a 3X3 matrix.
std::runtime_errorexception 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().

◆ GetFirstMaxIndex()

double IndexingUtils::GetFirstMaxIndex ( const double  magnitude_fft[],
size_t  N,
double  threshold 
)
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.

Parameters
magnitude_fftThe array containing the magnitude of the FFT values.
NThe size of the FFT array.
thresholdThe required threshold for the first peak. This must be positive.
Returns
The centroid (index) where the first maximum occurs, or -1 if no point in the FFT (beyond the DC term) equals or exceeds the required threshold.

Definition at line 1646 of file IndexingUtils.cpp.

Referenced by FFTScanFor_Directions().

◆ GetIndexedPeaks()

int IndexingUtils::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 
)
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.

Parameters
UBThe UB matrix that determines the indexing of the peaks.
q_vectorsList of V3D peaks in reciprocal space
required_toleranceThe maximum allowed error (as a faction of the corresponding Miller index) for a peak q_vector to be counted as indexed.
miller_indicesList of the Miller indices (h,k,l) of peaks that were indexed in all specified directions.
indexed_qsList of Qxyz value for the peaks that were indexed indexed in all specified directions.
fit_errorThe 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.
Returns
The number of q_vectors that are indexed to within the specified tolerance, in the specified direction.

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().

◆ GetIndexedPeaks_1D()

int IndexingUtils::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 
)
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.

Parameters
directionDirection 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_vectorsList of V3D peaks in reciprocal space
required_toleranceThe maximum allowed error (as a faction of the corresponding Miller index) for a peak q_vector to be counted as indexed.
index_valsList of the one-dimensional Miller indices peaks that were indexed in the specified direction.
indexed_qsList of Qxyz value for the peaks that were indexed indexed in the specified direction.
fit_errorThe 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.
Returns
The number of q_vectors that are indexed to within the specified tolerance, in the specified direction.

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().

◆ GetIndexedPeaks_3D()

int IndexingUtils::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 
)
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).

Parameters
direction_1Direction vector in the direction of the normal vector for the first family of parallel planes.
direction_2Direction vector in the direction of the normal vector for the second family of parallel planes.
direction_3Direction vector in the direction of the normal vector for the third family of parallel planes.
q_vectorsList of V3D peaks in reciprocal space
required_toleranceThe maximum allowed error (as a faction of the corresponding Miller index) for a peak q_vector to be counted as indexed.
miller_indicesList of the Miller indices (h,k,l) of peaks that were indexed in all specified directions.
indexed_qsList of Qxyz value for the peaks that were indexed indexed in all specified directions.
fit_errorThe 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.
Returns
The number of q_vectors that are indexed to within the specified tolerance, in the specified direction.

Definition at line 2464 of file IndexingUtils.cpp.

References fabs, Mantid::Kernel::V3D::scalar_prod(), and ValidIndex().

◆ GetLatticeParameters()

bool IndexingUtils::GetLatticeParameters ( const Kernel::DblMatrix UB,
std::vector< double > &  lattice_par 
)
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.

Parameters
UBA non-singular matrix representing an orientation matrix.
lattice_parstd::vector of doubles that will contain the lattice parameters and cell volume as it's first seven entries.
Returns
true if the lattice_par vector was filled with the lattice parameters and false if the matrix could not be inverted.

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().

◆ GetLatticeParameterString()

std::string IndexingUtils::GetLatticeParameterString ( const Kernel::DblMatrix 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.

Returns
a string listing the cell parameters.

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().

◆ GetMagFFT()

double IndexingUtils::GetMagFFT ( const std::vector< Kernel::V3D > &  q_vectors,
const Kernel::V3D current_dir,
const size_t  N,
double  projections[],
double  index_factor,
double  magnitude_fft[] 
)
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.

Parameters
q_vectorsThe list of Q vectors to project on the specified direction.
current_dirThe direction the Q vectors will be projected on.
NThe size of the projections[] array. This MUST BE a power of 2. The magnitude_fft[] array must be half the size.
projectionsArray 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_factorFactor that when multiplied by a projected Q vector will give a valid index into the projections array.
magnitude_fftArray that will be filled out with the magnitude of the FFT of the projections.
Returns
The largest value in the magnitude_fft, that is stored in position 5 or more.

Definition at line 1600 of file IndexingUtils.cpp.

References fabs, index, and Mantid::Kernel::V3D::scalar_prod().

Referenced by FFTScanFor_Directions().

◆ GetModulationVector()

bool IndexingUtils::GetModulationVector ( const Kernel::DblMatrix UB,
const Kernel::DblMatrix ModUB,
Kernel::V3D ModVec,
const int  j 
)
static

◆ GetModulationVectors()

int IndexingUtils::GetModulationVectors ( const Kernel::DblMatrix UB,
const Kernel::DblMatrix ModUB,
Kernel::V3D ModVec1,
Kernel::V3D ModVec2,
Kernel::V3D ModVec3 
)
static

◆ IndexingError()

double IndexingUtils::IndexingError ( const Kernel::DblMatrix UB,
const std::vector< Kernel::V3D > &  hkls,
const std::vector< Kernel::V3D > &  q_vectors 
)
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().

◆ makeCDir()

V3D IndexingUtils::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 
)
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.

Parameters
a_dirV3D object with length "a" in the direction of the rotated cell edge "a"
b_dirV3D object with length "b" in the direction of the rotated cell edge "b"
cThe length of the third cell edge, c.
cosAlphacos angle between edges b and c in radians.
cosBetacos angle between edges c and a in radians.
cosGammacos angle between edges a and b in radians.
sinGammasin angle between edges a and b in radians.
Returns
A new V3D object with length "c", in the direction of the third rotated unit cell edge, "c".

Definition at line 1884 of file IndexingUtils.cpp.

References Mantid::Kernel::toV3D(), and Mantid::Kernel::toVector3d().

Referenced by ScanFor_UB().

◆ MakeCircleDirections()

std::vector< V3D > IndexingUtils::MakeCircleDirections ( int  n_steps,
const Kernel::V3D axis,
double  angle_degrees 
)
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.

Parameters
n_stepsThe number of vectors to generate around the circle.
axisThe specified axis
angle_degreesThe specified angle
Returns
A std::vector containing direction vectors forming the same angle with the axis.
Exceptions
std::invalid_argumentexception 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().

◆ MakeHemisphereDirections()

std::vector< V3D > IndexingUtils::MakeHemisphereDirections ( int  n_steps)
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.

Parameters
n_stepsThe number of subdivisions in latitude in the upper hemisphere.
Returns
A std::vector containing directions distributed over the hemisphere with y-coordinate at least zero.
Exceptions
std::invalid_argumentexception 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().

◆ NumberIndexed()

int IndexingUtils::NumberIndexed ( const Kernel::DblMatrix UB,
const std::vector< Kernel::V3D > &  q_vectors,
double  tolerance 
)
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.

Parameters
UBA 3x3 matrix of doubles holding the UB matrix. The UB matrix must not be singular.
q_vectorsstd::vector of V3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors.
toleranceThe 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.
Returns
A non-negative integer giving the number of peaks indexed by UB.
Exceptions
std::invalid_argumentexception 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().

◆ NumberIndexed_1D()

int IndexingUtils::NumberIndexed_1D ( const Kernel::V3D direction,
const std::vector< Kernel::V3D > &  q_vectors,
double  tolerance 
)
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.

Parameters
directionA V3D specifying a possible edge vector in real space.
q_vectorsstd::vector of V3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors.
toleranceThe maximum allowed distance to an integer from the dot products of peaks with the specified direction.
Returns
A non-negative integer giving the number of q-vectors indexed in one direction by the specified direction vector.

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().

◆ NumberIndexed_3D()

int IndexingUtils::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 
)
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.

Parameters
a_dirA Vector3D representing unit cell edge vector a
b_dirA Vector3D representing unit cell edge vector b
c_dirA Vector3D representing unit cell edge vector c
q_vectorsVector of Vector3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors.
toleranceThe maximum allowed distance to an integer from the dot products of peaks with the specified direction.
Returns
A non-negative integer giving the number of peaks simultaneously indexed in all three directions by the specified direction vectors.

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().

◆ NumberOfValidIndexes()

int IndexingUtils::NumberOfValidIndexes ( const std::vector< Kernel::V3D > &  hkls,
double  tolerance,
double &  average_error 
)
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.

Parameters
hklsList of V3D objects containing hkl values
toleranceThe maximum acceptable deviation from integer values for the Miller indices.
average_errorThis 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().

◆ Optimize_6dUB() [1/2]

double IndexingUtils::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

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.

Parameters
UB3x3 matrix that will be set to the UB matrix
ModUB3x3 matrix that will be set to the ModUB matrix
hkl_vectorsstd::vector of V3D objects that contains the list of hkl values
mnp_vectorsstd::vector of V3D objects that contains the list of mnp values
ModDimint value from 1 to 3, defines the number of dimensions of modulation.
q_vectorsstd::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.
Returns
This will return the sum of the squares of the residual differences between the Q vectors provided and the UB*hkl values, in reciprocal space.
Exceptions
std::invalid_argumentexception 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_errorexception 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.

◆ Optimize_6dUB() [2/2]

double IndexingUtils::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

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.

Parameters
UB3x3 matrix that will be set to the UB matrix
ModUB3x3 matrix that will be set to the ModUB matrix
hkl_vectorsstd::vector of V3D objects that contains the list of hkl values
mnp_vectorsstd::vector of V3D objects that contains the
ModDimint value from 1 to 3, defines the number of dimensions of modulation.
q_vectorsstd::vector of V3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors.
sigabcerror 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.
sigqerror in the modulation vectors.
Returns
This will return the sum of the squares of the residual differences between the Q vectors provided and the UB*hkl values, in reciprocal space.
Exceptions
std::invalid_argumentexception 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_errorexception 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().

◆ Optimize_Direction()

double IndexingUtils::Optimize_Direction ( Kernel::V3D best_vec,
const std::vector< int > &  index_values,
const std::vector< Kernel::V3D > &  q_vectors 
)
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).

Parameters
best_vecV3D 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_valuesstd::vector of ints that contains the list of indices
q_vectorsstd::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.
Returns
This will return the sum of the squares of the residual errors.
Exceptions
std::invalid_argumentexception 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_errorexception 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().

◆ Optimize_UB() [1/2]

double IndexingUtils::Optimize_UB ( Kernel::DblMatrix UB,
const std::vector< Kernel::V3D > &  hkl_vectors,
const std::vector< Kernel::V3D > &  q_vectors 
)
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.

Parameters
UB3x3 matrix that will be set to the UB matrix
hkl_vectorsstd::vector of V3D objects that contains the list of hkl values
q_vectorsstd::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.
Returns
This will return the sum of the squares of the residual differences between the Q vectors provided and the UB*hkl values, in reciprocal space.
Exceptions
std::invalid_argumentexception 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_errorexception 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.

◆ Optimize_UB() [2/2]

double IndexingUtils::Optimize_UB ( Kernel::DblMatrix UB,
const std::vector< Kernel::V3D > &  hkl_vectors,
const std::vector< Kernel::V3D > &  q_vectors,
std::vector< double > &  sigabc 
)
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.

Parameters
UB3x3 matrix that will be set to the UB matrix
hkl_vectorsstd::vector of V3D objects that contains the list of hkl values
q_vectorsstd::vector of V3D objects that contains the list of q_vectors that are indexed by the corresponding hkl vectors.
sigabcerror 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.
Returns
This will return the sum of the squares of the residual differences between the Q vectors provided and the UB*hkl values, in reciprocal space.
Exceptions
std::invalid_argumentexception 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_errorexception 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().

◆ RoundHKL()

void IndexingUtils::RoundHKL ( Kernel::V3D hkl)
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.

Parameters
hklV3D object whose components will be rounded.

Definition at line 2005 of file IndexingUtils.cpp.

Referenced by RoundHKLs().

◆ RoundHKLs()

void IndexingUtils::RoundHKLs ( std::vector< Kernel::V3D > &  hkl_list)
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.

Parameters
hkl_listVector of V3D objects whose components will be rounded.

Definition at line 2020 of file IndexingUtils.cpp.

References RoundHKL().

◆ ScanFor_Directions()

size_t IndexingUtils::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 
)
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.

Parameters
directionsVector that will be filled with the directions that may correspond to unit cell edges.
q_vectorsVector of new Vector3D objects that contains the list of q_vectors that are to be indexed.
min_dLower 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_dUpper 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_toleranceThe maximum allowed deviation of Miller indices from integer values for a peak to be indexed.
degrees_per_stepThe 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().

◆ ScanFor_UB()

double IndexingUtils::ScanFor_UB ( Kernel::DblMatrix UB,
const std::vector< Kernel::V3D > &  q_vectors,
const UnitCell cell,
double  degrees_per_step,
double  required_tolerance 
)
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.

Parameters
UBThis will be set to the UB matrix that best indexes the supplied list of q_vectors.
q_vectorsList of locations of peaks in "Q".
cellUnit cell defining the parameters a,b,c and alpha, beta, gamma.
degrees_per_stepThe 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_toleranceThe maximum distance from an integer that the calculated h,k,l values can have if a peak is to be considered indexed.
Exceptions
std::invalid_argumentexception if the UB matrix is not a 3X3 matrix.
std::runtime_errorexception 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().

◆ SelectDirection()

int IndexingUtils::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 
)
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.

Parameters
best_directionThis 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_vectorsList of peak positions, specified according to the convention that |q| = 1/d. (i.e. Q/2PI)
direction_listList of possible directions for plane normals. Initially, this will be a long list of possible directions from MakeHemisphereDirections().
plane_spacingThe required spacing between planes in reciprocal space.
required_toleranceThe 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.
Returns
The number of peaks that lie within the specified tolerance of the family of planes with normal direction = best_direction and with spacing given by plane_spacing.
Exceptions
invalid_argumentexception 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().

◆ ValidIndex()

bool IndexingUtils::ValidIndex ( const Kernel::V3D hkl,
double  tolerance 
)
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).

Parameters
hklA V3D object containing what may be valid Miller indices for a peak.
toleranceThe maximum acceptable deviation from integer values for the Miller indices.
Returns
true if all components of the vector are within the tolerance of integer values (h,k,l) and (h,k,l) is NOT (0,0,0)

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().

Member Data Documentation

◆ g_Log

Mantid::Kernel::Logger& Mantid::Geometry::IndexingUtils::g_Log
staticprivate

Static reference to the logger class.

Definition at line 199 of file IndexingUtils.h.


The documentation for this class was generated from the following files: