From 65eb4d005b72948245c9bae9823930fb85bde708 Mon Sep 17 00:00:00 2001 From: Sam Hocevar Date: Sat, 21 Jun 2014 15:30:34 +0000 Subject: [PATCH] image: generate Gaussian kernels. --- src/image/filter/blur.cpp | 92 --------------------------------------- src/image/stock.cpp | 79 +++++++++++++++++++++++++++++++++ src/lol/image/image.h | 1 + 3 files changed, 80 insertions(+), 92 deletions(-) diff --git a/src/image/filter/blur.cpp b/src/image/filter/blur.cpp index 27e5d360..ba76b828 100644 --- a/src/image/filter/blur.cpp +++ b/src/image/filter/blur.cpp @@ -34,98 +34,6 @@ #define TEMPLATE_FILE "filter/blur.c" #include "pipi-template.h" -/* Any standard deviation below this value will be rounded up, in order - * to avoid ridiculously low values. exp(-1/(2*0.2*0.2)) is < 10^-5 so - * there is little chance that any value below 0.2 will be useful. */ -#define BLUR_EPSILON 0.2 - -pipi_image_t *pipi_gaussian_blur(pipi_image_t *src, float radius) -{ - return pipi_gaussian_blur_ext(src, radius, radius, 0.0, 0.0, 0.0); -} - -pipi_image_t *pipi_gaussian_blur_ext(pipi_image_t *src, float rx, float ry, - float angle, float dx, float dy) -{ - pipi_image_t *ret; - double *kernel; - double Kx, Ky, t = 0.0, sint, cost, bbx, bby; - int i, j, krx, kry, m, n; - - if(rx < BLUR_EPSILON) rx = BLUR_EPSILON; - if(ry < BLUR_EPSILON) ry = BLUR_EPSILON; - - sint = sin(angle * M_PI / 180.); - cost = cos(angle * M_PI / 180.); - - /* Compute the final ellipse's bounding box */ - bbx = sqrt(rx * rx * cost * cost + ry * ry * sint * sint); - bby = sqrt(ry * ry * cost * cost + rx * rx * sint * sint); - - /* FIXME: the kernel becomes far too big with large values of dx, because - * we grow both left and right. Fix the growing direction. */ - krx = (int)(3. * bbx + .99999 + ceil(abs(dx))); - m = 2 * krx + 1; - Kx = -1. / (2. * rx * rx); - - kry = (int)(3. * bby + .99999 + ceil(abs(dy))); - n = 2 * kry + 1; - Ky = -1. / (2. * ry * ry); - - kernel = malloc(m * n * sizeof(double)); - - for(j = -kry; j <= kry; j++) - { - for(i = -krx; i <= krx; i++) - { - /* FIXME: this level of interpolation sucks. We should - * interpolate on the full NxN grid for better quality. */ - static double const samples[] = - { - .0, .0, 1, - -.40, -.40, 0.8, - -.30, .0, 0.9, - -.40, .40, 0.8, - .0, .30, 0.9, - .40, .40, 0.8, - .30, .0, 0.9, - .40, -.40, 0.8, - .0, -.30, 0.9, - }; - double u, v, ex, ey, d = 0.; - unsigned int k; - - for(k = 0; k < sizeof(samples) / sizeof(*samples) / 3; k++) - { - u = ((double)i + samples[k * 3]) * cost - - ((double)j + samples[k * 3 + 1]) * sint + dx; - v = ((double)i + samples[k * 3]) * sint - + ((double)j + samples[k * 3 + 1]) * cost + dy; - ex = Kx * u * u; - ey = Ky * v * v; - d += samples[k * 3 + 2] * exp(ex + ey); - - /* Do not interpolate if this is a standard gaussian. */ - if(!dx && !dy && !angle) - break; - } - - kernel[(j + kry) * m + i + krx] = d; - t += d; - } - } - - for(j = 0; j < n; j++) - for(i = 0; i < m; i++) - kernel[j * m + i] /= t; - - ret = pipi_convolution(src, m, n, kernel); - - free(kernel); - - return ret; -} - pipi_image_t *pipi_box_blur(pipi_image_t *src, int size) { return pipi_box_blur_ext(src, size, size); diff --git a/src/image/stock.cpp b/src/image/stock.cpp index 78715dee..153a4b16 100644 --- a/src/image/stock.cpp +++ b/src/image/stock.cpp @@ -277,5 +277,84 @@ Array2D Image::EdiffKernel(EdiffAlgorithm algorithm) return ret; } +/* Any standard deviation below this value will be rounded up, in order + * to avoid ridiculously low values. exp(-1/(2*0.2*0.2)) is < 10^-5 so + * there is little chance that any value below 0.2 will be useful. */ +#define BLUR_EPSILON 0.2f + +Array2D Image::GaussianKernel(vec2 radius, float angle, vec2 delta) +{ + Array2D kernel; + + if (radius.x < BLUR_EPSILON) + radius.x = BLUR_EPSILON; + if (radius.y < BLUR_EPSILON) + radius.y = BLUR_EPSILON; + + float const sint = lol::sin(angle); + float const cost = lol::cos(angle); + + /* Compute the final ellipse's bounding box */ + float const bbx = lol::sqrt(sq(radius.x * cost) + sq(radius.y * sint)); + float const bby = lol::sqrt(sq(radius.y * cost) + sq(radius.x * sint)); + + /* FIXME: the kernel becomes far too big with large values of dx, because + * we grow both left and right. Fix the growing direction. */ + int const krx = (int)(3.f * bbx + .99999f + lol::ceil(lol::abs(delta.x))); + int const kry = (int)(3.f * bby + .99999f + lol::ceil(lol::abs(delta.y))); + ivec2 size(2 * krx + 1, 2 * kry + 1); + float const Kx = -1.f / (2.f * radius.x * radius.x); + float const Ky = -1.f / (2.f * radius.y * radius.y); + + kernel.SetSize(size); + + float t = 0.f; + + for (int j = -kry; j <= kry; j++) + { + for (int i = -krx; i <= krx; i++) + { + /* FIXME: this level of interpolation sucks. We should + * interpolate on the full NxN grid for better quality. */ + static vec3 const samples[] = + { + vec3( 0.0f, 0.0f, 1.0f), + vec3(-0.4f, -0.4f, 0.8f), + vec3(-0.3f, 0.0f, 0.9f), + vec3(-0.4f, 0.4f, 0.8f), + vec3( 0.0f, 0.3f, 0.9f), + vec3( 0.4f, 0.4f, 0.8f), + vec3( 0.3f, 0.0f, 0.9f), + vec3( 0.4f, -0.4f, 0.8f), + vec3( 0.0f, -0.3f, 0.9f), + }; + + float d = 0.f; + + for (auto p : samples) + { + float u = (i + p.x) * cost - (j + p.y) * sint + delta.x; + float v = (i + p.x) * sint + (j + p.y) * cost + delta.y; + float ex = Kx * u * u; + float ey = Ky * v * v; + d += p.z * lol::exp(ex + ey); + + /* Do not interpolate if this is a standard gaussian. */ + if (!delta.x && !delta.y && !angle) + break; + } + + kernel[i + krx][j + kry] = d; + t += d; + } + } + + for (int j = 0; j < size.y; j++) + for (int i = 0; i < size.x; i++) + kernel[i][j] *= (1.f / t); + + return kernel; +} + } /* namespace lol */ diff --git a/src/lol/image/image.h b/src/lol/image/image.h index 514dde55..75dcdb1d 100644 --- a/src/lol/image/image.h +++ b/src/lol/image/image.h @@ -95,6 +95,7 @@ public: static Array2D HalftoneKernel(ivec2 size); static Array2D EdiffKernel(EdiffAlgorithm algorithm); static Array2D NormalizeKernel(Array2D const &kernel); + static Array2D GaussianKernel(vec2 radius, float angle, vec2 delta); /* Rendering */ bool Stock(char const *desc);