PET Operator

[mehrhardt]

---

A SIRF Acquisition Model <=> A CIL Operator A.direct(x) A.adjoint(x)

[paskino]

https://github.com/CCPPETMR/Hackathon-SIRF/blob/48c182ab2a127469384dfff682682839def8405a/src/xSTIR/pSTIR/STIR.py#L1097-L1122

direct = forward backward = adjoint

[mehrhardt]

As discussed with @KrisThielemans we need to change direct slightly. It should contain only the linear part of the forward operator (important: needs to be documented!!!!).

[mehrhardt]

We also want to have subset functionality in some way. It would be good to write in python A[i].direct(x)to project the ith subset.

[paskino]

Specifically this needs to create a AcquisitionModel and not call

    acq_model = AcquisitionModelUsingRayTracingMatrix()
    ## avoid calling these to have only linear term
    acq_model.set_additive_term(add)
    acq_model.set_background_term(bck)

    print('projecting image...')
    # project the image to obtain simulated acquisition data
    # data from raw_data_file is used as a template
    acq_template = AcquisitionData(raw_data_file)
    acq_model.set_up(acq_template, image)
    simulated_data = acq_model.forward(image)

[mehrhardt]

For the subsets, one could "hack" this into python as

class AcquisitionModelSubset():
    def __init__(subset_num, num_subsets):
        self.AcquisitionModel = pet.AcquisitionModel(...)
        self.subset_num = subset_num
        self.num_subsets = num_subsets

    def direct(self, x):
        return self.AcquisitionModel.direct(x, subset_num=self.subset_num, num_subsets=self.num_subsets)

    def adjoint(self x):
        return self.AcquisitionModel.adjoint(x, subset_num=self.subset_num, num_subsets=self.num_subsets) 

then with A = [pet.AcquisitionModelSubset(subset_num=i, num_subsets=m) for i in range(m)] one can iterate over all subsets as

for i in range(m):
    A[i].direct(x)

This solution could be similarly adapted to the c code if people are interested.

[mehrhardt]

I suppose this is done? @paskino