Algorytmy przetwarzania sygnałów zdefiniowane w ciągłym czasie / przestrzeni / częstotliwości są zwykle implementowane przez próbkowanie sygnału na dyskretnej siatce i przekształcanie całek w sumy (i pochodne w różnice). Filtry przestrzenne są implementowane poprzez splot z jądrem splotu (tj. Ważoną sumą sąsiadów).
Istnieje ogromna wiedza na temat filtrowania próbkowanych sygnałów w dziedzinie czasu; filtry w dziedzinie czasu są implementowane jako filtry o skończonej odpowiedzi impulsowej , w których bieżąca próbka wyjściowa jest obliczana jako ważona suma poprzednich N próbek wejściowych; lub nieskończone filtry odpowiedzi impulsowej, w których prąd wyjściowy jest ważoną sumą poprzednich sygnałów wejściowych i poprzednich sygnałów wyjściowych . Formalnie, filtry czasu dyskretnego są opisywane przy użyciu transformaty z , która jest analogowym czasem dyskretnym do transformaty Laplace'a . Metoda tustina odwzorowuje jeden do drugiego ( c2d
i d2c
w Matlab).
Jak oceniasz funkcje w dowolnych punktach?
Gdy potrzebujesz wartości sygnału w punkcie, który nie leży bezpośrednio na siatce próbkowania, interpolujesz jego wartość z pobliskich punktów. Interpolacja może być tak prosta, jak wybranie najbliższej próbki, obliczenie średniej ważonej najbliższych próbek lub dopasowanie arbitralnie skomplikowanej funkcji analitycznej do próbkowanych danych i ocena tej funkcji przy wymaganych współrzędnych. Interpolacja na jednolitą drobniejszą siatkę to próbkowanie w górę . Jeśli twój oryginalny (ciągły) sygnał nie zawiera szczegółów (tj. Częstotliwości) drobniejszych niż połowa siatki próbkowania, to funkcję ciągłą można doskonale zrekonstruować na podstawie próbkowanej wersji ( twierdzenie Nyquista-Shannona o próbkowaniu ). Aby zobaczyć przykład interpolacji w 2D, zobaczinterpolacja dwuliniowa .
W Matlab możesz używać interp1
lub interp2
interpolować dane 2D lub regularnie próbkowane dane 2D (odpowiednio) lub griddata
interpolować dane 2D nieregularnie próbkowane.
Czy miałbyś pętlę for przechodzącą przez każdy woksel i obliczającą odpowiednią formułę?
Tak, dokładnie.
Matlab ratuje Cię przed koniecznością robienia tego za pomocą jawnych pętli for, ponieważ jest zaprojektowany do pracy na macierzach i wektorach (tj. Tablicach wielowymiarowych). W Matlabie nazywa się to „wektoryzacją”. Całka mogą być przybliżone sum
, cumsum
, trapz
, cumtrapz
, itd.
Przeczytałem książkę „Digital Image Processing” Gonzalez i Woods, ale wciąż jestem zagubiony. Czytałem również o serii książek z numerycznymi przepisami. Czy to byłby właściwy sposób?
Tak, przepisy numeryczne byłyby świetnym początkiem. Jest to bardzo praktyczne i obejmuje większość metod numerycznych, których będziesz potrzebować. (Przekonasz się, że Matlab już implementuje wszystko, czego potrzebujesz, ale Przepisy numeryczne zapewnią doskonałe tło).
Wziąłem klasę „algorytmy i struktury danych”, ale nie widzę związku między prezentowanym tam materiałem a wdrażaniem algorytmów naukowych.
Materiał omawiany na kursach „Algorytmy i struktury danych” ma tendencję do koncentrowania się na strukturach takich jak listy, tablice, drzewa i wykresy zawierające liczby całkowite lub łańcuchy oraz operacje, takie jak sortowanie i wybieranie: problemy, dla których zwykle występuje jeden poprawny wynik. Jeśli chodzi o algorytmy naukowe, to tylko połowa historii. Druga połowa dotyczy metod szacowania liczb rzeczywistych i funkcji analitycznych. Znajdziesz to na kursie „Metody numeryczne” (lub „Analiza numeryczna”; taki jak ten- przewiń w dół po slajdy): jak oszacować funkcje specjalne, jak oszacować całki i pochodne itp. Tutaj jednym z głównych zadań jest oszacowanie dokładności wyniku, a jednym wspólnym wzorem jest iteracja rutyny, która poprawia oszacuj, aż będzie wystarczająco dokładny. (Możesz zadać sobie pytanie, skąd Matlab wie, jak zrobić coś tak prostego, jak oszacowanie wartości sin(x)
dla niektórych x
).
Jako prosty przykład, oto krótki skrypt, który oblicza transformatę radonową obrazu w Matlabie. Transformacja radonu przenosi projekcje obrazu na zbiór kątów projekcji. Zamiast próbować obliczyć projekcję pod dowolnym kątem, zamiast tego obracam cały obraz za pomocą imrotate
, aby projekcja była zawsze pionowa. Następnie możemy po prostu użyć rzutu sum
, ponieważ sum
macierz zwraca wektor zawierający sumę z każdej kolumny.
Możesz napisać własny, imrotate
jeśli wolisz, używając interp2
.
%%# Home-made Radon Tranform
%# load a density map (image).
A = phantom;
n_pixels = size(A, 1); %# image width (assume square)
%# At what rotation angles do we want to take projections?
n_thetas = 101;
thetas = linspace(0, 180, n_thetas);
result = zeros(n_thetas, n_pixels);
%# Loop over angles
for ii=1:length(thetas)
theta = thetas(ii);
rotated_image = imrotate(A, theta, 'crop');
result(ii, :) = sum(rotated_image);
end
%# display the result
imagesc(thetas, 1:n_pixels, result.');
xlabel('projection angle [degrees]');
To, co kiedyś było całką gęstości wzdłuż promienia, jest teraz sumą nad kolumną dyskretnie próbkowanego obrazu, który z kolei został znaleziony przez interpolację oryginalnego obrazu nad przekształconym układem współrzędnych.