From b2fc6c70434674d74551c3a6c01ffb3233499312 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Mon, 1 Jul 2013 22:34:11 +0000 Subject: Update version to 1.3 --- samples/s009_projection_matrix.m | 45 ++++++++++++++++++++++++++++++++++++++++ 1 file changed, 45 insertions(+) create mode 100644 samples/s009_projection_matrix.m (limited to 'samples/s009_projection_matrix.m') diff --git a/samples/s009_projection_matrix.m b/samples/s009_projection_matrix.m new file mode 100644 index 0000000..c0df79e --- /dev/null +++ b/samples/s009_projection_matrix.m @@ -0,0 +1,45 @@ +%------------------------------------------------------------------------ +% This file is part of the +% All Scale Tomographic Reconstruction Antwerp Toolbox ("ASTRA-Toolbox") +% +% Copyright: iMinds-Vision Lab, University of Antwerp +% License: Open Source under GPLv3 +% Contact: mailto:astra@ua.ac.be +% Website: http://astra.ua.ac.be +%------------------------------------------------------------------------ + +vol_geom = astra_create_vol_geom(256, 256); +proj_geom = astra_create_proj_geom('parallel', 1.0, 384, linspace2(0,pi,180)); + +% For CPU-based algorithms, a "projector" object specifies the projection +% model used. In this case, we use the "strip" model. +proj_id = astra_create_projector('strip', proj_geom, vol_geom); + +% Generate the projection matrix for this projection model. +% This creates a matrix W where entry w_{i,j} corresponds to the +% contribution of volume element j to detector element i. +matrix_id = astra_mex_projector('matrix', proj_id); + +% Get the projection matrix as a Matlab sparse matrix. +W = astra_mex_matrix('get', matrix_id); + + +% Manually use this projection matrix to do a projection: +P = phantom(256)'; +s = W * P(:); +s = reshape(s, [proj_geom.DetectorCount size(proj_geom.ProjectionAngles, 2)])'; +figure(1), imshow(s,[]); + +% Because Matlab's matrices are stored transposed in memory compared to C++, +% reshaping them to a vector doesn't give the right ordering for multiplication +% with W. We have to take the transpose of the input and output to get the same +% results (up to numerical noise) as using the toolbox directly. + +% Each row of the projection matrix corresponds to a detector element. +% Detector t for angle p is for row 1 + t + p*proj_geom.DetectorCount. +% Each column corresponds to a volume pixel. +% Pixel (x,y) corresponds to column 1 + x + y*vol_geom.GridColCount. + + +astra_mex_projector('delete', proj_id); +astra_mex_matrix('delete', matrix_id); -- cgit v1.2.3