function pSino=genParalellSinogram(fanbeamsino, cameraParameters)
% first dimension: distance from center
% second dimension: number of angles

% field of camera parameters and the default values
% cameraParameters.numberOfSamples=100
% cameraParameters.rotationSteps=1.8
% cameraParameters.acquisitionTime=1000
% cameraParameters.verticalSteps
% cameraParameters.verticalRange
% cameraParameters.detectorWidth=51.2
% cameraParameters.detectorHeight=0.6
% cameraParameters.detectorPosition=[0 46 1]
% cameraParameters.XRayTubePosition=[0 -169 1]
% cameraParameters.centerOfRotation=[0 0 1]
% cameraParameters.extraRotationSteps
% cameraParameters.sliceNumber
% cameraParameters.numberOfDetectorPixels=128

% Size of a detector pixel in mm
detectorPixelSize=cameraParameters.detectorWidth/cameraParameters.numberOfDetectorPixels;

% Center of each detector pixel in mm
detectorPixelCentersX=[-cameraParameters.numberOfDetectorPixels/2+0.5:cameraParameters.numberOfDetectorPixels/2-0.5].*detectorPixelSize+cameraParameters.detectorPosition(1);

% Angular positions of the detector pixels related to the X-Ray tube
beamAngularPosition=atan((detectorPixelCentersX-cameraParameters.XRayTubePosition(1))/(cameraParameters.detectorPosition(2)-cameraParameters.XRayTubePosition(2)))*180/pi;

% Beam distance from the center of rotation
beamDistance=zeros(1,cameraParameters.numberOfDetectorPixels);
for i=1:cameraParameters.numberOfDetectorPixels
    
    detectorPixelCenter=[detectorPixelCentersX(i) cameraParameters.detectorPosition(2) 1];
    % Equation of the line along p1 and p2 equals to the cross produt of p1 and p2
    beamEquation=[  cameraParameters.XRayTubePosition(2)*detectorPixelCenter(3)-cameraParameters.XRayTubePosition(3)*detectorPixelCenter(2) ...
        cameraParameters.XRayTubePosition(3)*detectorPixelCenter(1)-cameraParameters.XRayTubePosition(1)*detectorPixelCenter(3) ...
        cameraParameters.XRayTubePosition(1)*detectorPixelCenter(2)-cameraParameters.XRayTubePosition(2)*detectorPixelCenter(1)];
    
    % distance between the line and a point is the
    % normalization of their scalar product
    beamDistance(i)=sum(cameraParameters.centerOfRotation.*beamEquation)/sqrt(beamEquation(1)^2+beamEquation(2)^2);
end

% The parallel beam sinogram
pSino=zeros(size(fanbeamsino,1),size(fanbeamsino,2)-cameraParameters.extraRotationSteps*2);

% Normalization matrix (if two different fan beam positions calculated to
% the same parallel position
normMatrix=zeros(size(pSino));


% Calculate parallel beam sinograms from fan beam data
for detpos=1:cameraParameters.numberOfDetectorPixels
    % Distance of the parallel beam to the center of rotation
    parallelDistance=round(beamDistance(detpos)/detectorPixelSize+cameraParameters.numberOfDetectorPixels/2);
    
    for angularPos=1:size(fanbeamsino,2)
        % Angular position of the parallel beam
        parallelAngle=beamAngularPosition(detpos)+(angularPos-1)*cameraParameters.rotationSteps-cameraParameters.extraRotationSteps*2*cameraParameters.rotationSteps;
        % Angular position in pixels
        parallelAngleIdx=round(parallelAngle/cameraParameters.rotationSteps)+1+cameraParameters.extraRotationSteps;
        % Underflow and Overflow and protection of the array indexes
        if (parallelAngleIdx<1 || parallelAngleIdx>size(pSino,2))
            continue;
        end
        
        % Add fan beam data to the parallel sinogram
        pSino(parallelDistance,parallelAngleIdx)=pSino(parallelDistance,parallelAngleIdx)+fanbeamsino(detpos,angularPos);
        % Increment in the normalization matrix
        normMatrix(parallelDistance,parallelAngleIdx)=normMatrix(parallelDistance,parallelAngleIdx)+1;
    end
end

% Normalize the sinogram
% Zero division protection
normMatrix(normMatrix==0)=1;
pSino=pSino./normMatrix;

