  We develop a positronium imaging method for the Jagiellonian PET (J-PET)
scanners based on the time-of-flight maximum likelihood expectation
maximisation (TOF MLEM). The system matrix elements are calculated on-the-fly
for the coincidences comprising two annihilation and one de-excitation photons
that originate from the ortho-positronium (o-Ps) decay. Using the Geant4
library, a Monte Carlo simulation was conducted for four cylindrical 22Na
sources of beta plus decay with diverse o-Ps mean lifetimes, placed
symmetrically inside the two JPET prototypes. The estimated time differences
between the annihilation and the positron emission were aggregated into
histograms (one per voxel), updated by the weights of the activities
reconstructed by TOF MLEM. The simulations were restricted to include only the
o-Ps decays into back-to-back photons, allowing a linear fitting model to be
employed for the estimation of the mean lifetime from each histogram built in
the log scale. To suppress the noise, the exclusion of voxels with activity
below 2-10 percent of the peak was studied. The estimated o-Ps mean lifetimes
were consistent with the simulation and distributed quasi-uniformly at high
MLEM iterations. The proposed positronium imaging technique can be further
upgraded to include various correction factors, as well as be modified
according to realistic o-Ps decay models.
