Skip to content

Commit

Permalink
WIP
Browse files Browse the repository at this point in the history
  • Loading branch information
kordejong committed Mar 13, 2024
1 parent 74f999b commit 2ef98de
Show file tree
Hide file tree
Showing 5 changed files with 90 additions and 19 deletions.
2 changes: 2 additions & 0 deletions source/data_model/gdal/include/lue/gdal/driver.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -17,4 +17,6 @@ namespace lue::gdal {

auto driver(std::string const& name) -> DriverPtr;

auto driver_name(std::string const& dataset_name) -> std::string;

} // namespace lue::gdal
3 changes: 2 additions & 1 deletion source/data_model/gdal/src/compare_rasters.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -305,7 +305,8 @@ namespace lue::gdal {
}
default:
{
differences.emplace_back(fmt::format("Unsupported GDAL data type: {}", data_type1));
differences.emplace_back(
fmt::format("Unsupported GDAL data type: {}", as_string(data_type1)));
}
}
}
Expand Down
7 changes: 7 additions & 0 deletions source/data_model/gdal/src/data_type.cpp
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
#include "lue/gdal/data_type.hpp"
#include <cassert>


namespace lue::gdal {
Expand Down Expand Up @@ -88,6 +89,12 @@ namespace lue::gdal {
result = "GDT_CFloat64";
break;
}
case GDT_TypeCount:
{
// This should never happen
assert(false);
break;
}
}

return result;
Expand Down
44 changes: 44 additions & 0 deletions source/data_model/gdal/src/driver.cpp
Original file line number Diff line number Diff line change
@@ -1,9 +1,29 @@
#include "lue/gdal/driver.hpp"
#include <fmt/format.h>
#include <algorithm>
#include <cctype>
#include <filesystem>
#include <map>
#include <stdexcept>


namespace lue::gdal {
namespace {

// GDAL driver name mapped by lower-cased file extension. Add to it whenever useful.
const std::map<std::string, std::string> driver_name_by_extension{
{".tif", "GTiff"},
{".hdf4", "HDF4"},
{".h4", "HDF4"},
{".hdf5", "HDF5"},
{".h5", "HDF5"},
{".nc", "netCDF"},
{".map", "PCRaster"},
{".csf", "PCRaster"},
{".img", "HFA"}, // Erdas Imagine
};

} // namespace

/*!
@brief Register all GDAL drivers
Expand Down Expand Up @@ -39,4 +59,28 @@ namespace lue::gdal {
return driver_ptr;
}


/*!
@brief Return the name of the driver that is (likely) able to handle the dataset
*/
auto driver_name(std::string const& dataset_name) -> std::string
{
std::string extension{std::filesystem::path(dataset_name).extension().string()};

std::transform(
extension.begin(),
extension.end(),
extension.begin(),
[](unsigned char const character) { return std::tolower(character); });

std::string result{"Unknown"};

if (auto it = driver_name_by_extension.find(extension); it != driver_name_by_extension.end())
{
result = (*it).second;
}

return result;
}

} // namespace lue::gdal
53 changes: 35 additions & 18 deletions source/data_model/translate/src/format/gdal.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -262,6 +262,9 @@ namespace lue::utility {
template<typename T>
void lue_to_gdal(data_model::Array const& array, gdal::Raster::Band& raster_band)
{
// TODO If a no-data value is set in the array, set it in the raster band


// It is assumed here that array only contains a 2D array
auto const blocks = utility::blocks(raster_band.shape(), raster_band.block_shape());
auto const& block_shape = blocks.block_shape();
Expand Down Expand Up @@ -361,8 +364,7 @@ namespace lue::utility {
return;
}

// TODO std::string const driver_name{gdal::driver_name(raster_name)};
std::string const driver_name{"GTiff"};
std::string const driver_name{gdal::driver_name(raster_name)};

// If the constant raster view finds a raster with the property name
// requested, export it to a single GDAL raster
Expand All @@ -388,21 +390,23 @@ namespace lue::utility {

RasterLayer layer{raster_view.layer(property_name)};
auto const& space_box{raster_view.space_box()};
gdal::Shape const raster_shape{hdf5_shape_to_gdal_shape(raster_view.grid_shape())};
auto const [nr_rows, nr_cols] = raster_shape;

gdal::Raster raster{gdal::create_dataset(
driver_name,
raster_name,
hdf5_shape_to_gdal_shape(raster_view.grid_shape()),
nr_bands,
memory_data_type_to_gdal_data_type(layer.memory_datatype()))};

// TODO
double const cell_size{0.000992063492063};
double const west{space_box[0]};
double const south{space_box[1]};
double const east{space_box[2]};
double const north{space_box[3]};
gdal::GeoTransform geo_transform{west, cell_size, 0, north, 0, -cell_size};

raster.set_geo_transform(geo_transform);
assert(east >= west);
assert(north >= south);
assert(nr_rows > 0);
assert(nr_cols > 0);

double const cell_width{(east - west) / nr_cols};
double const cell_height{(north - south) / nr_rows};

gdal::GeoTransform geo_transform{west, cell_width, 0, north, 0, -cell_height};

// TODO
// OGRSpatialReference oSRS;
Expand All @@ -415,6 +419,15 @@ namespace lue::utility {
// poDstDS->SetProjection( pszSRS_WKT );
// CPLFree( pszSRS_WKT );

gdal::Raster raster{gdal::create_dataset(
driver_name,
raster_name,
raster_shape,
nr_bands,
memory_data_type_to_gdal_data_type(layer.memory_datatype()))};

raster.set_geo_transform(geo_transform);

gdal::Raster::Band raster_band{raster.band(band_nr)};
lue_to_gdal(layer, raster_band);

Expand All @@ -439,13 +452,12 @@ namespace lue::utility {
else if (data_model::variable::contains_raster(dataset, phenomenon_name, property_set_name))
{
using RasterView = data_model::variable::RasterView<data_model::Dataset*>;
using RasterLayer = RasterView::Layer;

RasterView raster_view{&dataset, phenomenon_name, property_set_name};
// using RasterLayer = RasterView::Layer;

[[maybe_unused]] RasterView raster_view{&dataset, phenomenon_name, property_set_name};

// TODO

throw std::runtime_error("Importing stacks not supported yet. WIP.");

/// if (!raster_view.contains(property_name))
/// {
Expand Down Expand Up @@ -955,7 +967,12 @@ namespace lue::utility {
auto const [west, cell_width, row_rotation, north, col_rotation, cell_height] =
gdal_raster.geo_transform();
double const east = west + nr_cols * cell_width;
double const south = north - nr_rows * cell_height;
double const south = north - (nr_rows * std::abs(cell_height));

assert(east >= west);
assert(north >= south);
assert(nr_rows > 0);
assert(nr_cols > 0);

// Create / open raster view
using RasterView = data_model::constant::RasterView<data_model::Dataset*>;
Expand Down

0 comments on commit 2ef98de

Please sign in to comment.