PRAKTIKUM EL4027 PENGOLAHAN CITRA BIOMEDIKA EB7031 PENGOLAHAN CITRA BIOMEDIKA LANJUT MODUL 1 − IMAGE ENHANCEMENT
Dalam bidang pengolahan citra, istilah citra mengacu pada suatu fungsi intensitas dalam bidang dua dimensi. Pemrosesan citra menggunakan komputer membutuhkan citra digital sebagai masukannya. Oleh karena itu, kita mengenal beberapa macam format citra digital, yang masing-masing memiliki format penyimpanan dan pembacaan data yang berbeda-beda. Dengan menggunakan toolbox image processing di Matlab, proses pembacaan citra dapat dilakukan dengan mudah, sebagaimana diperlihatkan oleh beberapa contoh berikut:
•
Windows Bitmap (*.bmp) a=imread('taz_ref.bmp'); figure,imshow(a);
•
Joint Photographic Experts Group (*.jpg) b=imread('oranges.jpg'); figure,imshow(b);
•
Tagged Image File Format (*.tif) c=imread('saturn.tif'); figure,imshow(c);
Jumlah warna yang terdapat dalam suatu citra digital berkaitan erat dengan format data digital yang digunakan. Suatu citra 8 bit (misalnya pada format *.bmp 8-bit) dapat memiliki paling banyak 256 jenis warna di dalamnya. Sedangkan suatu citra 24 bit (misalnya pada format *.png 24-bit) memiliki 224 kombinasi data warna yang dapat ditampilkan. Suatu citra berwarna dapat dinyatakan dalam banyak cara. Salah satu metode yang paling populer adalah dengan menggunakan metode kanal RGB. Dalam format ini, suatu citra berwarna dinyatakan dalam bentuk gabungan dari tiga buah citra monochrome merah, hijau, dan biru yang berukuran sama. Warna dari setiap piksel citra digital bergantung pada kombinasi ketiga nilai intensitas piksel monochrome yang bersesuaian. Contoh program untuk mengekstraksi citra monochrome pada toolbox image processing di Matlab adalah sebagai berikut: d=imread('flowers.tif'); figure, imshow(d); dR=d(:,:,1); figure,imshow(dR),title('Kanal Merah'); dG=d(:,:,2); figure,imshow(dG),title('Kanal Hijau'); dB=d(:,:,3); figure,imshow(dB),title('Kanal Biru ');
Teknik image enhancement digunakan untuk meningkatkan kualitas suatu citra digital, baik dalam tujuan untuk menonjolkan suatu ciri tertentu dalam citra tersebut, maupun untuk memperbaiki aspek tampilan. Proses ini biasanya didasarkan pada prosedur yang bersifat eksperimental, subjektif, dan amat bergantung pada tujuan yang hendak dicapai.
1.
Operasi Titik
Operasi titik dalam image enhancement dilakukan dengan memodifikasi histogram citra masukan agar sesuai dengan karakteristik yang diharapkan. Histogram dari suatu citra adalah grafik yang menunjukkan distribusi frekuensi dari nilai intensitas piksel dalam citra tersebut. Contoh perintah untuk menampilkan histogram adalah: I=imread('rice.tif'); figure,imshow(I); figure,imhist(I);
Beberapa teknik image enhancement melalui operasi titik antara lain adalah intensity adjustment (termasuk brightening dan darkening), histogram equalization, dan thresholding. 1.1 Intensity Adjustment Intensity adjusment bekerja dengan cara melakukan pemetaan linear terhadap nilai intensitas pada histogram awal menjadi nilai intensitas pada histogram yang baru. Perintah umum untuk melakukan pemetaan linear tersebut adalah: J = imadjust(I,[low_in,high_in),[low_out,high_out])
dimana : low_in merupakan nilai intensitas yang akan dipetakan sebagai low_out high_in merupakan nilai intensitas yang akan dipetakan sebagai high_out
Contoh: Citra rice.tif di samping memiliki nilai kekontrasan yang rendah. Berdasarkan histogramnya, dapat diketahui bahwa citra ini tidak memiliki piksel dengan intensitas di bawah 40 dan di atas 225. Untuk memperbaikinya, kita dapat memetakan histogram secara linear sehingga diperoleh sebuah citra baru yang memiliki rentang histogram antara 0 hingga 255. Perintah untuk melakukan intensity adjustment: I=imread('rice.tif'); J=imadjust(I,[0.15 0.9],[0 1]); figure,imshow(I); figure,imhist(I); figure,imshow(J); figure,imhist(J);
1.2
Histogram Equalization
Teknik histogram equalization bertujuan untuk menghasilkan suatu citra keluaran yang memiliki nilai histogram yang relatif sama. Contoh perintah untuk melakukan histogram equalization: I=imread('rice.tif'); J=histeq(I); figure,imshow(I); figure,imhist(I); figure,imshow(J); figure,imhist(J);
1.3 Thresholding Thresholding merupakan proses pemisahan piksel-piksel berdasarkan derajat keabuan yang dimilikinya. Piksel yang memiliki derajat keabuan lebih kecil dari nilai batas yang ditentukan akan diberikan nilai 0, sementara piksel yang memiliki derajat keabuan yang lebih besar dari batas akan diubah menjadi bernilai 1 Contoh perintah untuk melakukan thresholding: I=imread('rice.tif'); J=im2bw(I,0.4); K=im2bw(I,0.5); figure,imshow(I); figure,imhist(I); figure,imshow(J); figure,imshow(K);
2.
Operasi Spasial
Operasi spasial dalam pengolahan citra digital dilakukan melalui penggunaan suatu kernel konvolusi 2-dimensi. Beberapa metode image enhancement yang termasuk dalam keluarga ini adalah neighborhood averaging, median filtering, dan high-pass filtering. 2.1 Neighborhood Averaging Pada prinsipnya, filter yang digunakan dalam neighborhood averaging merupakan salah satu jenis low-pass filter, yang bekerja dengan cara mengganti nilai suatu piksel pada citra asal dengan nilai rata-rata dari piksel tersebut dan lingkungan tetangganya. Contoh perintah untuk melakukan neighborhood averaging dengan kernel berukuran 3×3: I=imread('eight.tif'); kernel=[1 1 1;1 1 1;1 1 1]/9; J=uint8(conv2(double(I),kernel,'same')); figure,imshow(I); figure,imshow(J);
Salah satu persoalan pada penggunaan neighborhood averaging adalah apabila citra masukan telah terkontaminasi noise, sebagaimana diperlihatkan oleh contoh berikut: I=imread('eight.tif'); IN=imnoise(I,'salt & pepper',0.02); kernel=[1 1 1;1 1 1;1 1 1]/9; J=uint8(conv2(double(I),kernel,'same')); JN=uint8(conv2(double(IN),kernel,'same')); figure,imshow(I); figure,imshow(J); figure,imshow(IN); figure,imshow(JN);
2.2 Median Filtering Median filter merupakan salah satu jenis low-pass filter, yang bekerja dengan mengganti nilai suatu piksel pada citra asal dengan nilai median dari piksel tersebut dan lingkungan tetangganya. Dibandingkan dengan neighborhood averaging, filter ini lebih tidak sensitif terhadap perbedaan intensitas yang ekstrim.
Contoh perintah untuk melakukan median filtering dengan kernel berukuran 3×3: I=imread('eight.tif'); IN=imnoise(I,'salt & pepper',0.02); J=medfilt2(I,[3 3]); JN=medfilt2(IN,[3 3]); figure,imshow(I); figure,imshow(J); figure,imshow(IN); figure,imshow(JN);
2.3 High-pass Filtering Sebagaimana pada proses pengolahan sinyal satu dimensi, high-pass filter dua dimensi akan melewatkan komponen citra frekuensi tinggi dan meredam komponen citra frekuensi rendah. Contoh perintah untuk melakukan high-pass filtering: I=imread('saturn.tif'); hpf1=[ 1 -2 1;-2 5 -2; 1 -2 1]; hpf2=[ 0 -1 0;-1 5 -1; 0 -1 0]; hpf3=[-1 -1 -1;-1 9 -1;-1 -1 -1]; J1=uint8(conv2(double(I),hpf1,'same')); J2=uint8(conv2(double(I),hpf2,'same')); J3=uint8(conv2(double(I),hpf3,'same')); figure,imshow(I); figure,imshow(J1); figure,imshow(J2); figure,imshow(J3);
3.
Operasi Transformasi
Berbeda dengan beberapa metode yang telah dibahas sebelumnya, proses image enhancement berbasis transformasi citra dilakukan dengan: a. mentransformasi citra asal ke dalam domain yang sesuai bagi proses enhancement b. melakukan proses enhancement pada domain tersebut c. mengembalikan citra ke dalam domain spasial untuk ditampilkan/diproses lebih lanjut Salah satu metode transformasi yang paling populer dalam aplikasi pengolahan citra digital adalah Fast Fourier Transform (FFT). Transformasi ini memindahkan informasi citra dari domain spasial ke dalam domain frekuensi, yaitu dengan merepresentasikan citra spasial sebagai suatu penjumlahan eksponensial kompleks dari beragam frekuensi, magnituda, dan fasa. Contoh perintah untuk melakukan low-pass filtering dan high-pass filtering melalui FFT adalah: I=double(imread('cameraman.tif')); IF=fft2(I); mask_high=double(imread('maskpojok.bmp')); mask_low =not(mask_high); IFH=IF.*mask_high; IFL=IF.*mask_low; hasil_high=abs(ifft2(IFH)); hasil_low =abs(ifft2(IFL)); figure,imagesc(I) ,colormap('gray'),colorbar('vert'); figure,imagesc(log(abs(IF )),[0 17]),colormap('hot') ,colorbar('vert'); figure,imagesc(hasil_high) ,colormap('gray'),colorbar('vert'); figure,imagesc(log(abs(IFH)),[0 17]),colormap('hot') ,colorbar('vert'); figure,imagesc(hasil_low ) ,colormap('gray'),colorbar('vert'); figure,imagesc(log(abs(IFL)),[0 17]),colormap('hot') ,colorbar('vert');
Pada contoh program tersebut, proses filtering dilakukan melalui masking terhadap komponen frekuensi yang ditentukan. Agar tercipta karakteristik high-pass filter, maka komponen frekuensi rendah (koefisien frekuensi yang berada pada bagian pojok dari citra hasil FFT) dimasking menggunakan nilai 0. Demikian pula sebaliknya untuk memunculkan sifat low-pass filter, komponen frekuensi tinggi (koefisien frekuensi yang berada pada bagian tengah dari citra hasil FFT) dibuat menjadi 0 melalui perkalian dengan mask low-pass.
4.
Tugas
4.1 Operasi Titik Jalankan program berikut dan berikan analisis atas hasil yang diperoleh. I=imread('rice.tif'); J=imadjust(I,[0.15 0.9],[1 0]); figure,imshow(I); figure,imhist(I); figure,imshow(J); figure,imhist(J);
4.2 Operasi Spasial Jalankan program berikut dan berikan analisis atas hasil yang diperoleh. I=double(imread('saturn.tif')); kernel=[1 1 1;1 1 1;1 1 1]/9; J=conv2(I,kernel,'same'); K=I-J; figure,imagesc(I),colormap('gray'),colorbar('vert'); figure,imagesc(J),colormap('gray'),colorbar('vert'); figure,imagesc(K),colormap('gray'),colorbar('vert');
4.3 Operasi Transformasi Jalankan program berikut dan berikan analisis atas hasil yang diperoleh.¹ I =double(imread('segitiga.bmp')); F =fft2(I); IF=fliplr(I); FF=fft2(IF); IR=rot90 (I); FR=fft2(IR); figure,imagesc(I) ,colormap('gray'),colorbar('vert'); figure,imagesc(log(abs(F ))),colormap('hot') ,colorbar('vert'); figure,imagesc(IF) ,colormap('gray'),colorbar('vert'); figure,imagesc(log(abs(FF))),colormap('hot') ,colorbar('vert'); figure,imagesc(IR) ,colormap('gray'),colorbar('vert'); figure,imagesc(log(abs(FR))),colormap('hot') ,colorbar('vert');
1
Gunakan file segitiga.bmp sebagai citra masukan. Pastikan terlebih dahulu agar file tersebut telah berada di direktori kerja yang digunakan.
4.4 Studi Kasus
Image enhancement pada dasarnya merupakan proses yang bersifat subjektif, sehingga parameter keberhasilannya bersifat subjektif pula. Namun demikian, apabila terdapat suatu citra referensi yang merepresentasikan hasil image enhancement yang diinginkan, pengukuran keberhasilan prosedur image enhancement secara kuantitatif dapat dilakukan melalui perhitungan nilai parameter MSE (Mean Squared Error) sebagai berikut:
a.
Untuk citra yang terkontaminasi noise berikut (gambar sebelah kiri), lakukanlah prosedur image enhancement sehingga diperoleh hasil yang paling mendekati citra referensi yang bersesuaian (gambar sebelah kanan).² Prosedur image enhancement yang dilakukan dapat merupakan gabungan dari berbagai metode image enhancement yang tersedia.
b.
Tuliskan langkah-langkah image enhancement yang memberikan hasil optimum dalam sebuah m-file. Prosedur image enhancement yang optimum adalah langkahlangkah yang menghasilkan keluaran citra dengan nilai MSE paling rendah. Perhitungan nilai MSE antara citra yang telah di-enhance dengan citra referensi dapat dilakukan melalui program hitungmse.m berikut: %HITUNGMSE – Fungsi untuk menghitung nilai MSE dari dua buah citra function NilaiMSE=HitungMSE(Citra,CitraRef); [tinggi,lebar]=size(Citra); NilaiSSE=sum(sum((double(Citra)-double(CitraRef)).^2)); NilaiMSE=NilaiSSE/tinggi/lebar;
Dengan menggunakan file hitungmse.m tersebut, perhitungan MSE dapat dilakukan dari command window melalui perintah: nilai_mse = hitungmse(citra,citra_referensi)
2
Gunakan file taz_noise.bmp sebagai citra yang akan di-enhance, dan file taz_ref.bmp sebagai citra referensi.
PRAKTIKUM
EL4027 PENGOLAHAN CITRA BIOMEDIKA EB7031 PENGOLAHAN CITRA BIOMEDIKA LANJUT MODUL 2 − EDGE DETECTION
Penentuan tepian suatu objek dalam citra merupakan salah satu wilayah pengolahan citra digital yang paling awal dan paling banyak diteliti. Proses ini seringkali ditempatkan sebagai langkah pertama dalam aplikasi segmentasi citra, yang bertujuan untuk mengenali objekobjek yang terdapat dalam citra ataupun konteks citra secara keseluruhan. Deteksi tepi berfungsi untuk mengidentifikasi garis batas (boundary) dari suatu objek yang terdapat pada citra. Tepian dapat dipandang sebagai lokasi piksel dimana terdapat nilai perbedaan intensitas citra secara ekstrem. Sebuah edge detector bekerja dengan cara mengidentifikasi dan menonjolkan lokasi-lokasi piksel yang memiliki karakteristik tersebut.
1. Operator Gradien Pada citra digital f(x,y), turunan berarah sepanjang tepian objek akan bernilai maksimum pada arah normal dari kontur tepian yang bersesuaian. Sifat ini dipergunakan sebagai dasar pemanfaatan operator gradien sebagai edge detector. Operator gradien citra konvensional melakukan diferensiasi intensitas piksel pada arah baris dan kolom, mengikuti persamaan local intensity variation berikut :
Nilai magnitudo gradien |∇(x,y)| dari persamaan di atas dapat dinyatakan sebagai berikut:
Operator gradien dapat direpresentasikan oleh dua buah kernel konvolusi Gx dan Gy, yang masing-masing mendefinisikan operasi penghitungan gradien dalam arah sumbu x dan sumbu y yang saling tegak lurus. Dalam kasus penghitungan gradien dengan persamaan local intensity variation, maka kernel Gx dan Gy dapat dirumuskan seperti berikut:
Berikut adalah contoh fungsi Matlab untuk operasi penghitungan gradien orde satu:
I = double(imread('cameraman.tif')); % Gradien orde satu pada arah horizontal gx = [-1 1]; Ix = conv2(I,gx,'same'); % Gradien orde satu pada arah vertikal gy = [-1;1]; Iy = conv2(I,gy,'same'); % Magnitudo gradien J = sqrt((Ix.^2)+(Iy.^2)); % Gambar hasil figure,imagesc(I ),colormap('gray'),colorbar('vert'); figure,imagesc(Ix),colormap('gray'),colorbar('vert'); figure,imagesc(Iy),colormap('gray'),colorbar('vert'); figure,imagesc(J ),colormap('gray'),colorbar('vert');
Dari operator gradien konvensional di atas, dapat diturunkan berbagai operator gradien berikut: 1.1 Operator Selisih Terpusat Operator selisih terpusat juga dikenal sebagai Centered Difference Edge Detector Mask, dan dinyatakan sebagai kernel:
I = double(imread('cameraman.tif')); % Konvolusi dengan operator selisih terpusat d1x = [-1 0 1]; d1y = [-1;0;1]; Ix = conv2(I,d1x,'same'); Iy = conv2(I,d1y,'same'); J = sqrt((Ix.^2)+(Iy.^2)); % Gambar Hasil figure,imagesc(I ),colormap('gray'),colorbar('vert'); figure,imagesc(Ix),colormap('gray'),colorbar('vert'); figure,imagesc(Iy),colormap('gray'),colorbar('vert'); figure,imagesc(J ),colormap('gray'),colorbar('vert');
1.2 Operator Roberts Operator Roberts memiliki ukuran kernel sebesar 2×2, yang direpresentasikan sebagai:
Contoh perintah menggunakan operator Roberts: I = double(imread('cameraman.tif'));
% Konvolusi dengan operator Roberts robertshor = [0 1; -1 0]; robertsver = [1 0; 0 -1]; Ix = conv2(I,robertshor,'same'); Iy = conv2(I,robertsver,'same'); J = sqrt((Ix.^2)+(Iy.^2)); % Gambar Hasil figure,imagesc(I ),colormap('gray'),colorbar('vert'); figure,imagesc(Ix),colormap('gray'),colorbar('vert'); figure,imagesc(Iy),colormap('gray'),colorbar('vert'); figure,imagesc(J ),colormap('gray'),colorbar('vert');
Contoh perintah menggunakan operator Roberts (matlab toolbox) 1: I = imread('cameraman.tif'); J = edge(I,'roberts'); figure,imagesc(I),colormap('gray'),colorbar('vert'); figure,imagesc(J),colormap('gray'),colorbar('vert');
1.3 Operator Prewitt
Contoh perintah menggunakan operator Prewitt: I = double(imread('cameraman.tif')); %Konvolusi dengan operator Prewitt prewitthor = [-1 0 1; -1 0 1; -1 0 1]; prewittver = [-1 -1 -1; 0 0 0; 1 1 1]; Ix = conv2(I,prewitthor,'same'); Iy = conv2(I,prewittver,'same'); J = sqrt((Ix.^2)+(Iy.^2)); %Gambar Hasil figure,imagesc(I ),colormap('gray'),colorbar('vert'); figure,imagesc(Ix),colormap('gray'),colorbar('vert'); figure,imagesc(Iy),colormap('gray'),colorbar('vert'); figure,imagesc(J ),colormap('gray'),colorbar('vert');
Contoh perintah menggunakan operator Prewitt (matlab toolbox): I = imread('cameraman.tif'); J = edge(I,'prewitt'); figure,imagesc(I),colormap('gray'),colorbar('vert'); figure,imagesc(J),colormap('gray'),colorbar('vert'); 1
Perintah edge (pada toolbox Matlab) untuk mensimulasikan operator Prewitt, Roberts, Sobel, dan lainnya memiliki konsep dasar yang sama dengan operasi konvolusi kernel setiap operator secara manual. Bedanya, perintah edge menambahkan suatu skema thresholding secara otomatis, sehingga dihasilkan citra keluaran yang bersifat biner (bernilai 0 atau 1).
1.4 Operator Sobel
Contoh perintah menggunakan operator Sobel: I = double(imread('cameraman.tif')); %Konvolusi dengan operator Sobel sobelhor = [-1 0 1; -2 0 2; -1 0 1]; sobelver = [-1 -2 -1; 0 0 0; 1 2 1]; Ix = conv2(I,sobelhor,'same'); Iy = conv2(I,sobelver,'same'); J = sqrt((Ix.^2)+(Iy.^2)); %Gambar Hasil figure,imagesc(I ),colormap('gray'),colorbar('vert'); figure,imagesc(Ix),colormap('gray'),colorbar('vert'); figure,imagesc(Iy),colormap('gray'),colorbar('vert'); figure,imagesc(J ),colormap('gray'),colorbar('vert');
Contoh perintah menggunakan operator Sobel (matlab toolbox): I = imread('cameraman.tif'); J = edge(I,'sobel'); figure,imagesc(I),colormap('gray'),colorbar('vert'); figure,imagesc(J),colormap('gray'),colorbar('vert');
1.5 Operator Isotropic
Contoh perintah menggunakan operator isotropic: I = double(imread('cameraman.tif')); %Konvolusi dengan operator isotropic isohor = [-1 0 1; -sqrt(2) 0 sqrt(2); -1 0 1]; isover = [-1 -sqrt(2) -1; 0 0 0; 1 sqrt(2) 1]; Ix = conv2(I,isohor,'same'); Iy = conv2(I,isover,'same'); J = sqrt((Ix.^2)+(Iy.^2)); %Gambar Hasil figure,imagesc(I ),colormap('gray'),colorbar('vert'); figure,imagesc(Ix),colormap('gray'),colorbar('vert'); figure,imagesc(Iy),colormap('gray'),colorbar('vert'); figure,imagesc(J ),colormap('gray'),colorbar('vert');
1.6 Operator Compass
Operator Compass bekerja menggunakan pola empat arah mata angin:
1.7 Operator Kirsch Operator Kirsch bekerja menggunakan pola delapan arah mata angin:
2. Operator Laplacian Dalam kondisi transisi tepian yang lebih tidak ekstrem, penggunaan operator turunan kedua lebih dianjurkan.
Representasi turunan kedua dalam bentuk kernel operator Laplacian adalah sebagai berikut:
Dengan berbagai macam pembobotan, kernel Laplacian tersebut dapat dimodifikasi menjadi beberapa kernel konvolusi berikut :
Laplacian of Gaussian Turunan kedua memiliki sifat lebih sensitif terhadap noise, selain itu juga menghasilkan double edge. Oleh karena itu, operator Laplacian dalam deteksi tepi pada umumnya tidak dipergunakan secara langsung, namun dikombinasikan dengan suatu kernel Gaussian menjadi sebuah operator Laplacian of Gaussian. Fungsi transfer dari kernel Laplacian of Gaussian dapat dirumuskan sebagai berikut:
dimana σ merupakan standar deviasi dari kernel Gaussian. Contoh perintah untuk mendeteksi tepian dengan menggunakan operator Laplacian of Gaussian (Matlab toolbox) adalah: I = imread('cameraman.tif'); J = edge(I,'log'); figure,imagesc(I),colormap('gray'),colorbar('vert'); figure,imagesc(J),colormap('gray'),colorbar('vert');
3. Operator Canny Salah satu algoritma deteksi tepi modern adalah deteksi tepi dengan menggunakan metoda Canny. Berikut adalah diagram blok algoritma Canny :
Contoh perintah untuk mendeteksi tepian menggunakan metoda Canny (matlab toolbox) : I = imread('cameraman.tif'); J = edge(I,'canny'); figure,imagesc(I),colormap('gray'),colorbar('vert'); figure,imagesc(J),colormap('gray'),colorbar('vert');
4. Tugas A. Pada percobaan dengan kernel Sobel, kernel Isotropic, dan kernel Compass, jelaskan mengenai garis abu-abu yang ditemukan pada pinggiran citra hasil konvolusi. B. Untuk citra cameraman.tif, edge detector manakah yang memberikan hasil terbaik? Berikan analisis mengenai pendapat anda tersebut. C. Sensitivitas edge detector terhadap noise dapat diukur dengan menggunakan parameter error rate sebagai berikut:
dimana: nR : jumlah piksel yang dinyatakan sebagai edge pada citra referensi nN : jumlah piksel yang dinyatakan sebagai edge pada citra noisy Nilai P yang besar menyatakan sensitivitas edge detector yang tinggi terhadap noise. Dalam kasus citra cameraman.tif terkontaminasi oleh noise salt & pepper dengan distribusi 0.02, edge detector manakah yang paling terpengaruh performansinya? Simulasikan dengan Matlab. %Contoh perhitungan error rate untuk edge detector Sobel %(kerjakan pula untuk jenis edge detector lain pada Matlab toolbox) Ia = imread('cameraman.tif'); Ja = edge(Ia,'sobel'); In = imnoise(Ia,'salt & pepper',0.02); Jn = edge(In,'sobel'); nr = sum(sum(Ja)); nn = sum(sum(Jn)); P = abs(nn-nr)/nr
PRAKTIKUM EL4027 PENGOLAHAN CITRA BIOMEDIKA EB7031 PENGOLAHAN CITRA BIOMEDIKA LANJUT MODUL 3 − ANALISIS TEKSTUR
Tekstur merupakan karakteristik intrinsik dari suatu citra yang terkait dengan tingkat kekasaran (roughness), granularitas (granulation), dan keteraturan (regularity) susunan struktural piksel. Aspek tekstural dari sebuah citra dapat dimanfaatkan sebagai dasar dari segmentasi, klasifikasi, maupun interpretasi citra. Tekstur dapat didefinisikan sebagai fungsi dari variasi spasial intensitas piksel (nilai keabuan) dalam citra. Berdasarkan strukturnya, tekstur dapat diklasifikasikan dalam dua golongan : • Makrostruktur Tekstur makrostruktur memiliki perulangan pola lokal secara periodik pada suatu daerah citra, biasanya terdapat pada pola-pola buatan manusia dan cenderung mudah untuk direpresentasikan secara matematis. •
Mikrostruktur Pada tekstur mikrostruktur, pola-pola lokal dan perulangan tidak terjadi begitu jelas, sehingga tidak mudah untuk memberikan definisi tekstur yang komprehensif.
Contoh gambar berikut ini menunjukkan perbedaan tekstur makrostruktur dan mikrostruktur yang diambil dari album tekstur Brodatz.
Gambar 1 Contoh tekstur visual dari Album Tekstur Brodatz . Atas: makrostruktur Bawah: mikrostruktur
Analisis tekstur bekerja dengan mengamati pola ketetanggaan antar piksel dalam domain spasial. Dua persoalan yang seringkali berkaitan dengan analisis tekstur adalah: •
Ekstraksi ciri Ekstraksi ciri merupakan langkah awal dalam melakukan klasifikasi dan interpretasi citra. Proses ini berkaitan dengan kuantisasi karakteristik citra ke dalam sekelompok nilai ciri yang sesuai. Dalam praktikum ini kita akan mengamati metoda ekstraksi ciri statistik orde pertama dan kedua, serta mengenali performansi masing-masing skema dalam mengenali citra dengan karakteristik tekstural yang berlainan.
•
Segmentasi citra Segmentasi citra merupakan proses yang bertujuan untuk memisahkan suatu daerah pada citra dengan daerah lainnya. Berbeda dengan pada citra non-tekstural, segmentasi citra tekstural tidak dapat didasarkan pada intensitas piksel per piksel, tetapi perlu mempertimbangkan perulangan pola dalam suatu wilayah ketetanggaan lokal. Dalam praktikum ini kita akan mencoba menerapkan filter Gabor untuk melakukan segmentasi citra tekstural berdasarkan perulangan pola lokal pada orientasi dan frekuensi tertentu.
1. Ekstraksi Ciri Statistik Analisis tekstur lazim dimanfaatkan sebagai proses antara untuk melakukan klasifikasi dan interpretasi citra. Suatu proses klasifikasi citra berbasis analisis tekstur pada umumnya membutuhkan tahapan ekstraksi ciri, yang dapat terbagi dalam tiga macam metode berikut: •
Metode statistik Metode statistik menggunakan perhitungan statistik distribusi derajat keabuan (histogram) dengan mengukur tingkat kekontrasan, granularitas, dan kekasaran suatu daerah dari hubungan ketetanggaan antar piksel di dalam citra. Paradigma statistik ini penggunaannya tidak terbatas, sehingga sesuai untuk teksturtekstur alami yang tidak terstruktur dari sub pola dan himpunan aturan (mikrostruktur).
•
Metode spektral Metode spektral berdasarkan pada fungsi autokorelasi suatu daerah atau power distribution pada domain transformasi Fourier dalam mendeteksi periodisitas tekstur.
•
Metode struktural Analisis dengan metode ini menggunakan deskripsi primitif tekstur dan aturan sintaktik. Metode struktural banyak digunakan untuk pola-pola makrostruktur.
Bagian ini akan membahas metode ekstraksi ciri statistik orde pertama dan kedua. Ekstraksi ciri orde pertama dilakukan melalui histogram citra. Ekstraksi ciri statistik orde kedua dilakukan dengan matriks kookurensi, yaitu suatu matriks antara yang merepresentasikan hubungan ketetanggaan antar piksel dalam citra pada berbagai arah orientasi dan jarak spasial.
Gambar 2 Ilustrasi ekstraksi ciri statistik Kiri : Histogram citra sebagai fungsi probabilitas kemunculan nilai intensitas pada citra Kanan : Hubungan ketetanggaan antar piksel sebagai fungsi orientasi dan jarak spasial
1.1 Ekstraksi ciri orde pertama Ekstraksi ciri orde pertama merupakan metode pengambilan ciri yang didasarkan pada karakteristik histogram citra. Histogram menunjukkan probabilitas kemunculan nilai derajat keabuan piksel pada suatu citra. Dari nilai-nilai pada histogram yang dihasilkan, dapat dihitung beberapa parameter ciri orde pertama, antara lain adalah mean, skewness, variance, kurtosis, dan entropy. a.
Mean (μ) Menunjukkan ukuran dispersi dari suatu citra
dimana fn merupakan suatu nilai intensitas keabuan, sementara p(fn) menunjukkan nilai histogramnya (probabilitas kemunculan intensitas tersebut pada citra).
b.
Variance (σ2) Menunjukkan variasi elemen pada histogram dari suatu citra
c.
Skewness (α3) Menunjukkan tingkat kemencengan relatif kurva histogram dari suatu citra
d.
Kurtosis (α4) Menunjukkan tingkat keruncingan relatif kurva histogram dari suatu citra
e.
Entropy (H) Menunjukkan ukuran ketidakaturan bentuk dari suatu citra
Berikut adalah fungsi ciriordesatu yang dipergunakan untuk menghitung ciri orde satu dari citra:
1.2 Ekstraksi ciri orde kedua
Pada beberapa kasus, ciri orde pertama tidak lagi dapat digunakan untuk mengenali perbedaan antar citra. Pada kasus seperti ini, kita membutuhkan pengambilan ciri statistik orde dua. Salah satu teknik untuk memperoleh ciri statistik orde dua adalah dengan menghitung probabilitas hubungan ketetanggaan antara dua piksel pada jarak dan orientasi sudut tertentu. Pendekatan ini bekerja dengan membentuk sebuah matriks kookurensi dari data citra, dilanjutkan dengan menentukan ciri sebagai fungsi dari matriks antara tersebut. Kookurensi berarti kejadian bersama, yaitu jumlah kejadian satu level nilai piksel bertetangga dengan satu level nilai piksel lain dalam jarak (d) dan orientasi sudut (θ) tertentu. Jarak dinyatakan dalam piksel dan orientasi dinyatakan dalam derajat. Orientasi dibentuk dalam empat arah sudut dengan interval sudut 45°, yaitu 0°, 45°, 90°, dan 135°. Sedangkan jarak antar piksel biasanya ditetapkan sebesar 1 piksel. Matriks kookurensi merupakan matriks bujursangkar dengan jumlah elemen sebanyak kuadrat jumlah level intensitas piksel pada citra. Setiap titik (p,q) pada matriks kookurensi berorientasi θ berisi peluang kejadian piksel bernilai p bertetangga dengan piksel bernilai q pada jarak d serta orientasi θ dan (180−θ).
Gambar 3 Ilustrasi pembuatan matriks kookurensi
(a) Citra masukan (b) Nilai intensitas citra masukan (c) Hasil matriks kookurensi 0°
(d) Hasil matriks kookurensi 45° (e) Hasil matriks kookurensi 90° (f) Hasil matriks kookurensi 135°
Setelah memperoleh matriks kookurensi tersebut, kita dapat menghitung ciri statistik orde dua yang merepresentasikan citra yang diamati. Haralick et al mengusulkan berbagai jenis ciri tekstural yang dapat diekstraksi dari matriks kookurensi. Dalam modul ini dicontohkan
perhitungan 6 ciri statistik orde dua, yaitu Angular Second Moment, Contrast, Correlation, Variance, Inverse Difference Moment, dan Entropy. a. Angular Second Moment Menunjukkan ukuran sifat homogenitas citra.
dimana p(i,j) merupakan menyatakan nilai pada baris i dan kolom j pada matriks kookurensi. b.
Contrast Menunjukkan ukuran penyebaran (momen inersia) elemen-elemen matriks citra. Jika letaknya jauh dari diagonal utama, nilai kekontrasan besar. Secara visual, nilai kekontrasan adalah ukuran variasi antar derajat keabuan suatu daerah citra.
c.
Correlation Menunjukkan ukuran ketergantungan linear derajat keabuan citra sehingga dapat memberikan petunjuk adanya struktur linear dalam citra.
d.
Variance Menunjukkan variasi elemen-elemen matriks kookurensi. Citra dengan transisi derajat keabuan kecil akan memiliki variansi yang kecil pula.
e.
Inverse Different Moment Menunjukkan kehomogenan citra yang berderajat keabuan sejenis. Citra homogen akan memiliki harga IDM yang besar.
f.
Entropy Menunjukkan ukuran ketidakteraturan bentuk. Harga ENT besar untuk citra dengan transisi derajat keabuan merata dan bernilai kecil jika struktur citra tidak teratur (bervariasi).
Berikut adalah fungsi ciriordedua yang dipergunakan untuk menghitung ciri orde dua dari citra:
2. Filter Gabor Kemampuan sistem visual manusia dalam membedakan berbagai tekstur didasarkan atas kapabilitas dalam mengidentifikasikan berbagai frekuensi dan orientasi spasial dari tekstur yang diamati. Filter Gabor merupakan salah satu filter yang mampu mensimulasikan karakteristik sistem visual manusia dalam mengisolasi frekuensi dan orientasi tertentu dari citra. Karakteristik ini membuat filter Gabor sesuai untuk aplikasi pengenalan tekstur dalam computer vision. Secara spasial, sebuah fungsi Gabor merupakan sinusoida yang dimodulasi oleh fungsi Gauss. Respon impuls sebuah filter Gabor kompleks dua dimensi adalah :
dan dapat digambarkan sebagai berikut:
Gambar 4 Respon impuls filter Gabor dua dimensi.
Dalam domain frekuensi spasial, filter Gabor dapat direpresentasikan sebagai berikut:
Dalam domain frekuensi spasial, parameter-parameter filter Gabor dapat digambarkan sebagai:
Gambar 5 Parameter filter Gabor dalam domain frekuensi spasial
Tabel 1 Enam parameter filter Gabor
Ada enam parameter yang harus ditetapkan dalam implementasi filter Gabor. Keenam parameter tersebut adalah: F, θ, σx, σy, BF, dan Bθ •
Frekuensi (F) dan orientasi (θ) mendefinisikan lokasi pusat filter.
•
BF dan Bθ menyatakan konstanta lebar pita frekuensi dan jangkauan angular filter.
•
Variabel σx berkaitan dengan respon sebesar -6 dB untuk komponen frekuensi spasial.
•
Variabel σy berkaitan dengan respon sebesar -6dB untuk komponen angular.
•
•
•
Posisi (F, θ) dan lebar pita (σx, σy) dari filter Gabor dalam domain frekuensi harus ditetapkan dengan cermat agar dapat menangkap informasi tekstural dengan benar. Frekuensi tengah dari filter kanal harus terletak dekat dengan frekuensi karakteristik tekstur. Setelah mendapatkan ciri Gabor maka dapat dilakukan ekstraksi ciri. Salah satu ciri yang dapat dipilih adalah ciri energi, yang didefinisikan sebagai:
Dalam modul ini digunakan lebar pita frekuensi (BF), dan jarak frekuensi tengah (S F) sebesar satu oktaf, serta lebar pita angular (B θ) dan jarak angular (Sθ) sebesar 30° dan 45°. Pemilihan lebar pita angular sebesar 30° dan 45° adalah karena nilai ini dianggap mendekati karakteristik sistem visual manusia.
Berikut adalah fungsi gb dan en yang dipergunakan dalam percobaan filter Gabor ini:
2.1 Percobaan Mencari Frekuensi dan Orientasi Dominan pada Suatu Tekstur Berdasarkan program filter Gabor di atas, lakukan urutan kerja sebagai berikut : a.
Baca file citra straw.tif
1 2
b. Lakukan proses filtering terhadap citra tersebut menggunakan filter Gabor. Parameter frekuensi dan orientasi ditentukan sebagai:
c. Tampilan hasil citra yang telah difilter figure, imshow(G1);
d. Hitung nilai energi citra yang telah difilter E1 = en(G1);
e. Ulangi proses filtering pada beberapa frekuensi dan orientasi yang berbeda. Ambil kesimpulan berdasarkan besarnya energi dari masing-masing citra keluaran.
2.2 Percobaan Segmentasi Citra Tekstural Berdasarkan program filter Gabor di atas, lakukan urutan kerja sebagai berikut : a. Baca file citra tex2.tif
b. Analisis tekstur mengunakan filter Gabor
c. Lakukan thresholding dengan nilai threshold = 0.25.
d. Lakukan median filtering sebanyak dua kali
3. Tugas 3.1 Ciri Orde Satu a. Tampilkan histogram citra Taz1.bmp, Taz2.bmp, dan Taz3.bmp. Selanjutnya jalankan fungsi ciriordesatu terhadap masing-masing citra. b. Lakukan hal yang sama terhadap masing-masing citra Tekstur1.bmp, Tekstur2.bmp, dan Tekstur3.bmp. c. Berikan analisis mengenai proses yang telah dilakukan. 3.2 Ciri Orde Dua a. Jalankan fungsi ciriordedua terhadap masukan citra Tekstur1.bmp, Tekstur2.bmp, dan Tekstur3.bmp. b. Berikan analisis mengenai proses yang telah dilakukan. 3.3 Filter Gabor a. Jalankan program mikimos.m berikut. Berikan analisis mengenai proses yang dijalankan.
b. Gantilah baris perintah: X=imread('mikimos1.bmp'); dengan: X=imread('mikimos2.bmp');
Bandingkan hasilnya dengan (a), berikan analisis mengenai hasil yang diperoleh. c. Gantilah baris perintah: Y=gb(X,2,3);
dengan: Y=gb(X,2,6);
Bandingkan hasilnya dengan (a), berikan analisis mengenai hasil yang diperoleh.
PRAKTIKUM EL4027 PENGOLAHAN CITRA BIOMEDIKA EB7031 PENGOLAHAN CITRA BIOMEDIKA LANJUT MODUL 4 − 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 transaksial 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), -∞
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π ξ0s) /π2(d2 – 4s2) • Respons impuls diskrit : h(m) =2 / (π2d(1 - 4m2))
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 (x2 +y2)−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 1 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.
• • • •
• d.
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
e.
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 a. Perhatikan Gambar 3. Berikan analisis mengapa transformasi Radon menghasilkan gambar seperti di sebelah kanan untuk citra phantom yang diberikan di sebelah kiri! b. Lakukan kembali transformasi Radon untuk berbagai jangkauan nilai theta dan berbagai nilai step theta, kemudian bandingkan hasilnya. Jelaskan dimana dan bagaimana terdapat perbedaan! c. 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 a. Jelaskan parameter-parameter dari fungsi iradon dalam Matlab! b. Lakukan transformasi Radon invers terhadap citra phantom dengan input hasil transformasi Radon dalam berbagai jangkauan theta yang berbeda. Analisis hasilnya! c. 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.
2
PRAKTIKUM EL4027 PENGOLAHAN CITRA BIOMEDIKA EB7031 PENGOLAHAN CITRA BIOMEDIKA LANJUT MODUL 5 − TRANSFORMASI CITRA IMAGE TRANSFORM Ada banyak metode yang digunakan untuk melakukan tranformasi image. Dua diantaranya adalah DFT (Discrete Fourier Transform) dan DCT (Discrete Cosinus Transform. Dengan menggunakan MatLab hal ini dapat kita lakukan dengan mudah. DFT (Discrete Fourier Transform) Transformasi Fourier merupakan representasi dari sebuah citra sebagai penjumlahan eksponensial kompleks dari beragam magnitude, frekuensi dan fasa. Contoh perintah untuk melakukan transformasi Fourier pada suatu citra: f=checkerboard(8,8); f(5:24,13:17)=1; imshow(f,'notruesize') f=fft2(f); f2=log(abs(f)); figure,imshow(f2,[-1,5],'notruesize'); colormap(jet); colorbar
Output program diberikan sebagai berikut:
Gambar 1. Sebelum dikenai FFT
Gambar 2. Setelah dikenai FFT
DCT (Discrete Cosine Transform) Discrete Cosine Transform merepresentasikan sebuah citra dari penjumlahan sinusoida dari magnitude dan frekuensi yang berubah-ubah. Sifat dari DCT adalah mengubah informasi citra yang signifkan dikonsentrasikan hanya pada beberapa koefisien DCT. Oleh karena itu DCT sering digunakan untuk kompresi citra seperti pada JPEG. Contoh perintah untuk melakukan DCT pada suatu citra: i=checkerboard(8,8); j=dct2(i); imshow(i); figure,imshow(log(abs(j)),[]),colormap(jet(64)), colorbar j(abs(j)<10)=0; k=idct2(j)/255;
figure,imshow(k);
Output program diatas adalah sebagai berikut:
Transformasi Wavelet Wavelet merupakan klas dari suatu fungsi yang digunakan untuk melokalisasi suatu fungsi dalam ruang dan skala. Wavelet dapat dibentuk dari satu fungsi Ψ(x) , dikenal sebagai “mother wavelet” dalam suatu interval berhingga. “Daughter wavelet” Ψa,b (x) dibentuk oleh translasi (b) dan kontraksi (a). Wavelet sangat berguna untuk kompresi data image, karena transformasi wavelet mempunyai property yang ada dalam beberapa cara terhadap transformasi Fourier konvensional. Wavelet tunggal dapat didefinisikan dengan:
Selanjutnya
dan formula Calderon memberikan :
Wavelet yang sering digunakan didefinisikan dengan fungsi Haar sebagai berikut :
dan dengan j integer dan nonnegative, 0≤k≤2j-1 Transformasi dengan Image Processing Toolbox pada MatLab. Berikut ini diberikan contoh transformasi dengan Toolbox pada MatLab. Output memberikan diberikan gambar asli dan enam gambar hasil transformasi dengan enam metode yang berbeda yakni transformasi linear conformal, transformasi affine, transformasi proyeksi, transformasi polinomial, transformasi piecewise linear, transformasi barel, transformasi sinusoid dan transformasi pin cushion. clc;clear;
gbr = checkerboard(6,6); figure subplot(331) imshow(gbr) title('gb asli') % Transformasi linear dg 4 parameter skala = 1.2; % faktor skala sudut = 40*pi/180; % sudut rotasi tx = 0; % translasi x ty = 0; % translasi y sc = skala*cos(sudut); ss = skala*sin(sudut); T = [ sc -ss; ss sc; tx ty]; t_lc = maketform('affine',T); gbr_lc = imtransform(gbr,t_lc,'FillValues',.3); subplot(332) imshow(gbr_lc); title('transf. linear conformal') % Matrik T dg 6 elemen T = [1 0.1; 1 1; 0 0]; t_aff = maketform('affine',T); I_affine = imtransform(gbr,t_aff,'FillValues',.3); subplot(333) imshow(I_affine) title('affine') T = [1 0 0.008; 1 1 0.01; 0 0 1]; t_proj = maketform('projective',T); I_projective = imtransform(gbr,t_proj,'FillValues',.3); subplot(334) imshow(I_projective) title('transf. proyeksi') xybase = reshape(randn(12,1),6,2); t_poly = cp2tform(xybase,xybase,'polynomial',2); % matrik T dg 12 elemen T = [0 0; 1 0; 0 1; 0.001 0; 0.02 0; 0.01 0]; t_poly.tdata = T; I_polynomial = imtransform(gbr,t_poly,'FillValues',.3); subplot(335) imshow(I_polynomial) title('transf. polinomial') imid = round(size(gbr,2)/2); % mencari indeks tengah I_left = gbr(:,1:imid); stretch = 1.5; % faktor tarik size_right = [size(gbr,1) round(stretch*imid)]; I_right = gbr(:,imid+1:end); I_right_stretched = imresize(I_right,size_right); I_piecewiselinear = [I_left I_right_stretched]; subplot(336) imshow(I_piecewiselinear)
title('transf. piecewise linear') % variasi dg sinusoid [nrows,ncols] = size(gbr); [xi,yi] = meshgrid(1:ncols,1:nrows); a1 = 5; % amplitudo sinusoid a2 = 3; u = xi + a1*sin(pi*xi/imid); v = yi - a2*sin(pi*yi/imid); tmap_B = cat(3,u,v); resamp = makeresampler('linear','fill'); I_sinusoid = tformarray(gbr,[],resamp,[2 1],[1 2],[],tmap_B,.3); subplot(337) imshow(I_sinusoid) title('transf. sinusoid') % distorsi barel xt = xi(:) - imid; yt = yi(:) - imid; [theta,r] = cart2pol(xt,yt); a = .001; % amplitudo s = r + a*r.^3; [ut,vt] = pol2cart(theta,s); u = reshape(ut,size(xi)) + imid; v = reshape(vt,size(yi)) + imid; tmap_B = cat(3,u,v); I_barrel = tformarray(gbr,[],resamp,[2 1],[1 2],[],tmap_B,.3); subplot(338) imshow(I_barrel) title('transf. barrel') % distorsi radial pin cushion xt = xi(:) - imid; yt = yi(:) - imid; [theta,r] = cart2pol(xt,yt); a = -.0005; % amplitudo s = r + a*r.^3; [ut,vt] = pol2cart(theta,s); u = reshape(ut,size(xi)) + imid; v = reshape(vt,size(yi)) + imid; tmap_B = cat(3,u,v); I_pin = tformarray(gbr,[],resamp,[2 1],[1 2],[],tmap_B,.3); subplot(339) imshow(I_pin) title('transf. pin cushion')
Gambar 5. Output dengan enam model transformasi yang berbeda.
KOMPRESI CITRA GRAY SCALE DENGAN MODIFIKASI ALGORITMA KUANTISASI Krisnawati STMIK AMIKOM Yogyakarta Abstract Suatu file yang kapasitasnya besar dapat diperkecil dengan pemampatan (compression). Untuk file image metode kuantisasi bisa menjadi salah satu pilihan. Pada citra grayscale metode kuantisasi bekerja dengan mengurangi derajat keabuan, sehingga jumlah bit yang digunakan untuk merepresentasikan citra menjadi berkurang. Oleh karena jumlah bit berkurang maka ukuran file menjadi lebih kecil. Dari algoritma kuantisasi yang ada, dilakukan modifikasi pada saat konversi derajat keabuan lama ke derajat keabuan yang baru, dengan mengabaikan unsur penyebaran derajat keabuan. Hasilnya kapasitas citra hasil kompresi tidak berubah. Histogram citra sangat berpengaruh terhadap hasil akhir. Keyword : kompresi, citra grayscale, kuantisasi. 1.
Pendahuluan
Perkembangan media penyimpan berkapasitas besar mengakibatkan orang tidak lagi menemui masalah jika mempunyai file dengan ukuran yang besar. Lebih-lebih jika file yang kita punya merupakan file image. Walaupun demikian, adakalanya ukuran file yang besar tersebut terasa mengganggu jika kita harus memanage media penyimpan yang kita punya untuk bermacammacam data. Apalagi jika file tersebut akan akan kita kirim secara elektronik, tentunya kapasitas file menjadi masalah tersendiri. Kompresi citra (image compression) adalah proses untuk meminimalkan jumlah bit yang merepresentasikan suatu citra sehingga ukuran citra menjadi lebih kecil. Pada dasarnya teknik kompresi citra digunakan untuk proses transmisi data (data transmission) dan penyimpanan
data (storage). Kompresi citra banyak diaplikasikan pada penyiaran televisi, penginderaan jarak jauh (remote sensing), komunikasi militer, radar dan lain-lain. Ada banyak metode kompresi data, salah satunya adalah metode kuantisasi. Metode kuantisasi bekerja dengan mengurangi jumlah intensitas warna, sehingga jumlah bit yang digunakan untuk merepresentasikan citra menjadi berkurang. Oleh karena jumlah bit berkurang maka ukuran file menjadi lebih kecil. Dengan berkurangnya intensitas warna tentu saja ada informasi yang hilang dari citra asal. Oleh karena itu metode ini termasuk dalam loossy compression. Oleh karena itu citra yang sudah dikompresi sulit didekompresi kembali karena adanya informasi yang hilang. Proses kompresi tentunya akan berdampak kepada banyak hal. Yang pertama adalah ukuran citra hasil kompresi. Ukuran citra diharapkan lebih kecil dari citra asal. Kedua adalah kualitas citra untuk input terhadap proses berikutnya. Sampai berapa persenkah citra asli bisa dikompresi ini tentunya tergantung pada banyak factor. Faktor inilah yang ingin diketahui pula dalam penelitian ini. Ada beberapa teknik kompresi yang dapat dikategorikan ke dalam Lossless maupun Lossy Compression, antara lain: 1. Kompresi berbasis Statistik (Lossless) Merepresentasikan citra dengan frekuensi kemunculan nilai intensitas tertentu. 2. Kompresi berbasis Kuantisasi (Lossy) Mengurangi jumlah intensitas warna. 3. Kompresi berbasis Transformasi (Lossless/Lossy) Mengoptimalkan kinerja kompresi berbasis statistik dan kuantisasi dengan cara melakukan transformasi terlebih dahulu sebelum menerapkan salah satu teknik tersebut. Sehingga kompresi bersifat lossy atau lossles tergantung teknik mana yang digunakan setelah transformasi apakah itu statistik (lossless) atau kuantisasi (lossy). 4. Kompresi berbasis fraktal (Lossy) Fraktal merupakan bentuk rekursif yang merepresentasikan komponen dasar objek. Dalam konsep kompresi, data direpresentasikan sebagai pasangan antar elemen fraktal, pola umum konfigurasi yang membentuk objek secara keseluruhan, dan koefisien transformasi spasial (affine) untuk masing-masing fraktal sesuai dengan posisinya dalam konfigurasi pembentuk objek. 1.1 Proses Kompresi 1.1.1
Kompresi Berbasis Kuantisasi
Kompresi berbasis kuantisasi menggunakan metode pengurangan jumlah intensitas warna, sehingga dapat mengurangi jumlah bit yang digunakan untuk merepresentasikan citra. Kompresi ini bersifat lossy, karena intensitas warna berkurang, sehingga kualitas gambar hasil kompresi menjadi kurang baik. Secara algoritma teknik kompresi ini diberikan sebagai berikut: 1. Buat histogram citra asal. Histogram citra adalah grafik yang menunjukkan distribusidari intensitas citra. Histogram citra menyatakan frekuensi kemunculan berbagai derajat keabuan dalam citra. Teknik pemodelan histogram mengubah citra hingga memiliki histogram sesuai keinginan. Untuk meningkatkan kualitas citra salah satunya dapat dilakukan dengan ekualisasi histogram. Dengan ekualisasi histogram dapat diperoleh histogram citra dengan distribusi seragam. 2. Tentukan jumlah kelompok dalam histogram. Jumlah tersebut menunjukkan tingkat intensitas warna citra kompresi. 3. Buat kelompok intensitas warna baru dengan menghitungjumlah piksel/jumlah kelompok. 4. Atur pengelompokan. 5. Lakukan kuantisasi.
6. Kodekan nilai intensitas warna hasil kuantisasi ke dalam citra kompresi. Contoh: Diketahui citra array berukuran 6x6 piksel (8 derajat keabuan, 3 bit) sebagai berikut :
Citra diatas memiliki histogram sebagai berikut:
Misalnya citra akan dikompresi menjadi 4 derajat keabuan (0 – 3), sehingga setiap derajat keabuan direpresentasikan dalam 2 bit. Citra asli terdiri dari 36 piksel, sehingga jika dibagi menjadi 4 kelompok maka dibuat 36/4= 9 piksel tiap kelompok sehingga histogram baru menjadi sebagai berikut:
Dari histogram baru tersebut diperoleh citra baru (dengan 4 derajat keabuan, 2 bit) sebagai berikut:
Ukuran citra sebelum kompresi = 36 piksel x 3 bit = 108 bit Ukuran citra setelah kompresi = 36 piksel x 2 bit = 72 bit
1.1.2
Pengembangan Algoritma Kompresi
Pengembangan dilakukan dengan mengubah proses pengkodean warna baru. Pengkodean dilakukan dengan cara yang lebih sederhana dengan mengelompokkan intensitas warna yang ada ke intensitas warna kompresi yang diinginkan. Contoh: Citra asal terdiri dari 8 derajat keabuan akan dikompresi menjadi 4 derajat keabuan.
Histogram citra terkompresi menjadi sebagai berikut:
Dari histogram diatas didapat matrik citra terkompresi sebagai berikut:
Untuk ukuran citra kompresi tidak berubah dari algoritma asal yakni 36 piksel x 3 bit = 108 bit.
2.
Pembahasan
Allgoritma diatas diimplementasikan dengan MatLab 6.5. Gambar lama yang semula mempunyai 256 derajat keabuan (8 bit), dikompresi menjadi 128 derajat keabuan (7 bit). Langkah pertama adalah melakukan pembacaan terhadap file citra sebagai berikut: clc;clear;
j=imread('anakblm.jpg'); gblama=j;
Selanjutnya dicari ukuran matrik citra. Ukuran ini nantinya akan digunakan sebagai indek perulangan dalam proses selanjutnya. Selain itu juga harus diketahui derajat keabuan minimum dan maksimun untuk menentukan ukuran vektor histogram citra. [m,n]=size(gblama); % ukuran gblama array=double(gblama); % indeks perulangan figure,imshow(gblama); % tampilan gambar gblama %cari derajat keabuan min dan max abumax=max(max(array)); % indeks dari abu2 maximum abumin=min(min(array)); % indeks dari abu2 maximum
Histogram citra asal dapat diketahui dengan cara menghitung kemunculan untuk setiap derajat keabuan yang ada di dalam citra. %bikin histogram image yg pake fungsi imhist arraybaru=zeros(1,256); for i=1:256 for j=1:m for k=1:n if (array(j,k)==i) arraybaru(i)=arraybaru(i)+1; end end end end; %hitung jumlah pixel jpixel=m*n; %gambar histogram figure,plot(arraybaru); %cari rentang derajat keabuan derajatabu=(abumax-abumin)+1; %membuat histogram kompresi %256 drjt keabuan (8 bit) menjadi %128 derajat keabuan (7bit) hiskompresi=zeros(1,128); j=1; for i=1:128 hiskompresi(i)=arraybaru(j)+arraybaru(j+1); j=j+2; end %bikin array image baru dari hasil kompresi gbbaru=zeros(m,n); for i=1:m for j=1:n if mod(double(gblama(i,j)),2)==0 gbbaru(i,j)=((double(gblama(i,j))+1))/2; else gbbaru(i,j)=double(gblama(i,j))/2; end end end %konversi array gambar terkompresi %dari double ke uint8 gbbaru1=uint8(gbbaru); %Tampilkan histogram gambar hasil kompresi %dan gambar hasil kompresi figure, plot(hiskompresi);
figure, imshow(gbbaru1); %Simpan gambar ke media penyimpan imwrite(gbbaru1,'anakkomp.jpg');
Hasil pada saat program dijalankan dengan beberapa sampel gambar yang berbeda adalah sbb:
Rerata rasio kompresi didapat sebesar 30.33% 3.
Kesimpulan
Modifikasi yang dilakukan terhadap algoritma asal tidak mengubah kapasitas file dari citra hasil kompresi. Presentase kompresi dipengaruhi oleh histogram citra yang diberikan citra asal. Tingkat kemunculan derajat keabuan suatu citra, antara citra yang satu dengan citra yang lainpun berbeda-beda. Daftar Pustaka Departmen Teknik Elektro, Modul Praktikum Pengolahan Citra dan Pengenalan Pola, Institut Teknologi Bandung. Paul Wintz, 2000, Digital Image Processing, Prentice-Hall. MatLab 6 Help. William J Palm, 2004, Introduction to MatLab 6 for Engineers, The McGraw-Hill Companies, Inc. www.iprg.ee.itb.ac.id/lectures.html www.datacompression.info/Quantization.shtml www.cbloom.com/src/index_im.html
