《CT原理与图像重建.docx》由会员分享,可在线阅读,更多相关《CT原理与图像重建.docx(4页珍藏版)》请在三一办公上搜索。
1、CT原理与图像重建Radon变换入门matlab CT原理 简介 图像投影,就是说将图像在某一方向上做线性积分。如果将图像看成二维函数f(x, y),则其投影就是在特定方向上的线性积分,比如f(x, y)在垂直方向上的线性积分就是其在x轴上的投影;f(x, y)在水平方向上的线积分就是其在y轴上的投影。通过这些投影,可以获取图像在指定方向上的突出特性,这在图像模式识别等处理中可能会用到。 Radon变换,就是将数字图像矩阵在某一指定角度射线方向上做投影变换。这就是说可以沿着任意角度theta来做Radon变换。 实例 % 实验Radon变换 % By lyqmath % Dalian Univ
2、ersity of Technology % School of Mathematical Sciences clc; clear all; close all; I = zeros(256, 256); r, c = size(I); I(floor(1/5*r:4/5*r), floor(3/5*c:4/5*c) = 1; figure; subplot(2, 2, 1); imshow(I); title(原图像); R, xt = radon(I, 0 45 90); % 在0、45、90度方向做radon变换 subplot(2, 2, 2); plot(xt, R(:, 1); t
3、itle(水平方向的radon变换曲线); subplot(2, 2, 3); plot(xt, R(:, 2); title(45度方向的radon变换曲线); subplot(2, 2, 4); plot(xt, R(:, 3); title(垂直方向的radon变换曲线); 结果: 总结 由于radon变换将图像变换到按角度投影区域,和有名的hough类似,可以应用与检测直线。个人认为,通过将图像矩阵在多角度做积分投影,再对得到的数据做统计分析,可以确定出图像的一些基本性质。 =附= I = zeros(100,100); I(25:75, 25:75) = 1; imshow(I) R
4、,xp = radon(I,0 45); figure; plot(xp,R(:,1); title(R_0o (xprime) theta = 0:180; R,xp = radon(I,theta); imagesc(theta,xp,R); title(R_theta (Xprime); xlabel(theta (degrees); ylabel(Xprime); set(gca,XTick,0:20:180); colormap(hot); colorbar 结果 原图: 变换后: =逆变换= I = zeros(100,100); I(25:75, 25:75) = 1; P=I;
5、 imshow(P) theta1=0:10:170;R1,xp=radon(P,theta1); %存在18个角度投影 theta2=0:5:175;R2,xp=radon(P,theta2); %存在36个角度投影 theta3=0:2:178;R3,xp=radon(P,theta3); %存在90个角度投影 figure,imagesc(theta3,xp,R3);colormap(hot);colorbar; xlabel(theta);ylabel(xprime); % 定义坐标轴 I1=iradon(R1,10); I2=iradon(R2,5); I3=iradon(R3,2)
6、; figure,imshow(I1),title(I1); figure,imshow(I2),title(I2); figure,imshow(I3),title(I3); =反变换= I = imread(mm.JPG); I = rgb2gray(I); P=I; imshow(P) theta1=0:10:170;R1,xp=radon(P,theta1); %存在18个角度投影 theta2=0:5:175;R2,xp=radon(P,theta2); %存在36个角度投影 theta3=0:1:179;R3,xp=radon(P,theta3); %存在90个角度投影 figur
7、e,imagesc(theta3,xp,R3);colormap(hot);colorbar; xlabel(theta);ylabel(xprime); % 定义坐标轴 I1=iradon(R1,10); % 开始反变换,还原图像 I2=iradon(R2,5); I3=iradon(R3,1); figure,imshow(I1),title(I1); figure,imshow(I2),title(I2); figure,imshow(I3),title(I3); k=input(Your Thresh(0,1):) thres1=max(I3(:)*k %maybe 0.5 is perfect figure,imshow(I3thres1),title(I4); iradon函数的使用 iradon:图像处理命令 功能:进行反radon变换 语法: I=iradon(P,theta) I=iradon(P,theta, interp,filter,d,n) I,h=iradon() 举例: P=phantom(128); R=radon(P,0:179); I=iradon(R,0:179,nearest,Hann); imshow(P) figure,imshow(I) 结果: 相关命令:iradon phantom