Add support for origin/cell count instead of bounds.

Close #2415

doc/stages/writers.gdal.rst

@@ -16,14 +16,11 @@ The technique used to create the raster is a simple interpolation where
 each point that falls within a given radius_ of a raster cell center
 potentially contributes to the raster's value.  If no radius is provided,
 it is set to the product of the resolution_ and the square root of two.
-This is consistent with the original Points2Grid_ application from which
-this algorithm has its roots.  If a circle with the provided radius
+If a circle with the provided radius
 doesn't encompass the entire cell, it is possible that some points will
 not be considered at all, including those that may be within the bounds
 of the raster cell.
 
-.. _Points2Grid: http://www.opentopography.org/otsoftware/points2grid
-
 The GDAL writer creates rasters using the data specified in the dimension_
 option (defaults to `Z`). The writer creates up to six rasters based on
 different statistics in the output dataset.  The order of the layers in the
@@ -165,11 +162,22 @@ window_size
 dimension
   A dimension name to use for the interpolation. [Default: "Z"]
 
-.. _bounds:
-
 bounds
   The bounds of the data to be written.  Points not in bounds are discarded.
-  The format is ([minx, maxx],[miny,maxy]).
+  The format is ([minx, maxx],[miny,maxy]). [Deprecated]
+
+origin_x
+  X origin of the grid. [Default: None]
+
+origin_y
+  Y origin of the grid. [Default: None]
+
+width
+  Number of cells in the X direction. [Default: None]
+
+height
+  Number of cells in the Y direction. [Default: None]
 
 .. note::
-  The bounds_ option is required when a pipeline is run in streaming mode.
+    The bounds option is deprecated.  Use 'origin_x', 'origin_y', 'width'
+    and 'height' instead.

io/GDALWriter.cpp

@@ -80,8 +80,13 @@ void GDALWriter::addArgs(ProgramArgs& args)
     args.add("nodata", "No data value", m_noData,
         std::numeric_limits<double>::quiet_NaN());
     args.add("dimension", "Dimension to use", m_interpDimString, "Z");
-    args.add("bounds", "Bounds of data.  Required in streaming mode.",
-        m_bounds);
+    args.add("bounds", "Bounds of data. [deprecated]", m_bounds);
+    m_xOriginArg = &args.add("origin_x", "X origin for grid.", m_xOrigin);
+    m_yOriginArg = &args.add("origin_y", "Y origin for grid.", m_yOrigin);
+    m_widthArg = &args.add("width", "Number of cells in the X direction.",
+        m_width);
+    m_heightArg = &args.add("height", "Number of cells in the Y direction.",
+        m_height);
 }
 
 
@@ -111,6 +116,41 @@ void GDALWriter::initialize()
             throwError("Invalid output type: '" + ts + "'.");
     }
 
+    if (!m_radiusArg->set())
+        m_radius = m_edgeLength * sqrt(2.0);
+
+    int args = 0;
+    if (m_xOriginArg->set())
+        args |= 1;
+    if (m_yOriginArg->set())
+        args |= 2;
+    if (m_heightArg->set())
+        args |= 4;
+    if (m_widthArg->set())
+        args |= 8;
+    if (args != 0 && args != 15)
+        throwError("Must specify all or none of 'origin_x', 'origin_y', "
+            "'width' and 'height'.");
+    if (args == 15)
+    {
+        if (m_bounds.to2d().valid())
+            throwError("Specify either 'bounds' or 'origin_x'/'origin_y'/"
+                "'width'/'height' options -- not both");
+
+        // Subtracting .5 gets to the middle of the last cell.  This
+        // should get us back to the same place when figuring the
+        // cell count.
+        m_bounds = Bounds({m_xOrigin, m_yOrigin,
+            m_xOrigin + (m_edgeLength * (m_width - .5)),
+            m_yOrigin + (m_edgeLength * (m_height - .5))});
+    }
+
+    m_fixedGrid = m_bounds.to2d().valid();
+    // If we've specified a grid, we don't expand by point.  We also
+    // don't expand by point if we're running in standard mode.  That's
+    // set later in writeView.
+    m_expandByPoint = !m_fixedGrid;
+
     gdal::registerDrivers();
 }
 
@@ -121,13 +161,6 @@ void GDALWriter::prepared(PointTableRef table)
     if (m_interpDim == Dimension::Id::Unknown)
         throwError("Specified dimension '" + m_interpDimString +
             "' does not exist.");
-    if (!m_radiusArg->set())
-        m_radius = m_edgeLength * sqrt(2.0);
-    m_fixedGrid = m_bounds.to2d().valid();
-    // If we've specified a grid, we don't expand by point.  We also
-    // don't expand by point if we're running in standard mode.  That's
-    // set later in writeView.
-    m_expandByPoint = !m_fixedGrid;
 }
 
 

io/GDALWriter.hpp

@@ -86,6 +86,14 @@ class PDAL_DLL GDALWriter : public FlexWriter, public Streamable
     Position m_origin;
     double m_edgeLength;
     Arg *m_radiusArg;
+    double m_xOrigin;
+    double m_yOrigin;
+    size_t m_width;
+    size_t m_height;
+    Arg *m_xOriginArg;
+    Arg *m_yOriginArg;
+    Arg *m_heightArg;
+    Arg *m_widthArg;
     double m_radius;
     StringList m_options;
     StringList m_outputTypeString;

test/unit/io/GDALWriterTest.cpp

@@ -37,6 +37,7 @@
 #include <pdal/util/FileUtils.hpp>
 #include <filters/RangeFilter.hpp>
 #include <io/BufferReader.hpp>
+#include <io/FauxReader.hpp>
 #include <io/GDALWriter.hpp>
 #include <io/LasReader.hpp>
 #include <io/TextReader.hpp>
@@ -223,7 +224,11 @@ TEST(GDALWriterTest, min2)
     wo.add("filename", outfile);
 
     Options wo2 = wo;
-    wo2.add("bounds", "([-2, 4.7],[-2, 6.5])");
+//    wo2.add("bounds", "([-2, 4.7],[-2, 6.5])");
+    wo2.add("origin_x", -2);
+    wo2.add("origin_y", -2);
+    wo2.add("width", 7);
+    wo2.add("height", 9);
 
     const std::string output =
         "-9999.000 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00    -1.00"
@@ -629,7 +634,11 @@ TEST(GDALWriterTest, bounds)
     wo.add("resolution", 1);
     wo.add("radius", .7071);
     wo.add("filename", outfile);
-    wo.add("bounds", "([0, 4.5],[0, 4.5])");
+//    wo.add("bounds", "([0, 4.5],[0, 4.5])");
+    wo.add("origin_x", 0);
+    wo.add("origin_y", 0);
+    wo.add("height", 5);
+    wo.add("width", 5);
 
     const std::string output =
         "5.000 -9999.000     7.000     8.000     8.967 "
@@ -758,7 +767,11 @@ TEST(GDALWriterTest, issue_2545)
     wOpts.add("radius", 10);
     wOpts.add("output_type", "idw");
     wOpts.add("gdaldriver", "GTiff");
-    wOpts.add("bounds", BOX2D(.5, .5, 6.5, 6.5));
+    wOpts.add("origin_x", .5);
+    wOpts.add("origin_y", .5);
+    wOpts.add("width", 7);
+    wOpts.add("height", 7);
+//    wOpts.add("bounds", BOX2D(.5, .5, 6.5, 6.5));
     wOpts.add("filename", outfile);
 
     w.setOptions(wOpts);
@@ -790,4 +803,69 @@ TEST(GDALWriterTest, issue_2545)
         EXPECT_TRUE(data[i] <= 10.0 && data[i] >= 0);
 }
 
+// Test that using the alternate grid specification works
+TEST(GDALWriterTest, alternate_grid)
+{
+    FauxReader r;
+    Options opts;
+    opts.add("count", 1000);
+    opts.add("mode", "constant");
+    r.setOptions(opts);
+
+    std::string outfile(Support::temppath("test.tif"));
+
+    {
+        Options opts;
+        opts.add("origin_x", "10.0");
+        opts.add("filename", outfile);
+        opts.add("resolution", 1.1);
+
+        PointTable t;
+        GDALWriter w;
+        w.setOptions(opts);
+        w.setInput(r);
+        EXPECT_THROW(w.prepare(t), pdal_error);
+    }
+
+    {
+        Options opts;
+        opts.add("origin_x", "10.0");
+        opts.add("origin_y", "20.0");
+        opts.add("width", 10);
+        opts.add("height", 5);
+        opts.add("bounds", "([1, 10],[5, 25])");
+        opts.add("filename", outfile);
+        opts.add("resolution", 1.1);
+
+        PointTable t;
+        GDALWriter w;
+        w.setOptions(opts);
+        w.setInput(r);
+        EXPECT_THROW(w.prepare(t), pdal_error);
+    }
+
+    {
+        Options opts;
+        opts.add("origin_x", "10.0");
+        opts.add("origin_y", "20.0");
+        opts.add("width", 10);
+        opts.add("height", 5);
+        opts.add("filename", outfile);
+        opts.add("resolution", 1);
+
+        PointTable t;
+        GDALWriter w;
+        w.setOptions(opts);
+        w.setInput(r);
+        w.prepare(t);
+        w.execute(t);
+
+        gdal::registerDrivers();
+        gdal::Raster raster(outfile, "GTiff");
+        raster.open();
+        EXPECT_EQ(raster.width(), 10);
+        EXPECT_EQ(raster.height(), 5);
+    }
+}
+
 } // namespace pdal