Properly handle Bounds as both 2D and 3D in reprojection (#3117)

* Fix "Bounds" reprojection.

* Add test.

pdal/GDALUtils.cpp

@@ -163,9 +163,8 @@ GDALDataType toGdalType(Dimension::Type t)
   \param x  X coordinate of point to be reprojected in-place.
   \param y  Y coordinate of point to be reprojected in-place.
   \param z  Z coordinate of point to be reprojected in-place.
-  \param srcSrs  String in WKT or other suitable format of box coordinates.
-  \param dstSrs  String in WKT or other suitable format to which
-    coordinates should be projected.
+  \param srcSrs  Source SRS
+  \param dstSrs  Destination SRS
   \return  Whether the reprojection was successful or not.
 */
 bool reproject(double& x, double& y, double& z, const SpatialReference& srcSrs,
@@ -178,9 +177,8 @@ bool reproject(double& x, double& y, double& z, const SpatialReference& srcSrs,
 /**
   Reproject a bounds box from a source projection to a destination.
   \param box  Bounds box to be reprojected in-place.
-  \param srcSrs  String in WKT or other suitable format of box coordinates.
-  \param dstSrs  String in WKT or other suitable format to which
-    coordinates should be projected.
+  \param srcSrs  Source SRS.
+  \param dstSrs  Destination SRS.
   \return  Whether the reprojection was successful or not.
 */
 bool reprojectBounds(BOX3D& box, const SpatialReference& srcSrs,
@@ -195,24 +193,37 @@ bool reprojectBounds(BOX3D& box, const SpatialReference& srcSrs,
 }
 
 
+/**
+  Reproject a bounds box from a source projection to a destination.
+  \param box  2D or 3D bounds box to be reprojected.
+  \param srcSrs  Source SRS.
+  \param dstSrs  Destination SRS.
+  \return  Whether the reprojection was successful or not.
+*/
 bool reprojectBounds(Bounds& box, const SpatialReference& srcSrs,
     const SpatialReference& dstSrs)
 {
-    SrsTransform transform(srcSrs, dstSrs);
-
-    BOX3D b = box.to3d();
-    bool ok = transform.transform(b.minx, b.miny, b.minz);
-    if (ok)
-        ok = transform.transform(b.maxx, b.maxy, b.maxz);
+    bool ok = false;
+    if (box.is3d())
+    {
+        BOX3D b3 = box.to3d();
+        ok = reprojectBounds(b3, srcSrs, dstSrs);
+        box.reset(b3);
+    }
+    else
+    {
+        BOX2D b2 = box.to2d();
+        ok = reprojectBounds(b2, srcSrs, dstSrs);
+        box.reset(b2);
+    }
     return ok;
 }
 
 /**
   Reproject a bounds box from a source projection to a destination.
   \param box  2D Bounds box to be reprojected in-place.
-  \param srcSrs  String in WKT or other suitable format of box coordinates.
-  \param dstSrs  String in WKT or other suitable format to which
-    coordinates should be projected.
+  \param srcSrs  Source SRS.
+  \param dstSrs  Destination SRS.
   \return  Whether the reprojection was successful or not.
 */
 bool reprojectBounds(BOX2D& box, const SpatialReference& srcSrs,

pdal/util/Bounds.cpp

@@ -138,6 +138,23 @@ Bounds::Bounds(const BOX2D& box) : m_box(box)
     m_box.maxz = LOWEST;
 }
 
+void Bounds::reset(const BOX3D& box)
+{
+    m_box = box;
+}
+
+
+void Bounds::reset(const BOX2D& box)
+{
+    m_box.minx = box.minx;
+    m_box.maxx = box.maxx;
+    m_box.miny = box.miny;
+    m_box.maxy = box.maxy;
+    m_box.minz = HIGHEST;
+    m_box.maxz = LOWEST;
+}
+
+
 // We don't allow implicit conversion from a BOX2D to BOX3D.  Use the explicit
 // BOX3D ctor that takes a BOX2D if that's what you want.
 BOX3D Bounds::to3d() const

pdal/util/Bounds.hpp

@@ -630,6 +630,8 @@ class PDAL_DLL Bounds
     BOX3D to3d() const;
     BOX2D to2d() const;
     bool is3d() const;
+    void reset(const BOX3D& box);
+    void reset(const BOX2D& box);
     void grow(double x, double y);
     void grow(double x, double y, double z);
     void parse(const std::string& s, std::string::size_type& pos);

test/unit/filters/CropFilterTest.cpp

@@ -46,6 +46,10 @@
 #include <filters/StreamCallbackFilter.hpp>
 #include "Support.hpp"
 
+//ABELL
+#include <pdal/util/Bounds.hpp>
+#include <pdal/GDALUtils.hpp>
+
 using namespace pdal;
 
 TEST(CropFilterTest, create)
@@ -616,3 +620,54 @@ TEST(CropFilterTest, bounds_inside_outside)
     // Expect 1026 points when cropping to the outside of the bounds.
     EXPECT_EQ(nStreamPoints, 1026U);
 }
+
+// Make sure that transformed 2D and 3D bounds work.
+TEST(CropFilterTest, issue_3114)
+{
+    using namespace Dimension;
+
+    auto tst = [](const std::string& bounds, size_t count)
+    {
+        ColumnPointTable table;
+
+        table.layout()->registerDims({Id::X, Id::Y, Id::Z});
+        table.finalize();
+
+        PointViewPtr view(new PointView(table));
+
+        view->setField(Id::X, 0, 555000);
+        view->setField(Id::Y, 0, 4180000);
+        view->setField(Id::Z, 0, 100);
+
+        view->setField(Id::X, 1, 555000);
+        view->setField(Id::Y, 1, 4180000);
+        view->setField(Id::Z, 1, 1000);
+
+        view->setField(Id::X, 2, 565000);
+        view->setField(Id::Y, 2, 4180000);
+        view->setField(Id::Z, 2, 100);
+
+        BufferReader r;
+        r.addView(view);
+
+        Options ropts;
+        ropts.add("override_srs", "EPSG:26910");
+        r.setOptions(ropts);
+
+        CropFilter f;
+
+        Options copts;
+        copts.add("bounds", bounds);
+        copts.add("a_srs", "EPSG:4326");
+        f.setOptions(copts);
+        f.setInput(r);
+
+        f.prepare(table);
+        PointViewSet s = f.execute(table);
+        PointViewPtr v = *s.begin();
+        EXPECT_EQ(v->size(), count);
+    };
+
+    tst("([-122.530, -122.347], [37.695, 37.816])", 2);
+    tst("([-122.530, -122.347], [37.695, 37.816], [0,500])", 1);
+}