diff --git a/kernels/HausdorffKernel.cpp b/kernels/HausdorffKernel.cpp
index d4a23c80c1..50e6aee77b 100644
--- 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
index bbabb82332..d9c158de6e 100644
--- 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
index 62632a69d6..4c72255274 100644
--- 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
index 10756cc2d5..676af7d853 100644
--- 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);
 }