clear 
clc

load RawReconPt_Rotn36p0572.mat
RRR(RRR<0) = 0;
SSS = find(RRR==0);

Tsize = 241;
Thalf = (Tsize-1)/2;
Tcent = Thalf + 1;
PZ = Tcent;
Th = 3;

KKK = fftshift(fftn(RRR));

Rin = 45;
Rout = 75;
FFF = zeros(Tsize,Tsize,Tsize);
for a = 1:Tsize
    for b = 1:Tsize
        for c = 1:Tsize
            if (sqrt((a-Tcent)^2+(b-Tcent)^2+(c-Tcent)^2)>=Rin) & (sqrt((a-Tcent)^2+(b-Tcent)^2+(c-Tcent)^2)<=Rout)
                FFF(a,b,c) = 1;
            end
        end
    end
end

KKK = fftshift(fftn(RRR));
KKK = KKK.*FFF;
Amax = max(max(max(abs(KKK))))

GGG = abs(KKK);  
Cutoff = 0.1 * Amax
GGG(GGG<Cutoff) = 0;
GGG(GGG~=0) = 1;

ra = 1;
tc = ra+1;
box = zeros(2*ra+1,2*ra+1);
for a = 1:2*ra+1
    for b = 1:2*ra+1
        for c = 1:2*ra+1
            if sqrt((a-tc)^2+(b-tc)^2+(c-tc)^2)<=ra
                box(a,b,c) = 1;
            end
        end
    end
end
MMM = zeros(Tsize,Tsize,Tsize);
for a = 1:Tsize
    for b = 1:Tsize
        for c = 1:Tsize
            if GGG(a,b,c)==1;
               BBB = MMM(a-ra:a+ra,b-ra:b+ra,c-ra:c+ra);
               BBB(box==1) = 1; 
               MMM(a-ra:a+ra,b-ra:b+ra,c-ra:c+ra) = BBB;
            end
        end
    end
end
GGG = MMM;

KKK = KKK.*GGG;
PPPN = real(ifftn(ifftshift(KKK))); 
FFF = log10(abs(KKK));
PPPN(SSS) = 0;

Chunk = zeros(241,241,241);
Chunk(102:106,:,:) = PPPN(102:106,:,:);
save Chunk1_FourierFilter_10perCutoff.mat Chunk
Chunk = zeros(241,241,241);
Chunk(107:111,:,:) = PPPN(107:111,:,:);
save Chunk2_FourierFilter_10perCutoff.mat Chunk
