[X] Grundlegende Bildverarbeitungsalgorithmen
-
Korbinian schrieb:
ich warte mit dem ausbessern mal, bis noch mehr sachen zusammenkommen.
Was meinst du? Kommentare und technische Korrekturen oder noch eine Rechtschreibprüfung?
Mr. B
-
alles. möchte jetzt nur nicht mehr wegen jedem einzlenen fehler nachfummeln, sondern mach die paar sachen dann auf einmal. dass das dings dann mal fertig ist.
btw: gefällt dir die post-version oder das pdf besser?
-
Korbinian schrieb:
btw: gefällt dir die post-version oder das pdf besser?
die PDF version. allerdings wärs noch mal schicker, wenn du den code farbig (so wie in diesem forum hier) gestalten könntest.
Mr. B
-
lad dir nochmal die neue version des PDFs runter. ich hab ein farbiges syncol modul für tex gefunden
-
template<type T> void medianFilter(T *din, T *dout, unsigned sx, unsigned sy, unsigned ks)
Sollte wohl
template<typename T> void medianFilter(T *din, T *dout, unsigned sx, unsigned sy, unsigned ks)
oder
template<class T> void medianFilter(T *din, T *dout, unsigned sx, unsigned sy, unsigned ks)
sein
-
Korbinian schrieb:
lad dir nochmal die neue version des PDFs runter. ich hab ein farbiges syncol modul für tex gefunden
Ja, gut. Ich fänd das in den gewohnten Farben zwar besser, aber es ist ein Fortschritt!
Mr. B
-
Mr. B schrieb:
...in den gewohnten Farben zwar besser...
done, mann bin ich heut fleissig
ps: und für die forumleute die nicht dauernd quote wollen: www.korbinian-riedhammer.de/misc/BVAlgos_phpbb.txt
-
Hi
ich find die Idee ziemlich gut die Grundlagen zu erklaeren, was mir aber an Erklaerung noch fehlt:
- Warum werden die Tiefpassfilter (Mittelwert, Gauss) normiert? Das sollte wenigstens kurz erwaehnt werden.
- Vielleicht auch einen Hochpassfilter zeigen
- Randbetrachtung: Du schneidest den Rand einfach ab (was bei groesseren Filtern schnell zu nem deutlichen Bildverlist fuehrt), statt dessen gibt es ja auch die Moeglichkeit den Rand zu schmieren, also alle Werte < 0 auf 0 zu mappen und bei groesser breite/hoehe auf breite/hoehe (ich hoffe mal das kann man verstehen). Das liefert meiner Erfahrung nach gute Ergebnisse. Wenn du, weil einfacher zu implementieren, deine Loesung behalten willst dann weise vielleicht darauf hin, dass der Rand verloren geht.
- Medianfilter: Fuer binaeres Rauschen nimmt man moeglichst kleine, sprich 3x3, Filter. 7 oder mehr ist schon ein ziemlich starker Filter der das Bild u.U. sehr stark verwischt. Ausserdem waechst die Laufzeit ziemlich stark an je groesser der Filter wird.
- Zusaetzlicher Hinweis beim Histogrammausgleich, z.B. "funktioniert nur auf Grauwertbildern oder auf Farbbildern die einen Helligkeitskanal haben (YCrCb, HLS, ..) nicht auf Farbkanalbildern (RGB, CYMK)"
Willst du noch mehr erklaeren? Ich bin immer der Vertreter auch etwas von der Theorie dahinter reinzubringen und nicht nur die Anwendung zu erklaeren, aber ich weiss ja nicht wie du dir das Ziel gesetzt hast. Ich faende ein paar Erlaeuteringen zur Funktionsweise von Tief-/Hochpassfiltern und Gradientenberechnung gut, im Prinzip ja nicht viel anders als das was du zur Fouriertransformation geschrieben hast.
Gruss
pli
-
pli schrieb:
Vielleicht auch einen Hochpassfilter zeigen
ne das führt mir zu weit, wüsste auch nicht, wo der eingestzt wird.
Randbetrachtung...
jau ich bau noch ein sätzle ein
Medianfilter: Fuer binaeres Rauschen nimmt man moeglichst kleine, sprich 3x3, Filter. 7 oder mehr ist schon ein ziemlich starker Filter der das Bild u.U. sehr stark verwischt. Ausserdem waechst die Laufzeit ziemlich stark an je groesser der Filter wird.
das kommt ganz auf das einsatz gebiet drauf an... bei uns in der firma wurde das für gusseisenteilprüfung verwendet, da hatten die auf 512er bilder sogar ne 13er kernelgröße. es kommt ganz drauf an, was man will. der vorteil vom medianfilter ist ja gerade, dass er nicht blind einen arithmetisches mittel macht, sondern ausreissern keine chance gibt.
Zusaetzlicher Hinweis beim Histogrammausgleich, z.B. "funktioniert nur auf Grauwertbildern oder auf Farbbildern die einen Helligkeitskanal haben (YCrCb, HLS, ..) nicht auf Farbkanalbildern (RGB, CYMK)"
gute idee
Willst du noch mehr erklaeren? ...
vielleicht lass ich mich noch etwas mehr über die kantenerkennung aus. auf dem forumtreffen hat sich aber gezeigt, dass die leute, die ahnung von mathe haben (stichwort differenzenquotient) mit den paar worten klar kommen, aber leute die davon keinen schimmer haben, auch mit 1 seite erklärung kaum durchsteigen
ich setz mich heut abend nochmal hin, und besser das eine oder andere aus. habe jetzt auch ein testprojekt angelegt, wo die ganzen code-teile drin sind und noch ein paar andere sachen (musste noch ein paar sachen ausbessern), werde auch das noch mit dazu online stellen
-
so, pdf nochmal upgedated, und ich hab ein beispielprojekt hochgeladen, wo unter anderem die im artikel vorgestellten funktionen enthalten sind: www.korbinian-riedhammer.de/misc/ip-suite.tar.bz2
-
Korbinian schrieb:
jau ich bau noch ein sätzle ein
also ich wäre ja schon dafür, dass die rechtschreibprüfung wirklich am ende erfolgt
sonst ist alles für die katz. jetzt wo du noch ein paar sätzchen hinzufügen wolltest.Mr. B
-
macht ihr jetzt was ihr wollt mit dem artikel, ich hab fertig damit...
-
> Überarbeitete, eigentlich endgültige Version!
Für PDF-Fans: www.korbinian-riedhammer.de/misc/BVAlgos.pdf
Grundlegende Algorithmen in der Bildverarbeitung
Inhalt:
1 Einleitung 2 Rotation von Bildern 2.1 [kor]Octave-Implementierung[/kor] 2.2 [kor]C++-Implementierung[/kor] 3 Preprocessing von Bildern 3.1 Einfache Filter 3.2 Implentierung einfacher Filter in Octave 3.3 Implementierung einfacher Filter in c++ 3.4 Effektive Implementierung eines [kor]3x3-Gaußfilters[/kor] 4 Erweiterte Filter 4.1 Thresholding 4.2 Medianfilter 4.3 Histogram Equilisation 5 Defektpixelinterpolation 5.1 [kor]Octave-Implementierung[/kor]
Einleitung
In diesem Artikel möchte ich einige grundlegende Verfahren vorstellen, die in der Bildverarbeitung regen Einsatz finden. Ich werde jedes Verfahren zunächst beschreiben, dann anhand eines Octave Script veranschaulichen und schließlich (fragmentweise) in C++ implementieren. Diese Implementierungen sind i.d.R. nicht die schnellsten, da ich auf Anschaulichkeit Wert gelegt habe.
Warum benutze ich Octave, um die Algorithmen zu veranschaulichen, obwohl der Leser vermutlich nur an der [quote]C++-Variante[/kor] interessiert ist? Zum einen, weil Octave-Programmtext deutlich übersichtlicher und damit leichter nachvollziehbar ist. Meines Erachtens wichtiger ist aber, dass es gerade in der Bildverarbeitung praktisch ist, mit diesem Tool vertraut zu sein, da man oft ein Verfahren erst einmal ausprobieren möchte, bevor man es implementiert. Würde man sofort anfangen, die Idee in C++ zu implementieren, so hätte man erst einmal einen größeren Zeitaufwand, dann das größere Risiko, einen Fehler darin zu haben, und letztendlich verfängt man sich leicht in Startschwierigkeiten wie "Wie lese ich ein Bild ***ein?"***. All das wird durch die abstraktere Implementierung in Octave verhindert.
Neben den trivialen Implementierungen werde ich auch eine effektive am Beispiel eines Filters zeigen, die nicht nur wegen der Programmiertricks so gut ist, sondern auch, weil mathematische Umstände zu unserem Vorteil genutzt werden. Ich zeige das, um den Leser auch für die Theorie hinter den Verfahren zu interessieren. Oft entstehen neue schnelle Lösungen, weil die Probleme von einer anderen Seite betrachtet werden.
Alle hier vorgestellten Algorithmen sowie einige mehr finden sich in einem von mir zusammengestellten Beispielprojekt, somit entfällt das lästige abtippen
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Rotation von Bildern
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Oftmals sind die Bilder nicht in der Position, in der man sie gerne hätte. Das liegt entweder daran, dass die Kamera oder der Detektor nicht korrekt angebracht ist, oder ganz einfach daran, dass das Objekt nicht in der richtigen Lage ist, wie z.B. bei schiefen Fotos oder bei in der Schrifterkennung schiefen Zeilen. Rotationen im zweidimensionalen Raum lassen sich mathematisch recht einfach darstellen: Sei p ein zu rotierender und q der rotierte Punkt, so besteht folgender Zusammenhang mit der Rotationsmatrix R resultierend aus dem Rotationswinkel omega:
Möchte man nun eine Rotation eines ganzen Bildes durchführen, so berechnet man für jeden Punkt aus dem Urbild die Koordinaten im rotierten Bild. Bei der Implementierung ist auf ungültige Koordinaten zu achten: Das kann zum einen ein Nichtbeachten von nicht erreichbaren Koordinaten sein, zum anderen eine Expansion des Bildes (bei Winkeln ungleich 90, 180 und 270 Grad). Der Einfachheit halber werde ich hier ersteres zeigen. Punkte im Bild, denen kein Punkt aus dem Urbild zugewiesen wird, bleiben z.B. schwarz. Nimmt man nun die "Vorwärtsabbildung", so kann es sein, dass manche Pixel aufgrund von Rundungsfehlern nicht rotiert werden. Dem kann man entgegenkommen, indem man nicht vom Urbild zum Bild rechnet, sondern umgekehrt, und optional beim Zugriff interpoliert (z.B. bilinear). Damit wird jedem Pixel im Bild ein (interpoliertes) Pixel aus dem Urbild zugeordnet, was umgekehrt aufgrund der Diskretisierung nicht der Fall ist.Octave-Implemtierung
% Rotieren eines [kor]weißen[/kor] Vierecks um 30 Grad um % den Ursprung dim = 30; w = -(pi/6); % Bild und Urbild reservieren Bild = zeros(dim); Urbild = zeros(dim); % [kor]Weißen[/kor] Rand und Viereck zeichnen Urbild(10:20, 10:20) = 255; Urbild(:,2) = 255; Urbild(:,29) = 255; Urbild(2,:) = 255; Urbild(29,:) = 255; for i=1:dim for j=1:dim x = round(i*cos(w)-j*sin(w)); y = round(i*sin(w)+j*cos(w)); if (x > 0 && x <= dim && y > 0 && y <= dim) Bild(i,j) = Urbild(x,y); endif endfor endfor
Beispielbilder
Urbild
Bild, 30 Grad rotiertC++-Implementierung
template <typename T> void rotate(T* din, T* dout, unsigned sx, unsigned sy, double w) { for (unsigned i = 0; i < sy; ++i) { for (unsigned j = 0; j < sx; ++i) { signed rx = (signed)((double)j*cos(w) - (double)i*sin(w)); signed ry = (signed)((double)j*sin(w) + (double)i*cos(w)); if (rx < 0 || ry < 0 || rx >= sx || ry >= sy) continue; dout[i*sx + j] = din[ry*sx + rx]; } } return; }
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Preprocessing von Bildern
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Normalerweise werden Bilder vor ihrer Verarbeitung erst vorbereitet, um bessere Eingangsbedingungen für die eigentlichen Algorithmen zu haben. Unter einem Filter verstehe ich eine Operation, die ein Bild verändert. Ich möchte hier zuerst einfache Filter vorstellen, die auf einer Faltung zweier Bilder beruhen, und dann erweiterte Filter, die teils aus mathematischen Ideen resultieren und nichts mit Faltung zu tun haben.
Einfache Filter
Bei einfachen Filtern wird das Urbild mit einem Filterkernel gefaltet. Mathematisch ist die Faltung zweier Bilder im Diskreten für die Bildfunktionen f und g so beschrieben:
Da allerdings der Filterkernel normaleweise eine Dimension von 3x3 oder 5x5 hat, keinenfalls aber die Größe des Urbilds, wird das Urbild abschnittweise mit dem Kernel gefaltet. Praktisch gesehen ist das eine Multiplikation der jeweils überlappenden Matrixelemente (was sich im C++-Code gut sehen lässt). Kernelgrößen sind im Allgemeinen ungerade, damit beim Programmieren immer die Mitte des Kernelfensters eindeutig ist.
Bei der hier gezeigten c++ Implementierung wird der Rand (also die Pixel, für die die Filtermaske nicht vollständig gefüllt ist) nicht gefiltert. Ein besseres Verfahren ist es z.B. die fehlenden Werte mit 0 anzunehmen. Die Implementierung ist desweiteren nicht gerade die effektivste: Für jedes Pixel wird das Kernelfenster neu aufgestellt. Ist die Kernelgröße im Vorhinein bekannt, so kann man effektivere Implementierungen schreiben, indem man sich Pointerarithmetik zu Hilfe nimmt. Als Beispiel eines effektiv implementierten Filters zeige ich eine elegante Version des 3x3-Gaußfilters, der sich neben der Pointerarithmetik noch folgende Eigenschaft zunutze macht: Die zweifache Anwendung des 2x2-Mittelwertfilters entspricht der einfachen Anwendung eines 3x3-Gaußfilters.
Verknüpfung von Kernel- und Bildelementen beim [kor]3x3-Mittelwertfilter[kor]Beispiele für einfache Filter
- Mittelwert: Der Mittelwertfilter wird in verschiedenen Kernelgrößen verwendet und bildet jeweils den Mittelwert über die Umgebung. Kernelgrößen zwischen [kor]zwei[/kor] und fünf sind gebräuchlich.
- Gaußfilter: Der Gaußfilter wird in Kernelgrößen ab drei verwendet, die Elemente sind entsprechend der Gaußverteilung erstellt. Er erzielt eine bessere Glättung des Bildes als ein Mittelwertfilter, da die Pixel von außen nach innen immer mehr an Einfluss erhalten.
- Kantenerkennung, finite Differenz: Der einfache Kantenerkennungsfilter beruht auf der Idee der finiten Differenzen. Betrachtet man das Bild als zweidimensionale diskrete Funktion, so kann man mit Hilfe des Differenzenquotienten die Steigung für die x- und y-Richtung an jedem Punkt berechnen. Auf diese Weise erhält man kleine Werte, wenn sich der Farb-/Grauwert kaum ändert, und große Werte, wenn von einem Pixel zum nächsten ein starker An-/Abstieg ist. Die beiden gezeigten Kernel sind die sog. [kor]*forward difference-*Kernel[/kor] für die x- und y-Richtung.
Vereinfacht gesagt beruht die Kantenerkennung durch finite Differenzen auf den Unterschieden von einem zum nächsten Pixel in jede Bildrichtung: Ändert sich der Farbwert stark, so ist die Differenz groß und eine Kante liegt vor. Ändert er sich aber kaum, so wird ein homogenes Gebiet vorliegen. - Sobelfilter: Der Sobelfilter funktioniert ähnlich wie der einfache Kantenerkennungsfilter, nur dass er die Nachbarschaft mit in Betracht zieht:
Implementierung einfacher Filter in Octave
% Einfache Filter % [kor]Weißes[/kor] Viereck auf schwarzem Untergrund img = zeros(30); img(10:20, 10:20) = 1; M = (1/4)*[1 1; 1 1]; G = (1/16)*[1 2 1; 2 4 2; 1 2 1 ]; Dx = [-1 1]; S = (1/8)*[-1 0 1; -2 0 2; -1 0 1 ]; img_m = conv2(img, M, 'same'); img_g = conv2(img, G, 'same'); img_x = conv2(img, Dx, 'same'); img_s = conv2(img, S, 'same');
Originalbild
Mittelwertfilter
Gaußfilter
Kanten in x-Richtung
Sobel in x-RichtungImplementierung einfacher Filter in C++
void filter(double *din, double *dout, unsigned sx, unsigned sy, double *kernel, unsigned kdim) { for (unsigned y = (kdim/2); y < sy - (kdim/2); ++y) for (unsigned x = (kdim/2); x < sx - (kdim/2); ++x) { dout[y*sx+x] = 0.; for (unsigned i = 0; i < kdim; ++i) for (unsigned j = 0; j < kdim; ++j) dout[y*sx+x] += (kernel[i*kdim+j] * din[(y+j-kdim/2)*sx + x+i-kdim/2]; } }
-----------------------------------------------------------------
Den Rest gibts morgen. Ich kann leider nur zeitlich begrenzt ins Internet.
- Mittelwert: Der Mittelwertfilter wird in verschiedenen Kernelgrößen verwendet und bildet jeweils den Mittelwert über die Umgebung. Kernelgrößen zwischen [kor]zwei[/kor] und fünf sind gebräuchlich.
-
sei mir nicht böse, aber:
- ich hab alte rechtschreibung verwendet, daher z.b. auch "Effektive" groß, weil es anstelle des Substantives steht
- spielereien ob jetzt C++ oder c++ sind mir wayne gretzky, ich hab mich für die kleine variante entschieden, und die immer konsequent genommen. genauso wie X-Richtung statt x-Richtung.
- ich hab keine sonderzeichen in den octave files, weil ich auf dieser tastatur keine hatte. so einfach...
- über die anderen formulierungen kann man streiten, ich würd sie so lassen wie sie sind, aber wie vorher gepostet, macht was ihr wollt damit ich hab fertig
-
Kein Problem. Ich kann nur neue Rechtschreibung und passe den Text demnach auch entsprechend an. Ist ja jetzt seit sechs Tagen verbindlich :p X-Richtung hast du AFAIR nicht ganz konsequent verwendet. Ich schreib C++ statt c++, weil mans einfach häufiger antrifft, ist kein Rechtschreibfehler in dem Sinn.
So, muss weiterarbeiten, die Spachtelmasse wird hart.
-
Michael E. schrieb:
2.1 Octave-Implementierung
2.2 C++-Implementierung[...]
Michael E. schrieb:
Octave-Implemtierung
Heißt das jetzt "Implementierung" oder "Implemtierung"? Ich dachte, weil er das am Anfang "Implemtierung" geschrieben hatte, dass das auch richtig sei. Dachte, zwischen "Implemtierung" und "Implementierung" gäbe es einen Unterschied, wobei "Implemtierung" auch ein Fachbegriff ist, deshalb habe ich es auch nich korrigiert.
Mr. B
-
Es heißt Implementierung. Das Wort ist markiert, weil ich einen Bindestrich eingefügt habe. Hab dann wohl den Fehler übersehen.
-
eiei... die Implementierungen hab ich jetzt ausgebessert, war ja furchtbar
-
> Überarbeitete, eigentlich endgültige Version!
Für PDF-Fans: www.korbinian-riedhammer.de/misc/BVAlgos.pdf
Grundlegende Algorithmen in der Bildverarbeitung
Inhalt:
1 Einleitung 2 Rotation von Bildern 2.1 [kor]Octave-Implementierung[/kor] 2.2 [kor]C++-Implementierung[/kor] 3 Preprocessing von Bildern 3.1 Einfache Filter 3.2 [kor]Implementierung[/kor] einfacher Filter in Octave 3.3 Implementierung einfacher Filter in c++ 3.4 Effektive Implementierung eines [kor]3x3-Gaußfilters[/kor] 4 Erweiterte Filter 4.1 Thresholding 4.2 Medianfilter 4.3 Histogram [kor]Equalization[/kor] 5 Defektpixelinterpolation 5.1 [kor]Octave-Implementierung[/kor]
Einleitung
In diesem Artikel möchte ich einige grundlegende Verfahren vorstellen, die in der Bildverarbeitung regen Einsatz finden. Ich werde jedes Verfahren zunächst beschreiben, dann anhand eines Octave Script veranschaulichen und schließlich (fragmentweise) in C++ implementieren. Diese Implementierungen sind i.d.R. nicht die schnellsten, da ich auf Anschaulichkeit Wert gelegt habe.
Warum benutze ich Octave, um die Algorithmen zu veranschaulichen, obwohl der Leser vermutlich nur an der [quote]C++-Variante[/kor] interessiert ist? Zum einen, weil Octave-Programmtext deutlich übersichtlicher und damit leichter nachvollziehbar ist. Meines Erachtens wichtiger ist aber, dass es gerade in der Bildverarbeitung praktisch ist, mit diesem Tool vertraut zu sein, da man oft ein Verfahren erst einmal ausprobieren möchte, bevor man es implementiert. Würde man sofort anfangen, die Idee in C++ zu implementieren, so hätte man erst einmal einen größeren Zeitaufwand, dann das größere Risiko, einen Fehler darin zu haben, und letztendlich verfängt man sich leicht in Startschwierigkeiten wie "Wie lese ich ein Bild ***ein?"***. All das wird durch die abstraktere Implementierung in Octave verhindert.
Neben den trivialen Implementierungen werde ich auch eine effektive am Beispiel eines Filters zeigen, die nicht nur wegen der Programmiertricks so gut ist, sondern auch, weil mathematische Umstände zu unserem Vorteil genutzt werden. Ich zeige das, um den Leser auch für die Theorie hinter den Verfahren zu interessieren. Oft entstehen neue schnelle Lösungen, weil die Probleme von einer anderen Seite betrachtet werden.
Alle hier vorgestellten Algorithmen sowie einige mehr finden sich in einem von mir zusammengestellten Beispielprojekt, somit entfällt das lästige Abtippen
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Rotation von Bildern
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Oftmals sind die Bilder nicht in der Position, in der man sie gerne hätte. Das liegt entweder daran, dass die Kamera oder der Detektor nicht korrekt angebracht ist, oder ganz einfach daran, dass das Objekt nicht in der richtigen Lage ist, wie z.B. bei schiefen Fotos oder bei in der Schrifterkennung schiefen Zeilen. Rotationen im zweidimensionalen Raum lassen sich mathematisch recht einfach darstellen: Sei p ein zu rotierender und q der rotierte Punkt, so besteht folgender Zusammenhang mit der Rotationsmatrix R resultierend aus dem Rotationswinkel omega:
Möchte man nun eine Rotation eines ganzen Bildes durchführen, so berechnet man für jeden Punkt aus dem Urbild die Koordinaten im rotierten Bild. Bei der Implementierung ist auf ungültige Koordinaten zu achten: Das kann zum einen ein Nichtbeachten von nicht erreichbaren Koordinaten sein, zum anderen eine Expansion des Bildes (bei Winkeln ungleich 90, 180 und 270 Grad). Der Einfachheit halber werde ich hier ersteres zeigen. Punkte im Bild, denen kein Punkt aus dem Urbild zugewiesen wird, bleiben z.B. schwarz. Nimmt man nun die "Vorwärtsabbildung", so kann es sein, dass manche Pixel aufgrund von Rundungsfehlern nicht rotiert werden. Dem kann man entgegenkommen, indem man nicht vom Urbild zum Bild rechnet, sondern umgekehrt, und optional beim Zugriff interpoliert (z.B. bilinear). Damit wird jedem Pixel im Bild ein (interpoliertes) Pixel aus dem Urbild zugeordnet, was umgekehrt aufgrund der Diskretisierung nicht der Fall ist.Octave-Implementierung
% Rotieren eines [kor]weißen[/kor] Vierecks um 30 Grad um % den Ursprung dim = 30; w = -(pi/6); % Bild und Urbild reservieren Bild = zeros(dim); Urbild = zeros(dim); % [kor]Weißen[/kor] Rand und Viereck zeichnen Urbild(10:20, 10:20) = 255; Urbild(:,2) = 255; Urbild(:,29) = 255; Urbild(2,:) = 255; Urbild(29,:) = 255; for i=1:dim for j=1:dim x = round(i*cos(w)-j*sin(w)); y = round(i*sin(w)+j*cos(w)); if (x > 0 && x <= dim && y > 0 && y <= dim) Bild(i,j) = Urbild(x,y); endif endfor endfor
Beispielbilder
Urbild
Bild, 30 Grad rotiertC++-Implementierung
template <typename T> void rotate(T* din, T* dout, unsigned sx, unsigned sy, double w) { for (unsigned i = 0; i < sy; ++i) { for (unsigned j = 0; j < sx; ++i) { signed rx = (signed)((double)j*cos(w) - (double)i*sin(w)); signed ry = (signed)((double)j*sin(w) + (double)i*cos(w)); if (rx < 0 || ry < 0 || rx >= sx || ry >= sy) continue; dout[i*sx + j] = din[ry*sx + rx]; } } return; }
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Preprocessing von Bildern
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Normalerweise werden Bilder vor ihrer Verarbeitung erst vorbereitet, um bessere Eingangsbedingungen für die eigentlichen Algorithmen zu haben. Unter einem Filter verstehe ich eine Operation, die ein Bild verändert. Ich möchte hier zuerst einfache Filter vorstellen, die auf einer Faltung zweier Bilder beruhen, und dann erweiterte Filter, die teils aus mathematischen Ideen resultieren und nichts mit Faltung zu tun haben.
Einfache Filter
Bei einfachen Filtern wird das Urbild mit einem Filterkernel gefaltet. Mathematisch ist die Faltung zweier Bilder im Diskreten für die Bildfunktionen f und g so beschrieben:
Da allerdings der Filterkernel normaleweise eine Dimension von 3x3 oder 5x5 hat, keinenfalls aber die Größe des Urbilds, wird das Urbild abschnittweise mit dem Kernel gefaltet. Praktisch gesehen ist das eine Multiplikation der jeweils überlappenden Matrixelemente (was sich im C++-Code gut sehen lässt). Kernelgrößen sind im Allgemeinen ungerade, damit beim Programmieren immer die Mitte des Kernelfensters eindeutig ist.
Bei der hier gezeigten c++ Implementierung wird der Rand (also die Pixel, für die die Filtermaske nicht vollständig gefüllt ist) nicht gefiltert. Ein besseres Verfahren ist es z.B. die fehlenden Werte mit 0 anzunehmen. Die Implementierung ist des Weiteren nicht gerade die effektivste: Für jedes Pixel wird das Kernelfenster neu aufgestellt. Ist die Kernelgröße im Vorhinein bekannt, so kann man effektivere Implementierungen schreiben, indem man sich Pointerarithmetik zu Hilfe nimmt. Als Beispiel eines effektiv implementierten Filters zeige ich eine elegante Version des 3x3-Gaußfilters, der sich neben der Pointerarithmetik noch folgende Eigenschaft zunutze macht: Die zweifache Anwendung des 2x2-Mittelwertfilters entspricht der einfachen Anwendung eines 3x3-Gaußfilters.
Verknüpfung von Kernel- und Bildelementen beim [kor]3x3-Mittelwertfilter[kor]Beispiele für einfache Filter
- Mittelwert: Der Mittelwertfilter wird in verschiedenen Kernelgrößen verwendet und bildet jeweils den Mittelwert über die Umgebung. Kernelgrößen zwischen [kor]zwei[/kor] und fünf sind gebräuchlich.
- Gaußfilter: Der Gaußfilter wird in Kernelgrößen ab drei verwendet, die Elemente sind entsprechend der Gaußverteilung erstellt. Er erzielt eine bessere Glättung des Bildes als ein Mittelwertfilter, da die Pixel von außen nach innen immer mehr an Einfluss erhalten.
- Kantenerkennung, finite Differenz: Der einfache Kantenerkennungsfilter beruht auf der Idee der finiten Differenzen. Betrachtet man das Bild als zweidimensionale diskrete Funktion, so kann man mit Hilfe des Differenzenquotienten die Steigung für die x- und y-Richtung an jedem Punkt berechnen. Auf diese Weise erhält man kleine Werte, wenn sich der Farb-/Grauwert kaum ändert, und große Werte, wenn von einem Pixel zum nächsten ein starker An-/Abstieg ist. Die beiden gezeigten Kernel sind die sog. [kor]*forward difference-*Kernel[/kor] für die x- und y-Richtung.
Vereinfacht gesagt beruht die Kantenerkennung durch finite Differenzen auf den Unterschieden von einem zum nächsten Pixel in jede Bildrichtung: Ändert sich der Farbwert stark, so ist die Differenz groß und eine Kante liegt vor. Ändert er sich aber kaum, so wird ein homogenes Gebiet vorliegen. - Sobelfilter: Der Sobelfilter funktioniert ähnlich wie der einfache Kantenerkennungsfilter, nur dass er die Nachbarschaft mit in Betracht zieht:
Implementierung einfacher Filter in Octave
% Einfache Filter % [kor]Weißes[/kor] Viereck auf schwarzem Untergrund img = zeros(30); img(10:20, 10:20) = 1; M = (1/4)*[1 1; 1 1]; G = (1/16)*[1 2 1; 2 4 2; 1 2 1 ]; Dx = [-1 1]; S = (1/8)*[-1 0 1; -2 0 2; -1 0 1 ]; img_m = conv2(img, M, 'same'); img_g = conv2(img, G, 'same'); img_x = conv2(img, Dx, 'same'); img_s = conv2(img, S, 'same');
Originalbild
Mittelwertfilter
Gaußfilter
Kanten in x-Richtung
Sobel in x-RichtungImplementierung einfacher Filter in C++
void filter(double *din, double *dout, unsigned sx, unsigned sy, double *kernel, unsigned kdim) { for (unsigned y = (kdim/2); y < sy - (kdim/2); ++y) for (unsigned x = (kdim/2); x < sx - (kdim/2); ++x) { dout[y*sx+x] = 0.; for (unsigned i = 0; i < kdim; ++i) for (unsigned j = 0; j < kdim; ++j) dout[y*sx+x] += (kernel[i*kdim+j] * din[(y+j-kdim/2)*sx + x+i-kdim/2]; } }
Effektive Implementierung eines 3x3-Gaußfilters
Wie eingangs erwähnt, beruht diese Implementierung auf der Äquivalenz der zweifachen 2x2-Mittelwertfilterung zur 3x3-Gaußfilterung (wer's nicht glaubt, nehme ein Blatt Papier und rechne es nach ***:-) )***. Daraus ergeben sich zwei Vorteile: Zum einen ist das Fenster auf 2x2 geschrumpft, zum anderen sind die Gewichte für die Felder jetzt immer 0.25 statt vorher zwischen 0.0625 und 0.25, was die Berechnung erleichtert: Man kann im Voraus alle Werte mal 0.25 nehmen; so muss man beim Durchlaufen nur noch zusammenzählen. Es werden somit vier Pointer für das Kernelfenster und einer für das Output-Feld benutzt. Jetzt ist nur Vorsicht geboten: Wendet man den gleichen Filter zwei Mal (bin mir nicht sicher bei der Schreibweise) an, so stimmt das Ergebnis nicht. Das Output-Fenster muss einmal an (0,0) und einmal an (1,1) liegen.
Nimmt man all dies zusammen, so erhält man folgende Implementierung:void gaussian3filter(double *img, double *buf, unsigned sx, unsigned sy) { // this is a 3x3 gaussian filter realized by 2 2x2 mean value filters double *pass1 = buf; double *pass2 = img+sx+1; // first pass ptrs double *p1_1 = img; double *p1_2 = img+1; double *p1_3 = img+sx; double *p1_4 = img+sx+1; // second pass ptrs double *p2_1 = buf; double *p2_2 = buf+1; double *p2_3 = buf+sx; double *p2_4 = buf+sx+1; // perform 2 line of 1st pass for (unsigned i = 0; i < sx-1; ++i) *pass1++ = .25*((*p1_1++) + (*p1_2++) + (*p1_3++) + (*p1_4++)); // skip last field ++pass1; ++p1_1; ++p1_2; ++p1_3; ++p1_4; // 1 line more to go for (unsigned i = 0; i < sx-1; ++i) *pass1++ = .25*((*p1_1++) + (*p1_2++) + (*p1_3++) + (*p1_4++)); ++pass1; ++p1_1; ++p1_2; ++p1_3; ++p1_4; // now continue with the 1st pass, but work on 2nd in [kor]parallel[/kor] for (unsigned i = 2; i < sy-1; ++i) // we start in line 3 { for (unsigned j = 0; j < sx-1; ++j) { *pass1++ = .25*((*p1_1++) + (*p1_2++) + (*p1_3++) + (*p1_4++)); *pass2++ = .25*((*p2_1++) + (*p2_2++) + (*p2_3++) + (*p2_4++)); } ++pass1; ++p1_1; ++p1_2; ++p1_3; ++p1_4; ++pass2; ++p2_1; ++p2_2; ++p2_3; ++p2_4; } // 2 more lines for 2nd pass for (unsigned i = 0; i < sx-1; ++i) *pass2++ = .25*((*p2_1++) + (*p2_2++) + (*p2_3++) + (*p2_4++)); // skip last field ++pass2; ++p2_1; ++p2_2; ++p2_3; ++p2_4; // 1 line more to go for (unsigned i = 0; i < sx-1; ++i) *pass2++ = .25*((*p2_1++) + (*p2_2++) + (*p2_3++) + (*p2_4++)); }
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Erweiterte Filter
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Natürlich gibt es noch andere Filter, die nicht auf der Faltung von Urbild mit Kernel basieren. Sie machen sich meist andere Eigenschaften von realen Bildern oder unseres Perzeptionsmechanismus zunutze. Da die Implementierungen recht einfach sind und kaum Hilfsfunktionen benötigen, verzichte ich hier auf die Octave-Variante.
Thresholding
Beim Thresholding wird, wie der Name schon sagt, ein Schwellwert festgelegt und dementsprechend pixelweise für einen bestimmten Wert entschieden. Meistens benutzt man es, um s/w-Bilder zu erhalten, man kann das Verfahren allerdings auch auf mehrere Farben auslegen, z.B. bei der Farbraumumwandlung von 16bit auf 8bit bei Farb- oder Grauwertbildern. Da dieser Filter so trivial ist, werde ich keine Implementierung zeigen.Medianfilter
Dem Namen nach ähnelt er dem Mittelwertfilter aus dem vorigen Kapitel. Das dahinter liegende Verfahren ist jedoch ein komplett anderes: Die Farb-/Grauwerte werden über eine gewisse Fenstergröße (Kernelsize) in ein Array kopiert und sortiert, dann das Pixel auf den Wert gesetzt, der in der Mitte dieses Arrays steht. Sinn dieses Verfahrens ist z.B. die Feststellung von Inhomogenitäten in Kombination mit dem Differenzbildverfahren: Weicht ein Pixel zu stark von seiner Umgebung ab, so liegt vermutlich eine Unregelmäßigkeit vor (bei Gusseisenteilprüfung z.B. ein Gasbläschen).
Der Vorteil des Medianfilters gegenüber dem Mittelwertfilter ist, dass diese Methode sehr ausreißerstabil ist. Wenn sich unter der Filtermaske viele Pixel gleicher Helligkeit befinden und einzelne sehr verschiedene, so haben diese keinen Einfluss auf den resultierenden Helligkeitswert. Die empfohlene Filtergröße liegt bei sieben bis dreizehn, je nach Verwendung.template<type T> void medianFilter(T *din, T *dout, unsigned sx, unsigned sy, unsigned ks) { T *values = new T[ks*ks]; for (unsigned y = (ks/2); y < sy - (ks/2); ++y) { for (unsigned x = (ks/2); x < sx - (ks/2); ++x) { // Werte aufsammeln for (unsigned i = 0; i < ks; ++i) for (unsigned j = 0; j < ks; ++j) values[i*ks + j] = din[(y+j-ks/2)*sx + (x+i-ks/2)]; // Array sortieren std::sort(values, values+(ks*ks)); // Aktuelles Pixel auf Median setzen dout[y*sx + x] = (T)(*(values+(ks/2))); } } delete [] values; return; }
Histogram Equalization
Ziel der Histogrammnormalisierung ist eine Verbesserung des Kontrastes und Glättung des Histogramms. Letztere dient dazu, eigentlich gleiche, aber unter verschiedenen Lichtverhältnissen aufgenomme Bilder besser vergleichen zu können. Es wird zuerst ein Histogramm des Bildes erstellt und anschließend normalisiert, sodass das gesamte Spektrum des Farbraumes (Grauraumes) ausgenutzt wird. Wichtig ist jedoch, dass das Bildformat über einen Helligkeitskanal (YCrCb, HLS, ...) verfügt. Die hier vorgestellte Implementierung arbeitet mit [0;255]-Grauwertbildern.
Der Algorithmus:- Erstelle eine Tabelle mit Zuordnung Grauwert-Häufigkeit
- Erstelle die laufende Summe
- Normalisiere jedes Element der Summe durch Teilung durch die Gesamtsumme
- Multipliziere diese Werte mit dem höchsten Grauwert
- Ersetze die neue Grauwerttabelle mit der alten
- Aktualisiere das Bild
void [kor]histogramEqualization[/kor](unsigned char *din, unsigned char *dout, unsigned sx, unsigned sy) { int container[256]; std::memset(container, 0, sx*sy); // Container fuellen for (unsigned y = 0; y < sy; ++y) for (unsigned x = 0; x < sx; ++x) container[din[y*sx+x]]++; // Graymap erstellen std::map<unsigned char, int> graymap; for (unsigned i = 0; i < 256; ++i) if (container < 0) graymap[i] = container[i]; // Laufende Summe bilden unsigned char hgval = 0; double *sum = new double [graymap.size()] std::map<unsigned char, int>::iterator iter = graymap.begin(); sum[0] = (double)iter->second; for (unsigned i = 1, iter++; iter != graymap.end(); ++iter, ++i) { sum[i] = sum[i-1] + (double)iter->second; if (iter->first > hgval) hgval = iter->first; } for (unsigned i = 0; i < graymap.size(); ++i) sum[i] = (int)(sum[i] / sum[graymap.size()-1]) * (double)hgval); // Bild aktualisieren, dazu Grauwerttabelle anpassen for (unsigned i = 0; iter = graymap.begin(); iter != graymap.end(); ++iter, ++i) container[iter->first] = (int)sum[i]; // container [kor]enthält[/kor] jetzt die Zuordnung der neuen Grauwerte for (unsigned y = 0; y < sy; ++y) for (unsigned x = 0; x < sx; ++x) dout[y*sx+x] = container[din[y*sx+x]]; }
Originalbild
Bild nach Histogram Equalization%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Defektpixelinterpolation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Röntgendetektoren, CCD-Kameras oder Fingerprintleser haben oft defekte Pixel. An dieser Stelle liefern sie im Normalfall den Minimal- oder Maximalwert des Wertebereichs zurück. Bei Röntgendetektoren oder CCD-Kameras kann man diese defekten Pixel leicht identifizieren, indem man z.B. ein sog. Hellbild erstellt. Bei Röntgen ist das Vollbelichtung ohne durchstrahltes Objekt, bei CCD starke Überbelichtung ("Hellbild"). Das Ergebnis ist ein weißes Bild, das an den defekten Pixeln schwarz ist. Wiederholt man das mit einem Dunkelbild, so erhält man die entsprechend "anders" defekten Pixel. Hat man die defekten Pixel erst einmal identifiziert, kann man versuchen sie gezielt zu interpolieren. Oftmals genügt ein Mittelwertfilter im betroffenen Gebiet, möchte man allerdings mehr Information wiederherstellen oder ist das defekte Gebiet größer, so bietet sich das sog. [kor][i]band limitation*-Verfahren[/kor] an.
Wem die Fouriertransformation nicht geläufig ist, dem empfehle ich, zuerst eine gute Erklärung dazu z.B. auf www.wikipedia.org zu lesen
Die Idee: Das Bild wird als zweidimensionales Signal betrachtet. Man nimmt nun an, dass es in realen Bildern keine harten Sprünge in den Farb-/Grauwerten gibt, sondern immer mehr oder weniger weiche Übergänge. Wendet man nun die Fouriertransformation auf ein gewähltes Fenster des Bildes an, so bekommt man das Frequenzspektrum dieses Fensters, welches einem angibt, welche Sinus- oder Cosinusschwingungen man mit welcher Gewichtung überlagern müsste, um das Ausgangssignal zu erhalten. Harte Sprünge würden sich in hohen Frequenzen darstellen, da diese schnell oszillieren und man nur dadurch den geforderten schnellen Farbübergang erreichen kann. Beim [kor]band limitation-Verfahren[/kor] werden nun einfach alle Frequenzen oberhalb eines Grenzwertes auf 0 gesetzt, also entfernt. Rekonstruiert man nun das Bildfenster aus diesem modifizierten Frequenzspektrum, so erhält man eine Version des Bildes ohne die schnellen Frequenzen. Abschließend überträgt man die als defekt markierten Pixel von dem rekonstruierten Bild ins Originalbild. Wendet man dieses Verfahren nun iterativ an, so lässt sich die Qualität der Interpolation immer weiter verbessern. Will man das Verfahren noch optimieren, so kann man den Grenzwert, bei dem die Frequenzen "abgeschnitten" werden sollen, bei jeder Iteration neu bestimmen.
Zu diesem Verfahren biete ich nur eine Octave-Implementierung, da für die C++-Variante eine FFT-Bibliothek benötigt werden würde und dann der Code fast analog zur Octave-Variante wäre.
An den Bildern sieht man, dass die Frequenzmethode für größere Fehler nur bedingt nützlich ist. Man könnte das Ergebnis verbessern, indem man eine kleinere Fenstergröße nimmt. Im Beispiel wurde das ganze Bild genommen, was hier einer Fenstergröße von 128 entspricht. Fenstergrößen von 64 Pixel haben sich bewährt.Octave-Implementierung
% Bild einlesen org = imread("res/brain_org.png"); [sx, sy] = size(org); % Welche Frequenzen sollen abgeschnitten werden? m = 35; % Ein paar [kor]künstliche[/kor] Fehler... def = org; def(30, 1:sy) = 256; def(1:sx, 40) = 256; def(60:65, 60:65) = 256; def(90:95, 60:85) = 256; % Iterationsschleife for k = 1 : 10 defFT = fftshift(fft2(def)); defFT(1:sx, 1:m) = 0; % Frequenzen abschneiden defFT(1:sx, sy-m:sy) = 0; defFT(1:m, 1:sy) = 0; defFT(sx-m:sx, 1:sy) = 0; recon = ifft2(fftshift(defFT)); % inverse FFT def(30, 1:sy) = real(recon(30, 1:sy)); % Rekonstruierte Pixel [kor]übertragen[/kor] def(1:sx, 40) = real(recon(1:sx, 40)); def(60:65, 60:65) = real(recon(60:65, 60:65)); def(90:95, 60:85) = real(recon(90:95, 60:85)); endfor % Interpoliertes Bild in 'def'
Originalbild
Defektes Bild
Interpoliertes Bild--------------------------------------------------------------------
Edit: Fehler gefunden
Edit: kor-Tags eingefügt
Edit: gefundene Fehler von Mr. B verbessert
Edit: kor-Tags eingefügt
- Mittelwert: Der Mittelwertfilter wird in verschiedenen Kernelgrößen verwendet und bildet jeweils den Mittelwert über die Umgebung. Kernelgrößen zwischen [kor]zwei[/kor] und fünf sind gebräuchlich.
-
Bei Inhalt Punkt 3.2: Implentierung richtig?????
Außerdem gibts kein "desweiteren"; wenn schon, dann "des Weiteren".
Also sagt mir jetzt, ob Punkt 3.2 richtig ist oder nicht und wenn ja, dann was es denn nun wirklich gibt und was nicht.
Dann gucke ich mir die letzte Version noch einmal an, weil ich ja schon beim bloßen Drüberschauen Fehler gefunden habe...
Mr. B