Clear up some code in TileDataset

This commit is contained in:
Erik Broberg
2016-06-28 18:35:35 -04:00
parent dd8bdc69b7
commit cdc3daedca
3 changed files with 175 additions and 155 deletions

View File

@@ -51,8 +51,137 @@ namespace {
namespace openspace {
//////////////////////////////////////////////////////////////////////////////////
// Tile Data Layout //
//////////////////////////////////////////////////////////////////////////////////
TileDataLayout::TileDataLayout() {
}
TileDataLayout::TileDataLayout(GDALDataset* dataSet, GLuint _glType) {
// Assume all raster bands have the same data type
gdalType = _glType != 0 ? TileDataType::getGdalDataType(glType) : dataSet->GetRasterBand(1)->GetRasterDataType();
glType = TileDataType::getOpenGLDataType(gdalType);
numRasters = dataSet->GetRasterCount();
bytesPerDatum = TileDataType::numberOfBytes(gdalType);
bytesPerPixel = bytesPerDatum * numRasters;
textureFormat = TileDataType::getTextureFormat(numRasters, gdalType);
}
//////////////////////////////////////////////////////////////////////////////////
// GDAL Data Region //
//////////////////////////////////////////////////////////////////////////////////
GdalDataRegion::GdalDataRegion(GDALDataset * dataSet,
const ChunkIndex& chunkIndex, int tileLevelDifference) {
GDALRasterBand* firstBand = dataSet->GetRasterBand(1);
// Assume all raster bands have the same data type
// Level = overviewCount - overview (default, levels may be overridden)
int numOverviews = firstBand->GetOverviewCount();
// Generate a patch from the chunkIndex, extract the bounds which
// are used to calculated where in the GDAL data set to read data.
// pixelStart0 and pixelEnd0 defines the interval in the pixel space
// at overview 0
GeodeticPatch patch = GeodeticPatch(chunkIndex);
glm::uvec2 pixelStart0 = geodeticToPixel(dataSet, patch.getCorner(Quad::NORTH_WEST));
glm::uvec2 pixelEnd0 = geodeticToPixel(dataSet, patch.getCorner(Quad::SOUTH_EAST));
glm::uvec2 numPixels0 = pixelEnd0 - pixelStart0;
// Calculate a suitable overview to choose from the GDAL dataset
int minNumPixels0 = glm::min(numPixels0.x, numPixels0.y);
GDALRasterBand* maxOverview;
if (numOverviews <= 0) {
maxOverview = firstBand;
}
else {
maxOverview = firstBand->GetOverview(numOverviews - 1);
}
int sizeLevel0 = maxOverview->GetXSize();
// The dataset itself may not have overviews but even if it does not, an overview
// for the data region can be calculated and possibly be used to sample greater
// Regions of the original dataset.
int ov = std::log2(minNumPixels0) - std::log2(sizeLevel0 + 1) - tileLevelDifference;
if (numOverviews > 0) {
ov = glm::clamp(ov, 0, numOverviews - 1);
}
// Convert the interval [pixelStart0, pixelEnd0] to pixel space at
// the calculated suitable overview, ov. using a >> b = a / 2^b
int toShift = ov + 1;
// Set member variables
overview = ov;
pixelStart = glm::uvec2(pixelStart0.x >> toShift, pixelStart0.y >> toShift);
pixelEnd = glm::uvec2(pixelEnd0.x >> toShift, pixelEnd0.y >> toShift);
numPixels = pixelEnd - pixelStart;
}
glm::uvec2 GdalDataRegion::geodeticToPixel(GDALDataset* dataSet, const Geodetic2& geo) {
double padfTransform[6];
CPLErr err = dataSet->GetGeoTransform(padfTransform);
ghoul_assert(err != CE_Failure, "Failed to get transform");
Scalar Y = Angle<Scalar>::fromRadians(geo.lat).asDegrees();
Scalar X = Angle<Scalar>::fromRadians(geo.lon).asDegrees();
// convert from pixel and line to geodetic coordinates
// Xp = padfTransform[0] + P*padfTransform[1] + L*padfTransform[2];
// Yp = padfTransform[3] + P*padfTransform[4] + L*padfTransform[5];
// <=>
double* a = &(padfTransform[0]);
double* b = &(padfTransform[3]);
// Xp = a[0] + P*a[1] + L*a[2];
// Yp = b[0] + P*b[1] + L*b[2];
// <=>
double divisor = (a[2] * b[1] - a[1] * b[2]);
ghoul_assert(divisor != 0.0, "Division by zero!");
//ghoul_assert(a[2] != 0.0, "a2 must not be zero!");
double P = (a[0] * b[2] - a[2] * b[0] + a[2] * Y - b[2] * X) / divisor;
double L = (-a[0] * b[1] + a[1] * b[0] - a[1] * Y + b[1] * X) / divisor;
// ref: https://www.wolframalpha.com/input/?i=X+%3D+a0+%2B+a1P+%2B+a2L,+Y+%3D+b0+%2B+b1P+%2B+b2L,+solve+for+P+and+L
double Xp = a[0] + P*a[1] + L*a[2];
double Yp = b[0] + P*b[1] + L*b[2];
ghoul_assert(abs(X - Xp) < 1e-10, "inverse should yield X as before");
ghoul_assert(abs(Y - Yp) < 1e-10, "inverse should yield Y as before");
return glm::uvec2(glm::round(P), glm::round(L));
}
//////////////////////////////////////////////////////////////////////////////////
// Tile Dataset //
//////////////////////////////////////////////////////////////////////////////////
// INIT THIS TO FALSE AFTER REMOVED FROM TILEPROVIDER
bool TileDataset::GdalHasBeenInitialized = false;
TileDataset::TileDataset(const std::string& gdalDatasetDesc, int minimumPixelSize,
@@ -72,10 +201,11 @@ namespace openspace {
if (!_dataset) {
throw ghoul::RuntimeError("Failed to load dataset:\n" + gdalDatasetDesc);
}
_dataLayout = DataLayout(_dataset, dataType);
_dataLayout = TileDataLayout(_dataset, dataType);
_depthTransform = calculateTileDepthTransform();
_tileLevelDifference = calculateTileLevelDifference(_dataset, minimumPixelSize);
LDEBUG(gdalDatasetDesc << " - " << _tileLevelDifference);
_maxLevel = calculateMaxLevel(_tileLevelDifference);
}
@@ -211,7 +341,7 @@ namespace openspace {
}
char* TileDataset::getImageDataFlippedY(const GdalDataRegion& region,
const DataLayout& dataLayout, const char* imageData)
const TileDataLayout& dataLayout, const char* imageData)
{
size_t bytesPerLine = dataLayout.bytesPerPixel * region.numPixels.x;
size_t totalNumBytes = bytesPerLine * region.numPixels.y;
@@ -236,13 +366,13 @@ namespace openspace {
}
const TileDataset::DataLayout& TileDataset::getDataLayout() const {
const TileDataLayout& TileDataset::getDataLayout() const {
return _dataLayout;
}
std::shared_ptr<TilePreprocessData> TileDataset::preprocess(const char* imageData,
const GdalDataRegion& region, const DataLayout& dataLayout)
const GdalDataRegion& region, const TileDataLayout& dataLayout)
{
size_t bytesPerLine = dataLayout.bytesPerPixel * region.numPixels.x;
size_t totalNumBytes = bytesPerLine * region.numPixels.y;
@@ -283,113 +413,4 @@ namespace openspace {
glm::uvec2 TileDataset::geodeticToPixel(GDALDataset* dataSet, const Geodetic2& geo) {
double padfTransform[6];
CPLErr err = dataSet->GetGeoTransform(padfTransform);
ghoul_assert(err != CE_Failure, "Failed to get transform");
Scalar Y = Angle<Scalar>::fromRadians(geo.lat).asDegrees();
Scalar X = Angle<Scalar>::fromRadians(geo.lon).asDegrees();
// convert from pixel and line to geodetic coordinates
// Xp = padfTransform[0] + P*padfTransform[1] + L*padfTransform[2];
// Yp = padfTransform[3] + P*padfTransform[4] + L*padfTransform[5];
// <=>
double* a = &(padfTransform[0]);
double* b = &(padfTransform[3]);
// Xp = a[0] + P*a[1] + L*a[2];
// Yp = b[0] + P*b[1] + L*b[2];
// <=>
double divisor = (a[2]*b[1] - a[1]*b[2]);
ghoul_assert(divisor != 0.0, "Division by zero!");
//ghoul_assert(a[2] != 0.0, "a2 must not be zero!");
double P = (a[0]*b[2] - a[2]*b[0] + a[2]*Y - b[2]*X) / divisor;
double L = (-a[0]*b[1] + a[1]*b[0] - a[1]*Y + b[1]*X) / divisor;
// ref: https://www.wolframalpha.com/input/?i=X+%3D+a0+%2B+a1P+%2B+a2L,+Y+%3D+b0+%2B+b1P+%2B+b2L,+solve+for+P+and+L
double Xp = a[0] + P*a[1] + L*a[2];
double Yp = b[0] + P*b[1] + L*b[2];
ghoul_assert(abs(X - Xp) < 1e-10, "inverse should yield X as before");
ghoul_assert(abs(Y - Yp) < 1e-10, "inverse should yield Y as before");
return glm::uvec2(glm::round(P), glm::round(L));
}
TileDataset::GdalDataRegion::GdalDataRegion(GDALDataset * dataSet,
const ChunkIndex& chunkIndex, int tileLevelDifference) {
GDALRasterBand* firstBand = dataSet->GetRasterBand(1);
// Assume all raster bands have the same data type
// Level = overviewCount - overview (default, levels may be overridden)
int numOverviews = firstBand->GetOverviewCount();
// Generate a patch from the chunkIndex, extract the bounds which
// are used to calculated where in the GDAL data set to read data.
// pixelStart0 and pixelEnd0 defines the interval in the pixel space
// at overview 0
GeodeticPatch patch = GeodeticPatch(chunkIndex);
glm::uvec2 pixelStart0 = geodeticToPixel(dataSet, patch.getCorner(Quad::NORTH_WEST));
glm::uvec2 pixelEnd0 = geodeticToPixel(dataSet, patch.getCorner(Quad::SOUTH_EAST));
glm::uvec2 numPixels0 = pixelEnd0 - pixelStart0;
// Calculate a suitable overview to choose from the GDAL dataset
int minNumPixels0 = glm::min(numPixels0.x, numPixels0.y);
GDALRasterBand* maxOverview;
if (numOverviews <= 0) {
maxOverview = firstBand;
}
else {
maxOverview = firstBand->GetOverview(numOverviews - 1);
}
int sizeLevel0 = maxOverview->GetXSize();
// The dataset itself may not have overviews but even if it does not, an overview
// for the data region can be calculated and possibly be used to sample greater
// Regions of the original dataset.
int ov = std::log2(minNumPixels0) - std::log2(sizeLevel0 + 1) - tileLevelDifference;
if (numOverviews > 0) {
ov = glm::clamp(ov, 0, numOverviews - 1);
}
// Convert the interval [pixelStart0, pixelEnd0] to pixel space at
// the calculated suitable overview, ov. using a >> b = a / 2^b
int toShift = ov + 1;
// Set member variables
overview = ov;
pixelStart = glm::uvec2(pixelStart0.x >> toShift, pixelStart0.y >> toShift);
pixelEnd = glm::uvec2(pixelEnd0.x >> toShift, pixelEnd0.y >> toShift);
numPixels = pixelEnd - pixelStart;
}
TileDataset::DataLayout::DataLayout() {
}
TileDataset::DataLayout::DataLayout(GDALDataset* dataSet, GLuint _glType) {
// Assume all raster bands have the same data type
gdalType = _glType != 0 ? TileDataType::getGdalDataType(glType) : dataSet->GetRasterBand(1)->GetRasterDataType();
glType = TileDataType::getOpenGLDataType(gdalType);
numRasters = dataSet->GetRasterCount();
bytesPerDatum = TileDataType::numberOfBytes(gdalType);
bytesPerPixel = bytesPerDatum * numRasters;
textureFormat = TileDataType::getTextureFormat(numRasters, gdalType);
}
} // namespace openspace

View File

@@ -29,6 +29,7 @@
#include <set>
#include <queue>
#include <iostream>
#include <unordered_map>
#include "gdal_priv.h"
@@ -44,6 +45,38 @@ namespace openspace {
using namespace ghoul::opengl;
using namespace ghoul::filesystem;
struct TileDataLayout {
TileDataLayout();
TileDataLayout(GDALDataset* dataSet, GLuint glType);
GDALDataType gdalType;
GLuint glType;
size_t bytesPerDatum;
size_t numRasters;
size_t bytesPerPixel;
TextureFormat textureFormat;
};
struct GdalDataRegion {
GdalDataRegion(GDALDataset* dataSet, const ChunkIndex& chunkIndex, int tileLevelDifference);
glm::uvec2 pixelStart;
glm::uvec2 pixelEnd;
glm::uvec2 numPixels;
int overview;
static glm::uvec2 geodeticToPixel(GDALDataset* dataSet, const Geodetic2& geo);
};
class TileDataset {
public:
@@ -61,48 +94,18 @@ namespace openspace {
~TileDataset();
struct DataLayout {
DataLayout();
DataLayout(GDALDataset* dataSet, GLuint glType);
GDALDataType gdalType;
GLuint glType;
size_t bytesPerDatum;
size_t numRasters;
size_t bytesPerPixel;
TextureFormat textureFormat;
};
std::shared_ptr<TileIOResult> readTileData(ChunkIndex chunkIndex);
int getMaximumLevel() const;
TileDepthTransform getDepthTransform() const;
const DataLayout& getDataLayout() const;
const TileDataLayout& getDataLayout() const;
private:
struct GdalDataRegion {
GdalDataRegion(GDALDataset* dataSet, const ChunkIndex& chunkIndex, int tileLevelDifference);
glm::uvec2 pixelStart;
glm::uvec2 pixelEnd;
glm::uvec2 numPixels;
int overview;
};
//////////////////////////////////////////////////////////////////////////////////
// HELPER FUNCTIONS //
@@ -112,18 +115,14 @@ namespace openspace {
TileDepthTransform calculateTileDepthTransform();
static glm::uvec2 geodeticToPixel(GDALDataset* dataSet, const Geodetic2& geo);
std::shared_ptr<TilePreprocessData> preprocess(const char* imageData,
const GdalDataRegion& region, const DataLayout& dataLayout);
const GdalDataRegion& region, const TileDataLayout& dataLayout);
static int calculateTileLevelDifference(GDALDataset* dataset, int minimumPixelSize);
static char* getImageDataFlippedY(const GdalDataRegion& region,
const DataLayout& dataLayout, const char* imageData);
const TileDataLayout& dataLayout, const char* imageData);
//////////////////////////////////////////////////////////////////////////////////
@@ -139,7 +138,7 @@ namespace openspace {
TileDepthTransform _depthTransform;
GDALDataset* _dataset;
DataLayout _dataLayout;
TileDataLayout _dataLayout;
bool _doPreprocessing;
};

View File

@@ -146,7 +146,7 @@ namespace openspace {
void CachingTileProvider::initializeAndAddToCache(std::shared_ptr<TileIOResult> tileIOResult) {
ChunkHashKey key = tileIOResult->chunkIndex.hashKey();
TileDataset::DataLayout dataLayout = _asyncTextureDataProvider->getTextureDataProvider()->getDataLayout();
TileDataLayout dataLayout = _asyncTextureDataProvider->getTextureDataProvider()->getDataLayout();
Texture* texturePtr = new Texture(
tileIOResult->imageData,
tileIOResult->dimensions,