summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
Diffstat
-rwxr-xr-xWrappers/Python/ccpi/framework/BlockDataContainer.py300
-rwxr-xr-xWrappers/Python/ccpi/framework/__init__.py24
-rwxr-xr-x[-rw-r--r--]Wrappers/Python/ccpi/framework/framework.py (renamed from Wrappers/Python/ccpi/framework.py)65
-rwxr-xr-xWrappers/Python/ccpi/optimisation/operators/BlockOperator.py334
-rwxr-xr-xWrappers/Python/ccpi/optimisation/operators/LinearOperator.py19
-rwxr-xr-xWrappers/Python/ccpi/optimisation/operators/Operator.py25
-rwxr-xr-xWrappers/Python/ccpi/optimisation/operators/__init__.py10
-rw-r--r--Wrappers/Python/conda-recipe/conda_build_config.yaml2
-rw-r--r--Wrappers/Python/conda-recipe/meta.yaml7
-rw-r--r--Wrappers/Python/setup.py4
-rwxr-xr-xWrappers/Python/test/test_BlockDataContainer.py210
-rwxr-xr-xWrappers/Python/test/test_DataContainer.py4
12 files changed, 968 insertions, 36 deletions
diff --git a/Wrappers/Python/ccpi/framework/BlockDataContainer.py b/Wrappers/Python/ccpi/framework/BlockDataContainer.py
new file mode 100755
index 0000000..1bfc98c
--- /dev/null
+++ b/Wrappers/Python/ccpi/framework/BlockDataContainer.py
@@ -0,0 +1,300 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Tue Mar 5 16:04:45 2019
+
+@author: ofn77899
+"""
+from __future__ import absolute_import
+from __future__ import division
+from __future__ import print_function
+from __future__ import unicode_literals
+
+import numpy
+from numbers import Number
+import functools
+#from ccpi.framework import AcquisitionData, ImageData
+#from ccpi.optimisation.operators import Operator, LinearOperator
+
+class BlockDataContainer(object):
+ '''Class to hold DataContainers as blocks'''
+ __array_priority__ = 1
+ def __init__(self, *args, **kwargs):
+ '''containers must be passed row by row'''
+ self.containers = args
+ self.index = 0
+ shape = kwargs.get('shape', None)
+ if shape is None:
+ shape = (len(args),1)
+ self.shape = shape
+
+ n_elements = functools.reduce(lambda x,y: x*y, shape, 1)
+ if len(args) != n_elements:
+ raise ValueError(
+ 'Dimension and size do not match: expected {} got {}'
+ .format(n_elements,len(args)))
+
+
+ def __iter__(self):
+ '''BlockDataContainer is Iterable'''
+ return self
+ def next(self):
+ '''python2 backwards compatibility'''
+ return self.__next__()
+ def __next__(self):
+ try:
+ out = self[self.index]
+ except IndexError as ie:
+ raise StopIteration()
+ self.index+=1
+ return out
+
+ def is_compatible(self, other):
+ '''basic check if the size of the 2 objects fit'''
+ if isinstance(other, Number):
+ return True
+ elif isinstance(other, list):
+ # TODO look elements should be numbers
+ for ot in other:
+ if not isinstance(ot, (Number,\
+ numpy.int, numpy.int8, numpy.int16, numpy.int32, numpy.int64,\
+ numpy.float, numpy.float16, numpy.float32, numpy.float64, \
+ numpy.complex)):
+ raise ValueError('List/ numpy array can only contain numbers {}'\
+ .format(type(ot)))
+ return len(self.containers) == len(other)
+ elif isinstance(other, numpy.ndarray):
+ return self.shape == other.shape
+ return len(self.containers) == len(other.containers)
+ def get_item(self, row, col=0):
+ if row > self.shape[0]:
+ raise ValueError('Requested row {} > max {}'.format(row, self.shape[0]))
+ if col > self.shape[1]:
+ raise ValueError('Requested col {} > max {}'.format(col, self.shape[1]))
+
+ index = row*self.shape[1]+col
+ return self.containers[index]
+
+ def add(self, other, *args, **kwargs):
+ assert self.is_compatible(other)
+ out = kwargs.get('out', None)
+ if isinstance(other, Number):
+ return type(self)(*[ el.add(other, out, *args, **kwargs) for el in self.containers])
+ elif isinstance(other, list) or isinstance(other, numpy.ndarray):
+ return type(self)(*[ el.add(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other)])
+ return type(self)(*[ el.add(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other.containers)])
+
+ def subtract(self, other, *args, **kwargs):
+ assert self.is_compatible(other)
+ out = kwargs.get('out', None)
+ if isinstance(other, Number):
+ return type(self)(*[ el.subtract(other, out, *args, **kwargs) for el in self.containers])
+ elif isinstance(other, list) or isinstance(other, numpy.ndarray):
+ return type(self)(*[ el.subtract(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other)])
+ return type(self)(*[ el.subtract(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other.containers)])
+
+ def multiply(self, other, *args, **kwargs):
+ self.is_compatible(other)
+ out = kwargs.get('out', None)
+ if isinstance(other, Number):
+ return type(self)(*[ el.multiply(other, out, *args, **kwargs) for el in self.containers])
+ elif isinstance(other, list):
+ return type(self)(*[ el.multiply(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other)])
+ elif isinstance(other, numpy.ndarray):
+ return type(self)(*[ el.multiply(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other)])
+ return type(self)(*[ el.multiply(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other.containers)])
+
+ def divide(self, other, *args, **kwargs):
+ self.is_compatible(other)
+ out = kwargs.get('out', None)
+ if isinstance(other, Number):
+ return type(self)(*[ el.divide(other, out, *args, **kwargs) for el in self.containers])
+ elif isinstance(other, list) or isinstance(other, numpy.ndarray):
+ return type(self)(*[ el.divide(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other)])
+ return type(self)(*[ el.divide(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other.containers)])
+
+ def power(self, other, *args, **kwargs):
+ assert self.is_compatible(other)
+ out = kwargs.get('out', None)
+ if isinstance(other, Number):
+ return type(self)(*[ el.power(other, out, *args, **kwargs) for el in self.containers])
+ elif isinstance(other, list) or isinstance(other, numpy.ndarray):
+ return type(self)(*[ el.power(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other)])
+ return type(self)(*[ el.power(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other.containers)])
+
+ def maximum(self,other, *args, **kwargs):
+ assert self.is_compatible(other)
+ out = kwargs.get('out', None)
+ if isinstance(other, Number):
+ return type(self)(*[ el.maximum(other, out, *args, **kwargs) for el in self.containers])
+ elif isinstance(other, list) or isinstance(other, numpy.ndarray):
+ return type(self)(*[ el.maximum(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other)])
+ return type(self)(*[ el.maximum(ot, out, *args, **kwargs) for el,ot in zip(self.containers,other.containers)])
+
+ ## unary operations
+ def abs(self, *args, **kwargs):
+ out = kwargs.get('out', None)
+ return type(self)(*[ el.abs(out, *args, **kwargs) for el in self.containers])
+ def sign(self, *args, **kwargs):
+ out = kwargs.get('out', None)
+ return type(self)(*[ el.sign(out, *args, **kwargs) for el in self.containers])
+ def sqrt(self, *args, **kwargs):
+ out = kwargs.get('out', None)
+ return type(self)(*[ el.sqrt(out, *args, **kwargs) for el in self.containers])
+ def conjugate(self, out=None):
+ return type(self)(*[el.conjugate() for el in self.containers])
+
+ ## reductions
+ def sum(self, *args, **kwargs):
+ return numpy.asarray([ el.sum(*args, **kwargs) for el in self.containers])
+ def squared_norm(self):
+ y = numpy.asarray([el.squared_norm() for el in self.containers])
+ return y.sum()
+ def norm(self):
+ return numpy.sqrt(self.squared_norm())
+ def copy(self):
+ '''alias of clone'''
+ return self.clone()
+ def clone(self):
+ return type(self)(*[el.copy() for el in self.containers])
+
+ def __add__(self, other):
+ return self.add( other )
+ # __radd__
+
+ def __sub__(self, other):
+ return self.subtract( other )
+ # __rsub__
+
+ def __mul__(self, other):
+ return self.multiply(other)
+ # __rmul__
+
+ def __div__(self, other):
+ return self.divide(other)
+ # __rdiv__
+ def __truediv__(self, other):
+ return self.divide(other)
+
+ def __pow__(self, other):
+ return self.power(other)
+ # reverse operand
+ def __radd__(self, other):
+ '''Reverse addition
+
+ to make sure that this method is called rather than the __mul__ of a numpy array
+ the class constant __array_priority__ must be set > 0
+ https://docs.scipy.org/doc/numpy-1.15.1/reference/arrays.classes.html#numpy.class.__array_priority__
+ '''
+ return self + other
+ # __radd__
+
+ def __rsub__(self, other):
+ '''Reverse subtraction
+
+ to make sure that this method is called rather than the __mul__ of a numpy array
+ the class constant __array_priority__ must be set > 0
+ https://docs.scipy.org/doc/numpy-1.15.1/reference/arrays.classes.html#numpy.class.__array_priority__
+ '''
+ return (-1 * self) + other
+ # __rsub__
+
+ def __rmul__(self, other):
+ '''Reverse multiplication
+
+ to make sure that this method is called rather than the __mul__ of a numpy array
+ the class constant __array_priority__ must be set > 0
+ https://docs.scipy.org/doc/numpy-1.15.1/reference/arrays.classes.html#numpy.class.__array_priority__
+ '''
+ return self * other
+ # __rmul__
+
+ def __rdiv__(self, other):
+ '''Reverse division
+
+ to make sure that this method is called rather than the __mul__ of a numpy array
+ the class constant __array_priority__ must be set > 0
+ https://docs.scipy.org/doc/numpy-1.15.1/reference/arrays.classes.html#numpy.class.__array_priority__
+ '''
+ return pow(self / other, -1)
+ # __rdiv__
+ def __rtruediv__(self, other):
+ '''Reverse truedivision
+
+ to make sure that this method is called rather than the __mul__ of a numpy array
+ the class constant __array_priority__ must be set > 0
+ https://docs.scipy.org/doc/numpy-1.15.1/reference/arrays.classes.html#numpy.class.__array_priority__
+ '''
+ return self.__rdiv__(other)
+
+ def __rpow__(self, other):
+ '''Reverse power
+
+ to make sure that this method is called rather than the __mul__ of a numpy array
+ the class constant __array_priority__ must be set > 0
+ https://docs.scipy.org/doc/numpy-1.15.1/reference/arrays.classes.html#numpy.class.__array_priority__
+ '''
+ return other.power(self)
+
+ def __iadd__(self, other):
+ '''Inline addition'''
+ if isinstance (other, BlockDataContainer):
+ for el,ot in zip(self.containers, other.containers):
+ el += ot
+ elif isinstance(other, Number):
+ for el in self.containers:
+ el += other
+ elif isinstance(other, list) or isinstance(other, numpy.ndarray):
+ self.is_compatible(other)
+ for el,ot in zip(self.containers, other):
+ el += ot
+ return self
+ # __iadd__
+
+ def __isub__(self, other):
+ '''Inline subtraction'''
+ if isinstance (other, BlockDataContainer):
+ for el,ot in zip(self.containers, other.containers):
+ el -= ot
+ elif isinstance(other, Number):
+ for el in self.containers:
+ el -= other
+ elif isinstance(other, list) or isinstance(other, numpy.ndarray):
+ assert self.is_compatible(other)
+ for el,ot in zip(self.containers, other):
+ el -= ot
+ return self
+ # __isub__
+
+ def __imul__(self, other):
+ '''Inline multiplication'''
+ if isinstance (other, BlockDataContainer):
+ for el,ot in zip(self.containers, other.containers):
+ el *= ot
+ elif isinstance(other, Number):
+ for el in self.containers:
+ el *= other
+ elif isinstance(other, list) or isinstance(other, numpy.ndarray):
+ assert self.is_compatible(other)
+ for el,ot in zip(self.containers, other):
+ el *= ot
+ return self
+ # __imul__
+
+ def __idiv__(self, other):
+ '''Inline division'''
+ if isinstance (other, BlockDataContainer):
+ for el,ot in zip(self.containers, other.containers):
+ el /= ot
+ elif isinstance(other, Number):
+ for el in self.containers:
+ el /= other
+ elif isinstance(other, list) or isinstance(other, numpy.ndarray):
+ assert self.is_compatible(other)
+ for el,ot in zip(self.containers, other):
+ el /= ot
+ return self
+ # __rdiv__
+ def __itruediv__(self, other):
+ '''Inline truedivision'''
+ return self.__idiv__(other)
+
diff --git a/Wrappers/Python/ccpi/framework/__init__.py b/Wrappers/Python/ccpi/framework/__init__.py
new file mode 100755
index 0000000..083f547
--- /dev/null
+++ b/Wrappers/Python/ccpi/framework/__init__.py
@@ -0,0 +1,24 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Tue Mar 5 16:00:18 2019
+
+@author: ofn77899
+"""
+from __future__ import absolute_import
+from __future__ import division
+from __future__ import print_function
+from __future__ import unicode_literals
+
+import numpy
+import sys
+from datetime import timedelta, datetime
+import warnings
+from functools import reduce
+
+from .framework import DataContainer
+from .framework import ImageData, AcquisitionData
+from .framework import ImageGeometry, AcquisitionGeometry
+from .framework import find_key, message
+from .framework import DataProcessor
+from .framework import AX, PixelByPixelDataProcessor
+from .BlockDataContainer import BlockDataContainer
diff --git a/Wrappers/Python/ccpi/framework.py b/Wrappers/Python/ccpi/framework/framework.py
index 24f4ca6..3159cc7 100644..100755
--- a/Wrappers/Python/ccpi/framework.py
+++ b/Wrappers/Python/ccpi/framework/framework.py
@@ -382,7 +382,8 @@ class DataContainer(object):
return self.shape == other.shape
## algebra
- def __add__(self, other , out=None, *args, **kwargs):
+ def __add__(self, other, *args, **kwargs):
+ out = kwargs.get('out', None)
if issubclass(type(other), DataContainer):
if self.check_dimensions(other):
out = self.as_array() + other.as_array()
@@ -639,8 +640,8 @@ class DataContainer(object):
## binary operations
- def pixel_wise_binary(self,pwop, x2 , out=None, *args, **kwargs):
-
+ def pixel_wise_binary(self, pwop, x2, *args, **kwargs):
+ out = kwargs.get('out', None)
if out is None:
if isinstance(x2, (int, float, complex)):
out = pwop(self.as_array() , x2 , *args, **kwargs )
@@ -658,7 +659,8 @@ class DataContainer(object):
elif issubclass(type(out), DataContainer) and issubclass(type(x2), DataContainer):
if self.check_dimensions(out) and self.check_dimensions(x2):
- pwop(self.as_array(), x2.as_array(), out=out.as_array(), *args, **kwargs )
+ kwargs['out'] = out.as_array()
+ pwop(self.as_array(), x2.as_array(), *args, **kwargs )
#return type(self)(out.as_array(),
# deep_copy=False,
# dimension_labels=self.dimension_labels,
@@ -668,14 +670,15 @@ class DataContainer(object):
raise ValueError(message(type(self),"Wrong size for data memory: ", out.shape,self.shape))
elif issubclass(type(out), DataContainer) and isinstance(x2, (int,float,complex)):
if self.check_dimensions(out):
-
- pwop(self.as_array(), x2, out=out.as_array(), *args, **kwargs )
+ kwargs['out']=out.as_array()
+ pwop(self.as_array(), x2, *args, **kwargs )
return out
else:
raise ValueError(message(type(self),"Wrong size for data memory: ", out.shape,self.shape))
elif issubclass(type(out), numpy.ndarray):
if self.array.shape == out.shape and self.array.dtype == out.dtype:
- pwop(self.as_array(), x2 , out=out, *args, **kwargs)
+ kwargs['out'] = out
+ pwop(self.as_array(), x2, *args, **kwargs)
#return type(self)(out,
# deep_copy=False,
# dimension_labels=self.dimension_labels,
@@ -683,26 +686,27 @@ class DataContainer(object):
else:
raise ValueError (message(type(self), "incompatible class:" , pwop.__name__, type(out)))
- def add(self, other , out=None, *args, **kwargs):
- return self.pixel_wise_binary(numpy.add, other, out=out, *args, **kwargs)
+ def add(self, other, *args, **kwargs):
+ return self.pixel_wise_binary(numpy.add, other, *args, **kwargs)
- def subtract(self, other, out=None , *args, **kwargs):
- return self.pixel_wise_binary(numpy.subtract, other, out=out, *args, **kwargs)
+ def subtract(self, other, *args, **kwargs):
+ return self.pixel_wise_binary(numpy.subtract, other, *args, **kwargs)
- def multiply(self, other , out=None, *args, **kwargs):
- return self.pixel_wise_binary(numpy.multiply, other, out=out, *args, **kwargs)
+ def multiply(self, other, *args, **kwargs):
+ return self.pixel_wise_binary(numpy.multiply, other, *args, **kwargs)
- def divide(self, other , out=None ,*args, **kwargs):
- return self.pixel_wise_binary(numpy.divide, other, out=out, *args, **kwargs)
+ def divide(self, other, *args, **kwargs):
+ return self.pixel_wise_binary(numpy.divide, other, *args, **kwargs)
- def power(self, other , out=None, *args, **kwargs):
- return self.pixel_wise_binary(numpy.power, other, out=out, *args, **kwargs)
+ def power(self, other, *args, **kwargs):
+ return self.pixel_wise_binary(numpy.power, other, *args, **kwargs)
- def maximum(self,x2, out=None, *args, **kwargs):
- return self.pixel_wise_binary(numpy.maximum, x2=x2, out=out, *args, **kwargs)
+ def maximum(self, x2, *args, **kwargs):
+ return self.pixel_wise_binary(numpy.maximum, x2, *args, **kwargs)
## unary operations
- def pixel_wise_unary(self,pwop, out=None, *args, **kwargs):
+ def pixel_wise_unary(self, pwop, *args, **kwargs):
+ out = kwargs.get('out', None)
if out is None:
out = pwop(self.as_array() , *args, **kwargs )
return type(self)(out,
@@ -711,23 +715,25 @@ class DataContainer(object):
geometry=self.geometry)
elif issubclass(type(out), DataContainer):
if self.check_dimensions(out):
- pwop(self.as_array(), out=out.as_array(), *args, **kwargs )
+ kwargs['out'] = out.as_array()
+ pwop(self.as_array(), *args, **kwargs )
else:
raise ValueError(message(type(self),"Wrong size for data memory: ", out.shape,self.shape))
elif issubclass(type(out), numpy.ndarray):
if self.array.shape == out.shape and self.array.dtype == out.dtype:
- pwop(self.as_array(), out=out, *args, **kwargs)
+ kwargs['out'] = out
+ pwop(self.as_array(), *args, **kwargs)
else:
raise ValueError (message(type(self), "incompatible class:" , pwop.__name__, type(out)))
- def abs(self, out=None, *args, **kwargs):
- return self.pixel_wise_unary(numpy.abs, out=out, *args, **kwargs)
+ def abs(self, *args, **kwargs):
+ return self.pixel_wise_unary(numpy.abs, *args, **kwargs)
- def sign(self, out=None, *args, **kwargs):
- return self.pixel_wise_unary(numpy.sign , out=out, *args, **kwargs)
+ def sign(self, *args, **kwargs):
+ return self.pixel_wise_unary(numpy.sign, *args, **kwargs)
- def sqrt(self, out=None, *args, **kwargs):
- return self.pixel_wise_unary(numpy.sqrt, out=out, *args, **kwargs)
+ def sqrt(self, *args, **kwargs):
+ return self.pixel_wise_unary(numpy.sqrt, *args, **kwargs)
#def __abs__(self):
# operation = FM.OPERATION.ABS
@@ -735,7 +741,7 @@ class DataContainer(object):
# __abs__
## reductions
- def sum(self, out=None, *args, **kwargs):
+ def sum(self, *args, **kwargs):
return self.as_array().sum(*args, **kwargs)
def squared_norm(self):
'''return the squared euclidean norm of the DataContainer viewed as a vector'''
@@ -756,6 +762,7 @@ class DataContainer(object):
+
class ImageData(DataContainer):
'''DataContainer for holding 2D or 3D DataContainer'''
def __init__(self,
diff --git a/Wrappers/Python/ccpi/optimisation/operators/BlockOperator.py b/Wrappers/Python/ccpi/optimisation/operators/BlockOperator.py
new file mode 100755
index 0000000..b2af8fc
--- /dev/null
+++ b/Wrappers/Python/ccpi/optimisation/operators/BlockOperator.py
@@ -0,0 +1,334 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Thu Feb 14 12:36:40 2019
+
+@author: ofn77899
+"""
+#from ccpi.optimisation.ops import Operator
+import numpy
+from numbers import Number
+import functools
+from ccpi.framework import AcquisitionData, ImageData, BlockDataContainer
+from ccpi.optimisation.operators import Operator, LinearOperator
+
+
+
+class BlockOperator(Operator):
+ '''Class to hold a block operator
+
+ Class to hold a number of Operators in a block.
+ User may specify the shape of the block, by default is a row vector
+ '''
+ def __init__(self, *args, **kwargs):
+ '''
+ Class creator
+
+ Note:
+ Do not include the `self` parameter in the ``Args`` section.
+
+ Args:
+ vararg (Operator): Operators in the block.
+ shape (:obj:`tuple`, optional): If shape is passed the Operators in
+ vararg are considered input in a row-by-row fashion.
+ Shape and number of Operators must match.
+
+ Example:
+ BlockOperator(op0,op1) results in a row block
+ BlockOperator(op0,op1,shape=(1,2)) results in a column block
+ '''
+ self.operators = args
+ shape = kwargs.get('shape', None)
+ if shape is None:
+ shape = (len(args),1)
+ self.shape = shape
+ n_elements = functools.reduce(lambda x,y: x*y, shape, 1)
+ if len(args) != n_elements:
+ raise ValueError(
+ 'Dimension and size do not match: expected {} got {}'
+ .format(n_elements,len(args)))
+ def get_item(self, row, col):
+ if row > self.shape[0]:
+ raise ValueError('Requested row {} > max {}'.format(row, self.shape[0]))
+ if col > self.shape[1]:
+ raise ValueError('Requested col {} > max {}'.format(col, self.shape[1]))
+
+ index = row*self.shape[1]+col
+ return self.operators[index]
+
+ def norm(self):
+ norm = [op.norm() for op in self.operators]
+ b = []
+ for i in range(self.shape[0]):
+ b.append([])
+ for j in range(self.shape[1]):
+ b[-1].append(norm[i*self.shape[1]+j])
+ return numpy.asarray(b)
+
+ def direct(self, x, out=None):
+ shape = self.get_output_shape(x.shape)
+ res = []
+ for row in range(self.shape[0]):
+ for col in range(self.shape[1]):
+ if col == 0:
+ prod = self.get_item(row,col).direct(x.get_item(col))
+ else:
+ prod += self.get_item(row,col).direct(x.get_item(col))
+ res.append(prod)
+ return BlockDataContainer(*res, shape=shape)
+
+ def adjoint(self, x, out=None):
+ shape = self.get_output_shape(x.shape, adjoint=True)
+ res = []
+ for row in range(self.shape[1]):
+ for col in range(self.shape[0]):
+ if col == 0:
+ prod = self.get_item(col,row).adjoint(x.get_item(row))
+ else:
+ prod += self.get_item(col,row).adjoint(x.get_item(row))
+ res.append(prod)
+ return BlockDataContainer(*res, shape=shape)
+
+ def get_output_shape(self, xshape, adjoint=False):
+ sshape = self.shape[1]
+ oshape = self.shape[0]
+ if adjoint:
+ sshape = self.shape[0]
+ oshape = self.shape[1]
+ if sshape != xshape[0]:
+ raise ValueError('Incompatible shapes {} {}'.format(self.shape, xshape))
+ return (oshape, xshape[-1])
+
+
+
+class BlockLinearOperator(BlockOperator):
+ '''Class to hold a block operator
+
+ Class to hold a number of Operators in a block.
+ User may specify the shape of the block, by default is a row vector
+ '''
+ def __init__(self, *args, **kwargs):
+ '''
+ Class creator
+
+ Note:
+ Do not include the `self` parameter in the ``Args`` section.
+
+ Args:
+ vararg (Operator): LinearOperators in the block.
+ shape (:obj:`tuple`, optional): If shape is passed the Operators in
+ vararg are considered input in a row-by-row fashion.
+ Shape and number of Operators must match.
+
+ Example:
+ BlockLinearOperator(op0,op1) results in a row block
+ BlockLinearOperator(op0,op1,shape=(1,2)) results in a column block
+ '''
+ for i,op in enumerate(args):
+ if not op.is_linear():
+ raise ValueError('Operator {} must be LinearOperator'.format(i))
+ super(BlockLinearOperator, self).__init__(*args, **kwargs)
+
+ def adjoint(self, x, out=None):
+ '''Adjoint operation for the BlockOperator
+
+ only available on BlockLinearOperator
+ '''
+ shape = self.get_output_shape(x.shape, adjoint=True)
+ res = []
+ for row in range(self.shape[1]):
+ for col in range(self.shape[0]):
+ if col == 0:
+ prod = self.get_item(col,row).adjoint(x.get_item(row))
+ else:
+ prod += self.get_item(col,row).adjoint(x.get_item(row))
+ res.append(prod)
+ return BlockDataContainer(*res, shape=shape)
+
+
+
+
+
+
+
+
+if __name__ == '__main__':
+ #from ccpi.optimisation.Algorithms import GradientDescent
+ from ccpi.optimisation.algorithms import CGLS
+
+ from ccpi.plugins.ops import CCPiProjectorSimple
+ from ccpi.optimisation.ops import PowerMethodNonsquare
+ from ccpi.optimisation.ops import TomoIdentity
+ from ccpi.optimisation.funcs import Norm2sq, Norm1
+ from ccpi.framework import ImageGeometry, AcquisitionGeometry
+ from ccpi.optimisation.Algorithms import GradientDescent
+ #from ccpi.optimisation.Algorithms import CGLS
+ import matplotlib.pyplot as plt
+
+
+ # Set up phantom size N x N x vert by creating ImageGeometry, initialising the
+ # ImageData object with this geometry and empty array and finally put some
+ # data into its array, and display one slice as image.
+
+ # Image parameters
+ N = 128
+ vert = 4
+
+ # Set up image geometry
+ ig = ImageGeometry(voxel_num_x=N,
+ voxel_num_y=N,
+ voxel_num_z=vert)
+
+ # Set up empty image data
+ Phantom = ImageData(geometry=ig,
+ dimension_labels=['horizontal_x',
+ 'horizontal_y',
+ 'vertical'])
+ Phantom += 0.05
+ # Populate image data by looping over and filling slices
+ i = 0
+ while i < vert:
+ if vert > 1:
+ x = Phantom.subset(vertical=i).array
+ else:
+ x = Phantom.array
+ x[round(N/4):round(3*N/4),round(N/4):round(3*N/4)] = 0.5
+ x[round(N/8):round(7*N/8),round(3*N/8):round(5*N/8)] = 0.94
+ if vert > 1 :
+ Phantom.fill(x, vertical=i)
+ i += 1
+
+
+ perc = 0.02
+ # Set up empty image data
+ noise = ImageData(numpy.random.normal(loc = 0.04 ,
+ scale = perc ,
+ size = Phantom.shape), geometry=ig,
+ dimension_labels=['horizontal_x',
+ 'horizontal_y',
+ 'vertical'])
+ Phantom += noise
+
+ # Set up AcquisitionGeometry object to hold the parameters of the measurement
+ # setup geometry: # Number of angles, the actual angles from 0 to
+ # pi for parallel beam, set the width of a detector
+ # pixel relative to an object pixe and the number of detector pixels.
+ angles_num = 20
+ det_w = 1.0
+ det_num = N
+
+ angles = numpy.linspace(0,numpy.pi,angles_num,endpoint=False,dtype=numpy.float32)*\
+ 180/numpy.pi
+
+ # Inputs: Geometry, 2D or 3D, angles, horz detector pixel count,
+ # horz detector pixel size, vert detector pixel count,
+ # vert detector pixel size.
+ ag = AcquisitionGeometry('parallel',
+ '3D',
+ angles,
+ N,
+ det_w,
+ vert,
+ det_w)
+
+ # Set up Operator object combining the ImageGeometry and AcquisitionGeometry
+ # wrapping calls to CCPi projector.
+ A = CCPiProjectorSimple(ig, ag)
+
+ # Forward and backprojection are available as methods direct and adjoint. Here
+ # generate test data b and some noise
+
+ b = A.direct(Phantom)
+
+
+ #z = A.adjoint(b)
+
+
+ # Using the test data b, different reconstruction methods can now be set up as
+ # demonstrated in the rest of this file. In general all methods need an initial
+ # guess and some algorithm options to be set. Note that 100 iterations for
+ # some of the methods is a very low number and 1000 or 10000 iterations may be
+ # needed if one wants to obtain a converged solution.
+ x_init = ImageData(geometry=ig,
+ dimension_labels=['horizontal_x','horizontal_y','vertical'])
+ X_init = BlockDataContainer(x_init)
+ B = BlockDataContainer(b,
+ ImageData(geometry=ig, dimension_labels=['horizontal_x','horizontal_y','vertical']))
+
+ # setup a tomo identity
+ Ibig = 1e5 * TomoIdentity(geometry=ig)
+ Ismall = 1e-5 * TomoIdentity(geometry=ig)
+
+ # composite operator
+ Kbig = BlockOperator(A, Ibig, shape=(2,1))
+ Ksmall = BlockOperator(A, Ismall, shape=(2,1))
+
+ #out = K.direct(X_init)
+
+ f = Norm2sq(Kbig,B)
+ f.L = 0.00003
+
+ fsmall = Norm2sq(Ksmall,B)
+ f.L = 0.00003
+
+ simplef = Norm2sq(A, b)
+ simplef.L = 0.00003
+
+ gd = GradientDescent( x_init=x_init, objective_function=simplef,
+ rate=simplef.L)
+ gd.max_iteration = 10
+
+ cg = CGLS()
+ cg.set_up(X_init, Kbig, B )
+ cg.max_iteration = 1
+
+ cgsmall = CGLS()
+ cgsmall.set_up(X_init, Ksmall, B )
+ cgsmall.max_iteration = 1
+
+
+ cgs = CGLS()
+ cgs.set_up(x_init, A, b )
+ cgs.max_iteration = 6
+#
+ #out.__isub__(B)
+ #out2 = K.adjoint(out)
+
+ #(2.0*self.c)*self.A.adjoint( self.A.direct(x) - self.b )
+
+ for _ in gd:
+ print ("iteration {} {}".format(gd.iteration, gd.get_current_loss()))
+
+ cg.run(10, lambda it,val: print ("iteration {} objective {}".format(it,val)))
+
+ cgs.run(10, lambda it,val: print ("iteration {} objective {}".format(it,val)))
+
+ cgsmall.run(10, lambda it,val: print ("iteration {} objective {}".format(it,val)))
+ cgsmall.run(10, lambda it,val: print ("iteration {} objective {}".format(it,val)))
+# for _ in cg:
+# print ("iteration {} {}".format(cg.iteration, cg.get_current_loss()))
+#
+# fig = plt.figure()
+# plt.imshow(cg.get_output().get_item(0,0).subset(vertical=0).as_array())
+# plt.title('Composite CGLS')
+# plt.show()
+#
+# for _ in cgs:
+# print ("iteration {} {}".format(cgs.iteration, cgs.get_current_loss()))
+#
+ fig = plt.figure()
+ plt.subplot(1,5,1)
+ plt.imshow(Phantom.subset(vertical=0).as_array())
+ plt.title('Simulated Phantom')
+ plt.subplot(1,5,2)
+ plt.imshow(gd.get_output().subset(vertical=0).as_array())
+ plt.title('Simple Gradient Descent')
+ plt.subplot(1,5,3)
+ plt.imshow(cgs.get_output().subset(vertical=0).as_array())
+ plt.title('Simple CGLS')
+ plt.subplot(1,5,4)
+ plt.imshow(cg.get_output().get_item(0,0).subset(vertical=0).as_array())
+ plt.title('Composite CGLS\nbig lambda')
+ plt.subplot(1,5,5)
+ plt.imshow(cgsmall.get_output().get_item(0,0).subset(vertical=0).as_array())
+ plt.title('Composite CGLS\nsmall lambda')
+ plt.show()
diff --git a/Wrappers/Python/ccpi/optimisation/operators/LinearOperator.py b/Wrappers/Python/ccpi/optimisation/operators/LinearOperator.py
new file mode 100755
index 0000000..d0e7804
--- /dev/null
+++ b/Wrappers/Python/ccpi/optimisation/operators/LinearOperator.py
@@ -0,0 +1,19 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Tue Mar 5 15:57:52 2019
+
+@author: ofn77899
+"""
+
+from ccpi.optimisation.operators import Operator
+
+class LinearOperator(Operator):
+ '''Operator that maps from a space X -> Y'''
+ def is_linear(self):
+ '''Returns if the operator is linear'''
+ return True
+ def adjoint(self,x, out=None):
+ '''returns the adjoint/inverse operation
+
+ only available to linear operators'''
+ raise NotImplementedError
diff --git a/Wrappers/Python/ccpi/optimisation/operators/Operator.py b/Wrappers/Python/ccpi/optimisation/operators/Operator.py
new file mode 100755
index 0000000..ea08b30
--- /dev/null
+++ b/Wrappers/Python/ccpi/optimisation/operators/Operator.py
@@ -0,0 +1,25 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Tue Mar 5 15:55:56 2019
+
+@author: ofn77899
+"""
+
+class Operator(object):
+ '''Operator that maps from a space X -> Y'''
+ def __init__(self, **kwargs):
+ self.scalar = 1
+ def is_linear(self):
+ '''Returns if the operator is linear'''
+ return False
+ def direct(self,x, out=None):
+ raise NotImplementedError
+ def size(self):
+ # To be defined for specific class
+ raise NotImplementedError
+ def norm(self):
+ raise NotImplementedError
+ def range_geometry(self):
+ raise NotImplementedError
+ def domain_geometry(self):
+ raise NotImplementedError
diff --git a/Wrappers/Python/ccpi/optimisation/operators/__init__.py b/Wrappers/Python/ccpi/optimisation/operators/__init__.py
new file mode 100755
index 0000000..088f48c
--- /dev/null
+++ b/Wrappers/Python/ccpi/optimisation/operators/__init__.py
@@ -0,0 +1,10 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Tue Mar 5 15:56:27 2019
+
+@author: ofn77899
+"""
+
+from .Operator import Operator
+from .LinearOperator import LinearOperator
+from .BlockOperator import BlockOperator
diff --git a/Wrappers/Python/conda-recipe/conda_build_config.yaml b/Wrappers/Python/conda-recipe/conda_build_config.yaml
index 96a211f..30c8e9d 100644
--- a/Wrappers/Python/conda-recipe/conda_build_config.yaml
+++ b/Wrappers/Python/conda-recipe/conda_build_config.yaml
@@ -4,5 +4,5 @@ python:
- 3.6
numpy:
# TODO investigage, as it doesn't currently build with cvxp, requires >1.14
- #- 1.12
+ - 1.12
- 1.15
diff --git a/Wrappers/Python/conda-recipe/meta.yaml b/Wrappers/Python/conda-recipe/meta.yaml
index 8ded429..dd3238e 100644
--- a/Wrappers/Python/conda-recipe/meta.yaml
+++ b/Wrappers/Python/conda-recipe/meta.yaml
@@ -11,7 +11,7 @@ build:
test:
requires:
- python-wget
- - cvxpy # [not win]
+ - cvxpy # [ unix and py36 and np115 ]
source_files:
- ./test # [win]
@@ -24,8 +24,9 @@ test:
requirements:
build:
+ - {{ pin_compatible('numpy', max_pin='x.x') }}
- python
- - numpy {{ numpy }}
+ - numpy
- setuptools
run:
@@ -33,7 +34,7 @@ requirements:
- python
- numpy
- scipy
- - matplotlib
+ #- matplotlib
- h5py
about:
diff --git a/Wrappers/Python/setup.py b/Wrappers/Python/setup.py
index eaf124b..630e33e 100644
--- a/Wrappers/Python/setup.py
+++ b/Wrappers/Python/setup.py
@@ -31,7 +31,9 @@ if cil_version == '':
setup(
name="ccpi-framework",
version=cil_version,
- packages=['ccpi' , 'ccpi.io', 'ccpi.optimisation',
+ packages=['ccpi' , 'ccpi.io',
+ 'ccpi.framework', 'ccpi.optimisation',
+ 'ccpi.optimisation.operators',
'ccpi.optimisation.algorithms'],
# Project uses reStructuredText, so ensure that the docutils get
diff --git a/Wrappers/Python/test/test_BlockDataContainer.py b/Wrappers/Python/test/test_BlockDataContainer.py
new file mode 100755
index 0000000..824abf6
--- /dev/null
+++ b/Wrappers/Python/test/test_BlockDataContainer.py
@@ -0,0 +1,210 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Tue Mar 5 16:08:23 2019
+
+@author: ofn77899
+"""
+
+import unittest
+import numpy
+#from ccpi.plugins.ops import CCPiProjectorSimple
+from ccpi.optimisation.ops import PowerMethodNonsquare
+from ccpi.optimisation.ops import TomoIdentity
+from ccpi.optimisation.funcs import Norm2sq, Norm1
+from ccpi.framework import ImageGeometry, AcquisitionGeometry
+from ccpi.framework import ImageData, AcquisitionData
+#from ccpi.optimisation.algorithms import GradientDescent
+from ccpi.framework import BlockDataContainer
+#from ccpi.optimisation.Algorithms import CGLS
+
+class TestBlockDataContainer(unittest.TestCase):
+ def test_BlockDataContainer(self):
+ print ("test block data container")
+ ig0 = ImageGeometry(2,3,4)
+ ig1 = ImageGeometry(12,42,55,32)
+
+ data0 = ImageData(geometry=ig0)
+ data1 = ImageData(geometry=ig1) + 1
+
+ data2 = ImageData(geometry=ig0) + 2
+ data3 = ImageData(geometry=ig1) + 3
+
+ cp0 = BlockDataContainer(data0,data1)
+ cp1 = BlockDataContainer(data2,data3)
+ #
+ a = [ (el, ot) for el,ot in zip(cp0.containers,cp1.containers)]
+ print (a[0][0].shape)
+ #cp2 = BlockDataContainer(*a)
+ cp2 = cp0.add(cp1)
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == 2.)
+ assert (cp2.get_item(1,0).as_array()[0][0][0] == 4.)
+
+ cp2 = cp0 + cp1
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == 2.)
+ assert (cp2.get_item(1,0).as_array()[0][0][0] == 4.)
+ cp2 = cp0 + 1
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 1. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 2., decimal = 5)
+ cp2 = cp0 + [1 ,2]
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 1. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 3., decimal = 5)
+ cp2 += cp1
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , +3. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , +6., decimal = 5)
+
+ cp2 += 1
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , +4. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , +7., decimal = 5)
+
+ cp2 += [-2,-1]
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 2. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 6., decimal = 5)
+
+
+ cp2 = cp0.subtract(cp1)
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == -2.)
+ assert (cp2.get_item(1,0).as_array()[0][0][0] == -2.)
+ cp2 = cp0 - cp1
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == -2.)
+ assert (cp2.get_item(1,0).as_array()[0][0][0] == -2.)
+
+ cp2 = cp0 - 1
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , -1. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 0, decimal = 5)
+ cp2 = cp0 - [1 ,2]
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , -1. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , -1., decimal = 5)
+
+ cp2 -= cp1
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , -3. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , -4., decimal = 5)
+
+ cp2 -= 1
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , -4. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , -5., decimal = 5)
+
+ cp2 -= [-2,-1]
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , -2. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , -4., decimal = 5)
+
+
+ cp2 = cp0.multiply(cp1)
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == 0.)
+ assert (cp2.get_item(1,0).as_array()[0][0][0] == 3.)
+ cp2 = cp0 * cp1
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == 0.)
+ assert (cp2.get_item(1,0).as_array()[0][0][0] == 3.)
+
+ cp2 = cp0 * 2
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 2, decimal = 5)
+ cp2 = 2 * cp0
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 2, decimal = 5)
+ cp2 = cp0 * [3 ,2]
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 2., decimal = 5)
+ cp2 = cp0 * numpy.asarray([3 ,2])
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 2., decimal = 5)
+
+ cp2 = [3,2] * cp0
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 2., decimal = 5)
+ cp2 = numpy.asarray([3,2]) * cp0
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 2., decimal = 5)
+ cp2 = [3,2,3] * cp0
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 2., decimal = 5)
+
+ cp2 *= cp1
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0 , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , +6., decimal = 5)
+
+ cp2 *= 1
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , +6., decimal = 5)
+
+ cp2 *= [-2,-1]
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , -6., decimal = 5)
+
+
+ cp2 = cp0.divide(cp1)
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == 0.)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0], 1./3., decimal=4)
+ cp2 = cp0/cp1
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == 0.)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0], 1./3., decimal=4)
+
+ cp2 = cp0 / 2
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 0.5, decimal = 5)
+ cp2 = cp0 / [3 ,2]
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 0.5, decimal = 5)
+ cp2 = cp0 / numpy.asarray([3 ,2])
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 0.5, decimal = 5)
+ cp3 = numpy.asarray([3 ,2]) / (cp0+1)
+ numpy.testing.assert_almost_equal(cp3.get_item(0,0).as_array()[0][0][0] , 3. , decimal=5)
+ numpy.testing.assert_almost_equal(cp3.get_item(1,0).as_array()[0][0][0] , 1, decimal = 5)
+
+ cp2 += 1
+ cp2 /= cp1
+ # TODO fix inplace division
+
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 1./2 , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 1.5/3., decimal = 5)
+
+ cp2 /= 1
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0.5 , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 0.5, decimal = 5)
+
+ cp2 /= [-2,-1]
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , -0.5/2. , decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , -0.5, decimal = 5)
+ ####
+
+ cp2 = cp0.power(cp1)
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == 0.)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0], 1., decimal=4)
+ cp2 = cp0**cp1
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == 0.)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0], 1., decimal=4)
+
+ cp2 = cp0 ** 2
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0] , 0., decimal=5)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0] , 1., decimal = 5)
+
+ cp2 = cp0.maximum(cp1)
+ assert (cp2.get_item(0,0).as_array()[0][0][0] == cp1.get_item(0,0).as_array()[0][0][0])
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0], cp2.get_item(1,0).as_array()[0][0][0], decimal=4)
+
+
+ cp2 = cp0.abs()
+ numpy.testing.assert_almost_equal(cp2.get_item(0,0).as_array()[0][0][0], 0., decimal=4)
+ numpy.testing.assert_almost_equal(cp2.get_item(1,0).as_array()[0][0][0], 1., decimal=4)
+
+ cp2 = cp0.subtract(cp1)
+ s = cp2.sign()
+ numpy.testing.assert_almost_equal(s.get_item(0,0).as_array()[0][0][0], -1., decimal=4)
+ numpy.testing.assert_almost_equal(s.get_item(1,0).as_array()[0][0][0], -1., decimal=4)
+
+ cp2 = cp0.add(cp1)
+ s = cp2.sqrt()
+ numpy.testing.assert_almost_equal(s.get_item(0,0).as_array()[0][0][0], numpy.sqrt(2), decimal=4)
+ numpy.testing.assert_almost_equal(s.get_item(1,0).as_array()[0][0][0], numpy.sqrt(4), decimal=4)
+
+ s = cp0.sum()
+ numpy.testing.assert_almost_equal(s[0], 0, decimal=4)
+ s0 = 1
+ s1 = 1
+ for i in cp0.get_item(0,0).shape:
+ s0 *= i
+ for i in cp0.get_item(1,0).shape:
+ s1 *= i
+
+ numpy.testing.assert_almost_equal(s[1], cp0.get_item(0,0).as_array()[0][0][0]*s0 +cp0.get_item(1,0).as_array()[0][0][0]*s1, decimal=4)
+ \ No newline at end of file
diff --git a/Wrappers/Python/test/test_DataContainer.py b/Wrappers/Python/test/test_DataContainer.py
index f23179c..47feb95 100755
--- a/Wrappers/Python/test/test_DataContainer.py
+++ b/Wrappers/Python/test/test_DataContainer.py
@@ -174,7 +174,7 @@ class TestDataContainer(unittest.TestCase):
def binary_add(self):
print("Test binary add")
X, Y, Z = 512, 512, 512
- X, Y, Z = 256, 512, 512
+ X, Y, Z = 1024, 512, 512
steps = [timer()]
a = numpy.ones((X, Y, Z), dtype='float32')
steps.append(timer())
@@ -561,4 +561,4 @@ class TestDataContainer(unittest.TestCase):
if __name__ == '__main__':
unittest.main()
- \ No newline at end of file
+
generated by cgit v1.2.3 (git 2.25.1) at 2025-07-08 22:18:01 +0000