I am getting the following error message: Index in position 1 is invalid. Array indices must be positive integers or logical values. My code is: clc; close all; clear; % The radon transform of Circle Phantom PD=phantom('Modified Shepp-Logan'); % pad the image with zeros padimage = [2,2]; PD= padarray(PD,padimage); % Original image subplot(2,3,1) imagesc(PD); colormap('gray'); title('Circle Phantom') xlabel('X') ylabel('Y') % Theta is 0:180 with the increment of 0.5 degrees freq = 2; thetas = 0:freq:180; % compute sinogram / radon transformation gtheta = length(thetas); gl = size(PD,1); sinogram = zeros(gl,gtheta); % loop for the number of angles for i = 1:length(thetas) tmpImage = imrotate(PD,-thetas(i),'bilinear','crop'); sinogram(:,i) = sum(tmpImage); end subplot(2,3,2) imagesc(sinogram); title('Circle Sinogram'); xlabel('l'); ylabel('\theta'); % The inverse radon transform with only 1 back projection that is theta = 0 thetas=0; Fl = size(sinogram,1); Ftheta = length(thetas); % convert thetas to radians thetas = (pi/180)*thetas; % set up the backprojected image g0 = zeros(Fl,Fl); % find the middle index of the projections Fmid = ceil(Fl/2); % set up the coords of the image [x,y] = meshgrid(ceil(-Fl/2):ceil(Fl/2-1)); % loop over each projection for i = 1:Ftheta % Using the back projection formula rotCoords = Fmid+round(x*sin(thetas(i)) + y*cos(thetas(i))); % % check which coords are in bounds % indices = find((rotCoords > 0) & (rotCoords <= Fl)); % newCoords = rotCoords(indices); % summation rotCoords=floor(abs(rotCoords)); g0 = g0 + sinogram(rotCoords,i)./Ftheta; end subplot(2,3,3); imagesc(g0); title('Simple backprojection') xlabel('X') ylabel('Y') The error message is occurring for line 54: Index in position 1 is invalid. Array indices must be positive integers or logical values. Error in CT_HW3_BP (line 54) g0 = g0 + sinogram(rotCoords,i)./Ftheta; Can anyone tell me what the issue is here?
Kshitij Singh answered .
2025-11-20
You have
Fmid = ceil(Fl/2); [x,y] = meshgrid(ceil(-Fl/2):ceil(Fl/2-1)); rotCoords = Fmid+round(x*sin(thetas(i)) + y*cos(thetas(i)));