summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
-rw-r--r--NOTICE.txt15
-rwxr-xr-xWrappers/Python/ccpi/optimisation/algorithms/Algorithm.py30
-rw-r--r--Wrappers/Python/ccpi/optimisation/functions/FunctionOperatorComposition.py64
3 files changed, 82 insertions, 27 deletions
diff --git a/NOTICE.txt b/NOTICE.txt
new file mode 100644
index 0000000..c107329
--- /dev/null
+++ b/NOTICE.txt
@@ -0,0 +1,15 @@
+CCPi Core Imaging Library (CIL).
+Copyright 2017 Rutherford Appleton Laboratory STFC
+Copyright 2017 University of Manchester
+
+This software product is developed for the Collaborative Computational
+Project in Tomographic Imaging CCPi (http://www.ccpi.ac.uk/) at RAL STFC (http://www.stfc.ac.uk), University of Manchester
+and other contributing institutions.
+
+Main contributors in alphabetical order:
+Evelina Ametova (UoM)
+Jakob Jorgensen (UoM)
+Daniil Kazantsev (Diamond Light Source)
+Srikanth Nagella (STFC)
+Edoardo Pasca (STFC)
+Evangelos Papoutsellis (UoM)
diff --git a/Wrappers/Python/ccpi/optimisation/algorithms/Algorithm.py b/Wrappers/Python/ccpi/optimisation/algorithms/Algorithm.py
index a14378c..0ba2c55 100755
--- a/Wrappers/Python/ccpi/optimisation/algorithms/Algorithm.py
+++ b/Wrappers/Python/ccpi/optimisation/algorithms/Algorithm.py
@@ -16,7 +16,7 @@
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
-import time
+import time, functools
from numbers import Integral
class Algorithm(object):
@@ -148,14 +148,36 @@ class Algorithm(object):
print ("Stop cryterion has been reached.")
i = 0
+ if verbose:
+ print ("{:>9} {:>10} {:>11} {:>20}".format('Iter',
+ 'Max Iter',
+ 's/Iter',
+ 'Objective Value'))
for _ in self:
if (self.iteration -1) % self.update_objective_interval == 0:
if verbose:
- print ("Iteration {}/{}, = {}".format(self.iteration-1,
- self.max_iteration, self.get_last_objective()) )
+ print (self.verbose_output())
if callback is not None:
callback(self.iteration -1, self.get_last_objective(), self.x)
i += 1
if i == iterations:
break
-
+
+ def verbose_output(self):
+ '''Creates a nice tabulated output'''
+ timing = self.timing[-self.update_objective_interval-1:-1]
+ if len (timing) == 0:
+ t = 0
+ else:
+ t = sum(timing)/len(timing)
+ el = [ self.iteration-1,
+ self.max_iteration,
+ "{:.3f} s/it".format(t),
+ self.get_last_objective() ]
+
+ if type(el[-1] ) == list:
+ string = functools.reduce(lambda x,y: x+' {:>15.5e}'.format(y), el[-1],'')
+ out = "{:>9} {:>10} {:>11} {}".format(*el[:-1] , string)
+ else:
+ out = "{:>9} {:>10} {:>11} {:>20.5e}".format(*el)
+ return out
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)
+
+
+
+
+
+
- )