From 25e4bb3f6ca28f7c15e880530ea9b103005b9d7f Mon Sep 17 00:00:00 2001 From: Blake Thompson Date: Tue, 20 Mar 2018 17:01:32 -0400 Subject: [PATCH] A fix for two distinct issues associated with gdal featuresets, the first is overviews were not properly being utilized based on the resolution of the final image requested. The second is that allocation of far too much memory could be possible in GDAL to attempt to do resampling internally in mapnik. This now has a hard cap so that we allocate less memory in these situations but are still able to resample internally in mapnik. --- include/mapnik/image.hpp | 2 +- include/mapnik/image_impl.hpp | 6 +- plugins/input/gdal/gdal_featureset.cpp | 140 +++++++++++++++++++++++-- src/image.cpp | 2 +- test/data | 2 +- test/data-visual | 2 +- 6 files changed, 140 insertions(+), 14 deletions(-) diff --git a/include/mapnik/image.hpp b/include/mapnik/image.hpp index 6f7b6d730..2cc952681 100644 --- a/include/mapnik/image.hpp +++ b/include/mapnik/image.hpp @@ -74,7 +74,7 @@ public: static constexpr image_dtype dtype = T::id; static constexpr std::size_t pixel_size = sizeof(pixel_type); private: - detail::image_dimensions<65535> dimensions_; + detail::image_dimensions<4294836225> dimensions_; detail::buffer buffer_; double offset_; double scaling_; diff --git a/include/mapnik/image_impl.hpp b/include/mapnik/image_impl.hpp index 002cef08c..856bc6990 100644 --- a/include/mapnik/image_impl.hpp +++ b/include/mapnik/image_impl.hpp @@ -39,8 +39,10 @@ image_dimensions::image_dimensions(int width, int height) : width_(width), height_(height) { - if (width < 0 || static_cast(width) > max_size) throw std::runtime_error("Invalid width for image dimensions requested"); - if (height < 0 || static_cast(height) > max_size) throw std::runtime_error("Invalid height for image dimensions requested"); + int64_t area = (int64_t)width * (int64_t)height; + if (width < 0) throw std::runtime_error("Invalid width for image dimensions requested"); + if (height < 0) throw std::runtime_error("Invalid height for image dimensions requested"); + if (area > max_size) throw std::runtime_error("Image area too large based on image dimensions"); } template diff --git a/plugins/input/gdal/gdal_featureset.cpp b/plugins/input/gdal/gdal_featureset.cpp index eb3c78d06..4363e2f93 100644 --- a/plugins/input/gdal/gdal_featureset.cpp +++ b/plugins/input/gdal/gdal_featureset.cpp @@ -153,6 +153,8 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) //size of resized output pixel in source image domain double margin_x = 1.0 / (std::fabs(dx_) * std::get<0>(q.resolution())); double margin_y = 1.0 / (std::fabs(dy_) * std::get<1>(q.resolution())); + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: margin_x=" << margin_x; + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: margin_y=" << margin_y; if (margin_x < 1) { margin_x = 1.0; @@ -167,6 +169,10 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) int y_off = rint(box.miny() - margin_y); int end_x = rint(box.maxx() + margin_x); int end_y = rint(box.maxy() + margin_y); + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: x_off=" << x_off; + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: y_off=" << y_off; + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: end_x=" << end_x; + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: end_y=" << end_y; //clip to available data if (x_off < 0) @@ -185,21 +191,139 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) { end_y = raster_height_; } + + // width and height of the portion of the source image we are requesting int width = end_x - x_off; int height = end_y - y_off; + // In many cases we want GDAL to simply return the exact image so we + // can handle resampling internally in mapnik. In other cases such as + // when overviews exist or when the image allocated might be too large + // we want to utilize some resampling in GDAL instead. + int im_height = height; + int im_width = width; + double im_offset_x = x_off; + double im_offset_y = y_off; + int current_width = (int)raster_width_; + int current_height = (int)raster_height_; + + // loop through overviews -- snap up in resolution to closest overview if necessary + // we find an image size that most resembles the resolution of our output image. + double width_res = std::get<0>(q.resolution()); + double height_res = std::get<1>(q.resolution()); + int res_adjusted_raster_width = static_cast(std::floor(((double)raster_width_ * width_res) + .5)); + int res_adjusted_raster_height = static_cast(std::floor(((double)raster_height_ * height_res) + .5)); + if (band_ > 0 && band_ < nbands_) + { + GDALRasterBand * band = dataset_.GetRasterBand(band_); + int band_overviews = band->GetOverviewCount(); + if (band_overviews > 0) + { + for (int b = 0; b < band_overviews; b++) + { + GDALRasterBand * overview = band->GetOverview(b); + int overview_width = overview->GetXSize(); + int overview_height = overview->GetYSize(); + if ((overview_width < current_width || overview_height < current_height) && + res_adjusted_raster_width <= overview_width && + res_adjusted_raster_height <= overview_height) + { + current_width = overview_width; + current_height = overview_height; + } + } + } + } + else + { + for (int i = 0; i < nbands_; ++i) + { + GDALRasterBand * band = dataset_.GetRasterBand(i + 1); + int band_overviews = band->GetOverviewCount(); + if (band_overviews > 0) + { + for (int b = 0; b < band_overviews; b++) + { + GDALRasterBand * overview = band->GetOverview(b); + int overview_width = overview->GetXSize(); + int overview_height = overview->GetYSize(); + if ((overview_width < current_width || overview_height < current_height) && + res_adjusted_raster_width <= overview_width && + res_adjusted_raster_height <= overview_height) + { + current_width = overview_width; + current_height = overview_height; + } + } + } + } + } + if (current_width != (int)raster_width_ || current_height != (int)raster_height_) + { + if (current_width != (int)raster_width_) + { + double ratio = (double)current_width / (double)raster_width_; + int adjusted_width = static_cast(std::floor((ratio * im_width) + 0.5)); + double adjusted_ratio = (double)adjusted_width / (double)im_width; + im_offset_x = adjusted_ratio * im_offset_x; + im_width = adjusted_width; + } + if (current_height != (int)raster_height_) + { + double ratio = (double)current_height / (double)raster_height_; + int adjusted_height = static_cast(std::floor((ratio * im_height) + 0.5)); + double adjusted_ratio = (double)adjusted_height / (double)im_height; + im_offset_y = adjusted_ratio * im_offset_y; + im_height = adjusted_height; + } + } + + // Maximum memory limitation for image will be simply based on the maximum + // area we allow for an image. The true memory footprint therefore will vary based + // on the type of imagery that exists. + // max_im_area based on 50 mb limit for RGBA + constexpr int64_t max_im_area = (50*1024*1024) / 4; + int64_t im_area = (int64_t)im_width * (int64_t)im_height; + if (im_area > max_im_area) + { + int adjusted_width = static_cast(std::round(std::sqrt(max_im_area * ((double)im_width / (double)im_height)))); + int adjusted_height = static_cast(std::round(std::sqrt(max_im_area * ((double)im_height / (double)im_width)))); + if (adjusted_width < 1) + { + adjusted_width = 1; + } + if (adjusted_height < 1) + { + adjusted_height = 1; + } + double ratio_x = (double)adjusted_width / (double)im_width; + double ratio_y = (double)adjusted_height / (double)im_height; + im_offset_x = ratio_x * im_offset_x; + im_offset_y = ratio_y * im_offset_y; + im_width = adjusted_width; + im_height = adjusted_height; + current_width = static_cast(std::floor((ratio_x * current_width) + 0.5)); + current_height = static_cast(std::floor((ratio_y * current_height) + 0.5)); + } + //calculate actual box2d of returned raster - box2d feature_raster_extent(x_off, y_off, x_off + width, y_off + height); - feature_raster_extent = t.backward(feature_raster_extent); + view_transform t2(current_width, current_height, raster_extent_, 0, 0); + box2d feature_raster_extent(im_offset_x, im_offset_y, im_offset_x + im_width, im_offset_y + im_height); + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: Feature Raster extent=" << feature_raster_extent; + feature_raster_extent = t2.backward(feature_raster_extent); MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: Raster extent=" << raster_extent_; + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: Feature Raster extent=" << feature_raster_extent; MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: View extent=" << intersect; MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: Query resolution=" << std::get<0>(q.resolution()) << "," << std::get<1>(q.resolution()); MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: StartX=" << x_off << " StartY=" << y_off << " Width=" << width << " Height=" << height; + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: IM StartX=" << im_offset_x << " StartY=" << im_offset_y << " Width=" << im_width << " Height=" << im_height; if (width > 0 && height > 0) { - MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: Image Size=(" << width << "," << height << ")"; + + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: Requested Image Size=(" << width << "," << height << ")"; + MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: Image Size=(" << im_width << "," << im_height << ")"; MAPNIK_LOG_DEBUG(gdal) << "gdal_featureset: Reading band=" << band_; if (band_ > 0) // we are querying a single band { @@ -215,7 +339,7 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) { case GDT_Byte: { - mapnik::image_gray8 image(width, height); + mapnik::image_gray8 image(im_width, im_height); image.set(std::numeric_limits::max()); raster_nodata = band->GetNoDataValue(&raster_has_nodata); raster_io_error = band->RasterIO(GF_Read, x_off, y_off, width, height, @@ -235,7 +359,7 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) case GDT_Float64: case GDT_Float32: { - mapnik::image_gray32f image(width, height); + mapnik::image_gray32f image(im_width, im_height); image.set(std::numeric_limits::max()); raster_nodata = band->GetNoDataValue(&raster_has_nodata); raster_io_error = band->RasterIO(GF_Read, x_off, y_off, width, height, @@ -254,7 +378,7 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) } case GDT_UInt16: { - mapnik::image_gray16 image(width, height); + mapnik::image_gray16 image(im_width, im_height); image.set(std::numeric_limits::max()); raster_nodata = band->GetNoDataValue(&raster_has_nodata); raster_io_error = band->RasterIO(GF_Read, x_off, y_off, width, height, @@ -274,7 +398,7 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) default: case GDT_Int16: { - mapnik::image_gray16s image(width, height); + mapnik::image_gray16s image(im_width, im_height); image.set(std::numeric_limits::max()); raster_nodata = band->GetNoDataValue(&raster_has_nodata); raster_io_error = band->RasterIO(GF_Read, x_off, y_off, width, height, @@ -295,7 +419,7 @@ feature_ptr gdal_featureset::get_feature(mapnik::query const& q) } else // working with all bands { - mapnik::image_rgba8 image(width, height); + mapnik::image_rgba8 image(im_width, im_height); image.set(std::numeric_limits::max()); for (int i = 0; i < nbands_; ++i) { diff --git a/src/image.cpp b/src/image.cpp index 8eceb1da9..d1aca683f 100644 --- a/src/image.cpp +++ b/src/image.cpp @@ -84,7 +84,7 @@ void buffer::swap(buffer & rhs) std::swap(owns_, rhs.owns_); } -template struct MAPNIK_DECL image_dimensions<65535>; +template struct MAPNIK_DECL image_dimensions<4294836225>; } // end ns detail diff --git a/test/data b/test/data index 96d52a79d..3a636431a 160000 --- a/test/data +++ b/test/data @@ -1 +1 @@ -Subproject commit 96d52a79df7e20c5543321baa177811dc72bc93a +Subproject commit 3a636431ab189285378686c596be16a440c0808f diff --git a/test/data-visual b/test/data-visual index d374baac6..a6ec8c091 160000 --- a/test/data-visual +++ b/test/data-visual @@ -1 +1 @@ -Subproject commit d374baac6cd72e7de665a27ded9b129e92661e18 +Subproject commit a6ec8c0919eda8ac2f34ede5c4a53240d0bb275a