I am trying to use MATLAB to implement a CT (computed tomography) projection operator, A, which I think is also referred as "system matrix" often times.
Basically, for a N x N image M, the projection data, P, can be obtained by multiplication of the project operator to the image:
P = AM
and the backprojection procedure can be performed by multiplying the (conjugate) transpose of the projection operator to the projection data:
M = A'P
Anyone has any idea/example/sample code on how to implement matrix A (for example: Radon transform)? I would really like to start with a small size of matrix, say 8 x 8, or 16 x 16, if possible.
My question really is: how to implement the projection operator, such that by multiplying the operator with an image, I can get the projections, and by multiplying the (conjugate) transpose of the operator with the projections, I can get the original image back.
EDIT:
Particularly, I would like to implement distance-driven projector, in which case beam trajectory (parallel, fan, or etc) would not matter. Very simple example (MATLAB preferred) will be the best for me to start.
As far as I'm aware, there are no freely available implementations of the distance-driven projector/backprojector (it is patented). You can, however, code it yourself without too much difficulty.
Start by reading the papers and understanding what the projector is doing. There are only a few key parts that you need:
The first two are simple geometry. The overlap kernel is described in good detail (and mostly usable pseudocode) in the papers.
Note that you won't wind up with an actual matrix that does the projection. The system would be too large for all but the tiniest examples. Instead, you should write a function that implements the linear operator corresponding to distance-driven projection.