ĐÀM ĐỨc cưỜng ứng dụng kỹ thuật kết hợp tần số nhằm nâng cao chất lưỢng ảnh siêU Âm cắt lớP



tải về 362.16 Kb.
trang5/6
Chuyển đổi dữ liệu29.10.2017
Kích362.16 Kb.
1   2   3   4   5   6

KẾT LUẬN


Luận văn này đã thành công trong việc nâng cao chất lượng ảnh chụp siêu âm cắt lớp bằng cách sử dụng kết hợp 2 tần số f1 và f2. Ảnh khôi phục theo phương pháp đề xuất cho chất lượng tốt hơn ảnh theo phương pháp truyền thống.

Tác giả cũng đã tìm được số bước lặp tối ưu với f1 sao cho việc kết hợp f1 và f2 cho chất lượng tốt nhất.

Đánh giá được tham số chất lượng Q được trình bày ở phần 2.3. Từ đó kết luận được ảnh tái tạo bởi việc sử dụng kết hợp 2 tần số f1 và f2, cho kết quả đánh giá về mặt sai số toán học thông dụng hay có xét cả đến vấn đề VHS (visual human system) đều tốt hơn so với chỉ sử dụng một tần số đơn.

Như vậy việc sử dụng kết hợp 2 tần số trong việc cải thiện chất lượng ảnh đã thành công, tạo điều kiện áp dụng trong lĩnh vực Y – Sinh. Bước tiếp theo của đề xuất này là việc thử nghiệm đề xuất trong tạo ảnh với những dữ liệu thực tế để có thể áp dụng trong ngành chuẩn đoán y khoa.


TÀI LIỆU THAM KHẢO


[1] C. F. Schueler, H. Lee, and G. Wade, “Fundamentals of digital ultrasonic processing,” IEEE Transactions on Sonics and Ultrasonics, vol. 31, no. 4, pp. 195–217, July 1984.

[2] N. Duric, P. Littrup, A. Babkin, D. Chambers, S. Azevedo, A. Kalinin, R.Pevzner, M. Tokarev, E. Holsapple, O. Rama, and R. Duncan, “Development of ultrasound tomography for breast imaging: Technical assessment,” Medical Physics, vol. 32, no. 5, pp. 1375–1386, May 2005.

[3] J.-W. Jeong, T.-S. Kim, D. C. Shin, S. Do, M. Singh, and V. Z. Marmarelis, “Soft tissue differentiation using multiband signatures of high resolution ul-trasonic transmission tomography,” IEEE Transactions on Medical Imaging, vol. 24, no. 3, pp. 399–408, March 2005.

[4] S. A. Johnson, T. Abbott, R. Bell, M. Berggren, D. Borup, D. Robinson, J. Wiskin, S. Olsen, and B. Hanover, “Noninvasive breast tissue charac-terization using ultrasound speed and attenuation,” in Acoustical Imaging, vol. 28, 2007, pp. 147–154.

[5] C. Li, N. Duric, and L. Huang, “Breast imaging using transmission ultra-sound: Reconstructing tissue parameters of sound speed and attenuation,” in International Conference on BioMedical Engineering and Informatics, vol. 2, 2008, pp. 708–712.

[6] R. J. Lavarello and M. L. Oelze: Tomographic Reconstruction of Three-Dimensional Volumes Using the Distorted Born Iterative Method. IEEE Transactions on Medical Imaging, 28, 2009, pp. 1643-1653.

[7] Lavarello Robert: New Developments on Quantitative Imaging Using Ultrasonic Waves. University of Illinois at Urbana-Champaign, 2009.

[8] http://en.wikipedia.org/wiki/Nonlinear_conjugate_gradient_method

[9] M. T. Heath, Scientific Computing: An Introductory Survey. New York, NY: McGraw-Hill, 2002.

[10] Martin, R., Noise power spectral density estimation based on optimal smoothing and minimum statistics, IEEE Transactions on Speech and Audio Processing, Vol. 9, 2001, pp. 504 - 512.

[11] http://www-stat.stanford.edu/~susan/courses/s60/split/node60.html

[12] Tran Duc Tan, N. Linh-Trung, M. L. Oelze, M. N. Do, Application of L1 regularization for high-quality reconstruction of ultrasound tomography, International Federation for Medical and Biological Engineering (IFMBE), NXB SPRINGER, ISSN: 1680-0737, Volume 40, 2013, pp. 309-312.

[13] Tran Duc Tan, Nguyen Linh-Trung, Minh N. Do, Modified Distorted Born Iterative Method for Ultrasound Tomography by Random Sampling, The 12th International Symposium on Communications and Information Technologies (ISCIT 2012), Australia, 2012, pp. 1065-1068.

[14] Tran Duc Tan, Automated Regularization Parameter Selection in Born Iterative Method for Ultrasound Tomography, Vietnam Conference on Control and Automation (VCCA-2011), ISBN 978-604-911-020-7, 2011, pp.786-791.

[15] Tran Duc Tan, Gian Quoc Anh, Improvement of Distorted Born Iterative Method for Reconstructing of Sound Speed, Vietnam Conference on Control and Automation (VCCA-2011), ISBN 978-604-911-020-7, 2011, pp.798-803.

[16] Zhou Wang, Student Member : A Universal Image Quality Index, IEEE Signal Processing Letters, Vol. 9, No. 3, March 2002



PHỤ LỤC 1: CODE MATLAB DBIM


Để tính toán được sai số của phép khôi phục, cũng như để so sánh ảnh tạo được với vật thể ta cần có một giá trị tham chiếu của thí nghiệp gọi là hàm mục tiêu lý tưởng. Hàm này có thể được tạo bởi phương trình (2.1) như vậy ta có thể tạo hàm mục tiêu lý tưởng cho vật thể hình trụ với tần số đo f0 như sau:

Hàm mục tiêu lý tưởng

SC=zeros(N,N);

xo=(N+1)/2;yo=(N+1)/2;

for m=1:N

for n=1:N

dis=sqrt((xo-m)^2+(yo-n)^2);

if dis>3.5355*N/10

SC(m,n)=0;

else

SC(m,n)=(2*pi*f0)^2*(1/(c1^2)-1/(co^2));%da doi f thanh f0



end;

end;


end;
Sau khi đã có hàm mục tiêu lý tưởng ta tạo cấu hình hệ đo với việc bố trí các máy phát và máy thu xung quanh vật thể:
L=N*N;

phi=linspace(-pi,pi,L);

No=10*N; % at first N=11, No can be changed when changing N

% however it is affected to distance from tranceivers to object

% because distance=No*h= constant

K2=cos(phi)*(No+.5)+x0(2);

K1=sin(phi)*(1.5*No-.5) + x0(1);

KK2=K2;KK1=K1;

%noise_flag: option 0: noise init, 2: no noise

noise_flag=0;

transmiter=1:N:L;%may phat co the thay doi duoc

detector=1:2*N:L;%may thu co the thay doi duoc

plot(KK2(detector),KK1(detector),'s')

hold on


plot(K2(transmiter),K1(transmiter),'r*')

hold on;


mesh(abs(SC))

legend('detector','transmiter','scatter area')


Như vậy ta đã có hàm mục tiêu lý tưởng cần khôi phục và một hệ đo, áp dụng thuật toán 1: Lặp vi phân Born, để viết chương trình trên Matlab ta có các hàm con phải tính sau:

Đầu tiên


Tính , , ,và tương ứng sử dụng (2.6) và (2.7)

Tính tín hiệu của sóng tới viết trên Matlab ta có

pix=[];

k=1:N;


for i=1:N

pix=[pix;k];

end;

PINC=[];



for l=transmiter

pinc=besselj(0,ko*h*sqrt((K1(l)-pix').^2+(K2(l)-pix).^2));

%pinc=besselj(0,4.1e3*h*sqrt((K1(l)-pix').^2+(K2(l)-pix).^2));

PINC=[PINC ; pinc];

end;

save PINC_2D_matrix PINC



Tín hiệu trong thực tế có thể đo được bằng cách lấy hiệu số của tín hiệu tại máy thu khi có đối tượng và khi không có đối tượng. Còn trong mô phỏng thì lại có thể tính bằng phương trình (2.7) sử dụng hàm mục tiêu lý tưởng. Như vậy theo phương trình ta còn phải tính hai ma trận B và C, ma trận B và C tính là ma trận hệ số của hàm Green từ các pixel tới máy thu và hệ số Green giữa các pixel:
calculate_B_matrix_DBIM:
BB=[];

ko_SC=sqrt(ko*ko+abs(SC)); % matrix

ko_SC1=sqrt(ko*ko+abs(SC1)); % matrix

BB_SC=[];BB_SC1=[];

for l=detector

B=-.25*j*h*h*besselj(0,ko*h*sqrt((KK1(l)-pix').^2+(KK2(l)-pix).^2));% no update

B_SC=-.25*j*h*h*besselj(0,ko_SC.*h*sqrt((KK1(l)-pix').^2+(KK2(l)-pix).^2));

B_SC1=-.25*j*h*h*besselj(0,ko_SC1.*h*sqrt((KK1(l)-pix').^2+(KK2(l)-pix).^2));

BB=[BB ; B];

BB_SC=[BB_SC ; B_SC];BB_SC1=[BB_SC1 ; B_SC1];

end;
calculate_C_matrix_DBIM:
CC=[];

CC_SC=[];CC_SC1=[];

for l1=1:N

for l2=1:N

C=-.25*j*h*h*besselj(0,ko*h*sqrt((l1-pix').^2+(l2-pix).^2)); %no update

C_SC =-.25*j*h*h*besselj(0,ko_SC.*h*sqrt((l1-pix').^2+(l2-pix).^2));

C_SC1=-.25*j*h*h*besselj(0,ko_SC1.*h*sqrt((l1-pix').^2+(l2-pix).^2));

CC=[CC ; C];

CC_SC=[CC_SC ; C_SC];CC_SC1=[CC_SC1 ; C_SC1];

end;


end;

Sau khi đã tính những tham số và hàm con trên ta tính được từ giá trị đo được và giá trị tiên đoán và tính RRE tương ứng sử dụng công thức (2.11), Tính giá trị mới sử dụng (2.9). Việc tính sử dụng công thức (2.11) ta phải áp dụng phương pháp NCG như sau:



Áp Dụng NCG Để Tính
function[delta_sound]=test_NCG(Mt,delta_sc_t,ni,RRE,gama)

[n1,n2]=size(Mt);


b=Mt'*delta_sc_t;

x=b;


r=b;

delta_sound=zeros(n2,1);

%delta_sound=zz;

for i=1:ni

q=Mt*x;

%alpha=transpose(r)*r/(transpose(q)*q+gama*transpose(x)*x);



alpha=r'*r/(q'*q+gama*x'*x);

%s=transpose(Mt)*q;

s=Mt'*q;

r_update=r-alpha*(s+gama*x);

%beta=(transpose(r_update)*r_update)/(transpose(r)*r);

beta=(r_update'*r_update)/(r'*r);

delta_sound=delta_sound+alpha*x;

x=r_update+beta*x;

r=r_update;

%e=sum(abs(delta_sc_t))/sum(abs(p_sc_exact_t));

temp=delta_sc_t-Mt*delta_sound;

e=temp'*temp/(delta_sc_t'*delta_sc_t);

%tol=delta_sc_t'*delta_sc_t/(p_sc_exact_t'*p_sc_exact_t)

if e


e

fprintf('convergent at step %d \n',i)

break

end


end;
Như vậy ta có chương chính chính như sau với kết quả đầu ra là hàm khôi phục và các giá trị sai số:

Chương Trình Chính

clear,clc,close all;

f0=2e6;%

f=1e6;% tan so don co the thay doi duoc

N=22

Niter=8 %so lan lap


co=1540; % m/s

c1=co*(1-0.02); % Khoa hoc tre

contrast

landa=co/f; % m/s

lamda=co/f0;

ko=2*pi/landa; % increase when f is increased

D=4*lamda

x0=[(N+1)/2;(N+1)/2];

fprintf('Ratio landa/h = %f \n',landa/h);
%tao muc ham muc tieu ly tuong

SC=zeros(N,N);

SC1=SC;

xo=(N+1)/2;yo=(N+1)/2;



for m=1:N

for n=1:N

dis=sqrt((xo-m)^2+(yo-n)^2);

%if dis

if dis>3.5355*N/10

SC(m,n)=0;

else

SC(m,n)=(2*pi*f0)^2*(1/(c1^2)-1/(co^2));%da doi f thanh f0



%SC(m,n)=dis;

end;


% if dis<3.5355*N/20

% SC(m,n)=2*(2*pi*f0)^2*(1/(c1^2)-1/(co^2));

%SC(m,n)=dis;

% end;


end;

end;


figure;

[X,Y] = meshgrid(linspace(0,D/landa,length(SC)));

mesh(X,Y,SC);

xlabel('\lambda')

ylabel('\lambda')
RRE=2^(-4);

%RRE=2^(-6);

L=N*N;

phi=linspace(-pi,pi,L);



No=10*N; % at first N=11, No can be changed when changing N

% however it is affected to distance from tranceivers to object

% because distance=No*h= constant

K2=cos(phi)*(No+.5)+x0(2);

K1=sin(phi)*(1.5*No-.5) + x0(1);

KK2=K2;KK1=K1;

%noise_flag: option 0: noise init, 2: no noise

noise_flag=0;

transmiter=1:N:L;%may phat co the thay doi duoc

detector=1:2*N:L;%may thu co the thay doi duoc


plot(KK2(detector),KK1(detector),'s')

hold on


plot(K2(transmiter),K1(transmiter),'r*')

hold on;


mesh(abs(SC))

legend('detector','transmiter','scatter area')

calculate_PINC_matrix_cavichi;

figure(100)

subplot(211)

imagesc(abs(SC))

subplot(212)

mesh(abs(SC))

Niter;

err1=zeros(1,Niter);



MSE1=zeros(1,Niter);

PSNR1=zeros(1,Niter);

gama1=zeros(1,Niter);

Quality=zeros(1,Niter);

for iter=1:Niter

iter


Mt=[];

delta_sc_t=[];

X=[];

p_sc_exact_t=[];



p_sc_t=[];

calculate_B_matrix_DBIM;

calculate_C_matrix_DBIM;

% update C matrix for each iteration

uu=-1; % transmit index

for dec1=transmiter % for each transmiter and detector

uu=uu+1;

u=0; % detector index, reset for each detector

for dec2=detector %

% create data for a single transmision and receving

% otain the p_sc_exact : scatter field measured at detector; (1 value)

% p_sc : scatter field predicted at detector; (1 value)

% p: presure calculated (predict at each pixel);matrix (NxN)

B_SC=BB_SC(1+u*N:N+u*N,:);

B=BB(1+u*N:N+u*N,:);

p_inc=PINC(1+uu*N:N+uu*N,:);

u3=0; % index cho C, reset for each new pixel

for n1=1:N % di vao tung pixel, xac dinh boi 2D: n1,n2

for n2=1:N

C_SC=CC_SC(1+u3*N:N+u3*N,:);

C_SC1=CC_SC1(1+u3*N:N+u3*N,:);

pp(n1,n2)=sum(sum(C_SC.*SC.*p_inc)); % C va p_inc trong moi se khac

pp1(n1,n2)=sum(sum(C_SC1.*SC1.*p_inc));

u3=u3+1;


end;

end;


p =p_inc+pp; % pp<>0; wave equation; % IDEAL

p1=p_inc+pp1; % at the first step p_inc=p1 because SC=zeros, pp1=0; PRDEICT

p_sc_exact=sum(sum(B.*SC.*p)); % using p matrix, IDEAL,1 point, IDEAL

p_sc=sum(sum(B.*SC1.*p)); % using p matrix, Predict,1 point

u=u+1;

% B change depends on detector



delta_sc=p_sc_exact-p_sc; % 1 valuse

M=reshape(B.*p1,1,N*N);

Mt=[Mt;M]; % add more detector

p_sc_exact_t=[p_sc_exact_t ; p_sc_exact]; % NCG

p_sc_t=[p_sc_t p_sc];

delta_sc_t=[delta_sc_t;delta_sc]; % add more scatter field in detector

end; % end one cycle of transmit abd detect

end;


if noise_flag==0;

n=0.05*sqrt(var(delta_sc_t))*randn(size(delta_sc_t));

%n=0.05*sqrt(var(p_sc_t))*randn(size(delta_sc_t))

noise_flag=1;

end;

if noise_flag==2;



n=zeros(size(delta_sc_t));

end;


% add noise

if iter==1

delta_sc_t=delta_sc_t+n;

end;


%%%%%%%

[n1,n2]=size(Mt);

gama=1.8755e-021; % co dinh gama

%%%%%%%%
delta_sound = test_NCG(Mt,delta_sc_t,1e3,RRE,gama);


SC1=SC1+reshape(delta_sound,N,N); % update the sound contrast

figure;


subplot(211);imagesc(abs(SC1))

subplot(212);mesh(abs(SC1))

err1(iter)=sum(sum(abs(SC-abs(SC1))))/sum(sum(SC));

MSE1(iter)=(1/N^2)*sum(sum(abs(SC-abs(SC1))^2));

%tinh chi so Quality

[nn1,nn2]=size(SC1);

N=nn1

xx1=(1/N^2)*sum(sum(SC));



xx=xx1*ones(N,N);

yy1=(1/N^2)*sum(sum(abs(SC1)));

yy=xx1*ones(N,N);

detal_y=sum(sum((abs(SC1)-yy)^2/(N^2-1)));

detal_x=sum(sum((abs(SC)-xx)^2/(N^2-1)));

detal_xy=(1/(N^2-1))*sum(sum((abs(SC)-xx)*(abs(SC1)-yy)));

Q=4*detal_xy*xx1*yy1/((detal_x+detal_y)*(xx1^2+yy1^2));

Quality(iter)=Q;

%het

end; % interation for convergen of sound contrast



newSC1=SC1(floor(length(SC1)/2),1:length(SC1));

newSC=SC(floor(length(SC)/2),1:length(SC));

figure()

plot(newSC);

hold on;

plot(abs(newSC1),'r');

figure

plot(1:Niter,err1,'-ro');



xlabel('interation');ylabel('error');

figure


plot(1:Niter,MSE1,'-ro');

xlabel('interation');ylabel('MSE');





1   2   3   4   5   6


Cơ sở dữ liệu được bảo vệ bởi bản quyền ©tieuluan.info 2017
được sử dụng cho việc quản lý

    Quê hương