Der Harris Eckendetektor

1988 beschrieben Chris Harris und Mike Stephens einen verbesserten Eckendetektor als Erweiterung des Eckendetektors von Moravec.

Während Moravec die Summe der quadratischen Differenzen

\[f(\Delta x, \Delta y) = \sum_{(x_k, y_k)\in W} \left( I(x_k, y_k) - I(x_k + \Delta x, y_k + \Delta y) \right)^2\]

über einem Fenster \(W\) für vier grundlegende Richtungen

\[(\Delta x, \Delta y) \in \{ (1,0), (0, 1), (1,1), (1,-1) \}\]

betrachtet, erweiteren Harris und Stephens diese Idee. Sie approximieren \(I(x_k + \Delta x, y_k + \Delta y)\) über die s.g. Taylorapproximation ersten Grades, also

\[I(x_k + \Delta x, y_k + \Delta y) \approx I(x_k, y_k) + \Delta x \frac{\partial I}{\partial x} + \Delta y \frac{\partial I}{\partial y}\]

und erhalten damit die Approximation

\[f(\Delta x, \Delta y) \approx \sum_{(x_k, y_k)\in W} \left( \Delta x \frac{\partial I}{\partial x} + \Delta y \frac{\partial I}{\partial y} \right)^2\]

Diese Gleichung läßt sich mit

\[\begin{split}M = \sum_{(x_k,y_k)\in W} \begin{pmatrix} I_x^2 & I_x I_y\\ I_x I_y & I_y^2 \end{pmatrix}\end{split}\]

schreiben als

\[\begin{split}f(\Delta x, \Delta y) \approx \begin{pmatrix} \Delta x & \Delta y \end{pmatrix} \cdot M \cdot \begin{pmatrix} \Delta x \\ \Delta y \end{pmatrix}\end{split}\]

Dabei nennt man \(M\) Strukturtensor. Harris und Stephens berechnen dann die s.g. Eckenstärke als

\[R = det(M) - \kappa tr(M)^2\]

mit \(\kappa \approx 0.04\). In diesem Praktikum implementieren wir den Eckendetektor „from scatch“.

Der Code

In diesem Praktikum arbeiten Sie in der Datei

kanten.py

und implementieren die Funktion processImage

Schritt 1: Das Bild in Graustufen umwandeln

Wir müssen zunächst das Bild in Graustufen umwandeln. Konvertieren Sie das Bild mit der cv2.cvtColor Methode. Wandeln Sie das Bild danach über np.float32 in ein Float-Bild um. Normieren Sie die Grauwerte vorher indem Sie durch 255.0 teilen.

Lösung anzeigen

frame_gray = cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY)
frame_gray = np.float32(frame_gray / 255.0)

Schritt 2: Die Bildgradienten berechnen

Berechnen Sie nun mit der cv2.Sobel Methode wieder die Bildableitungen in X- und Y-Richtung. Verwenden Sie ksize=3 und ddepth=cv2.CV_32F.

Lösung anzeigen

gx = cv2.Sobel(frame, cv2.CV_32F, 1, 0, ksize=3) / 4.0
gy = cv2.Sobel(frame, cv2.CV_32F, 0, 1, ksize=3) / 4.0

Schritt 3: Den Strukturtensor \(M\) vorbereiten

Um den Strukturtensor berechnen zu können berechnen Sie zunächst

\[\begin{split}\begin{align} I_x^2 &=& \left(\frac{\partial I}{\partial x}\right)^2\\ I_y^2 &=& \left(\frac{\partial I}{\partial y}\right)^2\\ I_x\cdot I_y &=& \left(\frac{\partial I}{\partial x}\right)\cdot\left(\frac{\partial I}{\partial x}\right) \end{align}\end{split}\]

und speichern die Ergebnisse in drei neuen Bildern.

Lösung anzeigen

Ix2  = gx ** 2
IxIy = gx * gy
Iy2  = gy ** 2

Schritt 4: Den Strukturtensor berechnen

Der Strukturtensor \(M\) für jeden einzelnen Pixel ergibt sich durch Summation der Bilder \(I_x^2\), \(I_y^2\) sowie \(I_x I_y\). Um dies zu erreichen falten wir diese Bilder mit einem s.g. Block-Kernel

\[\begin{split}\begin{pmatrix} 1&\dots&1\\ \vdots&\ddots&\vdots\\ 1&\dots&1 \end{pmatrix}\end{split}\]

bestehend nur aus einsen. Importieren Sie das scipy Paket und verwenden Sie die convolve2d Methode um diese Faltung durchzuführen. Sie können einen geeigneten Kernel mit np.ones erzeugen. Tip: Dividieren Sie durch die Anzahl einsen in ihrem Kernel um effektiv der Mittelwert anstatt der Summe zu bilden.

Lösung anzeigen

N = 7
Ix2  = signal.convolve2d(Ix2, np.ones((N,N))) / (N**2)
Iy2  = signal.convolve2d(Iy2, np.ones((N,N))) / (N**2)
IxIy = signal.convolve2d(IxIy, np.ones((N,N))) / (N**2)

Schritt 5: Die Eckenstärke berechnen

Um die Eckenstärke zu berechnen betrachten wir zunächst die nochmal Gleichung

\[R = det(M) - tr(M)^2\]

mit

\[\begin{split}M = \sum_{(x_k,y_k)} \begin{pmatrix} I_x^2 & I_x I_y \\ I_x I_y & I_y^2 \end{pmatrix}\end{split}\]

Für die Determinante finden wir

\[det(M) = \left(\sum_{(x_k,y_k)} I_x^2\right) \left( \sum_{(x_k,y_k)} I_y^2 \right) - \left( \sum_{(x_k,y_k)} I_x I_y \right)^2\]

und für die Spur finden wir entsprechend

\[tr(M) = \sum_{(x_k,y_k)} I_x^2 + \sum_{(x_k,y_k)} I_y^2\]

Berechnen Sie nun die Determinante sowie die Spur im Python-Code und berechnen Sie dann die Eckenstärke \(R\) für \(kappa = 0.04\).

Lösung anzeigen

kappa = 0.04
det = Ix2 * Iy2 - IxIy ** 2
trace = Ix2 + Iy2
strength = det - kappa *  trace**2

Schritt 6: Die Eckenstärke anzeigen

Die Eckenstärke kann sowohl negativ als auch positiv sein. Wir sind jedoch nur an großen positiven Werten interessiert. Normieren Sie die Eckenstärke, indem Sie durch das globale Maximum dividieren (verwenden Sie np.max).

Nun binarisieren wir das Bild indem wir zunächst ein gleichgroßes Bild erzeugen (z.B. mit np.zeros_like) und überall dort, wo die Eckenstärke größer ist als ein Schwellwert \(T\) eine 1 setzen. Zeigen Sie dann beide Bilder mit cv2.imshow an.

Lösung anzeigen

strength /= np.max(strength)

corners = np.zeros_like(strength)
corners[strength > 0.1] = 1.0

cv2.imshow("Harris Corner Strength", strength)
cv2.imshow("Harris Corners", corners)

Musterlösung

Harris Eckendetektor - Musterlösung