2011-11-20 14 views
29

Moje pytanie jest bardzo zbliżony do tego pytania: How do I gaussian blur an image without using any in-built gaussian functions?wykonawcze Gaussian Blur - Jak obliczyć macierz splot (jądro)

Odpowiedź na to pytanie jest bardzo dobra, ale to nie daje przykład faktycznie obliczania prawdziwe jądro filtru Gaussa. Odpowiedź daje dowolne jądro i pokazuje, jak zastosować filtr przy użyciu tego jądra, ale nie w jaki sposób obliczyć prawdziwe jądro. Próbuję zaimplementować rozmycie Gaussa w C++ lub Matlab od zera, więc muszę wiedzieć, jak obliczyć jądro od zera.

Byłbym wdzięczny, gdyby ktoś mógł obliczyć prawdziwe jądro filtra Gaussa za pomocą dowolnej małej przykładowej matrycy obrazu.

+0

Czy to czytasz: http://en.wikipedia.org/wiki/Gaussian_function? –

+0

Albo nawet to: http://en.wikipedia.org/wiki/Gaussian_blur – Bart

+0

Tak, spędziłem dużo czasu próbując je zrozumieć. Potrzebuję przykładu "krok po kroku". Po jego zrozumieniu prawdopodobnie dodam przykład do strony Rozmycie gaussowskie. – gsingh2011

Odpowiedz

33

Można utworzyć jądra Gaussa od podstaw, jak wskazano w dokumentacji MATLAB z fspecial. Zapoznaj się z formułą tworzenia jądra Gaussa w części algorytmów na tej stronie i postępuj zgodnie z poniższym kodem. Kod jest stworzenie macierzy M przez N 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 

Zauważmy, że może to być zrobione przez wbudowany fspecial następująco:

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 

myślę jest prosty do wdrożenia tego w dowolnym języku, który lubisz.

EDYCJA: Pozwól mi również dodać wartości h1 i h2 dla danego przypadku, ponieważ możesz być niezaznajomiony z meshgrid, jeśli kodujesz w 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

Wpisałem [h1, h2] = meshgrid (- (m-1)/2: (m-1)/2, - (n-1)/2: (n-1)/2) i uzyskałem zakres h1 od -2 do 2, a nie od -1,5 do 1,5. Ten sam problem z h2. Ale mój wynik jest taki sam. Ponadto, dlaczego użyłeś wartości siatki siatki jako wartości we wzorze? Co to oznacza, jeśli obliczałeś to dla obrazu? – gsingh2011

+0

Masz rację! Zmieniłem 'm' i' n' na 4, aby sprawdzić, czy kod działa, a następnie skopiowałem wartości dla tego przypadku zamiast dawać je dla wartości 5. Poprawiłem je, dziękuję. – petrichor

+0

Wartości są obliczane na siatce, w której punktem środkowym jest początek, czyli h1 == 0 i h2 == 0 w naszym przypadku. Wszystkie pozostałe pary reprezentują inne współrzędne, gdy spojrzymy na wartości h1, h2 element po elemencie. Podczas filtrowania możesz myśleć, że ta siatka zostanie umieszczona na pikselu obrazu, w którym początek siatki pasuje dokładnie do piksela. Możesz przeczytać odpowiedź Goza w linku, który podałeś w swoim pytaniu dla szczegółów. – petrichor

22

Jest to tak proste, jak to brzmi:

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

To jest wadliwe: musisz również znormalizować jądro, lub obraz staje się ciemniejszy w zależności od W i sigma. po prostu: pobierz sumę wartości jądra i podziel każdą wartość jądra przez tę sumę. – Rookie

+2

@Rookie - Postanowiłem zmodyfikować ten post i dodać normalizację. Ma to umożliwić tym, którzy chcą rozwiązania C/C++, aby użyć tego wprost. Dobry połów! – rayryeng

+0

Wydaje się nieprawidłowe, gdy m, n są parzyste, w porównaniu do wyniku 'fspecial'. –

13

do wdrożenia gaussian blur wystarczy wziąć gaussian function i obliczyć jedną wartość dla każdego z elementów w jądrze.

Zazwyczaj do głównego elementu jądra przypisuje się maksymalną masę, a wartości bliskie zeru dla elementów na krawędziach jądra. Oznacza to, że jądro powinno mieć nieparzystą wysokość (lub szerokość), aby upewnić się, że rzeczywiście jest element centralny.

Aby obliczyć rzeczywiste elementy jądra, można przeskalować dzwon Gaussa do siatki jądra (wybierz dowolny, np. sigma = 1 i dowolny zakres, np. -2*sigma ... 2*sigma) i znormalizuj go, s.t. elementy sumują się do jednego. Aby to osiągnąć, jeśli chcesz wspierać dowolne rozmiary jądra, możesz dostosować sigmę do wymaganego rozmiaru jądra.

Oto przykład 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'; 
    } 
} 

Wyjście jest:

$ 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 

Jako uproszczenia nie trzeba używać 2D-jądro. Łatwiejsze do wdrożenia, a także bardziej wydajne do obliczenia jest użycie dwóch ortogonalnych rdzeni 1d. Jest to możliwe dzięki asocjatywności tego typu liniowego splotu (separacja liniowa). Możesz także chcieć zobaczyć this section odpowiedniego artykułu w Wikipedii.


Oto samo w Pythonie (w nadziei, że ktoś może go znaleźć przydatne):

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] 

produkuje jądro:

['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

Gdzie jest PI? – user2023370

+0

@ user2023370 Dlaczego uważasz, że powinien tam być PI? – moooeeeep

0

rozmycie Gaussa w Pythonie przy użyciu biblioteki obrazów PIL. Aby uzyskać więcej informacji przeczytaj to: 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

Dodaj komentarze do kodu, będzie to pomocne dla innych osób. – HDJEMAI