diff options
| -rw-r--r-- | Wrappers/Python/ccpi/reconstruction/astra_ops.py | 9 | ||||
| -rw-r--r-- | Wrappers/Python/test/simple_mc_demo.py | 38 |
2 files changed, 44 insertions, 3 deletions
diff --git a/Wrappers/Python/ccpi/reconstruction/astra_ops.py b/Wrappers/Python/ccpi/reconstruction/astra_ops.py index cf138b4..1346f22 100644 --- a/Wrappers/Python/ccpi/reconstruction/astra_ops.py +++ b/Wrappers/Python/ccpi/reconstruction/astra_ops.py @@ -87,7 +87,7 @@ class AstraProjectorMC(Operator): device=device) # Initialise empty for singular value. - self.s1 = None + self.s1 = 50 def direct(self, IM): self.fp.setInput(IM) @@ -103,8 +103,11 @@ class AstraProjectorMC(Operator): # astra.data2d.delete(self.proj_id) def get_max_sing_val(self): - self.s1, sall, svec = PowerMethodNonsquare(self,10) - return self.s1 + if self.s1 is None: + self.s1, sall, svec = PowerMethodNonsquare(self,10) + return self.s1 + else: + return self.s1 def size(self): # Only implemented for 2D diff --git a/Wrappers/Python/test/simple_mc_demo.py b/Wrappers/Python/test/simple_mc_demo.py index 1888740..a6959f9 100644 --- a/Wrappers/Python/test/simple_mc_demo.py +++ b/Wrappers/Python/test/simple_mc_demo.py @@ -94,3 +94,41 @@ for k in range(3): axarro[k].imshow(out2.as_array()[:,:,0,k],vmin=0,vmax=3500) plt.show() +# Create least squares object instance with projector and data. +f = Norm2sq(Aop,b,c=0.5) + +# Initial guess +x_init = VolumeData(np.zeros(x.shape),geometry=vg) + +# FISTA options +opt = {'tol': 1e-4, 'iter': 200} + +# Run FISTA for least squares without regularization +x_fista0, it0, timing0, criter0 = FISTA(x_init, f, None, opt) + + +ff0, axarrf0 = plt.subplots(1,numchannels) +for k in numpy.arange(3): + axarrf0[k].imshow(x_fista0.as_array()[:,:,0,k],vmin=0,vmax=2.5) +plt.show() + +plt.semilogy(criter0) +plt.title('Criterion vs iterations, least squares') +plt.show() + +# Now least squares plus 1-norm regularization +lam = 0.1 +g0 = Norm1(lam) + + +# Run FISTA for least squares plus 1-norm function. +x_fista1, it1, timing1, criter1 = FISTA(x_init, f, g0, opt) + +ff1, axarrf1 = plt.subplots(1,numchannels) +for k in numpy.arange(3): + axarrf1[k].imshow(x_fista1.as_array()[:,:,0,k],vmin=0,vmax=2.5) +plt.show() + +plt.semilogy(criter1) +plt.title('Criterion vs iterations, least squares plus 1-norm regu') +plt.show()
\ No newline at end of file |