2011-11-20 23 views
29

Mi pregunta es muy cercano a esta pregunta: How do I gaussian blur an image without using any in-built gaussian functions?Ejecución gaussiano - ¿Cómo calcular matriz de convolución (kernel)

La respuesta a esta pregunta es muy buena, pero no da un ejemplo de cálculo de una realidad núcleo de filtro gaussiano real. La respuesta proporciona un kernel arbitrario y muestra cómo aplicar el filtro utilizando ese núcleo, pero no cómo calcular un kernel real. Estoy tratando de implementar un desenfoque gaussiano en C++ o Matlab desde cero, así que necesito saber cómo calcular el kernel desde cero.

Agradecería que alguien pudiera calcular un núcleo de filtro Gaussiano real utilizando cualquier matriz de imagen de ejemplo pequeña.

+0

¿Ha leído esto: http://en.wikipedia.org/wiki/Gaussian_function? –

+0

O incluso esto: http://en.wikipedia.org/wiki/Gaussian_blur – Bart

+0

Sí, he pasado mucho tiempo intentando entenderlos. Lo que necesito es un ejemplo paso a paso. Después de que lo entiendo, probablemente agregue el ejemplo a la página Gaussian Blur. – gsingh2011

Respuesta

33

Se puede crear un kernel gaussiana desde cero como se indica en la documentación de MATLAB fspecial. Lea la fórmula de creación de kernel de Gauss en la parte de algoritmos de esa página y siga el código que se muestra a continuación. El código es crear una matriz m-por-n con sigma = 1.

m = 5; n = 5; 
sigma = 1; 
[h1, h2] = meshgrid(-(m-1)/2:(m-1)/2, -(n-1)/2:(n-1)/2); 
hg = exp(- (h1.^2+h2.^2)/(2*sigma^2)); 
h = hg ./ sum(hg(:)); 

h = 

    0.0030 0.0133 0.0219 0.0133 0.0030 
    0.0133 0.0596 0.0983 0.0596 0.0133 
    0.0219 0.0983 0.1621 0.0983 0.0219 
    0.0133 0.0596 0.0983 0.0596 0.0133 
    0.0030 0.0133 0.0219 0.0133 0.0030 

Obsérvese que esto se puede hacer mediante la incorporada en fspecial de la siguiente manera:

fspecial('gaussian', [m n], sigma) 
ans = 

    0.0030 0.0133 0.0219 0.0133 0.0030 
    0.0133 0.0596 0.0983 0.0596 0.0133 
    0.0219 0.0983 0.1621 0.0983 0.0219 
    0.0133 0.0596 0.0983 0.0596 0.0133 
    0.0030 0.0133 0.0219 0.0133 0.0030 

creo que es fácil de implementar en cualquier idioma que desee.

EDITAR: Permítame también agregar los valores de h1 y h2 para el caso dado, ya que puede que no esté familiarizado con meshgrid si codifica en C++.

h1 = 

    -2 -1  0  1  2 
    -2 -1  0  1  2 
    -2 -1  0  1  2 
    -2 -1  0  1  2 
    -2 -1  0  1  2 

h2 = 

    -2 -2 -2 -2 -2 
    -1 -1 -1 -1 -1 
    0  0  0  0  0 
    1  1  1  1  1 
    2  2  2  2  2 
+1

Escribí [h1, h2] = meshgrid (- (m-1)/2: (m-1)/2, - (n-1)/2: (n-1)/2) y obtuve un rango de h1 de -2 a 2, no de -1.5 a 1.5. Mismo problema con h2. Pero mi resultado es el mismo. Además, ¿por qué usaste los valores de malla de malla como los valores en la fórmula? ¿Qué representa esto si estuvieras calculando esto para una imagen? – gsingh2011

+0

¡Tienes razón! Cambié 'm' y' n' a 4 para ver si el código funciona y luego copié los valores para este caso en lugar de darles el valor 5. Lo arreglé, gracias. – petrichor

+0

Los valores se calculan en una cuadrícula donde el que está en el centro es el origen, que es h1 == 0 y h2 == 0 en nuestro caso. Todos los otros pares representan las otras coordenadas cuando mira el elemento de valores h1, h2 por elemento. Durante el filtrado, puede pensar que esta cuadrícula se colocará en un píxel de la imagen donde el origen de la cuadrícula se ajusta exactamente en el píxel. Puede leer la respuesta de Goz en el enlace que proporcionó en su pregunta para conocer los detalles. – petrichor

22

Es tan simple como suena:

double sigma = 1; 
int W = 5; 
double kernel[W][W]; 
double mean = W/2; 
double sum = 0.0; // For accumulating the kernel values 
for (int x = 0; x < W; ++x) 
    for (int y = 0; y < W; ++y) { 
     kernel[x][y] = exp(-0.5 * (pow((x-mean)/sigma, 2.0) + pow((y-mean)/sigma,2.0))) 
         /(2 * M_PI * sigma * sigma); 

     // Accumulate the kernel values 
     sum += kernel[x][y]; 
    } 

// Normalize the kernel 
for (int x = 0; x < W; ++x) 
    for (int y = 0; y < W; ++y) 
     kernel[x][y] /= sum; 
+11

Esto tiene defectos: también es necesario normalizar el kernel, o la imagen se oscurece dependiendo de W y sigma. en pocas palabras: obtenga la suma de los valores del kernel y divida cada valor del kernel por esa suma. – Rookie

+2

@Rookie: he decidido modificar esta publicación y agregar la normalización. Esto es para permitir que aquellos que quieren una solución C/C++ usen esto directamente. ¡Buena atrapada! – rayryeng

+0

Parece no ser correcto cuando m, n son pares, en comparación con el resultado de 'fspecial'. –

13

Para poner en práctica el gaussian blur usted simplemente toma el gaussian function y calcular un valor para cada uno de los elementos en su núcleo.

Por lo general, desea asignar el peso máximo al elemento central en su kernel y valores cercanos a cero para los elementos en los bordes del kernel. Esto implica que el núcleo debe tener una altura impar (o ancho) para garantizar que realmente exista un elemento central.

para calcular los elementos del núcleo reales es posible que la escala de la campana de Gauss a la red núcleo (elegir una arbitraria por ejemplo sigma = 1 y un rango arbitrario por ejemplo -2*sigma ... 2*sigma) y normalizarla, S.T. los elementos suman uno. Para lograr esto, si desea admitir tamaños de kernel arbitrarios, es posible que desee adaptar sigma al tamaño de kernel requerido.

Aquí está un ejemplo de C++:

#include <cmath> 
#include <vector> 
#include <iostream> 
#include <iomanip> 

double gaussian (double x, double mu, double sigma) { 
    return std::exp(-(((x-mu)/(sigma))*((x-mu)/(sigma)))/2.0); 
} 

typedef std::vector<double> kernel_row; 
typedef std::vector<kernel_row> kernel_type; 

kernel_type produce2dGaussianKernel (int kernelRadius) { 
    double sigma = kernelRadius/2.; 
    kernel_type kernel2d(2*kernelRadius+1, kernel_row(2*kernelRadius+1)); 
    double sum = 0; 
    // compute values 
    for (int row = 0; row < kernel2d.size(); row++) 
    for (int col = 0; col < kernel2d[row].size(); col++) { 
     double x = gaussian(row, kernelRadius, sigma) 
       * gaussian(col, kernelRadius, sigma); 
     kernel2d[row][col] = x; 
     sum += x; 
    } 
    // normalize 
    for (int row = 0; row < kernel2d.size(); row++) 
    for (int col = 0; col < kernel2d[row].size(); col++) 
     kernel2d[row][col] /= sum; 
    return kernel2d; 
} 

int main() { 
    kernel_type kernel2d = produce2dGaussianKernel(3); 
    std::cout << std::setprecision(5) << std::fixed; 
    for (int row = 0; row < kernel2d.size(); row++) { 
    for (int col = 0; col < kernel2d[row].size(); col++) 
     std::cout << kernel2d[row][col] << ' '; 
    std::cout << '\n'; 
    } 
} 

La salida es:

$ g++ test.cc && ./a.out 
0.00134 0.00408 0.00794 0.00992 0.00794 0.00408 0.00134 
0.00408 0..02412 0.03012 0.02412 0..00408 
0.00794 0.02412 0.04698 0.05867 0.04698 0.02412 0.00794 
0.00992 0.03012 0.05867 0.07327 0.05867 0.03012 0.00992 
0.00794 0.02412 0.04698 0.05867 0.04698 0.02412 0.00794 
0.00408 0..02412 0.03012 0.02412 0..00408 
0.00134 0.00408 0.00794 0.00992 0.00794 0.00408 0.00134 

Como una simplificación que no es necesario utilizar un 2d-núcleo. Más fácil de implementar y también más eficiente de calcular es usar dos núcleos 1d ortogonales. Esto es posible debido a la asociatividad de este tipo de convolución lineal (separabilidad lineal). Es posible que también desee ver this section del artículo de wikipedia correspondiente.


Aquí es lo mismo en Python (con la esperanza de que alguien podría encontrar útil):

from math import exp 

def gaussian(x, mu, sigma): 
    return exp(-(((x-mu)/(sigma))**2)/2.0) 

#kernel_height, kernel_width = 7, 7 
kernel_radius = 3 # for an 7x7 filter 
sigma = kernel_radius/2. # for [-2*sigma, 2*sigma] 

# compute the actual kernel elements 
hkernel = [gaussian(x, kernel_radius, sigma) for x in range(2*kernel_radius+1)] 
vkernel = [x for x in hkernel] 
kernel2d = [[xh*xv for xh in hkernel] for xv in vkernel] 

# normalize the kernel elements 
kernelsum = sum([sum(row) for row in kernel2d]) 
kernel2d = [[x/kernelsum for x in row] for row in kernel2d] 

for line in kernel2d: 
    print ["%.3f" % x for x in line] 

produce el núcleo:

['0.001', '0.004', '0.008', '0.010', '0.008', '0.004', '0.001'] 
['0.004', '0.012', '0.024', '0.030', '0.024', '0.012', '0.004'] 
['0.008', '0.024', '0.047', '0.059', '0.047', '0.024', '0.008'] 
['0.010', '0.030', '0.059', '0.073', '0.059', '0.030', '0.010'] 
['0.008', '0.024', '0.047', '0.059', '0.047', '0.024', '0.008'] 
['0.004', '0.012', '0.024', '0.030', '0.024', '0.012', '0.004'] 
['0.001', '0.004', '0.008', '0.010', '0.008', '0.004', '0.001'] 
+0

¿Dónde está PI aquí? – user2023370

+0

@ user2023370 ¿Por qué crees que debería haber un PI allí? – moooeeeep

0

desenfoque gaussiano en Python usando la biblioteca de imágenes PIL. Para más información leer esto: http://blog.ivank.net/fastest-gaussian-blur.html

from PIL import Image 
import math 

# img = Image.open('input.jpg').convert('L') 
# r = radiuss 
def gauss_blur(img, r): 
    imgData = list(img.getdata()) 

    bluredImg = Image.new(img.mode, img.size) 
    bluredImgData = list(bluredImg.getdata()) 

    rs = int(math.ceil(r * 2.57)) 

    for i in range(0, img.height): 
     for j in range(0, img.width): 
      val = 0 
      wsum = 0 
      for iy in range(i - rs, i + rs + 1): 
       for ix in range(j - rs, j + rs + 1): 
        x = min(img.width - 1, max(0, ix)) 
        y = min(img.height - 1, max(0, iy)) 
        dsq = (ix - j) * (ix - j) + (iy - i) * (iy - i) 
        weight = math.exp(-dsq/(2 * r * r))/(math.pi * 2 * r * r) 
        val += imgData[y * img.width + x] * weight 
        wsum += weight 
      bluredImgData[i * img.width + j] = round(val/wsum) 

    bluredImg.putdata(bluredImgData) 
    return bluredImg 
0
function kernel = gauss_kernel(m, n, sigma) 
% Generating Gauss Kernel 

x = -(m-1)/2 : (m-1)/2; 
y = -(n-1)/2 : (n-1)/2; 

for i = 1:m 
    for j = 1:n 
     xx(i,j) = x(i); 
     yy(i,j) = y(j); 
    end 
end 

kernel = exp(-(xx.*xx + yy.*yy)/(2*sigma*sigma)); 

% Normalize the kernel 
kernel = kernel/sum(kernel(:)); 

% Corresponding function in MATLAB 
% fspecial('gaussian', [m n], sigma) 
+0

Agregue algunos comentarios a su código, será útil para otras personas. – HDJEMAI

Cuestiones relacionadas