Skip to content

Instantly share code, notes, and snippets.

@connormanning
Created October 1, 2018 20:41
Show Gist options
  • Save connormanning/78ab3758eb34ae4584d0a3b0c3d8b368 to your computer and use it in GitHub Desktop.
Save connormanning/78ab3758eb34ae4584d0a3b0c3d8b368 to your computer and use it in GitHub Desktop.
PDAL Reader spatialreference option handling

Description

https://entwine.io/data/red-rocks.laz and https://entwine.io/data/red-rocks.bpf contain the same point cloud data and the same SRS info. The SRS is UTM zone 13N (EPSG:26913).

In this sample program, the following pipeline behaviors are run in order:

  1. reproject the BPF file to EPSG:3857
  2. reproject the LAZ file to EPSG:3857
  3. reproject the BPF file to EPSG:3857 with readers.bpf.spatialreference set to EPSG:26912
  4. reproject the LAZ file to EPSG:3857 with readers.laz.spatialreference set to EPSG:26912

For each of these pipelines, the pipeline itself and the output bounds are logged. The expected behavior would be that the extents of 1 and 2 are equivalent, and that 3 and 4 are equivalent to each other but not equivalent to 1 and 2.

Instead, the behavior of 3 is equivalent to 1 and 2, since the readers.bpf.spatialreference option is not respected by the BPF reader.

$ g++ -std=c++11 srs.cpp -lpdalcpp -ljsoncpp && ./a.out
{
	"pipeline" :
	[
		{
			"filename" : "https://entwine.io/data/red-rocks.bpf"
		},
		{
			"out_srs" : "EPSG:3857",
			"type" : "filters.reprojection"
		}
	]
}
Bounds: (-11711828.98, 4816845.99), (-11710915.45, 4817997.54)


{
	"pipeline" :
	[
		{
			"filename" : "https://entwine.io/data/red-rocks.laz"
		},
		{
			"out_srs" : "EPSG:3857",
			"type" : "filters.reprojection"
		}
	]
}
Bounds: (-11711828.98, 4816845.99), (-11710915.45, 4817997.54)


{
	"pipeline" :
	[
		{
			"filename" : "https://entwine.io/data/red-rocks.bpf",
			"spatialreference" : "EPSG:26912"
		},
		{
			"out_srs" : "EPSG:3857",
			"type" : "filters.reprojection"
		}
	]
}
Bounds: (-11711828.98, 4816845.99), (-11710915.45, 4817997.54)


{
	"pipeline" :
	[
		{
			"filename" : "https://entwine.io/data/red-rocks.laz",
			"spatialreference" : "EPSG:26912"
		},
		{
			"out_srs" : "EPSG:3857",
			"type" : "filters.reprojection"
		}
	]
}
Bounds: (-12379745.92, 4816845.99), (-12378832.39, 4817997.54)
#include <limits>
#include <pdal/PipelineManager.hpp>
#include <pdal/Stage.hpp>
#include <json/json.h>
// Native (EPSG:26913) -> EPSG:3857
// (-11711830, 4816845.6), (-11710915, 4817998.9)
//
// EPSG:26912 -> EPSG:3857
// (-12379747, 4816845.6), (-12378831, 4817998.9)
void analyze(const Json::Value& pipeline)
{
pdal::PipelineManager pm;
std::cout << pipeline << std::endl;
std::istringstream iss(pipeline.toStyledString());
pm.readPipeline(iss);
pdal::Stage* stage(pm.getStage());
if (!stage) throw std::runtime_error("No stage");
pdal::PointTable table;
stage->prepare(table);
const auto views(stage->execute(table));
double x, y, xmin, ymin, xmax, ymax;
xmin = ymin = std::numeric_limits<double>::max();
xmax = ymax = std::numeric_limits<double>::lowest();
for (auto view : views)
{
for (uint64_t i(0); i < view->size(); ++i)
{
x = view->getFieldAs<double>(pdal::Dimension::Id::X, i);
y = view->getFieldAs<double>(pdal::Dimension::Id::Y, i);
xmin = std::min(xmin, x);
ymin = std::min(ymin, y);
xmax = std::max(xmax, x);
ymax = std::max(ymax, y);
}
}
std::cout << std::fixed << std::setprecision(2) << "Bounds: " <<
"(" << xmin << ", " << ymin << "), " <<
"(" << xmax << ", " << ymax << ")\n\n" << std::endl;
}
int main()
{
Json::Value wrapper;
Json::Value& pipeline(wrapper["pipeline"]);
Json::Value& reader(pipeline.append(Json::objectValue));
Json::Value& filter(pipeline.append(Json::objectValue));
filter["type"] = "filters.reprojection";
filter["out_srs"] = "EPSG:3857";
// First we'll run both readers using the default spatial reference from
// the input files themselves. We expect the same output bounds from both
// of these since the LAZ and BPF files contain the same data and same SRS.
reader["filename"] = "https://entwine.io/data/red-rocks.bpf";
analyze(wrapper);
reader["filename"] = "https://entwine.io/data/red-rocks.laz";
analyze(wrapper);
// Now we'll set the "spatialreference" option on the reader and run both
// pipelines again. While the LAS reader respects the setting, the BPF
// reader does not.
reader["spatialreference"] = "EPSG:26912";
reader["filename"] = "https://entwine.io/data/red-rocks.bpf";
analyze(wrapper);
reader["filename"] = "https://entwine.io/data/red-rocks.laz";
analyze(wrapper);
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment