From c102b119a1dd9444fba0c244ebcfe260cd679a7f Mon Sep 17 00:00:00 2001
From: epapoutsellis <epapoutsellis@gmail.com>
Date: Thu, 4 Apr 2019 17:36:04 +0100
Subject: add methods for precond

---
 .../optimisation/operators/GradientOperator.py     | 74 +++++++++++++++++++---
 .../optimisation/operators/IdentityOperator.py     |  4 +-
 .../optimisation/operators/SparseFiniteDiff.py     |  2 +-
 3 files changed, 69 insertions(+), 11 deletions(-)

(limited to 'Wrappers/Python')

diff --git a/Wrappers/Python/ccpi/optimisation/operators/GradientOperator.py b/Wrappers/Python/ccpi/optimisation/operators/GradientOperator.py
index cd65ee4..e00de0c 100644
--- a/Wrappers/Python/ccpi/optimisation/operators/GradientOperator.py
+++ b/Wrappers/Python/ccpi/optimisation/operators/GradientOperator.py
@@ -107,22 +107,80 @@ class Gradient(LinearOperator):
        
 if __name__ == '__main__':
     
+    
+    from ccpi.optimisation.operators import Identity, BlockOperator
+    
     M, N = 2, 3
     ig = ImageGeometry(M, N)
     arr = ig.allocate('random_int' )
     
-    G_neum = Gradient(ig)
+    # check direct of Gradient and sparse matrix
+    G = Gradient(ig)
+    G_sp = G.matrix()
+    
+    res1 = G.direct(arr)
+    res1y = numpy.reshape(G_sp[0].toarray().dot(arr.as_array().flatten('F')), ig.shape, 'F')
+    
+    print(res1[0].as_array())
+    print(res1y)
+    
+    res1x = numpy.reshape(G_sp[1].toarray().dot(arr.as_array().flatten('F')), ig.shape, 'F')
+    
+    print(res1[1].as_array())
+    print(res1x)    
     
-    d = G_neum.matrix()
-    print(d[1])
+    #check sum abs row
+    conc_spmat = numpy.abs(numpy.concatenate( (G_sp[0].toarray(), G_sp[1].toarray() )))
+    print(numpy.reshape(conc_spmat.sum(axis=0), ig.shape, 'F'))    
+    print(G.sum_abs_row().as_array())
     
-    d1 = G_neum.sum_abs_row()
-    print(d1.as_array())
+    print(numpy.reshape(conc_spmat.sum(axis=1), ((2,) + ig.shape), 'F'))
     
-    d2 = G_neum.sum_abs_col()
-    print(d2)    
+    print(G.sum_abs_col()[0].as_array())
+    print(G.sum_abs_col()[1].as_array())   
     
-    d1 * d2
+    # Check Blockoperator sum abs col and row
     
+    op1 = Gradient(ig)
+    op2 = Identity(ig)
     
+    B = BlockOperator( op1, op2)
     
+    Brow = B.sum_abs_row()
+    Bcol = B.sum_abs_col()
+    
+    concB = numpy.concatenate( (numpy.abs(numpy.concatenate( (G_sp[0].toarray(), G_sp[1].toarray() ))), op2.matrix().toarray()))
+    
+    print(numpy.reshape(concB.sum(axis=0), ig.shape, 'F'))
+    print(Brow.as_array())
+    
+    print(numpy.reshape(concB.sum(axis=1)[0:12], ((2,) + ig.shape), 'F'))
+    print(Bcol[1].as_array())    
+    
+        
+#    print(numpy.concatene(G_sp[0].toarray()+ ))
+#    print(G_sp[1].toarray())
+#    
+#    d1 = G.sum_abs_row()
+#    print(d1.as_array())
+#    
+#    d2 = G_neum.sum_abs_col()
+##    print(d2)    
+#    
+#    
+#    ###########################################################
+    a = BlockDataContainer( BlockDataContainer(arr, arr), arr)
+    b = BlockDataContainer( BlockDataContainer(arr+5, arr+3), arr+2)
+    c = a/b
+    
+    print(c[0][0].as_array(), (arr/(arr+5)).as_array())
+    print(c[0][1].as_array(), (arr/(arr+3)).as_array())
+    print(c[1].as_array(), (arr/(arr+2)).as_array())
+    
+    
+    a1 = BlockDataContainer( arr, BlockDataContainer(arr, arr))
+#    
+#    c1 = arr + a
+#    c2 = arr + a
+#    c2 = a1 + arr
+#    
diff --git a/Wrappers/Python/ccpi/optimisation/operators/IdentityOperator.py b/Wrappers/Python/ccpi/optimisation/operators/IdentityOperator.py
index 52c7c3b..a58a296 100644
--- a/Wrappers/Python/ccpi/optimisation/operators/IdentityOperator.py
+++ b/Wrappers/Python/ccpi/optimisation/operators/IdentityOperator.py
@@ -50,11 +50,11 @@ class Identity(LinearOperator):
     
     def sum_abs_row(self):
         
-        return ImageData(np.array(np.reshape(abs(self.matrix()).sum(axis=0), self.gm_domain.shape, 'F')))
+        return self.gm_domain.allocate(1)#ImageData(np.array(np.reshape(abs(self.matrix()).sum(axis=0), self.gm_domain.shape, 'F')))
  
     def sum_abs_col(self):
         
-        return ImageData(np.array(np.reshape(abs(self.matrix()).sum(axis=1), self.gm_domain.shape, 'F')))
+        return self.gm_domain.allocate(1)#ImageData(np.array(np.reshape(abs(self.matrix()).sum(axis=1), self.gm_domain.shape, 'F')))
             
     
 if __name__ == '__main__':
diff --git a/Wrappers/Python/ccpi/optimisation/operators/SparseFiniteDiff.py b/Wrappers/Python/ccpi/optimisation/operators/SparseFiniteDiff.py
index 0fb5efb..0b5e85f 100644
--- a/Wrappers/Python/ccpi/optimisation/operators/SparseFiniteDiff.py
+++ b/Wrappers/Python/ccpi/optimisation/operators/SparseFiniteDiff.py
@@ -68,7 +68,7 @@ class SparseFiniteDiff():
     
     def sum_abs_col(self):
         
-        return ImageData(np.array(np.reshape(abs(self.matrix()).sum(axis=1), self.gm_domain.shape, 'F')))
+        return ImageData(np.array(np.reshape(abs(self.matrix()).sum(axis=1), self.gm_domain.shape, 'C')))
         
 if __name__ == '__main__':
     
-- 
cgit v1.2.3