Address PR comments, retain original computeHausdorff utility

PDAL · Feb 22, 2021 · 39131a8 · 39131a8
1 parent 5c81e32
commit 39131a8
Show file tree

Hide file tree

Showing 4 changed files with 44 additions and 7 deletions.
diff --git a/kernels/HausdorffKernel.cpp b/kernels/HausdorffKernel.cpp
@@ -86,7 +86,7 @@ int HausdorffKernel::execute()
     ColumnPointTable candTable;
     PointViewPtr candView = loadSet(m_candidateFile, candTable);
 
-    std::pair<double, double> result = Utils::computeHausdorff(srcView, candView);
+    std::pair<double, double> result = Utils::computeHausdorffPair(srcView, candView);
 
     MetadataNode root;
     root.add("filenames", m_sourceFile);

diff --git a/pdal/PDALUtils.cpp b/pdal/PDALUtils.cpp
@@ -364,8 +364,44 @@ bool fileExists(const std::string& path)
     return FileUtils::fileExists(path);
 }
 
-std::pair<double, double> computeHausdorff(PointViewPtr viewA,
-                                           PointViewPtr viewB)
+double computeHausdorff(PointViewPtr srcView, PointViewPtr candView)
+{
+    using namespace Dimension;
+
+    KD3Index &srcIndex = srcView->build3dIndex();
+    KD3Index &candIndex = candView->build3dIndex();
+
+    double maxDistSrcToCand = std::numeric_limits<double>::lowest();
+    double maxDistCandToSrc = std::numeric_limits<double>::lowest();
+
+    for (PointRef p : *srcView)
+    {
+        PointIdList indices(1);
+        std::vector<double> sqr_dists(1);
+        candIndex.knnSearch(p, 1, &indices, &sqr_dists);
+
+        if (sqr_dists[0] > maxDistSrcToCand)
+            maxDistSrcToCand = sqr_dists[0];
+    }
+
+    for (PointRef q : *candView)
+    {
+        PointIdList indices(1);
+        std::vector<double> sqr_dists(1);
+        srcIndex.knnSearch(q, 1, &indices, &sqr_dists);
+
+        if (sqr_dists[0] > maxDistCandToSrc)
+            maxDistCandToSrc = sqr_dists[0];
+    }
+
+    maxDistSrcToCand = std::sqrt(maxDistSrcToCand);
+    maxDistCandToSrc = std::sqrt(maxDistCandToSrc);
+
+    return (std::max)(maxDistSrcToCand, maxDistCandToSrc);
+}
+
+std::pair<double, double> computeHausdorffPair(PointViewPtr viewA,
+                                               PointViewPtr viewB)
 {
     // Computes both the max and mean of all nearest neighbor distances from
     // each point in the PointView to those in the KD3Index.

diff --git a/pdal/PDALUtils.hpp b/pdal/PDALUtils.hpp
@@ -273,7 +273,8 @@ std::string PDAL_DLL fetchRemote(const std::string& path);
 bool PDAL_DLL isRemote(const std::string& path);
 bool PDAL_DLL fileExists(const std::string& path);
 std::vector<std::string> PDAL_DLL maybeGlob(const std::string& path);
-std::pair<double, double> PDAL_DLL computeHausdorff(PointViewPtr srcView, PointViewPtr candView);
+double PDAL_DLL computeHausdorff(PointViewPtr srcView, PointViewPtr candView);
+std::pair<double, double> PDAL_DLL computeHausdorffPair(PointViewPtr srcView, PointViewPtr candView);
 double PDAL_DLL computeChamfer(PointViewPtr srcView, PointViewPtr candView);
 
 } // namespace Utils

diff --git a/test/unit/apps/HausdorffTest.cpp b/test/unit/apps/HausdorffTest.cpp
@@ -78,20 +78,20 @@ TEST(Hausdorff, distance)
     cand->setField(Dimension::Id::Y, 1, 2.0);
     cand->setField(Dimension::Id::Z, 1, 0.0);
 
-    std::pair<double, double> result = Utils::computeHausdorff(src, cand);
+    std::pair<double, double> result = Utils::computeHausdorffPair(src, cand);
     EXPECT_EQ(2.0, result.first);
 
     cand->setField(Dimension::Id::X, 1, 0.0);
     cand->setField(Dimension::Id::Y, 1, 0.0);
     cand->setField(Dimension::Id::Z, 1, 3.0);
 
-    result = Utils::computeHausdorff(src, cand);
+    result = Utils::computeHausdorffPair(src, cand);
     EXPECT_EQ(3.0, result.first);
 
     src->setField(Dimension::Id::X, 0, 1.0);
     src->setField(Dimension::Id::Y, 0, 1.0);
     src->setField(Dimension::Id::Z, 0, 1.0);
 
-    result = Utils::computeHausdorff(src, cand);
+    result = Utils::computeHausdorffPair(src, cand);
     EXPECT_EQ(std::sqrt(6.0), result.first);
 }