TI/Programowanie dla Fizyków Medycznych/Manipulacja obrazem
Spis treści
Maniplacja obrazem
Ćwiczenia z mainuplacji obrazem rozpocznijmy od wczytania pliku który będziemy przetwarzać. Dla fizyków medycznych naturalnie będzie to plik DICOM.
import dicom
import numpy as np
import pylab as py
plik=dicom.read_file('I00001.dcm')
pixel=plik.pixel_array
print pixel.shape
>>>(2964, 2364)
Stworzona w ten sposób tablica ma dosyć duże rozmiary. O ile wyświetlenie takiego pliku nie jest problemem dla znajdujących się w pracowni komputerów, to bardziej złożona modyfikacja takiej grafiki mogła by z czasem obliczeń znacznie wykraczać poza czas przewidziany na zajęcia. Wyłącznie dla celów dydaktycznych, by ułatwić i przyspieszyć pracę zmniejszymy rozmiar badanego obrazu.
pixel=pixel[::5,::5]
print min(pixel.flatten()),max(pixel.flatten())
>>> 0 3907
Wartości zapisane w tablicy odpowiadającej zmniejszonemu obrazowi mają wartości w przedziale [0,3907]. Metoda wyświetlając imshow najwyższej wartości czyli 3908 przypisze kolor biały, wartości 0 kolor czarny. Aby ułatwić sobie manipulacje kolorami przypiszmy im wartości z zakresu [0,1].
pixel=pixel*1.0/max(pixel.flatten())
Tak powstały obraz możemy dowolnie modyfikować. Aby lepiej widoczne były skutki manipulacji obrazem dodajmy u góry obrazka þlynne przejście od czerni do bieli
def dodajPasek(tablica, n):
tablica[:n,:]=np.linspace(0,1,tablica.shape[1])
dodajPasek(pixel,30)
py.imshow(pixel,cmap=py.cm.gray,interpolation='nearest')
py.show()
Jasność
def jasnosc(tablica, a):
def f(x,a):
return min(a*x,1.0)
wynik=tablica.copy()
for x,y in np.ndindex(tablica.shape):
wynik[x,y]=f(tablica[x,y],a)
dodajPasek(wynik,15)
return wynik
py.imshow(jasnosc(pixel,2),cmap=py.cm.gray,interpolation='nearest')
py.show()
py.imshow(jasnosc(pixel,0.5),cmap=py.cm.gray,interpolation='nearest')
py.show()
Gamma
def gamma(tablica, a):
def f(x,a):
return min(x**a,1.0)
wynik=tablica.copy()
for x,y in np.ndindex(tablica.shape):
wynik[x,y]=f(tablica[x,y],a)
dodajPasek(wynik,15)
return wynik
py.imshow(gamma(pixel,2),cmap=py.cm.gray,interpolation='nearest')
py.show()
py.imshow(gamma(pixel,0.5),cmap=py.cm.gray,interpolation='nearest')
py.show()
Próg
def prog(tablica, a):
def f(x,a):
return 0 if x<a else 1
wynik=tablica.copy()
for x,y in np.ndindex(tablica.shape):
wynik[x,y]=f(tablica[x,y],a)
dodajPasek(wynik,15)
return wynik
py.imshow(prog(pixel,0.5),cmap=py.cm.gray,interpolation='nearest')
py.show()
Schodek
def schodek(tablica, a,b):
def f(x,a,b):
return 1 if ((x<a) or (x>b)) else 0
wynik=tablica.copy()
for x,y in np.ndindex(tablica.shape):
wynik[x,y]=f(tablica[x,y],a,b)
dodajPasek(wynik,15)
return wynik
py.imshow(schodek(pixel,0.2,0.6),cmap=py.cm.gray,interpolation='nearest')
py.show()
Rozmycie i wyostrzenie
from scipy import ndimage
imag=ndimage.gaussian_filter(pixel,sigma=4)
print imag
WYOSTRZ
def wyostrz(tablica, sigma, a):
rozmyty=np.array(ndimage.gaussian_filter(tablica,sigma=sigma))
roznica=tablica-rozmyty
wynik=rozmyty+a*roznica
wynik/=max(wynik.flatten())
dodajPasek(wynik,15)
return wynik
py.imshow(wyostrz(pixel,4,0),cmap=py.cm.gray,interpolation='nearest')
py.show()
py.imshow(wyostrz(pixel,4,3),cmap=py.cm.gray,interpolation='nearest')
py.show()