From 97fa168083106c2053ab48f0771dceb18a55838a Mon Sep 17 00:00:00 2001
From: Owen Arnold <owen.arnold@stfc.ac.uk>
Date: Thu, 16 Apr 2015 16:47:56 +0100
Subject: [PATCH] refs #11573. Enhance MDBoxImplicitFunction

We need to be able to determine the fractional overlap between a rectilinear box represented by a set of vertexes and the box implicit function. This is a general feature.

I've added quite a few functional tests for this, but have not done anything about performance yet.
---
 .../MDGeometry/MDBoxImplicitFunction.h        |  19 +-
 .../MDGeometry/MDImplicitFunction.h           |   2 +-
 .../src/MDGeometry/MDBoxImplicitFunction.cpp  |  63 ++++-
 .../Geometry/test/MDBoxImplicitFunctionTest.h | 253 ++++++++++++++++++
 4 files changed, 326 insertions(+), 11 deletions(-)

diff --git a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/MDGeometry/MDBoxImplicitFunction.h b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/MDGeometry/MDBoxImplicitFunction.h
index b03d1f986bf..4f6bca39d81 100644
--- a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/MDGeometry/MDBoxImplicitFunction.h
+++ b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/MDGeometry/MDBoxImplicitFunction.h
@@ -39,7 +39,6 @@ namespace Geometry {
 */
 class DLLExport MDBoxImplicitFunction : public MDImplicitFunction {
 public:
-  MDBoxImplicitFunction();
 
   MDBoxImplicitFunction(const Mantid::Kernel::VMD &min,
                         const Mantid::Kernel::VMD &max);
@@ -47,10 +46,24 @@ public:
   MDBoxImplicitFunction(const std::vector<coord_t> &min,
                         const std::vector<coord_t> &max);
 
-  void construct(const Mantid::Kernel::VMD &min,
-                 const Mantid::Kernel::VMD &max);
+  double volume() const;
+
+  double fraction(coord_t* box, const size_t nVertexes) const;
 
   virtual ~MDBoxImplicitFunction();
+
+private:
+
+  void construct(const Mantid::Kernel::VMD &min,
+                         const Mantid::Kernel::VMD &max);
+
+  /// Maximum extents of MDBox
+  const Mantid::Kernel::VMD m_max;
+  /// Minimum extents of MDBox
+  const Mantid::Kernel::VMD m_min;
+  /// Box volume
+  double m_volume;
+
 };
 
 } // namespace Geometry
diff --git a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/MDGeometry/MDImplicitFunction.h b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/MDGeometry/MDImplicitFunction.h
index 1d710812293..d24d009dcaf 100644
--- a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/MDGeometry/MDImplicitFunction.h
+++ b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/MDGeometry/MDImplicitFunction.h
@@ -245,7 +245,7 @@ public:
    * @param numPoints :: number of vertexes in the array.
    * @return eContact enum value
    */
-  eContact boxContact(const coord_t *vertexes, const size_t numPoints) {
+  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;
diff --git a/Code/Mantid/Framework/Geometry/src/MDGeometry/MDBoxImplicitFunction.cpp b/Code/Mantid/Framework/Geometry/src/MDGeometry/MDBoxImplicitFunction.cpp
index 49ee9bcbc7e..8d5046a56c9 100644
--- a/Code/Mantid/Framework/Geometry/src/MDGeometry/MDBoxImplicitFunction.cpp
+++ b/Code/Mantid/Framework/Geometry/src/MDGeometry/MDBoxImplicitFunction.cpp
@@ -8,11 +8,6 @@ using Mantid::Kernel::VMD;
 namespace Mantid {
 namespace Geometry {
 
-//----------------------------------------------------------------------------------------------
-/** Default ctor
- */
-MDBoxImplicitFunction::MDBoxImplicitFunction() : MDImplicitFunction() {}
-
 //----------------------------------------------------------------------------------------------
 /** Constructor with min/max dimensions.
  *
@@ -24,7 +19,8 @@ MDBoxImplicitFunction::MDBoxImplicitFunction() : MDImplicitFunction() {}
  * @param max :: nd-sized vector of the maximum edge of the box
  */
 MDBoxImplicitFunction::MDBoxImplicitFunction(const std::vector<coord_t> &min,
-                                             const std::vector<coord_t> &max) {
+                                             const std::vector<coord_t> &max)
+    : m_max(max), m_min(min) {
   construct(VMD(min), VMD(max));
 }
 
@@ -39,7 +35,8 @@ MDBoxImplicitFunction::MDBoxImplicitFunction(const std::vector<coord_t> &min,
  * @param max :: nd-sized vector of the maximum edge of the box
  */
 MDBoxImplicitFunction::MDBoxImplicitFunction(const Mantid::Kernel::VMD &min,
-                                             const Mantid::Kernel::VMD &max) {
+                                             const Mantid::Kernel::VMD &max)
+    : m_max(max), m_min(min) {
   construct(min, max);
 }
 
@@ -59,7 +56,10 @@ void MDBoxImplicitFunction::construct(const Mantid::Kernel::VMD &min,
     throw std::invalid_argument(
         "MDBoxImplicitFunction::ctor(): Invalid number of dimensions!");
 
+  double volume = 1;
   for (size_t d = 0; d < nd; d++) {
+    volume *= (max[d] - min[d]);
+
     // Make two parallel planes per dimension
 
     // Normal on the min side, so it faces towards +X
@@ -84,6 +84,7 @@ void MDBoxImplicitFunction::construct(const Mantid::Kernel::VMD &min,
     MDPlane p_max(normal_max, origin_max);
     this->addPlane(p_max);
   }
+  m_volume = volume;
 }
 
 //----------------------------------------------------------------------------------------------
@@ -91,5 +92,53 @@ void MDBoxImplicitFunction::construct(const Mantid::Kernel::VMD &min,
  */
 MDBoxImplicitFunction::~MDBoxImplicitFunction() {}
 
+/**
+ * Calculate volume
+ * @return box volume
+ */
+double MDBoxImplicitFunction::volume() const { return m_volume; }
+
+/**
+ * Calculate the fraction of a box residing inside this implicit function
+ * @param boxVertexes to get fraction for
+ * @return fraction 0 to 1
+ */
+double MDBoxImplicitFunction::fraction(coord_t *boxVertexes,
+                                       const size_t nVertexes) const {
+
+  double fraction = 0;
+  const eContact contact = this->boxContact(boxVertexes, nVertexes);
+  if (contact == NOT_TOUCHING) {
+    fraction = 0; //
+  } else if (contact == CONTAINED) {
+    fraction = 1;
+  } else {
+
+    size_t nd = m_min.size();
+    std::vector<coord_t> dmins(nd, std::numeric_limits<coord_t>::max());
+    std::vector<coord_t> dmaxs(nd, std::numeric_limits<coord_t>::min());
+    coord_t frac = 1;
+
+    for (size_t d = 0; d < nd; ++d) {
+
+      for (size_t i = 0; i < nVertexes; ++i) {
+
+        dmins[d] = std::min(dmins[d], boxVertexes[i + nVertexes * d]);
+        dmaxs[d] = std::max(dmaxs[d], boxVertexes[i + nVertexes * d]);
+      }
+      float dBoxRange = (dmaxs[d] - dmins[d]);
+
+      const coord_t dInnerMin = std::max(m_min[d], dmins[d]);
+      const coord_t dInnerMax = std::min(m_max[d], dmaxs[d]);
+      const coord_t dOverlap = dInnerMax - dInnerMin;
+
+      frac *= dOverlap / dBoxRange;
+    }
+    fraction = frac;
+  }
+  return fraction;
+
+}
+
 } // namespace Mantid
 } // namespace Geometry
diff --git a/Code/Mantid/Framework/Geometry/test/MDBoxImplicitFunctionTest.h b/Code/Mantid/Framework/Geometry/test/MDBoxImplicitFunctionTest.h
index d39923045c3..cd4da9ad6c3 100644
--- a/Code/Mantid/Framework/Geometry/test/MDBoxImplicitFunctionTest.h
+++ b/Code/Mantid/Framework/Geometry/test/MDBoxImplicitFunctionTest.h
@@ -51,6 +51,259 @@ public:
     TS_ASSERT( !try2Dpoint(f, 1.5,2.1) );
   }
 
+
+  void test_volume()
+  {
+    std::vector<coord_t> min, max;
+    min.push_back(0);
+    min.push_back(0);
+    min.push_back(0);
+    max.push_back(1);
+    max.push_back(2);
+    max.push_back(3);
+    MDBoxImplicitFunction box(min, max);
+    TS_ASSERT_EQUALS(1*2*3, box.volume());
+  }
+
+  void test_fraction_when_not_contained()
+  {
+      // The implicit function box
+      const coord_t areaMin = 1;
+      const coord_t areaMax = 2;
+      std::vector<coord_t> min;
+      min.push_back(areaMin);
+      min.push_back(areaMin);
+      std::vector<coord_t> max;
+      max.push_back(areaMax);
+      max.push_back(areaMax);
+      MDBoxImplicitFunction f(min,max);
+
+      // The box to test.
+      const coord_t boxMin = 0;
+      const coord_t boxMax = 0.1;
+      std::vector<coord_t> boxVertexes;
+      // vertex
+      boxVertexes.push_back(boxMin);
+      boxVertexes.push_back(boxMin);
+      // vertex
+      boxVertexes.push_back(boxMin);
+      boxVertexes.push_back(boxMax);
+      // vertex
+      boxVertexes.push_back(boxMax);
+      boxVertexes.push_back(boxMin);
+      // vertex
+      boxVertexes.push_back(boxMax);
+      boxVertexes.push_back(boxMax);
+
+      TS_ASSERT_EQUALS(0.0, f.fraction(&boxVertexes[0], 4));
+  }
+
+  void test_fraction_when_fully_contained()
+  {
+      // The implicit function box
+      const coord_t areaMin = 1;
+      const coord_t areaMax = 2;
+      std::vector<coord_t> min;
+      min.push_back(areaMin);
+      min.push_back(areaMin);
+      std::vector<coord_t> max;
+      max.push_back(areaMax);
+      max.push_back(areaMax);
+      MDBoxImplicitFunction f(min,max);
+
+      // The box to test.
+      const coord_t boxMin = 1.1;
+      const coord_t boxMax = 1.9;
+      std::vector<coord_t> boxVertexes;
+      // vertex
+      boxVertexes.push_back(boxMin);
+      boxVertexes.push_back(boxMin);
+      // vertex
+      boxVertexes.push_back(boxMin);
+      boxVertexes.push_back(boxMax);
+      // vertex
+      boxVertexes.push_back(boxMax);
+      boxVertexes.push_back(boxMin);
+      // vertex
+      boxVertexes.push_back(boxMax);
+      boxVertexes.push_back(boxMax);
+
+      TS_ASSERT_EQUALS(1.0, f.fraction(&boxVertexes[0], 4));
+  }
+
+  void test_fraction_when_partially_contained_1D_simple()
+  {
+      // The implicit function box
+      const coord_t areaMin = 0.9;
+      const coord_t areaMax = 2;
+      std::vector<coord_t> min;
+      min.push_back(areaMin);
+      std::vector<coord_t> max;
+      max.push_back(areaMax);
+      MDBoxImplicitFunction f(min,max);
+
+      // The box to test.
+      const coord_t boxMin = 0;
+      const coord_t boxMax = 1;
+      std::vector<coord_t> boxVertexes;
+      // vertex
+      boxVertexes.push_back(boxMin);
+      // vertex
+      boxVertexes.push_back(boxMax);
+
+      /*
+
+                box to test
+      (x = 0) *------------------* (x = 1)
+
+                            implicit function 1D
+                    (x = 0.9) *--------------------------* (x = 2)
+
+      */
+
+      TSM_ASSERT_DELTA("Overlap fraction is incorrectly calculated", 0.1, f.fraction(&boxVertexes[0], boxVertexes.size()), 1e-3);
+  }
+
+  void test_fraction_when_partially_contained_1D_complex()
+  {
+      // The implicit function box
+      const coord_t areaMin = 0.25;
+      const coord_t areaMax = 0.75;
+      std::vector<coord_t> min;
+      min.push_back(areaMin);
+      std::vector<coord_t> max;
+      max.push_back(areaMax);
+      MDBoxImplicitFunction f(min,max);
+
+      // The box to test.
+      const coord_t boxMin = 0;
+      const coord_t boxMax = 1.0;
+      std::vector<coord_t> boxVertexes;
+      // vertex
+      boxVertexes.push_back(boxMin);
+      // vertex
+      boxVertexes.push_back(boxMax);
+
+      /*
+
+                                    box to test
+      (x = 0) *------------------------------------------------------* (x = 1)
+
+                                 implicit function 1D
+                  (x = 0.25) *--------------------------* (x = 0.75)
+
+      */
+
+      TSM_ASSERT_DELTA("Overlap fraction is incorrectly calculated", 0.5, f.fraction(&boxVertexes[0], boxVertexes.size()), 1e-3);
+  }
+
+
+  void test_fraction_when_partially_contained_2d_simple()
+  {
+
+      /*
+
+        1/4 overlap
+
+                ---------------
+                |             |
+                |             |
+        ---------------       |
+        |       |     |       |
+        |       |     |       |
+        |       ------|--------
+        |             |
+        |             |
+        ---------------
+
+
+
+      */
+
+      // The implicit function box
+      const coord_t areaMin = 0.5;
+      const coord_t areaMax = 1.5;
+      std::vector<coord_t> min;
+      min.push_back(areaMin);
+      min.push_back(areaMin);
+      std::vector<coord_t> max;
+      max.push_back(areaMax);
+      max.push_back(areaMax);
+      MDBoxImplicitFunction f(min,max);
+
+      // The box to test.
+      const coord_t boxMin = 0;
+      const coord_t boxMax = 1;
+      std::vector<coord_t> boxVertexes;
+      // vertex
+      boxVertexes.push_back(boxMin);
+      boxVertexes.push_back(boxMin);
+      // vertex
+      boxVertexes.push_back(boxMin);
+      boxVertexes.push_back(boxMax);
+      // vertex
+      boxVertexes.push_back(boxMax);
+      boxVertexes.push_back(boxMin);
+      // vertex
+      boxVertexes.push_back(boxMax);
+      boxVertexes.push_back(boxMax);
+
+      TSM_ASSERT_DELTA("2d overlap incorrectly calculated", 1.0/4, f.fraction(&boxVertexes[0], 4), 1e-3);
+  }
+
+  void test_fraction_when_partially_contained_2d_complex()
+  {
+
+      /*
+
+        1/8 overlap
+
+                ---------------
+                |  function   |
+                |             |
+        ---------------       |
+        |       |     |       |
+        |       ------|--------
+        |             |
+        |   box       |
+        |             |
+        ---------------
+
+
+
+      */
+
+      // The implicit function box
+      const coord_t areaMin = 0.5;
+      const coord_t areaMax = 1.5;
+      std::vector<coord_t> min;
+      min.push_back(areaMin); // xmin at 0.5
+      min.push_back(areaMin + (areaMin/2)); // ymin at 0.75
+      std::vector<coord_t> max;
+      max.push_back(areaMax);
+      max.push_back(areaMax + (areaMin/2)); // ymax at 0.75
+      MDBoxImplicitFunction f(min,max);
+
+      // The box to test.
+      const coord_t boxMin = 0;
+      const coord_t boxMax = 1;
+      std::vector<coord_t> boxVertexes;
+      // vertex
+      boxVertexes.push_back(boxMin);
+      boxVertexes.push_back(boxMin);
+      // vertex
+      boxVertexes.push_back(boxMin);
+      boxVertexes.push_back(boxMax);
+      // vertex
+      boxVertexes.push_back(boxMax);
+      boxVertexes.push_back(boxMin);
+      // vertex
+      boxVertexes.push_back(boxMax);
+      boxVertexes.push_back(boxMax);
+
+      TSM_ASSERT_DELTA("2d overlap incorrectly calculated", 1.0/8, f.fraction(&boxVertexes[0], 4), 1e-3);
+  }
+
 };
 
 
-- 
GitLab