Commit e7deea7e authored by Antti Soininen's avatar Antti Soininen
Browse files

Use the normalize() free function where applicable. Re #25373

parent 117de286
......@@ -64,7 +64,6 @@ void DetectorSearcher::createDetectorCache() {
const auto frame = m_instrument->getReferenceFrame();
auto beam = frame->vecPointingAlongBeam();
auto up = frame->vecPointingUp();
beam.normalize();
for (size_t pointNo = 0; pointNo < m_detInfo.size(); ++pointNo) {
if (m_detInfo.isMonitor(pointNo) || m_detInfo.isMasked(pointNo))
......@@ -72,8 +71,7 @@ void DetectorSearcher::createDetectorCache() {
// Calculate a unit Q vector for each detector
// This follows a method similar to that used in IntegrateEllipsoids
auto pos = m_detInfo.position(pointNo);
pos.normalize();
const auto pos = normalize(m_detInfo.position(pointNo));
auto E1 = (pos - beam) * -m_crystallography_convention;
const auto norm = E1.norm();
if (norm == 0.) {
......@@ -196,13 +194,11 @@ DetectorSearcher::DetectorSearchResult DetectorSearcher::handleTubeGap(
auto gapDir = V3D(0., 0., 0.);
gapDir[i] = gap;
auto beam1 = detectorDir + gapDir;
beam1.normalize();
auto beam1 = normalize(detectorDir + gapDir);
const auto result1 = checkInteceptWithNeighbours(beam1, neighbours);
const auto hit1 = std::get<0>(result1);
auto beam2 = detectorDir - gapDir;
beam2.normalize();
const auto beam2 = normalize(detectorDir - gapDir);
const auto result2 = checkInteceptWithNeighbours(beam2, neighbours);
const auto hit2 = std::get<0>(result2);
......
......@@ -331,8 +331,7 @@ void AbsorptionCorrection::calculateDistances(const IDetector &detector,
for (size_t i = 0; i < m_numVolumeElements; ++i) {
// Create track for distance in cylinder between scattering point and
// detector
V3D direction = detectorPos - m_elementPositions[i];
direction.normalize();
const V3D direction = normalize(detectorPos - m_elementPositions[i]);
Track outgoing(m_elementPositions[i], direction);
int temp = m_sampleObject->interceptSurface(outgoing);
......
......@@ -245,8 +245,7 @@ void DetectorEfficiencyCor::correctForEfficiency(
// now get the sin of the angle, it's the magnitude of the cross product of
// unit vector along the detector tube axis and a unit vector directed from
// the sample to the detector centre
V3D vectorFromSample = det_member.getPos() - m_samplePos;
vectorFromSample.normalize();
const V3D vectorFromSample = normalize(det_member.getPos() - m_samplePos);
Quat rot = det_member.getRotation();
// rotate the original cylinder object axis to get the detector axis in the
// actual instrument
......@@ -377,9 +376,7 @@ void DetectorEfficiencyCor::getDetectorGeometry(const Geometry::IDetector &det,
double DetectorEfficiencyCor::distToSurface(const V3D &start,
const IObject *shape) const {
// get a vector from the point that was passed to the origin
V3D direction = V3D(0.0, 0.0, 0.0) - start;
// it needs to be a unit vector
direction.normalize();
const V3D direction = normalize(-start);
// put the point and the vector (direction) together to get a line, here
// called a track
Track track(start, direction);
......
......@@ -226,8 +226,7 @@ He3TubeEfficiency::calculateExponential(std::size_t spectraIndex,
// now get the sin of the angle, it's the magnitude of the cross product of
// unit vector along the detector tube axis and a unit vector directed from
// the sample to the detector center
Kernel::V3D vectorFromSample = idet.getPos() - m_samplePos;
vectorFromSample.normalize();
const Kernel::V3D vectorFromSample = normalize(idet.getPos() - m_samplePos);
Kernel::Quat rot = idet.getRotation();
// rotate the original cylinder object axis to get the detector axis in the
// actual instrument
......@@ -324,9 +323,7 @@ void He3TubeEfficiency::getDetectorGeometry(const Geometry::IDetector &det,
double He3TubeEfficiency::distToSurface(const Kernel::V3D start,
const Geometry::IObject *shape) const {
// get a vector from the point that was passed to the origin
Kernel::V3D direction = Kernel::V3D(0.0, 0.0, 0.0) - start;
// it needs to be a unit vector
direction.normalize();
const Kernel::V3D direction = normalize(-start);
// put the point and the vector (direction) together to get a line,
// here called a track
Geometry::Track track(start, direction);
......
......@@ -101,8 +101,7 @@ double MCInteractionVolume::calculateAbsorption(
scatterPos = m_sample->generatePointInObject(rng, m_activeRegion,
m_maxScatterAttempts);
}
auto toStart = startPos - scatterPos;
toStart.normalize();
const auto toStart = normalize(startPos - scatterPos);
Track beforeScatter(scatterPos, toStart);
int nlinks = m_sample->interceptSurface(beforeScatter);
if (m_env) {
......@@ -130,8 +129,7 @@ double MCInteractionVolume::calculateAbsorption(
};
// Now track to final destination
V3D scatteredDirec = endPos - scatterPos;
scatteredDirec.normalize();
const V3D scatteredDirec = normalize(endPos - scatterPos);
Track afterScatter(scatterPos, scatteredDirec);
m_sample->interceptSurface(afterScatter);
if (m_env) {
......
......@@ -51,9 +51,9 @@ void SofQWCentre::exec() {
const auto &detectorInfo = inputWorkspace->detectorInfo();
const auto &spectrumInfo = inputWorkspace->spectrumInfo();
V3D beamDir = detectorInfo.samplePosition() - detectorInfo.sourcePosition();
beamDir.normalize();
double l1 = detectorInfo.l1();
const V3D beamDir =
normalize(detectorInfo.samplePosition() - detectorInfo.sourcePosition());
const double l1 = detectorInfo.l1();
g_log.debug() << "Source-sample distance: " << l1 << '\n';
// Loop over input workspace bins, reassigning data to correct bin in output
......@@ -85,9 +85,8 @@ void SofQWCentre::exec() {
try {
size_t idet = detectorInfo.indexOf(detID);
// Calculate kf vector direction and then Q for each energy bin
V3D scatterDir =
(detectorInfo.position(idet) - detectorInfo.samplePosition());
scatterDir.normalize();
const V3D scatterDir = normalize(detectorInfo.position(idet) -
detectorInfo.samplePosition());
for (size_t j = 0; j < numBins; ++j) {
if (X[j] < xAxis.front() || X[j + 1] > xAxis.back())
continue;
......
......@@ -102,11 +102,9 @@ SXPeak::SXPeak(double t, double phi, double intensity,
const auto sourcePos = spectrumInfo.sourcePosition();
const auto detPos = spectrumInfo.position(m_wsIndex);
// Normalized beam direction
auto beamDir = samplePos - sourcePos;
beamDir.normalize();
const auto beamDir = normalize(samplePos - sourcePos);
// Normalized detector direction
auto detDir = (detPos - samplePos);
detDir.normalize();
const auto detDir = normalize(detPos - samplePos);
m_unitWaveVector = beamDir - detDir;
m_qConvention = Kernel::ConfigService::Instance().getString("Q.convention");
}
......
......@@ -97,8 +97,7 @@ bool PeaksOnSurface::pointOutsideAnyExtents(const V3D &) const { return true; }
bool lineIntersectsSphere(const V3D &line, const V3D &lineStart,
const V3D &peakCenter, const double peakRadius) {
V3D peakToStart = peakCenter - lineStart;
V3D unitLine = line;
unitLine.normalize();
const V3D unitLine = normalize(line);
double proj = peakToStart.scalar_prod(unitLine); // All we are doing here is
// projecting the peak to
// segment start vector onto
......
......@@ -741,8 +741,7 @@ V3D VesuvioCalculateMS::generateDetectorPos(
nominalPos[m_upIdx] + (m_randgen->flat() - 0.5) * m_detHeight;
// Distance to exit the sample for this order
V3D scToDet = detPos - scatterPt;
scToDet.normalize();
const V3D scToDet = normalize(detPos - scatterPt);
Geometry::Track scatterToDet(scatterPt, scToDet);
if (m_sampleShape->interceptSurface(scatterToDet) > 0) {
scang = direcBeforeSc.angle(scToDet);
......
......@@ -669,8 +669,7 @@ bool Peak::findDetector() {
*/
bool Peak::findDetector(const InstrumentRayTracer &tracer) {
// Scattered beam direction
V3D beam = detPos - samplePos;
beam.normalize();
const V3D beam = normalize(detPos - samplePos);
return findDetector(beam, tracer);
}
......
......@@ -1834,8 +1834,7 @@ bool IndexingUtils::FormUB_From_abc_Vectors(DblMatrix &UB,
V3D c_dir;
bool c_found = false;
V3D perp = a_dir.cross_prod(b_dir);
perp.normalize();
const V3D perp = normalize(a_dir.cross_prod(b_dir));
double perp_ang;
double alpha;
double beta;
......@@ -2749,8 +2748,7 @@ std::vector<V3D> IndexingUtils::MakeCircleDirections(int n_steps,
V3D second_vec(0, 0, 0);
second_vec[min_index] = 1;
V3D perp_vec = second_vec.cross_prod(axis);
perp_vec.normalize();
const V3D perp_vec = normalize(second_vec.cross_prod(axis));
// next get a vector that is the specified
// number of degrees away from the axis
......
......@@ -427,9 +427,7 @@ IComponent_const_sptr Instrument::getSample() const {
* @returns A unit vector denoting the direction of the beam
*/
Kernel::V3D Instrument::getBeamDirection() const {
V3D retval = getSample()->getPos() - getSource()->getPos();
retval.normalize();
return retval;
return normalize(getSample()->getPos() - getSource()->getPos());
}
//------------------------------------------------------------------------------------------
......
......@@ -2034,9 +2034,10 @@ void InstrumentDefinitionParser::makeXYplaneFaceComponent(
// vector from facing object to component we want to rotate
Kernel::V3D facingDirection = pos - facingPoint;
if (facingDirection.norm() == 0.0)
const auto facingDirLength = facingDirection.norm();
if (facingDirLength == 0.0)
return;
facingDirection.normalize();
facingDirection /= facingDirLength;
// now aim to rotate shape such that the z-axis of of the object we want to
// rotate points in the direction of facingDirection. That way the XY plane
......@@ -2047,10 +2048,12 @@ void InstrumentDefinitionParser::makeXYplaneFaceComponent(
R.rotate(facingDirection);
Kernel::V3D normal = facingDirection.cross_prod(z);
if (normal.norm() == 0.) {
normal = -facingDirection;
const auto normalLength = normal.norm();
if (normalLength == 0.) {
normal = normalize(-facingDirection);
} else {
normal /= normalLength;
}
normal.normalize();
double theta = (180.0 / M_PI) * facingDirection.angle(z);
if (normal.norm() > 0.0)
......
......@@ -194,11 +194,8 @@ with the along beam vector.
@return result of whther the along beam and test vector are parallel.
*/
bool ReferenceFrame::isVectorPointingAlongBeam(const V3D &v) const {
V3D vec = v;
vec.normalize();
// Normalized (unit) parallel vectors should produce a scalar product of 1
return m_vecPointingAlongBeam.scalar_prod(vec) == 1;
return m_vecPointingAlongBeam.scalar_prod(normalize(v)) == 1;
}
} // namespace Geometry
......
......@@ -309,10 +309,10 @@ double cylinderSolidAngle(const V3D &observer, const V3D &centre,
// facing away from the observer gives a negative solid angle and is excluded
// For simplicity the triangulation points are constructed such that the cone
// axis points up the +Z axis and then rotated into their final position
const V3D axis_direction = normalize(axis);
// Required rotation
constexpr V3D initial_axis(0., 0., 1.0);
const Quat transform(initial_axis, axis_direction);
const Quat transform(initial_axis, axis);
// Do the base cap which is a point at the centre and nslices points around it
constexpr double angle_step =
......@@ -770,12 +770,10 @@ bool CSGObject::isOnSide(const Kernel::V3D &point) const {
return true;
}
}
Kernel::V3D NormPair;
for (auto xs = Snorms.begin(); xs != Snorms.end(); ++xs)
for (auto ys = std::next(xs); ys != Snorms.end(); ++ys) {
NormPair = (*ys) + (*xs);
try {
NormPair.normalize();
const V3D NormPair = normalize((*ys) + (*xs));
if (!checkSurfaceValid(point, NormPair))
return true;
} catch (std::runtime_error &) {
......
......@@ -418,8 +418,7 @@ ShapeFactory::parseCylinder(Poco::XML::Element *pElem,
Element *pElemRadius = getShapeElement(pElem, "radius");
Element *pElemHeight = getShapeElement(pElem, "height");
V3D normVec = parsePosition(pElemAxis);
normVec.normalize();
const V3D normVec = normalize(parsePosition(pElemAxis));
// getDoubleAttribute can throw - put the calls above any new
const double radius = getDoubleAttribute(pElemRadius, "val");
......@@ -475,8 +474,7 @@ std::string ShapeFactory::parseSegmentedCylinder(
Element *pElemRadius = getShapeElement(pElem, "radius");
Element *pElemHeight = getShapeElement(pElem, "height");
V3D normVec = parsePosition(pElemAxis);
normVec.normalize();
const V3D normVec = normalize(parsePosition(pElemAxis));
// getDoubleAttribute can throw - put the calls above any new
const double radius = getDoubleAttribute(pElemRadius, "val");
......@@ -534,8 +532,7 @@ std::string ShapeFactory::parseHollowCylinder(
Element *pElemOuterRadius = getShapeElement(pElem, "outer-radius");
Element *pElemHeight = getShapeElement(pElem, "height");
V3D normVec = parsePosition(pElemAxis);
normVec.normalize();
const V3D normVec = normalize(parsePosition(pElemAxis));
const double innerRadius = getDoubleAttribute(pElemInnerRadius, "val");
if (innerRadius <= 0.0) {
throw std::runtime_error(
......@@ -667,8 +664,7 @@ CuboidCorners ShapeFactory::parseCuboid(Poco::XML::Element *pElem) {
// Use a quarternion to do a rotation for us, with respect to the default
// axis. Our "Quat" implementation requires that the vectors passed to
// it be normalised.
V3D axis = parsePosition(pElem_axis);
axis.normalize();
const V3D axis = normalize(parsePosition(pElem_axis));
const Quat rotation(DEFAULT_AXIS, axis);
rotation.rotate(result.lfb);
......@@ -704,8 +700,7 @@ ShapeFactory::parseCuboid(Poco::XML::Element *pElem,
int &l_id) {
auto corners = parseCuboid(pElem);
V3D pointTowardBack = corners.lbb - corners.lfb;
pointTowardBack.normalize();
const V3D pointTowardBack = normalize(corners.lbb - corners.lfb);
// add front plane cutoff
auto pPlaneFrontCutoff = boost::make_shared<Plane>();
......@@ -731,8 +726,7 @@ ShapeFactory::parseCuboid(Poco::XML::Element *pElem,
retAlgebraMatch << "-" << l_id << " ";
l_id++;
V3D pointTowardRight = corners.rfb - corners.lfb;
pointTowardRight.normalize();
const V3D pointTowardRight = normalize(corners.rfb - corners.lfb);
// add left plane cutoff
auto pPlaneLeftCutoff = boost::make_shared<Plane>();
......@@ -756,8 +750,7 @@ ShapeFactory::parseCuboid(Poco::XML::Element *pElem,
retAlgebraMatch << "-" << l_id << " ";
l_id++;
V3D pointTowardTop = corners.lft - corners.lfb;
pointTowardTop.normalize();
const V3D pointTowardTop = normalize(corners.lft - corners.lfb);
// add bottom plane cutoff
auto pPlaneBottomCutoff = boost::make_shared<Plane>();
......@@ -803,8 +796,7 @@ ShapeFactory::parseInfiniteCone(Poco::XML::Element *pElem,
Element *pElemAxis = getShapeElement(pElem, "axis");
Element *pElemAngle = getShapeElement(pElem, "angle");
V3D normVec = parsePosition(pElemAxis);
normVec.normalize();
const V3D normVec = normalize(parsePosition(pElemAxis));
// getDoubleAttribute can throw - put the calls above any new
const double angle = getDoubleAttribute(pElemAngle, "val");
......@@ -850,8 +842,7 @@ ShapeFactory::parseCone(Poco::XML::Element *pElem,
Element *pElemAngle = getShapeElement(pElem, "angle");
Element *pElemHeight = getShapeElement(pElem, "height");
V3D normVec = parsePosition(pElemAxis);
normVec.normalize();
const V3D normVec = normalize(parsePosition(pElemAxis));
// getDoubleAttribute can throw - put the calls above any new
const double angle = getDoubleAttribute(pElemAngle, "val");
......@@ -896,14 +887,9 @@ ShapeFactory::parseCone(Poco::XML::Element *pElem,
std::string ShapeFactory::parseHexahedronFromStruct(
Hexahedron &hex, std::map<int, boost::shared_ptr<Surface>> &prim,
int &l_id) {
V3D pointTowardBack = hex.lbb - hex.lfb;
pointTowardBack.normalize();
V3D normal;
// add front face
auto pPlaneFrontCutoff = boost::make_shared<Plane>();
normal = (hex.rfb - hex.lfb).cross_prod(hex.lft - hex.lfb);
auto normal = (hex.rfb - hex.lfb).cross_prod(hex.lft - hex.lfb);
// V3D jjj = (normal*(rfb-rbb));
if (normal.scalar_prod(hex.rfb - hex.rbb) < 0)
......@@ -1070,8 +1056,9 @@ ShapeFactory::parseTaperedGuide(Poco::XML::Element *pElem,
// For centre and axis we allow defaults.
const V3D centre = pElemCentre ? parsePosition(pElemCentre) : DEFAULT_CENTRE;
V3D axis = pElemAxis ? parsePosition(pElemAxis) : DEFAULT_AXIS;
axis.normalize(); // Quat requires normalised axes.
// Quat requires normalised axes.
const V3D axis =
normalize(pElemAxis ? parsePosition(pElemAxis) : DEFAULT_AXIS);
const double apertureStartWidth =
getDoubleAttribute(pElemApertureStart, "width");
......@@ -1146,8 +1133,7 @@ ShapeFactory::parseTorus(Poco::XML::Element *pElem,
getShapeElement(pElem, "radius-from-centre-to-tube");
Element *pElemRadiusTube = getShapeElement(pElem, "radius-tube");
V3D normVec = parsePosition(pElemAxis);
normVec.normalize();
const V3D normVec = normalize(parsePosition(pElemAxis));
// getDoubleAttribute can throw - put the calls above any new
const double radiusCentre = getDoubleAttribute(pElemRadiusFromCentre, "val");
......@@ -1453,8 +1439,7 @@ void ShapeFactory::createGeometryHandler(Poco::XML::Element *pElem,
Element *pElemAxis = getShapeElement(pElem, "axis");
Element *pElemRadius = getShapeElement(pElem, "radius");
Element *pElemHeight = getShapeElement(pElem, "height");
V3D normVec = parsePosition(pElemAxis);
normVec.normalize();
const V3D normVec = normalize(parsePosition(pElemAxis));
shapeInfo.setCylinder(parsePosition(pElemCentre), normVec,
std::stod(pElemRadius->getAttribute("val")),
std::stod(pElemHeight->getAttribute("val")));
......@@ -1464,10 +1449,9 @@ void ShapeFactory::createGeometryHandler(Poco::XML::Element *pElem,
Element *pElemAngle = getShapeElement(pElem, "angle");
Element *pElemHeight = getShapeElement(pElem, "height");
V3D normVec = parsePosition(pElemAxis);
normVec.normalize();
double height = std::stod(pElemHeight->getAttribute("val"));
double radius =
const V3D normVec = normalize(parsePosition(pElemAxis));
const double height = std::stod(pElemHeight->getAttribute("val"));
const double radius =
height * tan(M_PI * std::stod(pElemAngle->getAttribute("val")) / 180.0);
shapeInfo.setCone(parsePosition(pElemTipPoint), normVec, radius, height);
}
......
......@@ -216,8 +216,7 @@ Raster calculateCylinder(const V3D &beamDirection, const CSGObject &shape,
// Assume that z' = axis. Then select whatever has the smallest dot product
// with axis to be the x' direction
V3D z_prime = params.axis;
z_prime.normalize();
const V3D z_prime = normalize(params.axis);
const V3D x_prime = createPerpendicular(z_prime);
const V3D y_prime = z_prime.cross_prod(x_prime);
......
......@@ -80,17 +80,15 @@ void render(detail::GeometryTriangulator &triangulator) {
const auto &faces = triangulator.getTriangleFaces();
const auto &points = triangulator.getTriangleVertices();
glBegin(GL_TRIANGLES);
V3D normal;
for (size_t i = 0; i < triangulator.numTriangleFaces(); i++) {
auto index2 = static_cast<size_t>(faces[i * 3 + 1] * 3);
auto index3 = static_cast<size_t>(faces[i * 3 + 2] * 3);
auto index1 = static_cast<size_t>(faces[i * 3] * 3);
// Calculate normal and normalize
V3D v1(points[index1], points[index1 + 1], points[index1 + 2]);
V3D v2(points[index2], points[index2 + 1], points[index2 + 2]);
V3D v3(points[index3], points[index3 + 1], points[index3 + 2]);
normal = (v1 - v2).cross_prod(v2 - v3);
normal.normalize();
const V3D v1(points[index1], points[index1 + 1], points[index1 + 2]);
const V3D v2(points[index2], points[index2 + 1], points[index2 + 2]);
const V3D v3(points[index3], points[index3 + 1], points[index3 + 2]);
const auto normal = normalize((v1 - v2).cross_prod(v2 - v3));
glNormal3d(normal[0], normal[1], normal[2]);
glVertex3dv(&points[index1]);
glVertex3dv(&points[index2]);
......@@ -177,13 +175,11 @@ void renderCuboid(const detail::ShapeInfo &shapeInfo) {
{0, 4, 5, 1}, // front
{4, 7, 6, 5}, // bottom
};
V3D normal;
// first face
glBegin(GL_QUADS);
for (auto &row : faceindex) {
normal = (vertex[row[0]] - vertex[row[1]])
.cross_prod((vertex[row[0]] - vertex[row[2]]));
normal.normalize();
const auto normal = normalize((vertex[row[0]] - vertex[row[1]])
.cross_prod((vertex[row[0]] - vertex[row[2]])));
glNormal3d(normal[0], normal[1], normal[2]);
for (const int ij : row) {
if (ij == 0)
......
......@@ -139,9 +139,9 @@ int Plane::setPlane(const Kernel::V3D &P, const Kernel::V3D &N)
@retval 0 :: success
*/
{
NormV = N;
NormV.normalize();
if (NormV.norm2() == 0.0) {
try {
NormV = normalize(N);
} catch (std::runtime_error &) {
throw std::invalid_argument("Attempt to create Plane with zero normal");
}
Dist = P.scalar_prod(NormV);
......
......@@ -204,9 +204,9 @@ void Torus::setNorm(const Kernel::V3D &A)
@param A :: New Normal direction
*/
{
if (A.norm() > Tolerance) {
Normal = A;
Normal.normalize();
const auto length = A.norm();
if (length > Tolerance) {
Normal = A / length;
}
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment