CT平行束滤波反投影重建程序

clear all;
close all;
%%
%获得正向投影
N=400;
p=phantom（256）;
figureimagesc（p）;
title（‘Original image‘）;
theta=0.9:0.9:360;
[Rxp]=radon（ptheta）;
figureimagesc（R）;
title（‘Radon transform image‘）;
%%
%滤波
H=size（R1）;
h=zeros（（H*2-1）1）;
for i=0:H-1
    if i==0
        h（H-i）=1/4;
    elseif rem（i2）==0
        h（H-i）=0;
        h（H+i）=0;
    else
        h（H-i）=-1/（i*pi）^2;
        h（H+i）=-1/（i*pi）^2;
    end
end
x=zeros（HN）;
for i=1:N
    s=R（:i