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

This is a low-level class to construct a map with lists of events near each peak Q-vector in the lab frame. More...

#include <IntegrateQLabEvents.h>

Public Member Functions

void addEvents (SlimEvents const &event_qs)
 distribute the events among the cells of the partitioned QLab space. More...
 
PeakShape_const_sptr ellipseIntegrateEvents (const std::vector< V3D > &E1Vec, V3D const &peak_q, bool specify_size, double peak_radius, double back_inner_radius, double back_outer_radius, std::vector< double > &axes_radii, double &inti, double &sigi, std::pair< double, double > &backi)
 Integrate the events around the specified peak QLab vector. More...
 
 IntegrateQLabEvents (const SlimEvents &peak_q_list, double radius, const bool useOnePercentBackgroundCorrection=true)
 Store events within a certain radius of the specified peak centers, and sum these events to estimate pixel intensities. More...
 
void populateCellsWithPeaks ()
 Assign events to each of the cells occupied by events. More...
 
void setRadius (const double &radius)
 Set peak integration radius. More...
 

Static Public Member Functions

static bool isOrigin (const V3D &q, const double &cellSize)
 Determine if an input Q-vector lies in the cell associated to the origin. More...
 

Private Member Functions

void addEvent (const SlimEvent event)
 assign an event to one cell of the partitioned QLab space. More...
 
double detectorQ (const std::vector< V3D > &E1Vec, const V3D &QLabFrame, const std::vector< double > &r)
 Calculate if this Q is on a detector. More...
 
PeakShapeEllipsoid_const_sptr ellipseIntegrateEvents (const std::vector< V3D > &E1Vec, V3D const &peak_q, SlimEvents const &ev_list, std::vector< V3D > const &directions, std::vector< double > const &sigmas, bool specify_size, double peak_radius, double back_inner_radius, double back_outer_radius, std::vector< double > &axes_radii, double &inti, double &sigi, std::pair< double, double > &backi)
 Integrate a list of events associated to one peak. More...
 

Static Private Member Functions

static void getEigenVectors (Kernel::DblMatrix const &cov_matrix, std::vector< V3D > &eigen_vectors, std::vector< double > &eigen_values)
 Eigen vectors of a 3x3 real symmetric matrix using the GSL. More...
 
static void makeCovarianceMatrix (SlimEvents const &events, Kernel::DblMatrix &matrix, double radius)
 3x3 covariance matrix of a list of SlimEvent objects More...
 
static std::pair< double, double > numInEllipsoid (SlimEvents const &events, std::vector< V3D > const &directions, std::vector< double > const &sizes)
 Number of events in an ellipsoid. More...
 
static std::pair< double, double > numInEllipsoidBkg (SlimEvents const &events, std::vector< V3D > const &directions, std::vector< double > const &sizes, std::vector< double > const &sizesIn, const bool useOnePercentBackgroundCorrection)
 Number of events in an ellipsoid with background correction. More...
 

Private Attributes

double m_cellSize
 size of the square cell unit, holding at most one single peak More...
 
std::unordered_map< size_t, SlimEventsm_cellsWithEvents
 list of cells occupied with events More...
 
std::unordered_map< size_t, OccupiedCellm_cellsWithPeaks
 list the occupied cells in an unordered map for fast searching More...
 
double m_radius
 
const bool m_useOnePercentBackgroundCorrection
 if one percent culling of the background should be performed. More...
 

Detailed Description

This is a low-level class to construct a map with lists of events near each peak Q-vector in the lab frame.

The Q-vector of each event is shifted by the Q-vector of the associated peak. A method is also provided to find the principal axes of such a list of events and to find the net integrated counts using ellipsoids with axis lengths determined from the standard deviations in the directions of the principal axes.

Author
Dennis Mikkelson
Date
2012-12-19

Definition at line 88 of file IntegrateQLabEvents.h.

Constructor & Destructor Documentation

◆ IntegrateQLabEvents()

Mantid::MDAlgorithms::IntegrateQLabEvents::IntegrateQLabEvents ( const SlimEvents peak_q_list,
double  radius,
const bool  useOnePercentBackgroundCorrection = true 
)

Store events within a certain radius of the specified peak centers, and sum these events to estimate pixel intensities.

Parameters
peak_q_list: List of Q-vectors for peak centers.
radius: The maximum distance from a peak's Q-vector, for an event to be stored in the list associated with that peak.
useOnePercentBackgroundCorrection: flag if one percent background correction should be used.

Definition at line 45 of file IntegrateQLabEvents.cpp.

References Mantid::MDAlgorithms::CellCoords::getHash(), Mantid::MDAlgorithms::CellCoords::isOrigin(), m_cellSize, m_cellsWithPeaks, and m_radius.

Member Function Documentation

◆ addEvent()

void Mantid::MDAlgorithms::IntegrateQLabEvents::addEvent ( const SlimEvent  event)
private

assign an event to one cell of the partitioned QLab space.

Parameters
event: SlimEvent to be assigned

Definition at line 232 of file IntegrateQLabEvents.cpp.

References Mantid::MDAlgorithms::CellCoords::getHash(), Mantid::MDAlgorithms::CellCoords::isOrigin(), m_cellSize, and m_cellsWithEvents.

Referenced by addEvents().

◆ addEvents()

void Mantid::MDAlgorithms::IntegrateQLabEvents::addEvents ( SlimEvents const &  event_qs)

distribute the events among the cells of the partitioned QLab space.

Given QLab partitioned into a cubic lattice with unit cell of certain size, assign each event to one particular cell depending on its QLab vector.

Parameters
event_qs: List of SlimEvent objects to be distributed

Definition at line 75 of file IntegrateQLabEvents.cpp.

References addEvent().

Referenced by Mantid::MDAlgorithms::IntegrateEllipsoidsV2::qListFromEventWS(), and Mantid::MDAlgorithms::IntegrateEllipsoidsV2::qListFromHistoWS().

◆ detectorQ()

double Mantid::MDAlgorithms::IntegrateQLabEvents::detectorQ ( const std::vector< V3D > &  E1Vec,
const V3D QLabFrame,
const std::vector< double > &  r 
)
private

Calculate if this Q is on a detector.

The distance from C to OE is given by dv=C-E*(C.scalar_prod(E)) If dv.norm<integration_radius, one of the detector trajectories on the edge is too close to the peak. This method is applied to all masked pixels. If there are masked pixels trajectories inside an integration volume, the peak must be rejected.

Parameters
E1Vec: Vector of values for calculating edge of detectors
QLabFrameThe Peak center.
rPeak radius.

Definition at line 362 of file IntegrateQLabEvents.cpp.

References Mantid::API::E1(), Mantid::Kernel::V3D::norm(), and Mantid::Kernel::V3D::scalar_prod().

Referenced by ellipseIntegrateEvents().

◆ ellipseIntegrateEvents() [1/2]

Mantid::Geometry::PeakShape_const_sptr Mantid::MDAlgorithms::IntegrateQLabEvents::ellipseIntegrateEvents ( const std::vector< V3D > &  E1Vec,
V3D const &  peak_q,
bool  specify_size,
double  peak_radius,
double  back_inner_radius,
double  back_outer_radius,
std::vector< double > &  axes_radii,
double &  inti,
double &  sigi,
std::pair< double, double > &  backi 
)

Integrate the events around the specified peak QLab vector.

The principal axes of the events near this Q-vector and the standard deviations in the directions of these principal axes determine ellipsoidal regions for integrating the peak and estimating the background. Alternatively, if peak and background radii are specified, then those will be used for half the major axis length of the ellipsoids, and the other axes of the ellipsoids will be set proportionally, based on the standard deviations.

Parameters
E1Vec: Vector of values for calculating edge of detectors
peak_q: The QLab-vector for the peak center.
specify_size: If true the integration will be done using the ellipsoids with major axes determined by the peak, back_inner and back_outer radii parameters. If false, the integration will be done using a peak region with major axis chosen so that it covers +- three standard deviations of the data in each direction. In this case, the background ellipsoidal shell is chosen to have the same VOLUME as the peak ellipsoid, and to use the peak ellipsoid for the inner radius.
peak_radius: Size of half the major axis of the ellipsoidal peak region.
back_inner_radius: Size of half the major axis of the INNER ellipsoidal boundary of the background region
back_outer_radius: Size of half the major axis of the OUTER ellipsoidal boundary of the background region
axes_radii: The radii used for integration in the directions of the three principal axes.
inti: (output) collects the net integrated intensity
sigi: (output) collects an estimate of the standard deviation of the net integrated intensity
backi: (output) collects background intensity subtracted from inti and sigi

Definition at line 82 of file IntegrateQLabEvents.cpp.

References ellipseIntegrateEvents(), Mantid::MDAlgorithms::OccupiedCell::events, getEigenVectors(), Mantid::MDAlgorithms::CellCoords::getHash(), m_cellSize, m_cellsWithPeaks, m_radius, makeCovarianceMatrix(), and sigma.

Referenced by ellipseIntegrateEvents(), Mantid::MDAlgorithms::IntegrateEllipsoidsV2::exec(), and Mantid::MDAlgorithms::IntegrateEllipsoidsV2::integratePeaksCutoffISigI().

◆ ellipseIntegrateEvents() [2/2]

PeakShapeEllipsoid_const_sptr Mantid::MDAlgorithms::IntegrateQLabEvents::ellipseIntegrateEvents ( const std::vector< V3D > &  E1Vec,
V3D const &  peak_q,
SlimEvents const &  ev_list,
std::vector< V3D > const &  directions,
std::vector< double > const &  sigmas,
bool  specify_size,
double  peak_radius,
double  back_inner_radius,
double  back_outer_radius,
std::vector< double > &  axes_radii,
double &  inti,
double &  sigi,
std::pair< double, double > &  backi 
)
private

Integrate a list of events associated to one peak.

The QLab vector of the events are shifted by the QLab vector of the peak. Spatial distribution of the events in QLab space is described with principal axes of the ellipsoid, as well as the standard deviations in the the directions of the principal axes.

Parameters
E1Vec: Vector of values for calculating edge of detectors
peak_q: The Q-vector for the peak center.
ev_list: List of events centered around the peak (here with Q=[0,0,0]).
directions: The three principal axes of the list of events
sigmas: The standard deviations of the events in the directions of the three principal axes.
specify_size: If true the integration will be done using the ellipsoids with major axes determined by the peak, back_inner and back_outer radii parameters. If false, the integration will be done using a peak region with major axis chosen so that it covers +- three standard deviations of the data in each direction. In this case, the background ellipsoidal shell is chosen to have the same VOLUME as the peak ellipsoid, and to use the peak ellipsoid for the inner radius.
peak_radius: Size of half the major axis of the ellipsoid
back_inner_radius: Size of half the major axis of the INNER ellipsoidal boundary of the background region
back_outer_radius: Size of half the major axis of the OUTER ellipsoidal boundary of the background region
axes_radii: The radii used for integration in the directions of the three principal axes.
inti: (output) net integrated intensity
sigi: (output) estimate of the standard deviation the intensity
backi: (output) background intensity subtracted from inti and sigi

Definition at line 248 of file IntegrateQLabEvents.cpp.

References Mantid::Kernel::Logger::debug(), detectorQ(), Mantid::API::g_log, m_radius, m_useOnePercentBackgroundCorrection, Mantid::Kernel::Logger::notice(), numInEllipsoid(), numInEllipsoidBkg(), and Mantid::Kernel::QLab.

◆ getEigenVectors()

void Mantid::MDAlgorithms::IntegrateQLabEvents::getEigenVectors ( Kernel::DblMatrix const &  cov_matrix,
std::vector< V3D > &  eigen_vectors,
std::vector< double > &  eigen_values 
)
staticprivate

Eigen vectors of a 3x3 real symmetric matrix using the GSL.

Parameters
cov_matrix: 3x3 real symmetric matrix.
eigen_vectors: (output) returned eigen vectors
eigen_values: (output) three eigenvalues

Definition at line 202 of file IntegrateQLabEvents.cpp.

Referenced by ellipseIntegrateEvents().

◆ isOrigin()

bool Mantid::MDAlgorithms::IntegrateQLabEvents::isOrigin ( const V3D q,
const double &  cellSize 
)
static

Determine if an input Q-vector lies in the cell associated to the origin.

Definition at line 68 of file IntegrateQLabEvents.cpp.

Referenced by Mantid::MDAlgorithms::IntegrateEllipsoidsV2::exec(), and Mantid::MDAlgorithms::IntegrateEllipsoidsV2::pairBraggSatellitePeaks().

◆ makeCovarianceMatrix()

void Mantid::MDAlgorithms::IntegrateQLabEvents::makeCovarianceMatrix ( SlimEvents const &  events,
Kernel::DblMatrix matrix,
double  radius 
)
staticprivate

3x3 covariance matrix of a list of SlimEvent objects

the purpose of the covariance matrix is to find the principal axes of the SlimeEvents, associated with a particular peak. Their QLab vectors are already shifted by the QLab vector of the peak. Only events within the specified distance from the peak (here at Q=[0,0,0]) will be used. The covariance matrix can be easily constructed. X, Y, Z of each peak position are the variables we wish to determine the covariance. The mean position in each dimension has already been calculated on subtracted, since this corresponds to the QLab vector peak. The expected values of each correlation test X,X X,Y X,Z e.t.c form the elements of this 3x3 matrix, but since the probabilities are equal, we can remove them from the sums of the expected values, and simply divide by the number of events for each matrix element. Note that the diagonal elements form the variance X,X, Y,Y, Z,Z

Parameters
events: SlimEvents associated to one peak
matrix: (output) 3x3 covariance matrix
radius: Only events within this distance radius of the peak (here at Q=[0,0,0]) are used for calculating the covariance matrix.

Definition at line 182 of file IntegrateQLabEvents.cpp.

References radius.

Referenced by ellipseIntegrateEvents().

◆ numInEllipsoid()

std::pair< double, double > Mantid::MDAlgorithms::IntegrateQLabEvents::numInEllipsoid ( SlimEvents const &  events,
std::vector< V3D > const &  directions,
std::vector< double > const &  sizes 
)
staticprivate

Number of events in an ellipsoid.

The ellipsoid is centered at 0,0,0 with the three specified axes and the three specified sizes in the direction of those axes. NOTE: The three axes must be mutually orthogonal unit vectors.

Parameters
events: List of SlimEvents centered at 0,0,0
directions: List of 3 orthonormal directions for the axes of the ellipsoid.
sizes: List of three values a,b,c giving half the length of the three axes of the ellisoid.
Returns
number of events and estimated error

Definition at line 121 of file IntegrateQLabEvents.cpp.

References count.

Referenced by ellipseIntegrateEvents().

◆ numInEllipsoidBkg()

std::pair< double, double > Mantid::MDAlgorithms::IntegrateQLabEvents::numInEllipsoidBkg ( SlimEvents const &  events,
std::vector< V3D > const &  directions,
std::vector< double > const &  sizes,
std::vector< double > const &  sizesIn,
const bool  useOnePercentBackgroundCorrection 
)
staticprivate

Number of events in an ellipsoid with background correction.

The ellipsoid is centered at 0,0,0 with the three specified axes and the three specified sizes in the direction of those axes. NOTE: The three axes must be mutually orthogonal unit vectors.

Parameters
events: List of 3D events centered at 0,0,0
directions: List of 3 orthonormal directions for the axes of the ellipsoid.
sizes: List of three values a,b,c giving half the length of the three axes of the ellisoid.
sizesIn: List of three values a,b,c giving half the length of the three inner axes of the ellisoid.
useOnePercentBackgroundCorrection: flag if one percent background correction should be used.
Returns
number of events and estimated error

Definition at line 139 of file IntegrateQLabEvents.cpp.

References count.

Referenced by ellipseIntegrateEvents().

◆ populateCellsWithPeaks()

void Mantid::MDAlgorithms::IntegrateQLabEvents::populateCellsWithPeaks ( )

Assign events to each of the cells occupied by events.

Iterate over each QLab cell containing a peak and accumulate the list of events for the cell and for the first-neighbor cells into a single list of events. The QLab vectors for this events are shifted by the QLab vector of the peak.

Definition at line 375 of file IntegrateQLabEvents.cpp.

References Mantid::MDAlgorithms::OccupiedCell::events, m_cellSize, m_cellsWithEvents, m_cellsWithPeaks, Mantid::MDAlgorithms::CellCoords::nearbyCellHashes(), and Mantid::MDAlgorithms::OccupiedCell::peakQ.

Referenced by Mantid::MDAlgorithms::IntegrateEllipsoidsV2::qListFromEventWS(), and Mantid::MDAlgorithms::IntegrateEllipsoidsV2::qListFromHistoWS().

◆ setRadius()

void Mantid::MDAlgorithms::IntegrateQLabEvents::setRadius ( const double &  radius)

Set peak integration radius.

Parameters
radius:: double integration radius. radius must be larger than 0.

Definition at line 66 of file IntegrateQLabEvents.cpp.

References m_radius, and radius.

Referenced by Mantid::MDAlgorithms::IntegrateEllipsoidsV2::exec(), and Mantid::MDAlgorithms::IntegrateEllipsoidsV2::integratePeaksCutoffISigI().

Member Data Documentation

◆ m_cellSize

double Mantid::MDAlgorithms::IntegrateQLabEvents::m_cellSize
private

size of the square cell unit, holding at most one single peak

Definition at line 285 of file IntegrateQLabEvents.h.

Referenced by addEvent(), ellipseIntegrateEvents(), IntegrateQLabEvents(), and populateCellsWithPeaks().

◆ m_cellsWithEvents

std::unordered_map<size_t, SlimEvents> Mantid::MDAlgorithms::IntegrateQLabEvents::m_cellsWithEvents
private

list of cells occupied with events

Definition at line 289 of file IntegrateQLabEvents.h.

Referenced by addEvent(), and populateCellsWithPeaks().

◆ m_cellsWithPeaks

std::unordered_map<size_t, OccupiedCell> Mantid::MDAlgorithms::IntegrateQLabEvents::m_cellsWithPeaks
private

list the occupied cells in an unordered map for fast searching

Definition at line 287 of file IntegrateQLabEvents.h.

Referenced by ellipseIntegrateEvents(), IntegrateQLabEvents(), and populateCellsWithPeaks().

◆ m_radius

double Mantid::MDAlgorithms::IntegrateQLabEvents::m_radius
private

Definition at line 281 of file IntegrateQLabEvents.h.

Referenced by ellipseIntegrateEvents(), IntegrateQLabEvents(), and setRadius().

◆ m_useOnePercentBackgroundCorrection

const bool Mantid::MDAlgorithms::IntegrateQLabEvents::m_useOnePercentBackgroundCorrection
private

if one percent culling of the background should be performed.

Definition at line 283 of file IntegrateQLabEvents.h.

Referenced by ellipseIntegrateEvents().


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