sij/mapnik
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

gdal plugin: calculate constant raster properties used in featureset up front at datasource creation (no real speed boost but less repeated code)

Browse source
This commit is contained in:
Dane Springmeyer 2010-09-18 19:19:27 +00:00
parent 985e54379e
commit 00f473de63
4 changed files with 63 additions and 34 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

35
plugins/input/gdal/gdal_datasource.cpp
View file

@ -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_));
}

3
plugins/input/gdal/gdal_datasource.hpp
View file

@ -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;

47
plugins/input/gdal/gdal_featureset.cpp
View file

@ -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<double> 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) {}
@ -85,33 +93,22 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q)
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];
std::clog << "dx_: " << dx_ << " dx: " << dx << " dy_: " << dy_ << "dy: " << dy << "\n";
*/
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<double> raster_extent(x0,y0,x1,y1);
CoordTransform t(raster_width,raster_height,raster_extent,0,0);
box2d<double> intersect = raster_extent.intersect(q.get_bbox());
CoordTransform t(raster_width_,raster_height_,raster_extent_,0,0);
box2d<double> intersect = raster_extent_.intersect(q.get_bbox());
box2d<double> 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);

10
plugins/input/gdal/gdal_featureset.hpp
View file

@ -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<double> 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<double> raster_extent_;
unsigned raster_width_;
unsigned raster_height_;
double dx_;
double dy_;
int nbands_;
double filter_factor_;
bool first_;
};