MDImplicitFunction.h

Go to the documentation of this file.

// Mantid Repository : https://github.com/mantidproject/mantid
//
// Copyright © 2011 ISIS Rutherford Appleton Laboratory UKRI,
//   NScD Oak Ridge National Laboratory, European Spallation Source,
//   Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
// SPDX - License - Identifier: GPL - 3.0 +
#pragma once
 
#include "MantidGeometry/MDGeometry/MDPlane.h"
#include "MantidGeometry/MDGeometry/MDTypes.h"
#include "MantidKernel/Exception.h"
#include <algorithm>
#include <memory>
#include <string>
#include <vector>
 
namespace Mantid {
namespace Geometry {
 
class MANTID_GEOMETRY_DLL MDImplicitFunction {
public:
  MDImplicitFunction();
  virtual ~MDImplicitFunction() = default;
 
  void addPlane(const MDPlane &plane);
 
  const MDPlane &getPlane(size_t index) const { return m_planes[index]; }
 
  size_t getNumDims() const { return m_nd; }
 
  size_t getNumPlanes() const { return m_planes.size(); }
 
  virtual std::string getName() const {
    throw std::runtime_error("Cannot call MDImplicitFunction does not implement getName()");
  }
 
  virtual std::string toXMLString() const {
    throw std::runtime_error("Cannot call MDImplicitFunction does not implement toXMLString()");
  }
 
  //----------------------------------------------------------------------------------------------
  enum eContact {
    NOT_TOUCHING = 0,
    TOUCHING = 1,
    CONTAINED = 2
  };
 
  // ==================== Methods that are inline for performance
  // ================================
  //----------------------------------------------------------------------------------------------
  virtual bool isPointContained(const coord_t *coords) {
    const auto it = std::find_if(m_planes.begin(), m_planes.end(),
                                 [coords](const MDPlane &plane) { return !plane.isPointBounded(coords); });
    return it == m_planes.end();
  }
 
  //----------------------------------------------------------------------------------------------
  virtual bool isPointContained(const Mantid::Kernel::VMD &coords) {
    const auto it = std::find_if(m_planes.begin(), m_planes.end(),
                                 [coords](const MDPlane &plane) { return !plane.isPointBounded(coords); });
    return it == m_planes.end();
  }
 
  //----------------------------------------------------------------------------------------------
  virtual bool isPointContained(const std::vector<coord_t> &coords) {
    const auto it = std::find_if(m_planes.begin(), m_planes.end(),
                                 [coords](const MDPlane &plane) { return !plane.isPointBounded(coords); });
    return it == m_planes.end();
  }
 
  //----------------------------------------------------------------------------------------------
  bool isBoxTouching(const std::vector<std::vector<coord_t>> &vertexes) {
    // The first plane with NO points inside it means the box does NOT touch.
    return std::all_of(m_planes.cbegin(), m_planes.cbegin() + m_numPlanes, [&vertexes](const MDPlane &plane) {
      return std::any_of(vertexes.cbegin(), vertexes.cend(),
                         [&plane](const std::vector<coord_t> &vertex) { return plane.isPointBounded(vertex); });
    });
  }
 
  //----------------------------------------------------------------------------------------------
  bool isBoxTouching(const coord_t *vertexes, const size_t numPoints) {
    // As the description states, the first plane with NO points inside it
    // means the box does NOT touch. So iterate by planes
    for (size_t i = 0; i < m_numPlanes; i++) {
      size_t numBounded = 0;
      for (size_t j = 0; j < numPoints; j++) {
        if (m_planes[i].isPointBounded(vertexes + j * m_nd)) {
          numBounded++;
          // No need to evaluate any more points.
          break;
        }
      }
      // Not a single point is in this plane
      if (numBounded == 0)
        // That means the box CANNOT touch the implicit function
        return false;
    }
    return true;
  }
 
  //----------------------------------------------------------------------------------------------
  eContact boxContact(const coord_t *vertexes, const size_t numPoints) const {
    // For speed, we can stop looking when we know the box CANNOT be fully
    // contained.
    bool lookForFullyContained = true;
 
    // As the description states, the first plane with NO points inside it
    // means the box does NOT touch. So iterate by planes
    for (size_t i = 0; i < m_numPlanes; i++) {
      size_t numBounded = 0;
      for (size_t j = 0; j < numPoints; j++) {
        if (m_planes[i].isPointBounded(vertexes + j * m_nd)) {
          numBounded++;
          // No need to evaluate any more points, unless we look for fully
          // contained
          if (!lookForFullyContained)
            break;
        } else
          // One of the vertexes is not contained by one of the planes.
          // This means that the box CANNOT be fully contained.
          lookForFullyContained = false;
      }
      // Not a single point is in this plane
      if (numBounded == 0)
        // That means the box CANNOT touch the implicit function
        return NOT_TOUCHING;
      // If all points were within this plane, then there is still a chance that
      // the box is fully contained
      if (numBounded != numPoints)
        lookForFullyContained = false;
    }
    // If nothing said that the box might not be fully contained, then it is
    // fully contained!
    if (lookForFullyContained)
      return CONTAINED;
    else
      return TOUCHING;
  }
 
protected:
  size_t m_nd;
 
  std::vector<MDPlane> m_planes;
 
  size_t m_numPlanes;
};
 
using MDImplicitFunction_sptr = std::shared_ptr<MDImplicitFunction>;
 
} // namespace Geometry
} // namespace Mantid