IMAGE RECONSTRUCTION DALAM PENGOLAHAN CITRA MENGGUNAKAN MATLAB

loading...

IMAGE RECONSTRUCTION

1. Pendahuluan

Istilah rekonstruksi citra umum dipergunakan dalam lingkup pencitraan tomografi, yaitu teknik pencitraan yang menghasilkan gambaran potongan lintang suatu objek melalui pengolahan terhadap sinyal proyeksi trans-aksial dari objek tersebut. Sinyal proyeksi trans-aksial diperoleh dengan cara memberikan radiasi terhadap objek dalam berbagai sudut orientasi. Modalitas pencitraan tomografi pada umumnya memiliki sistem mekanik yang mengatur rotasi posisi transduser pemancar radiasi dan transduser pendeteksi sinyal sehingga pengambilan proyeksi objek dapat dilakukan dari berbagai sudut orientasi.

Sinyal proyeksi yang menjadi input dari proses rekonstruksi citra dapat dihasilkan dari 3 macam teknik pencitraan tomografi, yaitu:

• Transmission Tomography

Pada transmission tomography, sinyal informasi proyeksi diperoleh dari intensitas sinyal radiasi yang ditangkap setelah sinyal menembus objek. Contoh modalitas pencitraan yang temasuk kategori ini adalah CT-Scanner. CT-Scanner mempergunakan radiasi sinar X dalam berbagai orientasi untuk menghasilkan citra irisan lintang objek pencitraan.

• Reflection Tomography

Teknik tomografi refleksi banyak dipergunakan dalam aplikasi radar. Dalam aplikasi ini, sinyal informasi proyeksi objek merupakan intensitas sinyal refleksi yang dipantulkan oleh permukaan luar objek, misalnya seperti pada aplikasi doppler radar-imaging.

• Emission Tomography

Tomografi emisi memanfaatkan sinyal informasi berupa intensitas peluruhan radioaktif dari komponen radionuklida yang dimasukkan ke dalam objek. Dalam aplikasi medis, contoh modalitas pencitraan yang menggunakan metode tomografi emisi antara lain adalah PET (Positron Emission Tomography).

2. Rekonstruksi Citra 2 Dimensi dari Sinyal 1 Dimensi

Transformasi Radon dan invers transformasi Radon merupakan pendekatan paling langsung terhadap permasalahan proyeksi dan rekonstruksi citra, sehingga teknik ini dipergunakan sebagai teknik dasar dalam rekonstruksi citra 2D dari sinyal-sinyal 1D.

2.1 Transformasi Radon

Transformasi Radon merepresentasikan citra sebagai kumpulan sinyal 1D hasil proyeksi citra tersebut pada berbagai sudut orientasi. Transformasi Radon dari suatu citra pada suatu sumbu proyeksi dinyatakan sebagai integral garis dari fungsi representasi citra f(x,y) pada arah yang tegak lurus dengan sumbu proyeksinya. Hasil transformasi Radon dalam domain transformasi dinyatakan sebagai sinyal g(s,θ), dimana s merupakan jarak titik yang ditransformasikan terhadap titik pusat koordinat citra asal, dan θ merupakan sudut orientasi proyeksi.

Gambar 1. Ilustrasi transformasi Radon dalam koordinat asal kartesian

Apabila citra dinyatakan dalam format koordinat asal polar f(r,θ), integral garis transformasi Radon sepanjang r menghasikan kurva sinusoidal s = r cos(φ−θ). Bentuk sinusoidal ini menghasilkan istilah umum “sinogram” untuk kumpulan hasil transformasi Radon pada berbagai orientasi proyeksi.

Gambar 2. Ilustrasi transformasi Radon dalam koordinat asal polar

Dalam image processing toolbox Matlab, transformasi Radon dapat dilakukan dengan memanggil fungsi berikut:

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

dimana R merupakan hasil transformasi Radon terhadap citra asal I, pada jangkauan orientasi arah theta dan jangkauan lebar daerah proyeksi pada koordinat transformasi xp.

Berikut adalah code lengkap dari simulasi transformasi Radon terhadap gambar phantom irisan lintang kepala menggunakan Matlab: 

Gambar 3. Hasil transformasi Radon pada citra irisan kepala

2.2 Transformasi Radon Invers

Invers transformasi Radon dari proyeksi g(s,θ) ≡ R(f), -∞<s<∞, 0θ≤<π dirumuskan sebagai berikut:

Transformasi Radon invers terdiri dari dua buah tahapan, yaitu:

·         Proyeksi balik (backprojection)

·         Unblur filtering

Gambar 4. Skema transformasi Radon invers

Agar menghasilkan rekonstruksi sempurna, proses unblur filtering membutuhkan filter unblur 2 dimensi. Akan tetapi pemfilteran 2 dimensi dengan kernel unblur merupakan operasi yang sangat exhaustive secara komputasional, terutama untuk resolusi citra yang besar. Strategi yang biasa digunakan dalam aplikasi praktis adalah melalui unblur filtering pada sinyal proyeksi 1 dimensi, dan dilanjutkan dengan operasi backprojection. Strategi ini diilustrasikan pada Gambar 5 berikut:

Gambar 5. Skema transformasi Radon invers, apabila unblur filtering dilakukan terhadap sinyal satu dimensi

Penukaran posisi blok unblur filtering dan blok proyeksi balik (antara Gambar 4 dan Gambar 5) diperbolehkan oleh adanya teorema (irisan) – proyeksi. Dengan berpegang pada teorema ini, proses unblur filtering dapat dikerjakan pada domain frekuensi, sehingga akan mengurangi beban komputasi.

Gambar 6. Teorema (Irisan) - Proyeksi

2.2.1 Operasi unblur filtering

Mekanisme unblur filtering dapat dilakukan dengan menggunakan berbagai skema filter unblur yang dapat dipilih. Filter unblur yang umum dipergunakan adalah filter Shepp-Logan dengan properti sebagai berikut:

·         Respons frekuensi : H(ξ) = |ξ| sinc(ξd) rect(ξd)

·         Respons impuls : h(s) = 2(1+sin2π ξ₀s) /π²(d² – 4s²)

·         Respons impuls diskrit : h(m) =2 / (π²d(1 - 4m²))

2.2.2 Operasi proyeksi balik Operasi backprojection dapat didefinisikan sebagai berikut:

·            Koordinat kartesian:

·            Koordinat polar:

Hasil operasi proyeksi balik menunjukkan jumlah berkas proyeksi yang melewati suatu titik (x,y) dalam citra, sebagaimana diilustrasikan dalam Gambar 7 berikut:

Gambar 7. Ilustrasi operasi backprojection

Hasil operasi backprojection dapat dinyatakan sebagai citra asal yang dikaburkan dengan fungsi blurring dengan point spread function berupa (x² +y²)^−1/2.

Dalam code Matlab, operasi transformasi Radon invers dilakukan dengan memanggil fungsi sebagai berikut: (ilustrasi code sebagai lanjutan dari code transformasi Radon sebelumnya)

Fungsi ini sekaligus menjalankan operasi unblurring dan backprojection. Matlab menyediakan berbagai skema unblurring, yang dapat dipilih.

Gambar 8. Hasil transformasi Radon invers (citra irisan kepala)

3. Rekonstruksi Citra 3 Dimensi dari Citra 2 Dimensi

Operasi rekonstruksi citra 3 dimensi dari citra 2 dimensi pada dasarnya merupakan bentuk lanjut dari operasi rekonstruksi citra 2 dimensi dari citra 2 dimensi yang berbeda sumbu iris. Proses ini memiliki spesifikasi input berupa kumpulan citra irisan lintang dari suatu objek dalam koordinat irisan yang terurut (image slices) pada sumbu irisan tertentu (misalnya: potongan transversal), dan spesifikasi output berupa kumpulan citra irisan lintang objek yang sama pada sumbu irisan yang berbeda (misalnya sagittal atau koronal).

Apabila rekonstruksi dilakukan untuk 3 sumbu irisan yang mewakili sistem koordinat ruang 3 dimensi, maka dapat dilakukan graphical rendering untuk menghasilkan citra rekonstruksi volumetrik 3 dimensi. Dalam pembahasan ini akan dijelaskan teknik rekonstruksi citra 2 dimensi dari citra 2 dimensi berbeda sumbu iris.

PERCOBAAN

Rekonstruksi irisan sagittal dan coronal dari citra MRI dengan sumbu iris transversal.

a. Loading & viewing movie citra MRI sumbu iris transversal

             

truesizewarning = iptgetpref('TruesizeWarning');

iptsetpref('TruesizeWarning','off'); % Turn warning off

load mri;

figure;

immovie(D,map);

montage(D,map);

title('Horizontal Slices');

Penjelasan:

·         Gambar MRI yang dipergunakan terdiri atas 27 slices gambar MRI cranium bertipe indexed dengan ukuran per slice 128x128 piksel, disimpan sebagai mri.mat.

·         Perintah load mri akan menambahkan 2 variabel ke dalam workspace, yaitu variabel D (128x128x1x27) dan variabel map (89x3).

·         Variabel D berisi data citra, sedangkan variabel map berisi data colormap yang digunakan sebagai acuan tipe gambar indexed. Dimensi variabel D disusun sebagai berikut:

Dimensi 1: array data MRI terurut dari arah kepala depan ke kepala belakang

Dimensi 2: array data MRI terurut dari arah kepala kiri ke kepala kanan

Dimensi 4: array data MRI terurut dari arah kepala bawah ke kepala atas

·         Pada saat akuisisi citra, dimensi 1 dan 2 di-sampling 2.5 kali lebih rapat dari dimensi 4.

·         Perintah immovie akan menghasilkan movie citra MRI transversal menurut dimensi ke-4.

·         Perintah montage akan menghasilkan tampilan citra MRI lengkap sebanyak 27 slices dengan urutan slices dari kiri ke kanan.

 

Gambar 9. Citra irisan kepala – transversal

b. Ekstraksi citra irisan sagittal dari citra irisan transversal dengan imtransform

M1 = D(:,64,:,:);

M2 = reshape(M1,[128 27]);

figure,imshow(M2,map);

title('Sagittal - Raw Data');

T0 = maketform('affine',[0 -2.5; 1 0; 0 0]);

R2 = makeresampler({'cubic','nearest'},'fill');

M3 = imtransform(M2,T0,R2);

figure,imshow(M3,map);

title('Sagittal - IMTRANSFORM')

Penjelasan:

·         Variabel M1 dari code di atas berisi data MRI pada koordinat dimensi kedua bernilai 64, berarti berarti array data irisan ke arah sagittal untuk ke-27 slices.

·         Perintah reshape menghasilkan ekstraksi gambar 2 dimensi dari array data 4 dimensi M1. Data hasil reshape merupakan ‘Sagittal-Raw Data’ dengan ukuran matriks 128x27.

·         Untuk memperbaiki bentuk irisan agar diperoleh orientasi citra yang tepat dan ukuran yang lebih representatif, dilakukan transformasi affine. 4 buah parameter pertama dalam transformasi affine merupakan elemen matriks 2x2 untuk fungsi rotasi dan dilatasi. Dalam hal ini, dilakukan rotasi citra sebesar 90 derajat dan dilatasi sebanyak 2.5 kali ke arah horizontal. 2 buah parameter transformasi affine terakhir merupakan opeartor translasi. Dalam hal ini tidak dilakukan translasi sehingga operator di-assign dengan nilai 0.

·         Selanjutnya untuk memperbaiki kualitas dilatasi, dilakukan interpolasi dengan fungsi resample. Citra diinterpolasi secara cubic untuk arah vertikal dan nearest neighbour untuk arah horizontal. Hasilnya merupakan gambar rekonstruksi ke arah sagittal.

Gambar 10. Citra hasil rekonstruksi – sagittal

c. Ekstraksi citra irisan sagittal dari citra irisan transversal dengan tformarray

T1 = maketform('affine',[-2.5 0; 0 1; 68.5 0]);

T2 = maketform('custom',3,2,[],@ipex003,64);

Tc = maketform('composite',T1,T2);

R3 = makeresampler({'cubic','nearest','nearest'},'fill');

M4 = tformarray(D,Tc,R3,[4 1 2],[1 2],[66 128],[],0);

figure, imshow(M4,map);

title('Sagittal - TFORMARRAY');

Penjelasan:

·         Langkah ini menghasilkan citra rekonstruksi sagittal yang sama dengan langkah kedua, namun tidak diperlukan citra antara M1 dan M2.

·         Variabel T1 mendefinisikan spesifikasi transformasi affine yang sama dengan langkah kedua, namun diberikan translasi sebanyak 68.5 ke arah vertikal agar nilai koordinat tetap positif.

·         Variabel T2 mendefinisikan spesifikasi transformasi custom untuk mengekstraksi informasi sagittal pada koordinat horizontal ke-64 dari gambar.

·         Variabel Tc mengkompositkan transformasi T1 dan T2. Fungsi tformarray melakukan transformasi Tc kepada data gambar D, disusul dengan interpolasi dengan resampler R3.

·         Parameter [4 1 2] melakukan pengurutan ulang terhadap dimensi dari data D

·         Parameter [1 2] menunjukkan urutan dimensi dari data hasil transformasi

·         Parameter [66 128] menunjukkan ukuran data hasil transformasi per slice.

·         Parameter [ ] dan 0 merupakan parameter translasi array (di-assign dengan nilai 0 karena tidak adanya perubahan terhadap urutan slices pada gambar).

·         Tformarray menyatakan perintah untuk melakukan transformasi 3-D ke 2-D.

a.       Loading & viewing movie citra MRI sumbu iris sagittal

T3 = maketform('affine',[-2.5 0 0; 0 1 0; 0 0 0.5; 68.5 0 -14]);

S = tformarray(D,T3,R3,[4 1 2],[1 2 4],[66 128 35],[],0);

figure;

immovie(S,map);

S2 = padarray(S,[6 0 0 0],0,'both');

montage(S2,map);

title('Sagittal Slices');

Penjelasan:

·         Langkah ini melakukan ekstraksi citra sagittal untuk keseluruhan koordinat irisan sagittal, tidak hanya untuk irisan pada koordinat horizontal ke-64.

·         Fungsi tformarray dipergunakan kembali, namun kini dipergunakan untuk transformasi dari 3-D ke 3-D.

·         Parameter transformasi affine menjadi parameter transformasi 3-D.

·         Diinginkan movie bergerak dari potongan kiri ke potongan kanan dengan urutan tampilan 30 slices kiri disusul 35 slice tengah-genap dan diakhiri dengan 30 slices kanan, sementara montage hanya menampilkan 35 citra mid-sagittal. Oleh karena itu, ditambahkan nilai 0.5 dan –14 untuk menormalisasi urutan slices 30-98 pada dimensi keempat menjadi nomor slice 1-35.

·         Fungsi padarray memberikan jeda di antara urutan citra sagittal agar movie lebih mudah diamati.

Gambar 11. Citra hasil rekonstruksi – sagittal

b.       Loading & viewing movie citra MRI sumbu iris coronal

T4 = maketform('affine',[-2.5 0 0; 0 1 0; 0 0 -0.5; 68.5 0 61]);

C = tformarray(D,T4,R3,[4 2 1],[1 2 4],[66 128 45],[],0);

figure;

immovie(C,map);

C2 = padarray(C,[6 0 0 0],0,'both');

montage(C2,map);

title('Coronal Slices');

% Restore preference for TruesizeWarning

iptsetpref('TruesizeWarning',truesizewarning);

Penjelasan:

·         Langkah ini melakukan hal yang sama dengan langkah sebelumnya, namun mengekstraksi gambar ke arah irisan coronal. Perhatikan perubahan urutan dimensi pada parameter transformasi affine dan tformarray.

·         Diinginkan movie bergerak dari potongan depan ke potongan belakang dengan urutan tampilan 8 slices depan disusul 57 slice tengah-genap dan diakhiri dengan 8 slices belakang, sementara montage hanya menampilkan 45 citra mid-coronal. Untuk itu, ditambahkan nilai -0.5 dan 61 untuk menormalisasi urutan slices 8-120 pada dimensi keempat menjadi nomor slice 57-1, dan parameter [66 128 45] pada tformarray untuk menampilkan 45 slices penyusun montage.

·         Fungsi padarray memberikan jeda di antara urutan citra coronal agar movie lebih mudah diamati.

Gambar 12. Citra hasil rekonstruksi – coronal

Dari hasil data citra irisan transversal, sagittal, dan coronal, dapat dihasilkan rekonstruksi volumetrik 3 dimensi. Dalam aplikasi medis, hal ini sangat menguntungkan karena dengan melakukan satu kali pencitraan tomografi, dapat diekstraksi informasi citra dari irisan yang lain. Selain efisiensi biaya, kemudahan semacam ini juga memperkecil risiko prosedur pencitraan terhadap pasien. 

4. Tugas

4.1 Transformasi Radon

1. Perhatikan Gambar 3. Berikan analisis mengapa transformasi Radon menghasilkan gambar seperti di sebelah kanan untuk citra phantom yang diberikan di sebelah kiri!
2. Lakukan kembali transformasi Radon untuk berbagai jangkauan nilai theta dan berbagai nilai step theta, kemudian bandingkan hasilnya. Jelaskan dimana dan bagaimana terdapat perbedaan!
3. Citra phantom yang diberikan pada Gambar 3 merupakan gambar phantom ideal yang disediakan oleh Matlab. Simulasikan transformasi Radon terhadap gambar phantom di atas pada kondisi tidak ideal (terdapat noise) dan analisis hasilnya untuk berbagai jenis dan nilai distribusi noise!

Hint: refer ke Matlab Help>Image Processing Tolbox>phantom, radon, dan fungsi terkait, dan modul praktikum sebelumnya (image enhancement) tentang pemberian noise pada citra.

4.2 Transformasi Radon Invers

1. Jelaskan parameter-parameter dari fungsi iradon dalam Matlab!
2. Lakukan transformasi Radon invers terhadap citra phantom dengan input hasil transformasi Radon dalam berbagai jangkauan theta yang berbeda. Analisis hasilnya!
3. Matlab menyediakan berbagai macam filter unblur yang dapat dipergunakan dalam proses rekonstruksi citra hasil transformasi Radon. Analisis hasil rekonstruksi citra phantom dengan menggunakan berbagai model filter unblur!

Hint: refer ke Matlab Help>Image Processing Tolbox>irad.

loading...