From 00f473de63df7f1606c69769916d7b9d7c64d71c Mon Sep 17 00:00:00 2001 From: Dane Springmeyer Date: Sat, 18 Sep 2010 19:19:27 +0000 Subject: [PATCH] gdal plugin: calculate constant raster properties used in featureset up front at datasource creation (no real speed boost but less repeated code) --- plugins/input/gdal/gdal_datasource.cpp | 35 ++++++++++++++---- plugins/input/gdal/gdal_datasource.hpp | 3 ++ plugins/input/gdal/gdal_featureset.cpp | 49 ++++++++++++-------------- plugins/input/gdal/gdal_featureset.hpp | 10 +++++- 4 files changed, 63 insertions(+), 34 deletions(-) diff --git a/plugins/input/gdal/gdal_datasource.cpp b/plugins/input/gdal/gdal_datasource.cpp index 471cc8a94..b3d7a8c0e 100644 --- a/plugins/input/gdal/gdal_datasource.cpp +++ b/plugins/input/gdal/gdal_datasource.cpp @@ -48,6 +48,11 @@ using mapnik::datasource_exception; */ inline GDALDataset *gdal_datasource::open_dataset() const { + +#ifdef MAPNIK_DEBUG + std::clog << "GDAL Plugin: opening: " << dataset_name_ << "\n"; +#endif + GDALDataset *dataset; #if GDAL_VERSION_NUM >= 1600 if (shared_dataset_) @@ -88,18 +93,34 @@ gdal_datasource::gdal_datasource(parameters const& params) GDALDataset *dataset = open_dataset(); - // TODO: Make more class attributes from geotransform... + nbands_ = dataset->GetRasterCount(); width_ = dataset->GetRasterXSize(); height_ = dataset->GetRasterYSize(); + double tr[6]; dataset->GetGeoTransform(tr); - double dx = tr[1]; - double dy = tr[5]; + if (tr[2] !=0 || tr[4] != 0) + { + throw datasource_exception("GDAL Plugin: only 'north up' images are supported"); + } + + dx_ = tr[1]; + dy_ = tr[5]; + double x0 = tr[0]; double y0 = tr[3]; - double x1 = tr[0] + width_ * dx + height_ *tr[2]; - double y1 = tr[3] + width_ *tr[4] + height_ * dy; + double x1 = tr[0] + width_ * dx_ + height_ *tr[2]; + double y1 = tr[3] + width_ *tr[4] + height_ * dy_; + + /* + double x0 = tr[0] + (height_) * tr[2]; // minx + double y0 = tr[3] + (height_) * tr[5]; // miny + + double x1 = tr[0] + (width_) * tr[1]; // maxx + double y1 = tr[3] + (width_) * tr[4]; // maxy + */ + extent_.init(x0,y0,x1,y1); GDALClose(dataset); @@ -135,11 +156,11 @@ layer_descriptor gdal_datasource::get_descriptor() const featureset_ptr gdal_datasource::features(query const& q) const { gdal_query gq = q; - return featureset_ptr(new gdal_featureset(*open_dataset(), band_, gq, filter_factor_)); + return featureset_ptr(new gdal_featureset(*open_dataset(), band_, gq, extent_, width_, height_, nbands_, dx_, dy_, filter_factor_)); } featureset_ptr gdal_datasource::features_at_point(coord2d const& pt) const { gdal_query gq = pt; - return featureset_ptr(new gdal_featureset(*open_dataset(), band_, gq, filter_factor_)); + return featureset_ptr(new gdal_featureset(*open_dataset(), band_, gq, extent_, width_, height_, nbands_, dx_, dy_, filter_factor_)); } diff --git a/plugins/input/gdal/gdal_datasource.hpp b/plugins/input/gdal/gdal_datasource.hpp index 0f89f6c98..f7ad752ac 100644 --- a/plugins/input/gdal/gdal_datasource.hpp +++ b/plugins/input/gdal/gdal_datasource.hpp @@ -46,6 +46,9 @@ private: mapnik::layer_descriptor desc_; unsigned width_; unsigned height_; + double dx_; + double dy_; + int nbands_; bool shared_dataset_; double filter_factor_; inline GDALDataset *open_dataset() const; diff --git a/plugins/input/gdal/gdal_featureset.cpp b/plugins/input/gdal/gdal_featureset.cpp index 7cd7b1729..dbb914c13 100644 --- a/plugins/input/gdal/gdal_featureset.cpp +++ b/plugins/input/gdal/gdal_featureset.cpp @@ -37,10 +37,18 @@ using mapnik::point_impl; using mapnik::geometry2d; -gdal_featureset::gdal_featureset(GDALDataset & dataset, int band, gdal_query q, double filter_factor) +gdal_featureset::gdal_featureset(GDALDataset & dataset, int band, gdal_query q, + mapnik::box2d extent, double width, double height, int nbands, + double dx, double dy, double filter_factor) : dataset_(dataset), band_(band), gquery_(q), + raster_extent_(extent), + raster_width_(width), + raster_height_(height), + dx_(dx), + dy_(dy), + nbands_(nbands), filter_factor_(filter_factor), first_(true) {} @@ -84,34 +92,23 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) GDALRasterBand * blue = 0; GDALRasterBand * alpha = 0; GDALRasterBand * grey = 0; - - int nbands = dataset_.GetRasterCount(); - unsigned raster_width = dataset_.GetRasterXSize(); - unsigned raster_height = dataset_.GetRasterYSize(); - - // TODO - pull from class attributes... + /* double tr[6]; dataset_.GetGeoTransform(tr); double dx = tr[1]; double dy = tr[5]; - - double x0 = tr[0] + (raster_height) * tr[2]; // minx - double y0 = tr[3] + (raster_height) * tr[5]; // miny - - double x1 = tr[0] + (raster_width) * tr[1]; // maxx - double y1 = tr[3] + (raster_width) * tr[4]; // maxy - - box2d raster_extent(x0,y0,x1,y1); - - CoordTransform t(raster_width,raster_height,raster_extent,0,0); - box2d intersect = raster_extent.intersect(q.get_bbox()); + std::clog << "dx_: " << dx_ << " dx: " << dx << " dy_: " << dy_ << "dy: " << dy << "\n"; + */ + + CoordTransform t(raster_width_,raster_height_,raster_extent_,0,0); + box2d intersect = raster_extent_.intersect(q.get_bbox()); box2d box = t.forward(intersect); //size of resized output pixel in source image domain - double margin_x = 1.0/(fabs(dx)*boost::get<0>(q.resolution())); - double margin_y = 1.0/(fabs(dy)*boost::get<1>(q.resolution())); + double margin_x = 1.0/(fabs(dx_)*boost::get<0>(q.resolution())); + double margin_y = 1.0/(fabs(dy_)*boost::get<1>(q.resolution())); if (margin_x < 1) margin_x = 1.0; if (margin_y < 1) @@ -126,10 +123,10 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) x_off = 0; if (y_off < 0) y_off = 0; - if (end_x > (int)raster_width) - end_x = raster_width; - if (end_y > (int)raster_height) - end_y = raster_height; + if (end_x > (int)raster_width_) + end_x = raster_width_; + if (end_y > (int)raster_height_) + end_y = raster_height_; int width = end_x - x_off; int height = end_y - y_off; // don't process almost invisible data @@ -142,7 +139,7 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) intersect = t.backward(feature_raster_extent); #ifdef MAPNIK_DEBUG - std::clog << "GDAL Plugin: Raster extent=" << raster_extent << "\n"; + std::clog << "GDAL Plugin: Raster extent=" << raster_extent_ << "\n"; std::clog << "GDAL Plugin: View extent=" << intersect << "\n"; std::clog << "GDAL Plugin: Query resolution=" << boost::get<0>(q.resolution()) << "," << boost::get<1>(q.resolution()) << "\n"; std::clog << boost::format("GDAL Plugin: StartX=%d StartY=%d Width=%d Height=%d \n") % x_off % y_off % width % height; @@ -201,7 +198,7 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) else // working with all bands { - for (int i = 0; i < nbands; ++i) + for (int i = 0; i < nbands_; ++i) { GDALRasterBand * band = dataset_.GetRasterBand(i+1); diff --git a/plugins/input/gdal/gdal_featureset.hpp b/plugins/input/gdal/gdal_featureset.hpp index 248797492..d8c204940 100644 --- a/plugins/input/gdal/gdal_featureset.hpp +++ b/plugins/input/gdal/gdal_featureset.hpp @@ -36,7 +36,9 @@ class gdal_featureset : public mapnik::Featureset { public: - gdal_featureset(GDALDataset & dataset, int band, gdal_query q, double filter_factor); + gdal_featureset(GDALDataset & dataset, int band, gdal_query q, + mapnik::box2d extent, double width, double height, int nbands, + double dx, double dy, double filter_factor); virtual ~gdal_featureset(); mapnik::feature_ptr next(); private: @@ -46,6 +48,12 @@ private: GDALDataset & dataset_; int band_; gdal_query gquery_; + mapnik::box2d raster_extent_; + unsigned raster_width_; + unsigned raster_height_; + double dx_; + double dy_; + int nbands_; double filter_factor_; bool first_; };