Abstract: An image is produced of a target region from tomographic measurements through statistical inversion. The target region is mapped with a grid of voxels. Unknowns in the voxels are the components of a vector-valued variable. A next version of the vector-valued variable is repeatedly generated by drawing for each component a value from a conditional distribution. Elements of a theory matrix indicate, what is the contribution of each voxel in the grid to a detection result. The non-zero elements of the theory matrix are stored in a compressed form. A proper subset of the elements of the theory matrix is decompressed for either updating the quantity Ax and/or calculating the quantity AjT(Ax), where A is the theory matrix, x is the vector-valued variable and AjT is the transpose of the j:th column of the theory matrix.