Megabangetz’s Weblog

Rekontruksi 3D

Listing program stereo.m

% STEREO

%

% Usage: pt3D = stereo(im1, im2, c1, c2)

%

% Where: im1 and im2 are the two stereo images

% C1 and C2 are the calibration matrices

% for these two images respectively

%

% The function will prompt you to digitise some points in the first image

% (finishing by clicking the right button). The function then

% prompts you to digitise the equivalent points (which you must digitise

% in exactly the same sequence) in the second image.

% The function then solves the stereo equations

% and returns the 3D coordinates of the points in pt3D.

function [XYZ, uv1, uv2] = stereo(im1, im2, C1, C2)

fprintf(1, ‘Digitise some points in figure 1\n’);

figure(1)

imshow(im1);

[u1,v1] = digipts;

uv1 = [u1,v1]’;

fprintf(1, ‘Digitise some points in figure 2\n’);

figure(2)

imshow(im2);

[u2,v2] = digipts;

uv2 = [u2,v2]’;

% check if same number of points are selected

if length(u1) ~= length(u2)

fprintf(1, ‘Same number of points not selected\n’);

end

for i = 1:length(u1)

a = [C1(1:2,1:3) – [u1(i)*C1(3,1:3); v1(i)*C1(3,1:3)];

C2(1:2,1:3) – [u2(i)*C2(3,1:3); v2(i)*C2(3,1:3)]];

c = [u1(i) – C1(1,4);

v1(i) – C1(2,4);

u2(i) – C2(1,4);

v2(i) – C2(2,4)];

b(:, i) = a \ c;

end

XYZ = b’;

 

Image ’stereo1.jpg’

clip_image002

 

Image ”stereo2.jpg’

clip_image004

 

Dengan menjalankan program stereo.m, maka didapatkan koordinat 3D dari titik-titik sudut ketiga bangun ruang pada gambar diatas sebagai berikut:

pt3D =

-287.7542 156.7465 143.2885

-176.6269 162.0598 144.3627

-175.6786 19.8716 141.7267

-290.0593 13.4079 144.3193

-287.5655 157.5390 10.9533

-286.3338 13.6666 9.4036

-172.4321 17.1203 3.2873

-83.8519 -66.6201 1.3714

-66.1821 -141.0900 122.3533

20.2680 -179.4214 -6.0719

-134.6197 -204.2349 -2.2669

136.6957 -93.6632 63.9256

207.1254 -91.0740 63.7992

201.4685 -162.9965 60.2827

138.5328 -155.1594 59.2890

132.6638 -94.1806 -4.4145

134.7488 -160.4601 -9.8006

206.4264 -163.0088 -8.9430

 

Panjang sisi balok

slengths =

143.3608 111.2594 142.2157 114.5926

114.5926 138.5048 114.1181 134.9674

134.9674 143.8860 132.3377 143.3608

114.7735 132.3377 111.2594 134.9674

134.9674 142.2157 138.5048 145.3295

145.3295 114.1181 143.8860 114.7735

nface =

4 1 2 3 4

4 3 7 6 4

4 6 5 1 4

8 5 1 2 8

8 2 3 7 8

8 7 6 5 8

 

Panjang sisi limas

slengths =

143.1594 159.4866 153.6897

143.1594 155.5672 146.7257

153.6897 156.9089 146.7257

159.4866 156.9089 155.5672

nface =

1 2 3 1

1 2 4 1

1 3 4 1

2 3 4 2

 

Panjang sisi kubus

slengths =

61.6981 70.4773 72.2303 63.4295

63.4295 69.4031 71.7281 69.3959

69.3959 66.5307 68.4608 61.6981

70.7170 68.4608 70.4773 69.3959

69.3959 72.2303 69.4031 66.8054

66.8054 71.7281 66.5307 70.7170

nface =

4 1 2 3 4

4 3 7 6 4

4 6 5 1 4

8 5 1 2 8

8 2 3 7 8

8 7 6 5 8

 

Untuk menampilkan dua view yang berbeda dari rekonstruksi 3D kubus, balok dan limas dengan koordinat titik sudut estimasi dari kubus dan balok yang tidak kelihatan, juga sumbu koordinat yang mengindikasikan bidang dasar, jalankan fungsi rekon.m

 

Listing fungsi rekon.m

function rekon()

im1 = imread(‘stereo1.jpg’ );

im2 = imread(‘stereo2.jpg’ );

C1 = [0.6596 -0.7391 -0.0615 363.4235;

-0.1851 -0.1387 -0.9437 342.7417;

0.0005 0.0003 -0.0003 1.0000];

C2 = [0.9234 -0.2221 -0.0257 347.7796;

-0.0741 -0.2278 -0.9168 339.8960;

0.0002 0.0004 -0.0002 1.0000];

pt3D = stereo(im1, im2, C1, C2)

figure(3)

cube(pt3D(1:7,:));

tetrahedron(pt3D(8:11,:));

cube(pt3D(12:18,:));

% label coordinate axes

text(100,0,0,’x’);

text(0,100,0,’y’);

text(0,0,100,’z’);

% draw in a set of coordinate axes

axislength = 100*eye(3);

for i=1:3

line([0, axislength(i,1)], [0, axislength(i,2)], [0, axislength(i,3)]);

end

axis equal; box on; rotate3D on; grid on;

end

function cube(cubepts3D)

% determine hidden vertex

cubepts3D(8,: ) = – cubepts3D(4,: ) + cubepts3D(6,: ) + cubepts3D(2,: );

% define faces from standard numbering

cubefaces = [4 1 2 3

4 3 7 6

4 6 5 1

8 5 1 2

8 2 3 7

8 7 6 5];

% draw ‘patcheds’ from vertice and face matrix

patch(‘Faces’,cubefaces,’Vertices’,cubepts3D, ‘FaceColor’, ‘none’)

fprintf(1, ‘Length matrix of face sides\n’);

[slengths, nface] = sidelengths(cubepts3D, cubefaces)

end

function tetrahedron(tetrahedronpts3D)

% define faces from standard numbering

tetrahedronfaces = [1 2 3

1 2 4

1 3 4

2 3 4];

% draw ‘patcheds’ from vertice and face matrix

patch(‘Faces’,tetrahedronfaces,’Vertices’,tetrahedronpts3D, ‘FaceColor’, ‘none’)

fprintf(1, ‘Length matrix of face sides\n’);

[slengths, nface] = sidelengths(tetrahedronpts3D, tetrahedronfaces)

end

function [slengths, nface] = sidelengths(pt3D, face)

[rows, cols] = size(face);

nface = [face face(:,1)];

for i=1:cols

for j=1:rows

slengths(j,i) = norm(pt3D(nface(j,i),:)-pt3D(nface(j,i+1),:));

end

end

end

 

Hasilnya (dilihat dari 2 sisi yang berbeda):

clip_image006

clip_image008

Kalibrasi Kamera

Listing fungsi calibrate.m

% CALIBRATE

%

% Function to perform camera calibration

%

% Usage: C = calibrate(im, XYZ, uv)

%

% Where: im – is the image of the calibration target.

% XYZ – is a n x 3 array of XYZ coordinates

% of the calibration target points.

% uv – is a 2 x n array of the image coordinates

% of the calibration target points.

% C – is the 3 x 4 camera calibration matrix.

%

% This function plots the uv coordinates onto the image of

% the calibration target. It also projects the XYZ coordinates

% back into image coordinates using the calibration matrix

% and plots these points too as a visual check on the accuracy of

% the calibration process.

% Lines from the origin to the vanishing points in the X, Y and

% Z directions are overlaid on the image.

% The mean squared error between the positions of the uv coodinates

% and the projected XYZ coordinates is also reported.

%

% The function should also report the error in satisfying the

% camera calibration matrix constraint – the magnitude of

% (q1 x q3).(q2 x q3)

%

function C = calibrate(im, XYZ, uv)

% obtain rows so arbitrary number of points can be used

[rows, cols] = size(XYZ);

XYZ1 = [XYZ, ones(rows, 1)]; % makes it easier to work with

% build B matrix

for n = 1:rows

B(2*n-1, : )= [XYZ1(n, : ) 0 0 0 0 -uv(1, n)*XYZ(n, : )];

B(2*n, : ) = [0 0 0 0 XYZ1(n, : ) -uv(2, n)*XYZ(n, : )];

end

c = B \ uv(:);

c(12) = 1;

C = reshape(c,4,3)’

XYZ1 = XYZ1′;

for i = 1:rows

suv(:,i) = C*XYZ1(:,i);

suv(:,i) = suv(:,i)/suv(3,i);

end

% calculate the mean squared error between the positions of the uv coodinates and suv.

mse = mean(mean((uv – suv(1:2,:)).^2));

fprintf(1, ‘mean squared error is %d\n’, mse);

% calculate the error in satisfying the camera calibration matrix constraint

q1 = C(1,1:3)’;

q2 = C(2,1:3)’;

q3 = C(3,1:3)’;

error = abs(dot(cross(q1, q3), cross(q2, q3)));

fprintf(1, ‘error in satisfying the camera calibration matrix is %d\n’, error);

% annotate image

imshow(im);

hold;

% plot uv coordinates

plot(uv(1, :), uv(2, :),’r+’)

% plot XYZ coordinates

plot(suv(1,:), suv(2,:),’bx’)

% draw vanishing lines

for i = 1:3

line([C(1,4), C(1,i)/C(3,i)],[C(2,4), C(2,i)/C(3,i)])

end

hold;

Sebelum menjalankan fungsi calibrate.m, fungsi stereo.m ( ada di portofolio 8 ) harus dijalankan terlebih dahulu untuk mendapatkan XYZ dan uv.

Dari fungsi stereo.m didapatkan:

XYZ =

-7.8964 124.8445 228.2591

-7.8892 48.7874 227.7442

-2.3221 130.4471 147.1380

-0.7141 50.5747 146.7381

-5.0965 127.0813 68.1338

-5.3585 45.8466 71.0790

51.4819 7.7006 220.6782

134.5497 12.8033 222.6472

50.3066 4.4372 146.1820

132.0424 8.0994 143.9290

51.7956 7.9634 57.3379

137.3195 9.2072 60.4677

 

uv1 =

Columns 1 through 10

261 327 258 326 258 325 393 427 391 425

115 130 187 202 257 274 127 103 196 178

Columns 11 through 12

384 421

277 245

 

uv2 =

Columns 1 through 10

306 333 306 335 302 330 400 469 396 464

103 124 172 197 241 263 138 127 207 195

Columns 11 through 12

392 462

279 268

 

Kemudian dijalankan fungsi calibrate.m, dari ’stereo1.jpg’ dihasilkan:

clip_image002

 

Kemudian dijalankan fungsi calibrate.m, dari ’stereo2.jpg’ dihasilkan:

clip_image004

 

Nilai matriks kalibrasi untuk ’stereo1.jpg’:

C =

0.6590 -0.7035 -0.0556 363.9159

-0.1924 -0.1134 -0.9588 345.3708

0.0005 0.0005 -0.0003 1.0000

 

Nilai matriks kalibrasi untuk ’stereo2.jpg’:

C =

0.9040 -0.2829 -0.0331 348.0122

-0.0735 -0.2575 -0.8987 336.8385

0.0002 0.0002 -0.0002 1.0000

 

Nilai error reproyeksi untuk ’stereo1.jpg’:

mean squared error is 5.546909e-001

error in satisfying the camera calibration matrix is 4.894755e-010

 

Nilai error reproyeksi untuk ’stereo2.jpg’:

mean squared error is 5.742843e-001

error in satisfying the camera calibration matrix is 8.831331e-009

Edge Detection

% mendapatkan image

im = imread(‘polka.jpg’);

% image grayscale

im = rgb2gray(im);

% smoothing image dengan filter gausian dengan sigma 2.5

imgaus = filter2(fspecial(‘gaussian’, size(im), 2.5), im); % makin gede makin blur

figure,imagesc(imgaus),colormap gray;

% mendapatkan image dg gradien pada horizontal direction

imhoris = filter2([-1 0 1], imgaus);

figure,imagesc(imhoris),colormap gray;

% mendapatkan image dg gradien pada vertikal direction

imvert = filter2([-1; 0; 1], imgaus);

figure,imagesc(imvert),colormap gray;

% mendapatkan edge detection dari image dengan canny option

imed = EDGE(imgaus,’canny’,[0.035 0.140]); % makin gede mkin dikit egdenya

figure,imagesc(imed),colormap gray;

 

untuk sigma 0.5

smoothing image dengan filter gausian:

clip_image002

image dg gradien pada horizontal direction:

clip_image004

image dg gradien pada vertikal direction:

clip_image006

edge detection dari image dengan canny option:

clip_image008

 

Untuk sigma 2.5

smoothing image dengan filter gausian:

clip_image010

image dg gradien pada horizontal direction:

clip_image012

image dg gradien pada vertikal direction:

clip_image014

edge detection dari image dengan canny option:

clip_image016

 

Untuk sigma 10

smoothing image dengan filter gausian:

clip_image018

image dg gradien pada horizontal direction:

clip_image020

image dg gradien pada vertikal direction:

clip_image022

edge detection dari image dengan canny option:

clip_image024

Threshold yang digunakan pada fungsi edge MATLAB:

[0.035 0.140]

 

Kesimpulan:

Semakin besar threshold yang digunakan, maka semakin sedikit edge yang terbentuk.

Semakin besar sigma yang digunakan pada filter gausian, maka hasil image akan semakin tampak blur.

Image Enhancement

% mendapatkan image

im = imread(‘polka.jpg’);

% grayscale image

im = rgb2gray(im);

% lowpass filter di domain spasial, cutoff 0.04, n 60

lowpassffthigh = lowpassfilter(size(im), 0.04, 60);

% highpass filter di domain spasial, cutoff 0.05, n 3

highpassffthigh = highpassfilter(size(im), 0.05, 3);

% highboost filter di domain spasial, cutoff 0.1, n 10, boost 150

highboostfft = highboostfilter(size(im), 0.1, 10, 150);

% lowpass filter di domain frekuensi, cutoff 0.2, n 50

lowpassfftlow = lowpassfilter(size(im), 0.2, 50);

% highpass filter di domain frekuensi, cutoff 0.1, n 35

highpassfftlow = highpassfilter(size(im), 0.1, 35);

%FFT image

imfft = fft2(im);

% menampilkan plot permukaan lowpass filter di domain frekuensi tinggi

surfl(fftshift(lowpassffthigh)), shading interp;

print -dpng lowpassffthighfreq.png

% menampilkan plot permukaan highpass filter di domain frekuensi tinggi

surfl(fftshift(highpassffthigh)), shading interp;

print -dpng highpassffthighfreq.png

% menampilkan plot permukaan highboost filter di domain frekuensi

surfl(fftshift(highboostfft)), shading interp;

print -dpng highboostfftfreq.png

% menampilkan plot permukaan lowpass filter di domain frekuensi rendah

surfl(fftshift(lowpassfftlow)), shading interp;

print -dpng lowpassfftlowfreq.png

% menampilkan plot permukaan highpass filter di domain frekuensi rendah

surfl(fftshift(highpassfftlow)), shading interp;

print -dpng highpassfftlowfreq.png

% menampilkan plot permukaan lowpass filter di domain spasial tinggi

surfl(fftshift(real(ifft2(lowpassffthigh)))), shading interp;

print -dpng lowpassffthighspac.png

% menampilkan plot permukaan highpass filter di domain spasial tinggi

surfl(fftshift(real(ifft2(highpassffthigh)))), shading interp;

print -dpng highpassffthighspac.png

% menampilkan plot permukaan highboost filter di domain spasial

surfl(fftshift(real(ifft2(highboostfft)))), shading interp;

print -dpng highboostfftspac.png

% menampilkan plot permukaan lowpass filter di domain spasial rendah

surfl(fftshift(real(ifft2(lowpassfftlow)))), shading interp;

print -dpng lowpassfftlowspac.png

% menampilkan plot permukaan highpass filter di domain spasial rendah

surfl(fftshift(real(ifft2(highpassfftlow)))), shading interp;

print -dpng highpassfftlowspac.png

% perkalian filter dengan fft image

imglowpassffthigh = lowpassffthigh.*imfft;

imghighpassffthigh = highpassffthigh.*imfft;

imghighboostfft = highboostfft.*imfft;

imglowpassfftlow = lowpassfftlow.*imfft;

imghighpassfftlow = highpassfftlow.*imfft;

% menampilkan image hasil perkalian filter dengan fft image

figure,imagesc((real(ifft2(imglowpassffthigh)))),colormap gray;

figure,imagesc((real(ifft2(imghighpassffthigh)))),colormap gray;

figure,imagesc((real(ifft2(imghighboostfft)))),colormap gray;

figure,imagesc((real(ifft2(imglowpassfftlow)))),colormap gray;

figure,imagesc((real(ifft2(imghighpassfftlow)))),colormap gray;

%menampilkan image grayscale

figure,imagesc(im),colormap gray;

program lain yang digunakan: lowpassfilter.m, highpassfilter.m, dan highboostfilter.m

Output

plot permukaan lowpass filter di domain frekuensi tinggi:

clip_image002

plot permukaan highpass filter di domain frekuensi tinggi:

clip_image004

plot permukaan lowpass filter di domain frekuensi rendah:

clip_image006

plot permukaan highpass filter di domain frekuensi rendah:

clip_image008

plot permukaan highboost filter di domain frekuensi:

clip_image010

plot permukaan lowpass filter di domain spasial tinggi:

clip_image012

plot permukaan highpass filter di domain spasial tinggi:

clip_image014

plot permukaan lowpass filter di domain spasial rendah:

clip_image016

plot permukaan highpass filter di domain spasial rendah:

clip_image018

plot permukaan highboost filter di domain spasial:

clip_image020

Image dengan high order pada lowpass filter:

clip_image022

Image dengan high order pada highpass filter:

clip_image024

Image dengan high order pada high-boost filter :

clip_image026

Image dengan low order pada lowpass filter:

clip_image028

Image dengan low order pada highpass filter:

clip_image030

Parameter yang digunakan

Pada lowpass filter domain spasial:

cutoff : 0.04

n : 60

Pada highpass filter domain spasial:

cutoff : 0.05

n : 3

Pada lowpass filter domain frekuensi:

cutoff : 0.2

n : 50

Pada highpass filter domain frekuensi:

cutoff : 0.1

n : 35

Pada highboost filter:

cutoff : 0.1

n : 10

boost : 150

Kesimpulan: apabila nilai boost yang digunakan pada high-boost filtering sangat besar, maka efek blur pada image juga akansemakin besar.

Merekonstruksi foto diri dengan menggunakan freqcomp.m

Listing program freqcomp.m

function recon = freqcomp(im, Npts, delay)

if ndims(im) == 3

% mengubah image jadi grayscale

im = rgb2gray(im);

warning(‘converting colour image to greyscale’);

end

if nargin < 2

Npts = 50;

end

[rows,cols] = size(im);

% If necessary crop one row and/or column so that there are an even

% number of rows and columns, this makes things easier later.

if mod(rows,2) % odd

rows = rows-1;

end

if mod(cols,2) % odd

cols = cols-1;

end

im = im(1:rows,1:cols);

rc = fix(rows/2)+1; % Centre point

cc = fix(cols/2)+1;

% The following code constructs two arrays of coordinates within the FFT

% of the image that correspond to complex conjugate pairs of points that

% spiral out from the centre visiting every frequency component on

% the way.

p1 = zeros(Npts,2); % Path 1

p2 = zeros(Npts,2); % Path 2

m = zeros(rows,cols); % Matrix for marking visited points

m(rc,cc) = 1;

m(rc,cc-1) = 1;

m(rc,cc+1) = 1;

p1(1,:) = [rc cc-1];

p2(1,:) = [rc cc+1];

d1 = [0 -1]; % initial directions of the paths

d2 = [0 1];

% Mark out two symmetric spiral paths out from the centre (I wish I

% could think of a neater way of doing this)

for n = 2:Npts

l1 = [-d1(2) d1(1)]; % left direction

l2 = [-d2(2) d2(1)];

lp1 = p1(n-1,:) + l1; % coords of point in left direction

lp2 = p2(n-1,:) + l2;

if ~m(lp1(1), lp1(2)) % go left

p1(n,:) = lp1;

d1 = l1;

m(p1(n,1), p1(n,2)) = 1; % mark point as visited

else % go sraight ahead

p1(n,:) = p1(n-1,:) + d1;

m(p1(n,1), p1(n,2)) = 1; % mark point as visited

end

if ~m(lp2(1), lp2(2)) % go left

p2(n,:) = lp2;

d2 = l2;

m(p2(n,1), p2(n,2)) = 1; % mark point as visited

else % go sraight ahead

p2(n,:) = p2(n-1,:) + d2;

m(p2(n,1), p2(n,2)) = 1; % mark point as visited

end

end

% Having constructed the path of frequency components to be visited

% we take the FFT of the image and then enter a loop that

% incrementally reconstructs the image from its components.

IM = fftshift(fft2(im));

recon = zeros(rows,cols); % Initialise reconstruction matrix

if max(rows,cols) < 150

fontsze = 7;

else

fontsze = 10;

end

figure(1), clf

subplot(2,2,1),imagesc(im),colormap gray, axis image, axis off

title(‘Original Image’,’FontSize’,fontsze);

subplot(2,2,2),imagesc(log(abs(IM))),colormap gray, axis image

axis off, title(‘Fourier Transform + frequency component pair’,’FontSize’,fontsze);

warning off % Turn off warnings that might arise if the images cannot be

% displayed full size

truesize(1)

for n = 1:Npts

% Extract the pair of Fourier components

F = zeros(rows,cols);

F(p1(n,1), p1(n,2)) = IM(p1(n,1), p1(n,2));

F(p2(n,1), p2(n,2)) = IM(p2(n,1), p2(n,2));

% Invert and add to reconstruction

f = real(ifft2(fftshift(F)));

recon = recon+f;

% Display results

subplot(2,2,2),imagesc(log(abs(IM))),colormap gray, axis image

axis off, title(‘Fourier Transform + frequency component pair’,’FontSize’,fontsze);

hold on, plot([p1(n,2), p2(n,2)], [p1(n,1), p2(n,1)],’r.’); hold off

subplot(2,2,3),imagesc(recon),colormap gray, axis image, axis off, title(‘Reconstruction’,’FontSize’,fontsze);

subplot(2,2,4),imagesc(f),colormap gray, axis image, axis off

title(‘Basis function corresponding to frequency component pair’,’FontSize’,fontsze);

if nargin == 3

pause(delay);

else

fprintf(‘Hit any key to continue \n’); pause

end

end

Dengan jumlah komponen frekuensi 30, hasilnya seperti dibawah ini:

clip_image002

Dengan jumlah komponen frekuensi 50, hasilnya seperti dibawah ini:

clip_image004

Dengan jumlah komponen frekuensi 120, hasilnya seperti dibawah ini:

clip_image006

Filter Gausian

function tugasgausian = Gausian(images);

% mendapatkan image

i = imread(images);

% image grayscale

im = rgb2gray(i);

% filter gausian

gausianfilter = fspecial(‘gaussian’, size(im), 5);

% FFT image

imagefft = fft2(im);

% FFT dari filter gaussian

gausianfilterfft = fft2(gausianfilter);

% perkalian FFT image dengan FFT filter gaussian

result = gausianfilterfft.*imagefft;

figure,imshow(im);

clip_image008

figure,imagesc(gausianfilter),colormap gray;

clip_image010

figure,imagesc(fftshift(log(abs(imagefft)+eps))),colormap gray;

clip_image012

figure,imagesc(fftshift(log(abs(gausianfilterfft)+eps))),colormap gray;

clip_image014

figure,imagesc(fftshift(log(abs(result)+eps))),colormap gray;

clip_image016

% menampilkan hasil smoothing

figure,imagesc(fftshift(real(ifft2(result)))),colormap gray;

clip_image018

Filter Average

function tugasaverage = average(images);

% mendapatkan image

i = imread(images);

% image grayscale

im = rgb2gray(i);

[row col] = size(im);

averageFilter = zeros(row,col);

if (mod(row,2)==0 && mod(col,2)==0)

pusatx = row/2;

pusaty = col/2;

averageFilter(pusatx-10:pusatx+11, pusaty-10:pusaty+11) = 1;

elseif (mod(row,2)==0 && mod(col,2)==1)

pusatx = row/2;

pusaty = (col-1)/2;

averageFilter(pusatx-10:pusatx+11, pusaty-10:pusaty+10) = 1;

elseif (mod(row,2)==1 && mod(col,2)==0)

pusatx = (row-1)/2;

pusaty = col/2;

averageFilter(pusatx-10:pusatx+10, pusaty-10:pusaty+11) = 1;

elseif (mod(row,2)==1 && mod(col,2)==1)

pusatx = (row-1)/2;

pusaty = (col-1)/2;

averageFilter(pusatx-10:pusatx+10, pusaty-10:pusaty+10) = 1;

end

% filter average

averageFilter = averageFilter/(sum(sum(averageFilter)));

% FFT image

imagefft = fft2(im);

% FFT dari filter average

averagefilterfft = fft2(averageFilter);

% perkalian FFT image dengan FFT filter average

result = averagefilterfft.*imagefft;

figure,imshow(im);

clip_image019

figure,imagesc(averageFilter),colormap gray;

clip_image021

figure,imagesc(fftshift(log(abs(imagefft)+eps))),colormap gray;

clip_image012[1]

figure,imagesc(fftshift(log(abs(averagefilterfft)+eps))),colormap gray;

clip_image023

figure,imagesc(fftshift(log(abs(result)+eps))),colormap gray;

clip_image025

% menampilkan hasil smoothing

figure,imagesc(fftshift(real(ifft2(result)))),colormap gray;

clip_image027

Pertukaran Nilai Fase

Image1

clip_image029

Image2

clip_image031

informasi amplitudo dari image1 dan informasi fase dari image2, hasilnya lebih mirip image2

clip_image033

Informasi amplitudo dari image2 dan informasi fase dari image1, hasilnya lebih mirip image1

clip_image035

Fase dari image1

clip_image037

Fase dari image2

clip_image039

Listing program pertukaran nilai fase

% mendapatkan 2 gambar yg berbeda

image1 = imread(‘polka.jpg’);

image2 = imread(‘sok imut1.jpg’);

% image1 dirubah ke bentuk grayscale

image1 = rgb2gray(image1);

% image2 dirubah ke bentuk grayscale

image2 = rgb2gray(image2);

% FFT image1

image1fft = fft2(image1);

% FFT image2

image2fft = fft2(image2);

% mendapatkan informasi fase dari image1

image1fase = angle(image1fft);

% mendapatkan informasi fase dari image2

image2fase = angle(image2fft);

% mendapatkan amplitudo dari image1

image1amp = abs(image1fft);

% mendapatkan amplitudo dari image2

image2amp = abs(image2fft);

% informasi amplitudo dari image1, informasi fase dari image2

newimfft1 = image1amp.*(cos(image2fase) + i*sin(image2fase));

% informasi amplitudo dari image2, informasi fase dari image1

newimfft2 = image2amp.*(cos(image1fase) + i*sin(image1fase));

figure,imagesc(real(ifft2(newimfft1))),colormap gray;

figure,imagesc(real(ifft2(newimfft2))),colormap gray;

% mendapatkan fase dari image1 saja

image1faseaja = image1fft./image1amp;

% mendapatkan fase dari image2 saja

image2faseaja = image2fft./image2amp;

figure,imagesc(real(ifft2(image1faseaja))),colormap gray;

figure,imagesc(real(ifft2(image2faseaja))),colormap gray;

figure,imagesc(image1),colormap gray;

figure,imagesc(image2),colormap gray;

Tugas 3

Menggunakan fungsi opening dan closing

Foto asli:

% membuka image

im = imread (‘polka.jpg’);

imshow (im)

clip_image002

% mengubah image ke bentuk grayscale

grayscale = rgb2gray(im);

clip_image004

% mengubah image ke bentuk biner

biner = im2bw(grayscale,0.5);

clip_image006

% operasi closing

se = strel(‘disk’,3);

IMC = imclose(biner,se);

imshow(IMC), title(‘Closing’)

clip_image008

% operasi opening

IMO = imopen(biner,se);

imshow(IMO), title(‘Opening’)

clip_image010

% operasi open-close untuk menghilangkan noise

openclose = imclose(IMO,se);

imshow(openclose), title (‘open-close’)

clip_image012

% operasi close-open

closeopen = imopen (IMC,se);

imshow(closeopen), title (‘close-open’)

clip_image014

Tugas 4

Menentukan lokasi empat kotak di pojok form menggunakan program locatemarks.m

Image ’surf01.png’

clip_image016

Outputnya

clip_image018

Image ’surf02.png’

clip_image020

Outputnya

clip_image022

Image ’surf05.png’

clip_image024

Outputnya

clip_image026

Listing program locatemarks.m

% LOCATELANDMARKS – locates landmarks on SURF form

%

% Usage: [tl, tr, bl, br] = locatelandmarks(im)

%

% Argument: im – citra biner

%

% Returns: tl, tr, bl, br

% – Koordinat titik pusat dari landmarks di kiri-atas

% (top-left),kanan-atas (top-right), kiri-bawah

% (bottom-left), dan kanan-bawah (bottom-right).

% Koordinat ini berbentuk vektor kolom [row; col]

% untuk setiap landmark.

%

% Fungsi ini juga menampilkan image dengan titik pusat (centroids)

% dari landmarks overlayed.

function [tl, tr, bl, br] = locatelandmarks(im)

% membaca input image

gambar = imread(im)

[rows,cols] = size(gambar);

% koordinat landmarks

tl=[rows;cols];

tr=[rows;0 ];

bl=[0 ;cols];

br=[0 ;0 ];

% citra biner dengan background hitam dan objek putih

bw = ~gambar;

% faktor skala

a = 8/800;

SE = ones(a*rows, a*cols);

% operasi close-open

closeopen = imopen(imclose(bw,SE),SE);

x = ones(rows,1)*[1:cols];

y = [1:rows]’*ones(1,cols);

% bwlabel untuk menandai setiap blob dan menghasilkan jumlah blob

[L, num] = bwlabel(closeopen);

% mencari blob untuk centroid

for i = 1:num

img = L==i;

area(i) = sum(sum(img));

meanx = sum(sum(double(img).*x))/area(i);

meany = sum(sum(double(img).*y))/area(i);

% koordinat top-left

if tl(1)+tl(2) > meanx + meany

tl = [meanx; meany];

tlindex = i;

end

% koordinat top-right

if tr(1)-tr(2) > meanx – meany

tr = [meanx; meany];

trindex = i;

end

% koordinat bottom-left

if bl(1)-bl(2) < meanx – meany

bl = [meanx; meany];

blindex = i;

end

% koordinat bottom-right

if br(1)+br(2) < meanx + meany

br = [meanx; meany];

brindex = i;

end

end

imshow(gambar)

hold;

landmarks = [tl, tr, br, bl];

plot(landmarks(1,:), landmarks(2,:),’xr’,’MarkerSize’,12,’LineWidth’,1);

plot(landmarks(1,:), landmarks(2,:),’+r’,’MarkerSize’,12,’LineWidth’,1);

hold;

end

Tugas 1

Fungsi yang dapat mengubah citra lego berwarna menjadi citra lego biner

dengan nilai threshold tertentu.

Nama : Elvira Megasari S.

NRP : 5105100020

% memanggil gambar lego

im = imread(‘lego1.png’);

% memunculkan gambar lego

imshow(im);

% mengubah citra lego berwarna menjadi citra lego biner dengan nilai

% threshold 0.64

BW = im2bw(im,0.64);

% Citra biner yang terbentuk akan mempunyai warna hitam untuk obyek, dan

% warna putih untuk background

imshow(BW);

% diubah warna hitam untuk background, dan warna putih untuk obyek

imshow(~BW);

% untuk pemisahan dan pelebelan citra biner

BWL = bwlabel(~BW);

imagesc(BWL), colormap(gray);

% memilih salah satu objek

im = BWL==9;

imshow(im);

Hasilnya sebagai berikut:

lego1.png

clip_image002

Citra lego biner dengan nilai threshold 0.64

clip_image004

di negasi

clip_image006

clip_image008

Memilih objek 9

clip_image010

Tugas 2

% MOMENTS

%

% Fungsi menghitung momen sebuah citra biner dan mengembalikan

% area, centroid, sudut dari sumbu inertia minimum, dan sebuah

% ukuran ‘kebulatan’. Fungsi ini mengasumsikan bahwa hanya ada

% satu obyek pada citra biner.

%

% function [area, centroid, theta, roundness] = moments(im)

%

% Argument: im – citra biner berisi nilai 0 atau 1

%

% Returns: area – luasan citra biner

% centroid – vektor dengan 2 elemen

% theta – sudut sumbu kelembaman minimum (radian)

% roundness – ratio minimum inertia/maximum inertia.

%

% Fungsi ini juga menampilkan citra dan letak dari

% centroid dan sumbu inertia minimum. Area yang dihitung

% harus ditampilkan sebagai title di citra.

function [area, centroid, theta, roundness] = moments(im)

% mendapatkan ukuran matriks image

[rows,cols] = size(im);

% mendapatkan koordinat x

x = ones(rows,1)*[1:cols];

% mendapatkan koordinat y

y = [1:rows]’*ones(1,cols);

% menghitung luasan image/area

area = sum(sum(im));

area

% menghitung mean x

meanx = sum(sum(double(im).* x))/area;

% menghitung mean y

meany = sum(sum(double(im).*y))/area;

%mendapatkan vektor centroid

centroid = [meanx meany];

% mendapatkan invers x dan y

x = x – meanx;

y = y – meany;

a = sum(sum(double(im).*(x.^2)));

b = 2*(sum(sum(double(im).*x.*y)));

c = sum(sum((y.^2).*double(im)));

% mendapatkan nilai theta

denominator = sqrt(b^2 + (a-c)^2);

y1 = b/denominator;

x1 = (a-c)/denominator;

minTheta = 0.5 * atan2(y1,x1);

maxTheta = 0.5 * atan2(-y1,-x1);

theta = minTheta;

% hitung roundness

Imin = 0.5 * (c + a) – (0.5 * (a – c) * cos(2 * minTheta)) – (0.5 * b * sin(2 * minTheta));

Imax = 0.5 * (c + a) – (0.5 * (a – c) * cos(2 * maxTheta)) – (0.5 * b * sin(2 * maxTheta));

roundness = Imin/Imax;

xmin = [meanx-50 : meanx + 50];

ymin = ((xmin – meanx) * tan(minTheta)) + meany;

figure, imshow(im),

title([‘area : ‘, num2str(area),’ centroid: [‘, num2str(centroid), ‘]’, ‘theta: ‘, num2str(theta), ’roundness : ‘, num2str(roundness)]),

line(xmin, ymin, ‘Color’, ‘b’, ‘LineWidth’,1),

text(meanx-3, meany, ‘X’);

pixval on;

Hasilnya sebagai berikut:

Objek : 9

area = 5393

centroid = 339.1081 227.2199

theta = -0.8402

roundness = 0.1136

clip_image012

Objek : 6

area = 3158

centroid = 305.5741 342.6491

theta = 0.0336

roundness = 0.0319

clip_image014

Objek : 3

area = 1741

centroid = 137.0270 352.6588

theta = 0.1400

roundness = 0.8212

clip_image016


  • None
  • No comments yet

Archives