Mantid
Loading...
Searching...
No Matches
MDPlane.h
Go to the documentation of this file.
1// Mantid Repository : https://github.com/mantidproject/mantid
2//
3// Copyright © 2011 ISIS Rutherford Appleton Laboratory UKRI,
4// NScD Oak Ridge National Laboratory, European Spallation Source,
5// Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
6// SPDX - License - Identifier: GPL - 3.0 +
7#pragma once
8
9#include "MantidGeometry/DllConfig.h"
11#include "MantidKernel/VMD.h"
12#include <vector>
13
14namespace Mantid {
15namespace Geometry {
16
41class MANTID_GEOMETRY_DLL MDPlane {
42public:
43 MDPlane(const Mantid::Kernel::VMD &normal, const Mantid::Kernel::VMD &point);
44 MDPlane(const std::vector<coord_t> &normal, const std::vector<coord_t> &point);
45 MDPlane(const size_t nd, const float *normal, const float *point);
46 MDPlane(const size_t nd, const double *normal, const double *point);
47 MDPlane(const std::vector<Mantid::Kernel::VMD> &vectors, const Mantid::Kernel::VMD &origin,
48 const Mantid::Kernel::VMD &insidePoint);
49 MDPlane(const MDPlane &other);
50 MDPlane &operator=(const MDPlane &other);
51 ~MDPlane();
52
54 size_t getNumDims() const { return m_nd; }
55
57 const coord_t *getNormal() const { return m_normal; }
58
60 coord_t getInequality() const { return m_inequality; }
61
62 // ==================== Methods that are inline for performance
63 // ================================
64 //----------------------------------------------------------------------------------------------
71 inline bool isPointBounded(const coord_t *coords) const {
72 coord_t total = 0;
73 for (size_t d = 0; d < m_nd; ++d) {
74 total += m_normal[d] * coords[d];
75 }
76 return (total >= m_inequality);
77 }
78
79 //----------------------------------------------------------------------------------------------
86 inline bool isPointBounded(const Mantid::Kernel::VMD &coords) const {
87 coord_t total = 0;
88 for (size_t d = 0; d < m_nd; ++d) {
89 total += m_normal[d] * static_cast<coord_t>(coords[d]);
90 }
91 return (total >= m_inequality);
92 }
93
94 //----------------------------------------------------------------------------------------------
101 inline bool isPointBounded(const std::vector<coord_t> &coords) const {
102 coord_t total = 0;
103 for (size_t d = 0; d < m_nd; ++d) {
104 total += m_normal[d] * coords[d];
105 }
106 return (total >= m_inequality);
107 }
108
109 //----------------------------------------------------------------------------------------------
119 inline bool isPointInside(const coord_t *coords) const {
120 coord_t total = 0;
121 for (size_t d = 0; d < m_nd; ++d) {
122 total += m_normal[d] * coords[d];
123 }
124 return (total > m_inequality);
125 }
126
127 //----------------------------------------------------------------------------------------------
137 inline bool isPointInside(const std::vector<coord_t> &coords) const {
138 coord_t total = 0;
139 for (size_t d = 0; d < m_nd; ++d) {
140 total += m_normal[d] * coords[d];
141 }
142 return (total > m_inequality);
143 }
144
145 //----------------------------------------------------------------------------------------------
154 bool doesLineIntersect(const coord_t *pointA, const coord_t *pointB) const {
155 bool AisBounded = isPointBounded(pointA);
156 bool BisBounded = isPointBounded(pointB);
157 // The line crosses the plane if one point is bounded and not the other.
158 // Simple! :)
159 return (AisBounded != BisBounded);
160 }
161
162private:
163 //----------------------------------------------------------------------------------------------
170 template <typename IterableType> void construct(IterableType normal, IterableType point) {
171 m_normal = new coord_t[m_nd];
172 m_inequality = 0;
173 for (size_t d = 0; d < m_nd; ++d) {
174 m_normal[d] = static_cast<coord_t>(normal[d]);
175 m_inequality += static_cast<coord_t>(point[d]) * m_normal[d];
176 }
177 }
178
179protected:
181 size_t m_nd;
182
187
190};
191
192} // namespace Geometry
193} // namespace Mantid
A generalized description of a N-dimensional hyperplane.
Definition: MDPlane.h:41
bool isPointBounded(const Mantid::Kernel::VMD &coords) const
Is a point in MDimensions bounded by this hyperplane, that is, is (a1*x1 + a2*x2 + ....
Definition: MDPlane.h:86
coord_t getInequality() const
Definition: MDPlane.h:60
coord_t m_inequality
Right-hand side of the linear equation. aka b in : a1*x1 + a2*x2 + ... < b.
Definition: MDPlane.h:189
bool isPointBounded(const std::vector< coord_t > &coords) const
Is a point in MDimensions bounded by this hyperplane, that is, is (a1*x1 + a2*x2 + ....
Definition: MDPlane.h:101
bool doesLineIntersect(const coord_t *pointA, const coord_t *pointB) const
Given two points defining the start and end point of a line, is there an intersection between the hyp...
Definition: MDPlane.h:154
bool isPointInside(const std::vector< coord_t > &coords) const
Is a point in MDimensions bounded by this hyperplane, that is, is (a1*x1 + a2*x2 + ....
Definition: MDPlane.h:137
bool isPointInside(const coord_t *coords) const
Is a point in MDimensions bounded by this hyperplane, that is, is (a1*x1 + a2*x2 + ....
Definition: MDPlane.h:119
void construct(IterableType normal, IterableType point)
Templated method for constructing the MDPlane.
Definition: MDPlane.h:170
size_t m_nd
Number of dimensions in the MDEventWorkspace.
Definition: MDPlane.h:181
bool isPointBounded(const coord_t *coords) const
Is a point in MDimensions bounded by this hyperplane, that is, is (a1*x1 + a2*x2 + ....
Definition: MDPlane.h:71
coord_t * m_normal
Coefficients to multiply each coordinate, sized m_nd.
Definition: MDPlane.h:186
const coord_t * getNormal() const
Definition: MDPlane.h:57
size_t getNumDims() const
Definition: MDPlane.h:54
Helper class which provides the Collimation Length for SANS instruments.
float coord_t
Typedef for the data type to use for coordinate axes in MD objects such as MDBox, MDEventWorkspace,...
Definition: MDTypes.h:27