#!/usr/bin/env python
import time
import numpy as np
import tomograpy
import lo
# object
obj = tomograpy.centered_cubic_map(10, 64)
obj[:] = tomograpy.phantom.shepp_logan(obj.shape)
# data 
radius = 200
a = tomograpy.fov(obj, radius)
data = tomograpy.centered_stack(a, 128, n_images=60, radius=radius, max_lon=np.pi)
# model
kwargs = {"pb":"pb", "obj_rmin":1.5, "data_rmin":1.5}
P, D, obj_mask, data_mask = tomograpy.models.thomson(data, obj, u=.5, **kwargs)
# projection
t = time.time()
data[:] = (P * obj.ravel()).reshape(data.shape)
print("projection time : " + str(time.time() - t))
# data
# backprojection
t = time.time()
x0 = P.T * data.ravel()
bpj = x0.reshape(obj.shape)
print("backprojection time : " + str(time.time() - t))
# inversion using scipy.sparse.linalg
t = time.time()
sol = lo.acg(P, data.ravel(), D, 1e-3 * np.ones(3),  maxiter=100, tol=1e-8)
sol = sol.reshape(obj.shape)
print("inversion time : " + str(time.time() - t))

[(7, 'tomograpy.centered_cubic_map', 'tomograpy.centered_cubic_map', (['(10)', '(64)'], {}), False, 'import tomograpy\n'), (8, 'tomograpy.phantom.shepp_logan', 'tomograpy.phantom.shepp_logan', (['obj.shape'], {}), False, 'import tomograpy\n'), (11, 'tomograpy.fov', 'tomograpy.fov', (['obj', 'radius'], {}), False, 'import tomograpy\n'), (12, 'tomograpy.centered_stack', 'tomograpy.centered_stack', (['a', '(128)'], {'n_images': '(60)', 'radius': 'radius', 'max_lon': 'np.pi'}), False, 'import tomograpy\n'), (15, 'tomograpy.models.thomson', 'tomograpy.models.thomson', (['data', 'obj'], {'u': '(0.5)'}), False, 'import tomograpy\n'), (17, 'time.time', 'time.time', ([], {}), False, 'import time\n'), (22, 'time.time', 'time.time', ([], {}), False, 'import time\n'), (27, 'time.time', 'time.time', ([], {}), False, 'import time\n'), (28, 'numpy.ones', 'np.ones', (['(3)'], {}), True, 'import numpy as np\n'), (19, 'time.time', 'time.time', ([], {}), False, 'import time\n'), (25, 'time.time', 'time.time', ([], {}), False, 'import time\n'), (30, 'time.time', 'time.time', ([], {}), False, 'import time\n')]