From fab42230c06ea7b213042e0dea358a2f6a0dcd48 Mon Sep 17 00:00:00 2001
From: epapoutsellis <epapoutsellis@gmail.com>
Date: Thu, 23 May 2019 17:18:42 +0100
Subject: check accuracy
---
.../functions/FunctionOperatorComposition.py | 64 ++++++++++++++--------
1 file changed, 41 insertions(+), 23 deletions(-)
diff --git a/Wrappers/Python/ccpi/optimisation/functions/FunctionOperatorComposition.py b/Wrappers/Python/ccpi/optimisation/functions/FunctionOperatorComposition.py
index 8895f3d..e788d5a 100644
--- a/Wrappers/Python/ccpi/optimisation/functions/FunctionOperatorComposition.py
+++ b/Wrappers/Python/ccpi/optimisation/functions/FunctionOperatorComposition.py
@@ -50,41 +50,59 @@ class FunctionOperatorComposition(Function):
self.operator.adjoint(tmp, out=out)
-
-
- #TODO do not know if we need it
- #def call_adjoint(self, x):
- #
- # return self.function(self.operator.adjoint(x))
-
-
- #def convex_conjugate(self, x):
- #
- # ''' convex_conjugate does not take into account the Operator'''
- # return self.function.convex_conjugate(x)
-
-
-
if __name__ == '__main__':
- from ccpi.framework import ImageGeometry
+ from ccpi.framework import ImageGeometry, AcquisitionGeometry
from ccpi.optimisation.operators import Gradient
from ccpi.optimisation.functions import L2NormSquared
+ from ccpi.astra.ops import AstraProjectorSimple
+ import numpy as np
- M, N, K = 2,3
+ M, N= 50, 50
ig = ImageGeometry(voxel_num_x=M, voxel_num_y = N)
- G = Gradient(ig)
+ detectors = N
+ angles_num = N
+ det_w = 1.0
+
+ angles = np.linspace(0, np.pi, angles_num, endpoint=False)
+ ag = AcquisitionGeometry('parallel',
+ '2D',
+ angles,
+ detectors,det_w)
+
+
+ Aop = AstraProjectorSimple(ig, ag, 'cpu')
+
+ u = ig.allocate('random_int')
+ u1 = ig.allocate('random_int')
+ b = Aop.direct(u1)
+
+
+# G = Gradient(ig)
alpha = 0.5
- f = L2NormSquared()
- f_comp = FunctionOperatorComposition(G, alpha * f)
- x = ig.allocate('random_int')
- print(f_comp.gradient(x).shape
+ f1 = alpha * L2NormSquared(b=b)
+
+ f_comp = FunctionOperatorComposition(f1, Aop)
+
+ print(f_comp(u))
+
+
+ z1 = Aop.direct(u)
+ tmp = 0.5 * ((z1 - b)**2).sum()
+
+
+ print(tmp)
+
+
+
+
+
+
- )
--
cgit v1.2.3