Skip to content

Commit

Permalink
Radius expansion for grid. Tests broken.
Browse files Browse the repository at this point in the history
  • Loading branch information
@abellgithub
abellgithub committed Jul 20, 2018
1 parent 7271403 commit 5ef5b00
Show file tree
Hide file tree
Showing 6 changed files with 157 additions and 80 deletions.
8 changes: 8 additions & 0 deletions io/GDALGrid.cpp
View file
Expand Up @@ -43,6 +43,14 @@
namespace pdal
{

//ABELL - In the beginning this data needed to be contiguous, as it was passed
// directly to GDAL to write. Since we started supporting various data types
// in GDAL output, we end up copying/casting data to a block for output,
// so there's no reason that we must have contiguous data -- we just need
// an iterator that allows traversal of the data in row-major order. So,
// this should probably be re-implemented in some way that doesn't require
// moving data around every time the grid is resized.

GDALGrid::GDALGrid(size_t width, size_t height, double edgeLength,
double radius, int outputTypes, size_t windowSize) :
m_width(width), m_height(height), m_windowSize(windowSize),
Expand Down
67 changes: 41 additions & 26 deletions io/GDALWriter.cpp
View file
Expand Up @@ -120,13 +120,7 @@ void GDALWriter::prepared(PointTableRef table)
"' does not exist.");
if (!m_radiusArg->set())
m_radius = m_edgeLength * sqrt(2.0);
}


void GDALWriter::readyTable(PointTableRef table)
{
if (m_bounds.to2d().empty() && !table.supportsView())
throwError("Option 'bounds' required in streaming mode.");
m_fixedGrid = m_bounds.to2d().valid();
}


Expand All @@ -136,8 +130,11 @@ void GDALWriter::readyFile(const std::string& filename,
m_outputFilename = filename;
m_srs = srs;
m_grid.reset();
if (m_bounds.to2d().valid())
if (m_fixedGrid)
{
createGrid(m_bounds.to2d());
m_expandByPoint = false;
}
}


Expand All @@ -149,7 +146,7 @@ void GDALWriter::createGrid(BOX2D bounds)
try
{
m_grid.reset(new GDALGrid(width, height, m_edgeLength, m_radius,
m_outputTypes, m_windowSize));
m_outputTypes, m_windowSize));
}
catch (GDALGrid::error& err)
{
Expand Down Expand Up @@ -185,16 +182,21 @@ void GDALWriter::expandGrid(BOX2D bounds)

void GDALWriter::writeView(const PointViewPtr view)
{
BOX2D bounds;
if (m_bounds.to2d().valid())
bounds = m_bounds.to2d();
else
view->calculateBounds(bounds);
m_expandByPoint = false;

if (!m_grid)
createGrid(bounds);
else
expandGrid(bounds);
// When we're running in standard mode, it's better to get the bounds and
// expand once, rather than have to do this for every point, since an
// expansion causes data to move.
if (!m_fixedGrid)
{
BOX2D bounds;
view->calculateBounds(bounds);
bounds.grow(m_radius);
if (!m_grid)
createGrid(bounds);
else
expandGrid(bounds);
}

PointRef point(*view, 0);
for (PointId idx = 0; idx < view->size(); ++idx)
Expand All @@ -207,23 +209,36 @@ void GDALWriter::writeView(const PointViewPtr view)

bool GDALWriter::processOne(PointRef& point)
{
double x = point.getFieldAs<double>(Dimension::Id::X) -
m_curBounds.minx;
double y = point.getFieldAs<double>(Dimension::Id::Y) -
m_curBounds.miny;
double x = point.getFieldAs<double>(Dimension::Id::X);
double y = point.getFieldAs<double>(Dimension::Id::Y);
double z = point.getFieldAs<double>(m_interpDim);

if (m_expandByPoint)
{
if (!m_grid)
createGrid(BOX2D(x, y, x, y).grow(m_radius));
else if (!m_curBounds.contains(x, y))
expandGrid(BOX2D(x, y, x, y).grow(m_radius));
}
x -= m_curBounds.minx;
y -= m_curBounds.miny;
if (!m_fixedGrid)
{
x += m_radius;
y += m_radius;
}

m_grid->addPoint(x, y, z);
return true;
}


void GDALWriter::doneFile()
{
if (!m_grid) {
throw pdal_error("Unable to write GDAL data, grid is uninitialized. You "
"might have provided the GDALWriter zero points.");
}
if (!m_grid)
throw pdal_error("Unable to write GDAL data with no points "
"for output.");

std::array<double, 6> pixelToPos;

pixelToPos[0] = m_curBounds.minx;
Expand Down
10 changes: 8 additions & 2 deletions io/GDALWriter.hpp
View file
Expand Up @@ -49,14 +49,18 @@ class PDAL_DLL GDALWriter : public FlexWriter, public Streamable
public:
std::string getName() const;

GDALWriter() : m_outputTypes(0)
GDALWriter() : m_outputTypes(0), m_expandByPoint(true)
{}

size_t gridWidth() const
{ return m_grid->width(); }
size_t gridHeight() const
{ return m_grid->height(); }

private:
virtual void addArgs(ProgramArgs& args);
virtual void initialize();
virtual void prepared(PointTableRef table);
virtual void readyTable(PointTableRef table);
virtual void readyFile(const std::string& filename,
const SpatialReference& srs);
virtual void writeView(const PointViewPtr view);
Expand All @@ -82,6 +86,8 @@ class PDAL_DLL GDALWriter : public FlexWriter, public Streamable
Dimension::Id m_interpDim;
std::string m_interpDimString;
Dimension::Type m_dataType;
bool m_expandByPoint;
bool m_fixedGrid;
};

}
18 changes: 16 additions & 2 deletions pdal/util/Bounds.cpp
View file
Expand Up @@ -127,20 +127,34 @@ bool BOX3D::valid() const
return !empty();
}

void BOX2D::grow(double x, double y)

BOX2D& BOX2D::grow(double dist)
{
assert(valid());
minx -= dist;
maxx += dist;
miny -= dist;
maxy += dist;
return *this;
}


BOX2D& BOX2D::grow(double x, double y)
{
if (x < minx) minx = x;
if (x > maxx) maxx = x;

if (y < miny) miny = y;
if (y > maxy) maxy = y;
return *this;
}

void BOX3D::grow(double x, double y, double z)
BOX3D& BOX3D::grow(double x, double y, double z)
{
BOX2D::grow(x, y);
if (z < minz) minz = z;
if (z > maxz) maxz = z;
return *this;
}

const BOX2D& BOX2D::getDefaultSpatialExtent()
Expand Down
39 changes: 29 additions & 10 deletions pdal/util/Bounds.hpp
View file
Expand Up @@ -91,15 +91,19 @@ class PDAL_DLL BOX2D
void clear();

/**
Expand the bounds of the box if a value is less than the current
minimum or greater than the current maximum. If the bounds box is
currently empty, both minimum and maximum box bounds will be set to
the provided value.
Expand the bounds of the box to include the specified point.
\param x X dimension value.
\param y Y dimension value.
\param x X point location.
\param y Y point location.
*/
BOX2D& grow(double x, double y);

/**
Expand the bounds of the box in all directions by a specified amount.
\param dist Distance by which to expand the box.
*/
void grow(double x, double y);
BOX2D& grow(double dist);

/**
Determine if a bounds box contains a point.
Expand Down Expand Up @@ -157,13 +161,14 @@ class PDAL_DLL BOX2D
\param other Box that this box should contain.
*/
void grow(const BOX2D& other)
BOX2D& grow(const BOX2D& other)
{
if (other.minx < minx) minx = other.minx;
if (other.maxx > maxx) maxx = other.maxx;

if (other.miny < miny) miny = other.miny;
if (other.maxy > maxy) maxy = other.maxy;
return *this;
}

/**
Expand Down Expand Up @@ -347,7 +352,7 @@ class PDAL_DLL BOX3D : private BOX2D
\param y Y dimension value.
\param z Z dimension value.
*/
void grow(double x, double y, double z);
BOX3D& grow(double x, double y, double z);

/**
Clear the bounds box to an empty state.
Expand Down Expand Up @@ -428,11 +433,25 @@ class PDAL_DLL BOX3D : private BOX2D
\param other Box that this box should contain.
*/
void grow(const BOX3D& other)
BOX3D& grow(const BOX3D& other)
{
BOX2D::grow(other);
if (other.minz < minz) minz = other.minz;
if (other.maxz > maxz) maxz = other.maxz;
return *this;
}

/**
Expand this box by a specified amount.
\param dist Distance by which box should be expanded.
*/
BOX3D& grow(double dist)
{
BOX2D::grow(dist);
minz -= dist;
maxz += dist;
return *this;
}

/**
Expand Down
95 changes: 55 additions & 40 deletions test/unit/io/GDALWriterTest.cpp
View file
Expand Up @@ -54,58 +54,73 @@ namespace
void runGdalWriter(const Options& wo, const std::string& infile,
const std::string& outfile, const std::string& values)
{
FileUtils::deleteFile(outfile);
auto run = [=](bool streamMode)
{
FileUtils::deleteFile(outfile);

Options ro;
ro.add("filename", infile);
Options ro;
ro.add("filename", infile);

TextReader r;
r.setOptions(ro);
TextReader r;
r.setOptions(ro);

GDALWriter w;
w.setOptions(wo);
w.setInput(r);
GDALWriter w;
w.setOptions(wo);
w.setInput(r);

PointTable t;
if (streamMode)
{
FixedPointTable t(100);

w.prepare(t);
w.execute(t);
w.prepare(t);
w.execute(t);
}
else
{
PointTable t;

using namespace gdal;
w.prepare(t);
w.execute(t);
}

std::istringstream iss(values);
using namespace gdal;

std::vector<double> arr;
while (true)
{
double d;
iss >> d;
if (!iss)
break;
arr.push_back(d);
}
std::istringstream iss(values);

registerDrivers();
Raster raster(outfile, "GTiff");
if (raster.open() != GDALError::None)
{
throw pdal_error(raster.errorMsg());
}
std::vector<double> data;
raster.readBand(data, 1);
std::vector<double> arr;
while (true)
{
double d;
iss >> d;
if (!iss)
break;
arr.push_back(d);
}

int row = 0;
int col = 0;
for (size_t i = 0; i < arr.size(); ++i)
{
EXPECT_NEAR(arr[i], data[i], .001) << "Mismatch for row/col = " <<
row << "/" << col;
if (++col == raster.width())
registerDrivers();
Raster raster(outfile, "GTiff");
if (raster.open() != GDALError::None)
{
col = 0;
row++;
throw pdal_error(raster.errorMsg());
}
}
std::vector<double> data;
raster.readBand(data, 1);
int row = 0;
int col = 0;
for (size_t i = 0; i < arr.size(); ++i)
{
EXPECT_NEAR(arr[i], data[i], .001) << "Error row/col = " <<
row << "/" << col << std::endl;
if (++col == raster.width())
{
col = 0;
row++;
}
}
};

run(false);
run(true);
}

void runGdalWriter2(const Options& wo, const std::string& outfile,
Expand Down

0 comments on commit 5ef5b00

Please sign in to comment.