1 /* 2 * Copyright © 2008 Kristian Høgsberg 3 * Copyright © 2009 Chris Wilson 4 * 5 * Permission to use, copy, modify, distribute, and sell this software and its 6 * documentation for any purpose is hereby granted without fee, provided that 7 * the above copyright notice appear in all copies and that both that copyright 8 * notice and this permission notice appear in supporting documentation, and 9 * that the name of the copyright holders not be used in advertising or 10 * publicity pertaining to distribution of the software without specific, 11 * written prior permission. The copyright holders make no representations 12 * about the suitability of this software for any purpose. It is provided "as 13 * is" without express or implied warranty. 14 * 15 * THE COPYRIGHT HOLDERS DISCLAIM ALL WARRANTIES WITH REGARD TO THIS SOFTWARE, 16 * INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS, IN NO 17 * EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY SPECIAL, INDIRECT OR 18 * CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, 19 * DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER 20 * TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE 21 * OF THIS SOFTWARE. 22 */ 23 24 #include <math.h> 25 #include "blur.h" 26 /* Performs a simple 2D Gaussian blur of standard devation @sigma surface @surface. */ 27 void 28 blur_image_surface (cairo_surface_t *surface, int sigma) 29 { 30 cairo_surface_t *tmp; 31 int width, height; 32 uint32_t *src, *dst; 33 34 if (cairo_surface_status (surface)) 35 return; 36 37 width = cairo_image_surface_get_width (surface); 38 height = cairo_image_surface_get_height (surface); 39 40 switch (cairo_image_surface_get_format (surface)) { 41 case CAIRO_FORMAT_A1: 42 default: 43 /* Don't even think about it! */ 44 return; 45 46 case CAIRO_FORMAT_A8: 47 /* Handle a8 surfaces by effectively unrolling the loops by a 48 * factor of 4 - this is safe since we know that stride has to be a 49 * multiple of uint32_t. */ 50 width /= 4; 51 break; 52 53 case CAIRO_FORMAT_RGB24: 54 case CAIRO_FORMAT_ARGB32: 55 break; 56 } 57 58 tmp = cairo_image_surface_create (CAIRO_FORMAT_ARGB32, width, height); 59 if (cairo_surface_status (tmp)) 60 return; 61 62 src = (uint32_t*)cairo_image_surface_get_data (surface); 63 dst = (uint32_t*)cairo_image_surface_get_data (tmp); 64 65 // according to a paper by Peter Kovesi [1], box filter of width w, equals to Gaussian blur of following sigma: 66 // σ_av = sqrt((w*w-1)/12) 67 // for our 7x7 filter we have σ_av = 2.0. 68 // applying the same Gaussian filter n times results in σ_n = sqrt(n*σ_av*σ_av) [2] 69 // after some trivial math, we arrive at n = ((σ_d)/(σ_av))^2 70 // since it's a box blur filter, n >= 3 71 // 72 // [1]: http://www.peterkovesi.com/papers/FastGaussianSmoothing.pdf 73 // [2]: https://en.wikipedia.org/wiki/Gaussian_blur#Mathematics 74 75 int n = lrintf((sigma*sigma)/(SIGMA_AV*SIGMA_AV)); 76 if (n < 3) n = 3; 77 78 for (int i = 0; i < n; i++) 79 { 80 // horizontal pass includes image transposition: 81 // instead of writing pixel src[x] to dst[x], 82 // we write it to transposed location. 83 // (to be exact: dst[height * current_column + current_row]) 84 #ifdef __SSE2__ 85 blur_impl_horizontal_pass_sse2(src, dst, width, height); 86 blur_impl_horizontal_pass_sse2(dst, src, height, width); 87 #else 88 blur_impl_horizontal_pass_generic(src, dst, width, height); 89 blur_impl_horizontal_pass_generic(dst, src, height, width); 90 #endif 91 } 92 93 cairo_surface_destroy (tmp); 94 cairo_surface_flush (surface); 95 cairo_surface_mark_dirty (surface); 96 } 97 98 void blur_impl_horizontal_pass_generic(uint32_t *src, uint32_t *dst, int width, int height) { 99 uint32_t *o_src = src; 100 for (int row = 0; row < height; row++) { 101 for (int column = 0; column < width; column++, src++) { 102 uint32_t rgbaIn[KERNEL_SIZE + 1]; 103 104 // handle borders 105 int leftBorder = column < HALF_KERNEL; 106 int rightBorder = column > width - HALF_KERNEL; 107 int i = 0; 108 if (leftBorder) { 109 // for kernel size 7x7 and column == 0, we have: 110 // x x x P0 P1 P2 P3 111 // first loop mirrors P{0..3} to fill x's, 112 // second one loads P{0..3} 113 for (; i < HALF_KERNEL - column; i++) 114 rgbaIn[i] = *(src + (HALF_KERNEL - i)); 115 for (; i < KERNEL_SIZE; i++) 116 rgbaIn[i] = *(src - (HALF_KERNEL - i)); 117 } else if (rightBorder) { 118 for (; i < width - column; i++) 119 rgbaIn[i] = *(src + i); 120 for (int k = 0; i < KERNEL_SIZE; i++, k++) 121 rgbaIn[i] = *(src - k); 122 } else { 123 for (; i < KERNEL_SIZE; i++) { 124 if ((uintptr_t) ((src + 4*i - HALF_KERNEL) + 1) 125 > (uintptr_t) (o_src + (height * width))) 126 break; 127 rgbaIn[i] = *(src + i - HALF_KERNEL); 128 } 129 } 130 131 uint32_t acc[4] = {0}; 132 133 for (i = 0; i < KERNEL_SIZE; i++) { 134 acc[0] += (rgbaIn[i] & 0xFF000000) >> 24; 135 acc[1] += (rgbaIn[i] & 0x00FF0000) >> 16; 136 acc[2] += (rgbaIn[i] & 0x0000FF00) >> 8; 137 acc[3] += (rgbaIn[i] & 0x000000FF) >> 0; 138 } 139 140 for(i = 0; i < 4; i++) 141 acc[i] *= 1.0/KERNEL_SIZE; 142 143 *(dst + height * column + row) = (acc[0] << 24) | 144 (acc[1] << 16) | 145 (acc[2] << 8 ) | 146 (acc[3] << 0); 147 } 148 } 149 }
