Image Processing Toolbox



LUCRAREA NR 11

GENERAREA IMAGINII UNUI CAP FANTOMĂ - phantom

Sintaxa

P = phantom(def, n)

P = phantom(E, n)

[P, E] = phantom(...)

Descriere

P = phantom(def, n) generează o imagine a unui cap fantomă care poate fi folosită pentru a testa acurateţea numerică a unor algoritmi de reconstrucţie bidimensională.

P este o imagine grayscale care constă dintr-o elipsă mare (reprezentând creierul) care conţine mai multe elipse mai mici (reprezentând caracteristici ale acestuia).

def este un şir care specifică tipul de cap fantoma care va fi generat:

• 'Shepp-Logan' — Imagine test larg utilizată de către cercetătorii din domeniul tomografiei

• 'Modified Shepp-Logan' (predefinită) – Variantă a fantomei Shepp-Logan la care contrastul este îmbunătăţit pentru o mai bună percepţie vizuală

• n este un scalar care specifică numărul de linii şi coloane din P. Valoarea lui predefinită este 256.

P = phantom(E, n) generează o fantomă definită de câtre utilizator, în care fiecare linie a matricei E specifică o elipsă în imagine. E are şase coloane, fiecare dintre aceastea conţinând un parametru diferit al elipsei.

[P, E] = phantom(...) returnează matricea E utilizată pentru a genera fantoma.

Exemplu:

P = phantom('Modified Shepp-Logan',200);

imshow(P)

[pic]

Fig. nr. 1

TRANSFORMAREA RADON - radon

Sintaxa

R = radon(I, theta)

[R,xp] = radon(...)

Descriere

R = radon(I, theta) returnează transformarea Radon R a unei imagini I pentru un unghi de theta grade. Această transformare reprezintă proiecţia unei imagini de-a lungul unei linii radiale orientate sub un anumit unghi. Dacă theta este un scalar, R este un vector coloană care conţine transformarea Radon pentru theta grade. Dacă theta este un vector, R este o matrice în care fiecare coloană reprezintă transformarea Radon pentru unul dintre unghiurile din theta. Dacă argumentul theta este omis, valorile lui predefite variază între 0:179.

[R,xp] = radon(...) returnează un vector xp care conţine coordonatele radiale corespunzând fiecărei linii a R. Coordonatele radiale returnate în xp reprezintă valorile de-a lungul axei x', care este orientată în sens antiorar cu theta faţă de axa. Originea ambelor axe se află în pixelul central al imaginii, defint mai jos:

floor((size(I)+1)/2)

Exemplu (fig. nr. 2)

iptsetpref('ImshowAxesVisible','on')

I = zeros(100,100);

I(25:75, 25:75) = 1;

theta = 0:180;

[R,xp] = radon(I,theta);

imshow(R,[],'Xdata',theta,'Ydata',xp,...

'InitialMagnification','fit')

xlabel('\theta (degrees)')

ylabel('x''')

colormap(hot), colorbar

[pic]

Fig. nr. 2

Algoritm

Transformarea Radon a unei imagini reprezintă suma transformărilor Radon a fiecărui pixel individual. Algoritmul divide prima oară pixelii din imagine în patru subpixeli şi proiectează fiecare subpixel separat, ca în fig. nr. 3:

[pic]

Fig. nr. 3

Contribuţia fiecărui subpixel este împărţită proporţional între două fascicule învecinate, conform cu distanţa dintre locaţia proiectată şi centrul fascicolelor. Dacă proiecţia subpixelului corespunde punctului central al unui fascicul, fascicolul şi axele primesc valoarea întreagă a unui subpixel, sau o pătrime din valoarea pixelului. Dacă proiecţia subpixelului corespunde limitei dintre două fascicule, valoarea subpixelului este împărţită în mod egal între fascicule.

Funcţia radon calculează proiecţia unei matrice imagine de-a lungul unei direcţii specificate.

O proiecţie a unei funcţii bidimensionale f(x,y) este un set de integrale liniare. Funcţia radon calculează integrale liniare din mai multe surse, de-a lungul unor căi paralele, sau fascicule, intr-o anumită direcţie. Fasciculele sunt distanţate între ele la intervale de 1 pixel. Pentru a reprezenta o imagine, funcţia radon preia mai multe proiecţii de fascicule paralele, la diferite unghiuri, obţinute prin rotirea sursei în jurul centrului imaginii. Fig. nr. 4 prezintă o singură proiecţie la un anumit unghi de rotaţie.

[pic]

Fig. nr. 4

Integrala liniară a funcţiei f(x,y) pe direcţie verticală este proiecţia f(x,y) pe axa x; integrala liniară pe direcţie orizontală este proiecţia f(x,y) pe axa y. Figura nr. 5 prezintă proiecţiile verticale ale unei funcţii simple bidimensionale..

[pic]

Fig. nr. 5

Proiecţiile pot fi calculate de-a lungul oricărui unghi [[THETA]]. În general, transformarea Radon a funcţiei f(x,y) este integrala liniară a f, paralel cu axa y´ (fig. nr. 6)

[pic]

Unde:

[pic]

[pic]

Fig. nr. 6

Reprezentarea grafică a transformării Radon

Transformarea Radon a unei imagini I, pentru unghiurile specificate în vectorul theta , utilizând funcţia radon , poate fi calculată cu ajutorul sintaxei:

[R,xp] = radon(I,theta);

Coloana R conţine transformarea Radon pentru fiecare unghi din theta. Vectorul xp conţine coordonatele corespunzătoare de-a lungul axei x. Pixelul central al I este definit ca floor((size(I)+1)/2); acesta este pixelul de pe axa x´corespunzând la [pic].

Instrucţiunile de mai jos calculează şi reprezintă grafic transformarea Radon la 0° şi 45° a unei imagini care conţine un singur obiect de formă pătrată (Fig nr 7, 8, 9). xp este aceeaşi pentru toate unghiurile de proiecţie.

I = zeros(100,100);

I(25:75, 25:75) = 1;

imshow(I)

[R,xp] = radon(I,[0 45]);

figure; plot(xp,R(:,1)); title('R_{0^o} (x\prime)')

[pic]

Fig. nr. 7

[pic]

Fig. nr. 8 Transformarea Radon la 0 grade

figure; plot(xp,R(:,2)); title('R_{45^o} (x\prime)')

[pic]

Fig. nr. 9 Transformarea Radon la 45 grade

Vizualizarea transformării Radon sub forma unei imagini

Transformarea Radon pentru un număr mare de unghiuri este afişata de multe ori sub forma unei imagini. În exemplul următor, (fig. nr. 10), transformarea Radon a unei imagini pătrate este calculată pentru unghiuri variind între 0° şi 180°, cu incrementul de 1°.

theta = 0:180;

[R,xp] = radon(I,theta);

imagesc(theta,xp,R);

title('R_{\theta} (X\prime)');

xlabel('\theta (degrees)');

ylabel('X\prime');

set(gca,'XTick',0:20:180);

colormap(hot);

colorbar

[pic]

Fig. nr. 10

Utilizarea transformării Radon pentru detectarea liniilor

Transformarea Radon transform este strâns înrudită cu un operator uzual, cunoscut sub numele de transformarea Hough. Funcţia radon poate fi folosită pentru a implementa o formă a transformării Hough utilizată pentru detectarea liniilor drepte. Paşii necesari sunt prezentaţi mai jos:

Calculul unei imagini contur binare (fig. nr. 11), utilizând funcţia edge

I = fitsread('solarspectra.fts');

I = mat2gray(I);

BW = edge(I);

imshow(I), figure, imshow(BW)

[pic]

Fig. nr. 11

Calculul transformării Radon asupra imaginii contur (fig. nr. 12).

theta = 0:179;

[R,xp] = radon(BW,theta);

figure, imagesc(theta, xp, R); colormap(hot);

xlabel('\theta (degrees)'); ylabel('x\prime');

title('R_{\theta} (x\prime)');

colorbar

[pic]

Fig. nr. 12

Identificarea locaţiilor vârfurilor tari din matricea de transformare. Locaţiile acestor vârfuri corespund locaţiilor liniilor drepte din imaginea originală. În fig. nr. 13, cele mai intense vârfuri din R corespund pentru [pic] şi [pic]. Linia perpendiculară pe acel unghi şi localizată la [pic] este reprezentată mai jos, supraimprimată cu roşu în imaginea originală. Transformarea Radon este reprezentată în negru.

[pic]

Fig. nr. 13

TRANSFORMAREA RADON INVERSĂ - iradon

Sintaxa

I = iradon(R, theta)

I = iradon(P, theta, interp, filter, frequency_scaling, output_size)

[I,H] = iradon(...)

Descriere

I = iradon(R, theta) reconstruieşte imaginea I din datele stocate în matricea bidimensională R. În coloana R sunt datele de proiecţie ale fasciculelor paralele, definite ca ceil(size(R,1)/2).

theta descrie unghiurile (in grade) sub care au fost realizate proiecţiile. El poate fi fie un vector conţinând unghiurile sau un scalar care specifică D_theta, incrementul unghiului dintre proiecţii. Dacă theta este un vector, el trebuie să conţină unghiuri egal distanţate între ele. Dacă theta este un scalar care specifică incrementul D_theta, proiecţiile au fost realizate la unghiuri theta = m*D_theta, unde m = 0,1,2,...,size(R,2)-1. Dacă mărimea de intrare este o matrice goală ([]), D_theta capătă valoarea predefinită 180/size(R,2).

I = iradon(P, theta, interp, filter, frequency_scaling, output_size) specifică parametrii care se utilizează în transformarea Radon inversă. iradon utilizează valorile predefinite pentru argumentele omise.

interp specifică tipul de interpolare care poate fi folosită: 'nearest','linear','spline','cubic'.

filter specifică filtrele utilizate pentru filtrarea în domeniul frecvenţelor: 'Ram-Lak''Shepp-Logan''Cosine''Hamming''Hann''None'

frequency_scaling este un scalar cuprins în intervalul (0,1] care modifică parametrii de filtrare. Valoarea lui predefinită este 1.

output_size este un a scalar care specifică numărul de linii şi coloane ale imaginii reconstruite. Dacă output_size nu este specificată, dimensiunea ei este determinată din valorile proiecţiilor.

output_size = 2*floor(size(R,1)/(2*sqrt(2)))

Dacă se specifică output_size, iradon reconstruieşte o porţiune mai mică sau mai mare a imaginii, dar nu modifică scara acesteia. Dacă proiecţiile au fost calculate cu ajutorul funcţiei radon, imaginea reconstruită s-ar putea să nu aibă aceleaşi dimensiuni cu imaginea originală.

[I,H] = iradon(...) returnează răspunsul în frecvenţă al filtrului, în vectorul H.

Exemplu:

Compararea proiecţiilor filtrate şi nefiltrate (fig. nr. 14) backprojection.

P = phantom(128);

R = radon(P,0:179);

I1 = iradon(R,0:179);

I2 = iradon(R,0:179,'linear','none');

subplot(1,3,1), imshow(P), title('Original')

subplot(1,3,2), imshow(I1), title('Filtered backprojection')

subplot(1,3,3), imshow(I2,[]), title('Unfiltered backprojection')

[pic]

Fig. nr. 14

Calculează proiecţia inversă a unui singur vector proiecţie. Sintaxa funcţiei iradon nu ne permite să facem acest lucru în mod direct, deoarece dacă theta este un scalar, el este tratat ca un increment. .

P = phantom(128);

R = radon(P,0:179);

r45 = R(:,46);

I = iradon([r45 r45], [45 45])/2;

imshow(I, [])

title('Backprojection from the 45-degree projection')

Definiţia transformării Radon inverse

Funcţia iradon inversează transformarea Radon, prin urmare poate fi folosită pentru reconstrucţia imaginilor din datele de proiecţie. Transformarea Radon inversă este folosită în mod curent în aplicaţii tomografice.

Pentru o imagine dată I şi un set de unghiuri theta, funcţia radon poate fi folosită pentru calculul transformării Radon:

R = radon(I,theta);

Funcţia iradon poate fi prin urmare utilizată pentru reconstrucţia imaginii I.

IR = iradon(R,theta);

În exemplul următor, proiecţiile sunt calculate din imaginea originală I. Pentru cele mai multe aplicaţii, nu există nici o imagine originală din care sunt formate proiecţiile.De exemplu, în cazul tomografiei computerizate, proiecţiile se formează prin măsurarea atenuării fasciculelor de radiaţii care trec printr-un corp fizic, sub diferite unghiuri. Imaginea originală poate fi considerată drept o secţiune transversală prin corp, unde valorile intensităţii reprezintă densitatea acestuia. Proiecţiile sunt colectate cu ajutorul unor dispozitive special concepute în acest scop, după care imaginea obiectului este reconstituită cu ajutorul funcţiei iradon.

iradon reconstruieşte imaginea din proiecţiile fasciculelor paralele. In geometria fasciculelor paralele, fiecare proiecţie este formată prin combinarea unui set de integrale liniare.

Figura nr. 15 ilustrează modul în care geometria fasciculelor paralele se aplică în cazul tomografiei computerizate. Se remarcă faptul că există un număr egal de n emiţători şi n senzori. Fiecare senzor măsoară radiaţia emisă de emiţătorul care îi corespunde. Această aranjare corespunde integralei liniare care este calculată prin transformarea Radon.

[pic]

Fig. nr. 15

Îmbunătăţirea rezultatelor

iradon utilizează algoritmul filtered backprojection pentru a calcula transformarea Radon inversă. Acest algoritm formează o aproximare a imaginii I bazată pe proiecţiile din coloana R. Un rezultat mai bun poate fi obţinut prin utilizarea mai multor proiecţii pentru realizarea reconstrucţiei. Cu cât este mai mare numărul de proiecţii (lungimea lui theta), cu atât mai bine imaginea reconstituită IR aproximează imaginea originală I. Vectorul theta trebuie să conţină valori angulare care cresc monoton, cu un unghi incremental constant. Δ[[THETA]]. Atunci când valoarea acestei creşteri este cunoscută, ea poate fi transferată funcţiei iradon in locul matricei valorilor θ, ca mai jos:

IR = iradon(R,Dtheta);

Algoritmul filtrează proiecţiile din R după care reconstruieşte imaginea utilizând proiecţiile filtrate. În unele cazuri, zgomotul poate fi prezent în proiecţii. Pentru îndepărtarea zgomotelor de frecvenţă înaltă, se aplică o fereastră de filtrare care atenuează zgomotul. Funcţia iradon dispune de mai multe astfel de ferestre de filtrare. Exemplul următor apelează funcţia iradon şi aplică o fereastră de filtrare Hamming. Pentru a obţine date nefiltrate, în locul parametrului filtrului se specifică 'none'.

IR = iradon(R,theta,'Hamming');

iradon mai permite şi specificarea unei frecvente normalizate, D, deasupra căreia filtrul are un răspuns zero. D trebuie să fie un scalar cuprins în intervalul [0,1]. Cu această opţiune, axa frecvenţelor este rescalată, astfel încât întregul filtru este comprimat pentru a putea fi cuprins în domeniul de frecvenţe [0,D]. Această opţiune poate fi folositoare în cazurile în care proiecţiile conţin un volum mic de informaţii la frecvenţe înalte şi mult zgomot de înaltă frecvenţă. În această situaţie zgomotul poate fi complet înlăturat, fără a compromite reconstrucţia. Apelarea instrucţiunii iradon ca în exemplul de mai jos setează valoarea frecvenţei normalizate la 0.85.

IR = iradon(R,theta,0.85);

Exemple:

Instrucţiunile de mai jos ilustrează modul în care se realizează reconstrucţia unei imagini cu ajutorul unui set de date obţinute din proiecţii paralele. Imaginea de test este capul fantomă Shepp-Logan, care poate fi generat de funcţia phantom a Image Processing Toolbox.

Generarea imaginii Shepp-Logan (fig nr 16):

P = phantom(256);

imshow(P)

[pic]

Fig. nr. 16

Calculul transformării Radon Asupra capului fantomă pentru trei seturi diferite de valori ale lui teta: R1 are 18 proiecţii, R2 are 36 proiecţii şi R3 are 90 de proiecţii.

theta1 = 0:10:170; [R1,xp] = radon(P,theta1);

theta2 = 0:5:175; [R2,xp] = radon(P,theta2);

theta3 = 0:2:178; [R3,xp] = radon(P,theta3);

Afişarea unei reprezentări grafice a unei transformări Radon asupra capului fantomă Shepp-Logan. În figura nr. 17 este reprezentată R3, transformarea cu 90 de proiecţii.

figure, imagesc(theta3,xp,R3); colormap(hot); colorbar

xlabel('\theta'); ylabel('x\prime');

[pic]

Fig. nr. 17

Reconstrucţia imaginii fantomă din datele proiecţiilor create în mai sus şi afişarea rezultatelor.

I1 = iradon(R1,10);

I2 = iradon(R2,5);

I3 = iradon(R3,2);

imshow(I1)

figure, imshow(I2)

figure, imshow(I3)

Figura nr. 18 prezintă rezultatele celor trei reconstrucţii. De remarcat faptul că imaginea I1, care a fost reconstruită din doar 18 proiecţii, este reconstrucţia cu cea mai mică acurateţe. Imaginea I2, care a fost reconstruită din 36 proiecţii, are o calitate mai bună, dar nu suficientă pentru a distinge în mod clar micile elipse din porţiunea inferioară a imaginii. Imaginea I3, reconstruită cu ajutorul a 90 de proiecţii, se apropie cel mai mult de imaginea originală. De remarcat faptul că dacă numărul proiecţiilor este relativ mic (ca în cazul I1 şi I2), reconstrucţia poate să includă şi anumite artefacte.

[pic]

Fig. nr. 18

................
................

In order to avoid copyright disputes, this page is only a partial summary.

Google Online Preview   Download