CT图像三维重建附源码.doc
- 文档编号:1226548
- 上传时间:2023-04-30
- 格式:DOC
- 页数:8
- 大小:40.50KB
CT图像三维重建附源码.doc
《CT图像三维重建附源码.doc》由会员分享,可在线阅读,更多相关《CT图像三维重建附源码.doc(8页珍藏版)》请在冰点文库上搜索。
结束
计算重建图像和原图的性噪比,并进行输出
用函数ifanbeam根据扇束投影数据重建图像,并展示
用函数fanbeam进行映射,得到扇束的数据,并展示
对图片信息进行预处理,并进行展示
输入图片数字选择
生成128的图片信息
程序流图:
MATLAB源码:
clc;
clearall;
closeall;
%loadmri%载入mri数据,是matlab自带库
%ph=squeeze(D);%压缩载入的数据D,并赋值给ph
ph=phantom3d(128);
prompt={'EnterthePiecenum(1to128):
'};%提示信息“输入1到27的片的数字”
name='Inputnumber';%弹出框名称
defaultanswer={'1'};%默认数字
numInput=inputdlg(prompt,name,1,defaultanswer)%弹出框,并得到用户的输入信息
P=squeeze(ph(:
:
str2num(cell2mat(numInput))));%将用户的输入信息转换成数字,并从ph中得到相应的片信息P
imshow(P)%展示图片P
D=250;%将D赋值为250,是从扇束顶点到旋转中心的像素距离。
dsensor1=2;%正实数指定扇束传感器的间距2
F1=fanbeam(P,D,'FanSensorSpacing',dsensor1);%通过P,D等计算扇束的数据值
dsensor2=1;%正实数指定扇束传感器的间距1
F2=fanbeam(P,D,'FanSensorSpacing',dsensor2);%通过P,D等计算扇束的数据值
dsensor3=0.25%正实数指定扇束传感器的间距0.25
[F3,sensor_pos3,fan_rot_angles3]=fanbeam(P,D,...
'FanSensorSpacing',dsensor3);%通过P,D等计算扇束的数据值,并得到扇束传感器的位置sensor_pos3和旋转角度fan_rot_angles3
figure,%创建窗口
imagesc(fan_rot_angles3,sensor_pos3,F3)%根据计算出的位置和角度展示F3的图片
colormap(hot);%设置色图为hot
colorbar;%显示色栏
xlabel('FanRotationAngle(degrees)')%定义x坐标轴
ylabel('FanSensorPosition(degrees)')%定义y坐标轴
output_size=max(size(P));%得到P维数的最大值,并赋值给output_size
Ifan1=ifanbeam(F1,D,...
'FanSensorSpacing',dsensor1,'OutputSize',output_size);
%根据扇束投影数据F1及D等数据重建图像
figure,imshow(Ifan1)%创建窗口,并展示图片Ifan1
title('图一');
disp('图一和原图的性噪比为:
');
result=psnr1(Ifan1,P);
Ifan2=ifanbeam(F2,D,...
'FanSensorSpacing',dsensor2,'OutputSize',output_size);
%根据扇束投影数据F2及D等数据重建图像
figure,imshow(Ifan2)%创建窗口,并展示图片Ifan2
disp('图二和原图的性噪比为:
');
result=psnr1(Ifan2,P);
title('图二');
Ifan3=ifanbeam(F3,D,...
'FanSensorSpacing',dsensor3,'OutputSize',output_size);
%根据扇束投影数据F3及D等数据重建图像
figure,imshow(Ifan3)%创建窗口,并展示图片Ifan3
title('图三');
disp('图三和原图的性噪比为:
');
result=psnr1(Ifan3,P);
function[p,ellipse]=phantom3d(varargin)
%PHANTOM3DThree-dimensionalanalogueofMATLABShepp-Loganphantom
%P=PHANTOM3D(DEF,N)generatesa3Dheadphantomthatcan
%beusedtotest3-Dreconstructionalgorithms.
%
%DEFisastringthatspecifiesthetypeofheadphantomtogenerate.
%Validvaluesare:
%
%'Shepp-Logan'Atestimageusedwidelybyresearchersin
%tomography
%'ModifiedShepp-Logan'(default)AvariantoftheShepp-Loganphantom
%inwhichthecontrastisimprovedforbetter
%visualperception.
%
%NisascalarthatspecifiesthegridsizeofP.
%Ifyouomittheargument,Ndefaultsto64.
%
%P=PHANTOM3D(E,N)generatesauser-definedphantom,whereeachrow
%ofthematrixEspecifiesanellipsoidintheimage.Ehastencolumns,
%witheachcolumncontainingadifferentparameterfortheellipsoids:
%
%Column1:
Atheadditiveintensityvalueoftheellipsoid
%Column2:
athelengthofthexsemi-axisoftheellipsoid
%Column3:
bthelengthoftheysemi-axisoftheellipsoid
%Column4:
cthelengthofthezsemi-axisoftheellipsoid
%Column5:
x0thex-coordinateofthecenteroftheellipsoid
%Column6:
y0they-coordinateofthecenteroftheellipsoid
%Column7:
z0thez-coordinateofthecenteroftheellipsoid
%Column8:
phiphiEulerangle(indegrees)(rotationaboutz-axis)
%Column9:
thetathetaEulerangle(indegrees)(rotationaboutx-axis)
%Column10:
psipsiEulerangle(indegrees)(rotationaboutz-axis)
%
%Forpurposesofgeneratingthephantom,thedomainsforthex-,y-,and
%z-axesspan[-1,1].Columns2through7mustbespecifiedinterms
%ofthisrange.
%
%[P,E]=PHANTOM3D(...)returnsthematrixEusedtogeneratethephantom.
%
%ClassSupport
%-------------
%Allinputsmustbeofclassdouble.Alloutputsareofclassdouble.
%
%Remarks
%-------
%Foranygivenvoxelintheoutputimage,thevoxel'svalueisequaltothe
%sumoftheadditiveintensityvaluesofallellipsoidsthatthevoxelisa
%partof.Ifavoxelisnotpartofanyellipsoid,itsvalueis0.
%
%TheadditiveintensityvalueAforanellipsoidcanbepositiveornegative;
%ifitisnegative,theellipsoidwillbedarkerthanthesurroundingpixels.
%Notethat,dependingonthevaluesofA,somevoxelsmayhavevaluesoutside
%therange[0,1].
%
%Example
%-------
%ph=phantom3d(128);
%figure,imshow(squeeze(ph(64,:
:
)))
%
%Copyright2005MatthiasChristianSchabel(matthias@stanfordalumni.org)
%UniversityofUtahDepartmentofRadiology
%UtahCenterforAdvancedImagingResearch
%729ArapeenDrive
%SaltLakeCity,UT84108-1218
%
%ThiscodeisreleasedundertheGnuPublicLicense(GPL).Formoreinformation,
%see:
http:
//www.gnu.org/copyleft/gpl.html
%
%Portionsofthiscodearebasedonphantom.m,copyrightedbytheMathworks
%
[ellipse,n]=parse_inputs(varargin{:
});
p=zeros([nnn]);
rng=((0:
n-1)-(n-1)/2)/((n-1)/2);
[x,y,z]=meshgrid(rng,rng,rng);
coord=[flatten(x);flatten(y);flatten(z)];
p=flatten(p);
fork=1:
size(ellipse,1)
A=ellipse(k,1);%Amplitudechangeforthisellipsoid
asq=ellipse(k,2)^2;%a^2
bsq=ellipse(k,3)^2;%b^2
csq=ellipse(k,4)^2;%c^2
x0=ellipse(k,5);%xoffset
y0=ellipse(k,6);%yoffset
z0=ellipse(k,7);%zoffset
phi=ellipse(k,8)*pi/180;%firstEulerangleinradians
theta=ellipse(k,9)*pi/180;%secondEulerangleinradians
psi=ellipse(k,10)*pi/180;%thirdEulerangleinradians
cphi=cos(phi);
sphi=sin(phi);
ctheta=cos(theta);
stheta=sin(theta);
cpsi=cos(psi);
spsi=sin(psi);
%Eulerrotationmatrix
alpha=[cpsi*cphi-ctheta*sphi*spsicpsi*sphi+ctheta*cphi*spsispsi*stheta;
-spsi*cphi-ctheta*sphi*cpsi-spsi*sphi+ctheta*cphi*cpsicpsi*stheta;
stheta*sphi-stheta*cphictheta];
%rotatedellipsoidcoordinates
coordp=alpha*coord;
idx=find((coordp(1,:
)-x0).^2./asq+(coordp(2,:
)-y0).^2./bsq+(coordp(3,:
)-z0).^2./csq<=1);
p(idx)=p(idx)+A;
end
p=reshape(p,[nnn]);
return;
functionout=flatten(in)
out=reshape(in,[1prod(size(in))]);
return;
function[e,n]=parse_inputs(varargin)
%eisthem-by-10arraywhichdefinesellipsoids
%nisthesizeofthephantombrainimage
n=128;%Thedefaultsize
e=[];
defaults={'shepp-logan','modifiedshepp-logan','yu-ye-wang'};
fori=1:
nargin
ifischar(varargin{i})%Lookforadefaultphantom
def=lower(varargin{i});
idx=strmatch(def,defaults);
ifisempty(idx)
eid=sprintf('Images:
%s:
unknownPhantom',mfilename);
msg='Unknowndefaultphantomselected.';
error(eid,'%s',msg);
end
switchdefaults{idx}
case'shepp-logan'
e=shepp_logan;
case'modifiedshepp-logan'
e=modified_shepp_logan;
case'yu-ye-wang'
e=yu_ye_wang;
end
elseifnumel(varargin{i})==1
n=varargin{i};%ascalaristheimagesize
elseifndims(varargin{i})==2&&size(varargin{i},2)==10
e=varargin{i};%userspecifiedphantom
else
eid=sprintf('Images:
%s:
invalidInputArgs',mfilename);
msg='Invalidinputarguments.';
error(eid,'%s',msg);
end
end
%ellipseisnotyetdefined
ifisempty(e)
e=modified_shepp_logan;
end
return;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Defaultheadphantoms:
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
functione=shepp_logan
e=modified_shepp_logan;
e(:
1)=[1-.98-.02-.02.01.01.01.01.01.01];
return;
functione=modified_shepp_logan
%
%ThisheadphantomisthesameastheShepp-Loganexcept
%theintensitiesarechangedtoyieldhighercontrastin
%theimage.TakenfromToft,199-200.
%
%Aabcx0y0z0phithetapsi
%-----------------------------------------------------------------
e=[1.6900.920.810000000
-.8.6624.874.7800-.01840000
-.2.1100.310.220.2200-18010
-.2.1600.410.280-.220018010
.1.2100.250.4100.35-.15000
.1.0460.046.0500.1.25000
.1.0460.046.0500-.1.25000
.1.0460.023.050-.08-.6050000
.1.0230.023.0200-.6060000
.1.0230.046.020
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- CT 图像 三维重建 源码
![提示](https://static.bingdoc.com/images/bang_tan.gif)