OpenCV & Python – La binarizzazione di Otsu

la binarizzazione di Otsu OpenCV Python

Introduzione

In questo articolo vedremo come applicare con la libreria OpenCV un’altra importante tecnica: la binarizzazione di Otsu. Questa tecnica è molto importante nell’analisi delle immagini, specialmente nei casi in cui vogliamo applicare una soglia nelle tecniche di thresholding in maniera efficiente.

the Otsu Binarization OpenCV Python main

Le immagini bimodali

Per applicare al meglio questa tecnica, le immagini da utilizzare dovrebbero essere bimodali. Ma che cosa significa immagine bimodale?

Qualsiasi immagine è un insieme di puntini definiti come pixel, e ciascuno di questi pixel può essere rappresentato mediante una tripletta di valori RGB nel caso di una immagine a colori, e un valore singolo se l’immagine è nella scala dei grigi. Questi valori assumono tutti i valori compresi tra 0 e 255.  Se prendiamo tutti i pixel di una immagine e contiamo quanti di questi hanno valore 0, quanti 1, quanti 2, e così via fino a 255, otteniamo quello che si chiama un istogramma.

image histogram

Come si può vedere dalla figura un istogramma non è altro che un modo di rappresentare la distribuzione del grado di colore presente in una immagine. Dalla figura per esempio possiamo vedere che abbiamo un massimo nella parte centrale. Quindi in questo caso l’immagine ha una sola moda.

Invece un’immagine bimodale, una volta rappresentata sotto forma di istogramma, dovrà presentare due massimi distinti tra di loro (mode). Prendiamo per esempio questa immagine a colori che io ho prodotto aggiungendo un po’ di rumore di fondo.

Adesso vedremo come effettuare un’analisi tramite openCV per ottenere l’istogramma dell’immagine e capire se l’immagine è bimodale. Scrivete il seguente programma e salvatelo come otsu01.py.

import cv2
import numpy as np
from matplotlib import pyplot as plt

img = cv2.imread('noisy_leaf.jpg',0)
plt.subplot(2,1,1), plt.imshow(img,cmap = 'gray')
plt.title('Original Noisy Image'), plt.xticks([]), plt.yticks([])
plt.subplot(2,1,2), plt.hist(img.ravel(), 256)
plt.title('Histogram'), plt.xticks([]), plt.yticks([])

plt.show()

Eseguendo il programma, l’immagine dapprima verrà caricata in scala di grigi e poi verrà generato l’istogramma dell’immagine. Visualizzandolo subito dopo l’immagine tramite matplotlib è possibile vedere che l’immagine da noi scelta è perfettamente bimodale.

Come potete vedere presenta due mode ben distinte tra di loro, dividendo in due parti la distribuzione dei pixel dell’immagine.

La binarizzazione di Otsu

Adesso che abbiamo compreso insieme cosa è un’immagine bimodale e gli strumenti per poterla riconoscere come tale, andiamo a vedere in che cosa consiste la binarizzazoine di Otsu.

Ritornando all’immagine precedente abbiamo visto che l’immagine bimodale si presenta sull’istogramma in due distribuzioni ben distinte, quasi separabili tra di loro. Infatti tra le due mode ci sta un punto di minimo, in cui si potrebbe considerare l’eventualità di separare l’istogramma in due parti. Bene la binarizzazione di Otsu ci aiuta ad ottenere automaticamente quel valore.

Inoltre una volta acquisita questa tecnica non ci sarà più bisogno di visualizzare e studiare l’istogramma per ricavarci tale punto, ma il tutto verrà fatto automaticamente.

Il thresholding

Una delle tecniche più utilizzate per l’analisi delle immagini è quella del thresholding, cioè l’applicazione di una soglia lungo una particolare scala di valori, per filtrare in qualche modo un’immagine.

Una di queste tecniche è per esempio quella che converte una qualsiasi immagine in scala di grigi (o a colori) in una totalmente in bianco e nero. Spesso questo è molto utile per riconoscere delle forme regolari, dei contorni o delle figure all’interno di un’immagine, o anche per delimitarne e suddividere zone all’interno, per poi essere utilizzate in maniera differente nelle elaborazioni successive.

Quindi applicata ad un istogramma, si andrà a scegliere un valore in cui tutti i valori sottostanti verranno convertiti a 0 (bianco) e tutti quelli sovrastanti a 255 (nero), convertendo un’immagine da scala di grigi in bianco e nero.

In OpenCV per effettuare il thresholding si utilizza la funzione cv2.threshold()

Prendiamo il caso dell’immagine della foglia precedente. Facciamo caso di aver bisogno di riconoscere la forma della foglia ma di non poter (saper o voler) utilizzare un’istogramma. Come primo approccio si proverà ad applicare un threshold (una soglia) a caso, e poi dopo diversi tentativi trovare un valore ottimale.

Il primo valore sicuramente da provare è quello di 127 (che nella scala da 0 – 255) sta perfettamente nel mezzo. Applichiamo quindi questo valore alla funzione cv2.threshold() scrivendo il seguente programma e salvatelo come otsu02.py.

import cv2
import numpy as np
from matplotlib import pyplot as plt

img = cv2.imread('noisy_leaf.jpg',0)

#ret1, th1 = cv2.threshold(img, 127, 255, cv2.THRESH_BINARY)

ret, imgf = cv2.threshold(img, 0, 255, cv2.THRESH_BINARY+cv2.THRESH_OTSU)

#blur = cv2.GaussianBlur(img, (5,5), 0)
#ret3, th3 = cv2.threshold(blur, 0, 255, cv2.THRESH_BINARY+cv2.THRESH_OTSU)

plt.subplot(3,1,1), plt.imshow(img,cmap = 'gray')
plt.title('Original Noisy Image'), plt.xticks([]), plt.yticks([])
plt.subplot(3,1,2), plt.hist(img.ravel(), 256)
plt.axvline(x=ret, color='r', linestyle='dashed', linewidth=2)
plt.title('Histogram'), plt.xticks([]), plt.yticks([])
plt.subplot(3,1,3), plt.imshow(imgf,cmap = 'gray')
plt.title('Otsu thresholding'), plt.xticks([]), plt.yticks([])
plt.show()

Eseguitelo

 

Come possiamo vedere dalla figura, in realtà il valore scelto da noi a caso è completamente fuori dal valore desiderato, e l’unica immagine che riusciamo a visualizzare sono le venature più grosse della foglia.

A questo punto potremmo scegliere dei valori di soglia facendo diversi tentativi, ma impiegheremmo tempo. Meglio sarebbe se potessimo visualizzare l’istogramma e scegliere un valore minimo tra le due mode. In realtà anche in questo caso il valore da noi scelto non sarebbe quello ottimale (infatti non è il minimo il valore corretto di soglia), ma quello che tiene conto dei pesi delle distribuzioni. Inoltre spesso la situazione degli istogrammi non è proprio così netta…

La binarizzazione di Otsu

Ecco qui che entra in gioco la binarizzazione di Otsu. Questo algoritmo ci permetterà di ottenere in maniera veloce ed automatica, il corretto valore di soglia da scegliere tra due mode nell’istogramma, potendo così applicare il thresholding in maniera ottimale.

In OpenCV, l’applicazione della binarizzazione di Otsu è molto semplice. Sarà sufficiente aggiungere una flag in più all’interno della funzione cv2.threshold(), chiamata

cv2.THRESH_OTSY

Vediamo insieme il codice. Scrivi il codice seguente come otsu03.py.

import cv2
import numpy as np
from matplotlib import pyplot as plt

img = cv2.imread('noisy_leaf.jpg',0)

ret, imgf = cv2.threshold(img, 0, 255, cv2.THRESH_BINARY+cv2.THRESH_OTSU)

plt.subplot(3,1,1), plt.imshow(img,cmap = 'gray')
plt.title('Original Noisy Image'), plt.xticks([]), plt.yticks([])
plt.subplot(3,1,2), plt.hist(img.ravel(), 256)
plt.axvline(x=ret, color='r', linestyle='dashed', linewidth=2)
plt.title('Histogram'), plt.xticks([]), plt.yticks([])
plt.subplot(3,1,3), plt.imshow(imgf,cmap = 'gray')
plt.title('Otsu thresholding'), plt.xticks([]), plt.yticks([])
plt.show()

Eseguendo il codice questa volta otterremo dei risultati ottimali. Inoltre è possibile vedere quale è il valore di soglia ottimale trovato dalla binarizzazione di Otsu all’interno dell’istogramma. Come possiamo vedere è tra le due mode, ma non si trova perfettamente nel punto di minimo.

OpenCV Otsu Binary threshold of noisy_leaf image 2

L’algoritmo: Come lavora la binarizzazione di Otsu

Dato che stiamo lavorando con immagini bimodali, l’algoritmo di Otsu cercherà di trovare il valore di soglia (threshold) t che minimizzi la varianza intra-classe (weighted within-class variance), definita come somma pesata delle varianze delle due classi:

 \sigma_w^2(t) = q_1(t)\sigma_1^2(t)+q_2(t)\sigma_2^2(t)

dove

 q_1(t) = \sum_{i=1}^{t} P(i) \quad \& \quad q_1(t) = \sum_{i=t+1}^{I} P(i)

e

 \mu_1(t) = \sum_{i=1}^{t} \frac{iP(i)}{q_1(t)} \quad \& \quad \mu_2(t) = \sum_{i=t+1}^{I} \frac{iP(i)}{q_2(t)}

e

 \sigma_1^2(t) = \sum_{i=1}^{t} [i-\mu_1(t)]^2 \frac{P(i)}{q_1(t)} \quad \& \quad \sigma_2^2(t) = \sum_{i=t+1}^{I} [i-\mu_1(t)]^2 \frac{P(i)}{q_2(t)}

In Internet ho trovato il seguente codice (vedi qui) che implementa le espressioni matematiche riportate sopra nell’algoritmo ed effettua tutti i calcoli necessari. Comunque ho dovuto apporre alcune modifiche dato che generava alcuni errori, come per esempio la divisione per zero. Scrivetelo e salvatelo come otsu04.py.

import cv2
import numpy as np

img = cv2.imread('noisy_leaf.jpg',0)
#blur = cv2.GaussianBlur(img,(5,5),0)

# find normalized_histogram, and its cumulative distribution functio
hist = cv2.calcHist([img],[0],None,[256],[0,256])
hist_norm = hist.ravel()/hist.max()
Q = hist_norm.cumsum()
bins = np.arange(256)
fn_min = np.inf
thresh = -1
for i in xrange(1,256):
 p1,p2 = np.hsplit(hist_norm,[i]) # probabilities
 q1,q2 = Q[i],Q[255]-Q[i] # cum sum of classes
 if q1 == 0:
    q1 = 0.00000001
 if q2 == 0:
    q2 = 0.00000001
 b1,b2 = np.hsplit(bins,[i]) # weights
 # finding means and variances
 m1,m2 = np.sum(p1*b1)/q1, np.sum(p2*b2)/q2
 v1,v2 = np.sum(((b1-m1)**2)*p1)/q1,np.sum(((b2-m2)**2)*p2)/q2
 # calculates the minimization function
 fn = v1*q1 + v2*q2
 if fn < fn_min:
 fn_min = fn
 thresh = i
# find otsu's threshold value with OpenCV function
ret, otsu = cv2.threshold(img,0,255,cv2.THRESH_BINARY+cv2.THRESH_OTSU)
print thresh,ret

Anche se eseguendo i risultati dei due threshold non mi vengono perfettamente uguali

204 202.0

Conclusioni

In questo articolo abbiamo visto come poter applicare al miglior modo la tecnica del thresholding nel caso delle immagini bimodali, e questo grazie alla binarizzazione di Otsu. In ulteriori articoli a venire, verranno approfondite altre tecniche riguardo al thresholding e all’analisi delle immagini in generale utilizzando la libreria OpenCV su Python.

Lascia un commento