diff options
-rw-r--r-- | Wrappers/Python/ccpi/astra/astra_ops.py | 6 | ||||
-rw-r--r-- | Wrappers/Python/ccpi/astra/astra_processors.py | 2 | ||||
-rwxr-xr-x | Wrappers/Python/ccpi/astra/ops.py | 132 | ||||
-rwxr-xr-x | Wrappers/Python/ccpi/astra/processors.py | 205 | ||||
-rwxr-xr-x[-rw-r--r--] | Wrappers/Python/ccpi/astra/utils.py (renamed from Wrappers/Python/ccpi/astra/astra_utils.py) | 0 | ||||
-rwxr-xr-x | Wrappers/Python/conda-recipe/meta.yaml | 2 | ||||
-rwxr-xr-x | Wrappers/Python/wip/demo_sophiabeads.py | 4 | ||||
-rwxr-xr-x | Wrappers/Python/wip/simple_demo_astra.py | 6 | ||||
-rwxr-xr-x | Wrappers/Python/wip/simple_mc_demo.py | 6 |
9 files changed, 350 insertions, 13 deletions
diff --git a/Wrappers/Python/ccpi/astra/astra_ops.py b/Wrappers/Python/ccpi/astra/astra_ops.py index 45b8238..cd0ef9e 100644 --- a/Wrappers/Python/ccpi/astra/astra_ops.py +++ b/Wrappers/Python/ccpi/astra/astra_ops.py @@ -15,12 +15,12 @@ # You should have received a copy of the GNU General Public License # along with this program. If not, see <http://www.gnu.org/licenses/>. -from ccpi.reconstruction.ops import Operator +from ccpi.optimisation.ops import Operator import numpy import astra from ccpi.framework import AcquisitionData, ImageData, DataContainer -from ccpi.reconstruction.ops import PowerMethodNonsquare -from ccpi.astra.astra_processors import AstraForwardProjector, AstraBackProjector, \ +from ccpi.optimisation.ops import PowerMethodNonsquare +from ccpi.astra.processors import AstraForwardProjector, AstraBackProjector, \ AstraForwardProjectorMC, AstraBackProjectorMC class AstraProjectorSimple(Operator): diff --git a/Wrappers/Python/ccpi/astra/astra_processors.py b/Wrappers/Python/ccpi/astra/astra_processors.py index 3de43bf..16c1f78 100644 --- a/Wrappers/Python/ccpi/astra/astra_processors.py +++ b/Wrappers/Python/ccpi/astra/astra_processors.py @@ -1,5 +1,5 @@ from ccpi.framework import DataSetProcessor, ImageData, AcquisitionData -from ccpi.astra.astra_utils import convert_geometry_to_astra +from ccpi.astra.utils import convert_geometry_to_astra import astra diff --git a/Wrappers/Python/ccpi/astra/ops.py b/Wrappers/Python/ccpi/astra/ops.py new file mode 100755 index 0000000..cd0ef9e --- /dev/null +++ b/Wrappers/Python/ccpi/astra/ops.py @@ -0,0 +1,132 @@ +# -*- coding: utf-8 -*- +# This work is independent part of the Core Imaging Library developed by +# Visual Analytics and Imaging System Group of the Science Technology +# Facilities Council, STFC +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see <http://www.gnu.org/licenses/>. + +from ccpi.optimisation.ops import Operator +import numpy +import astra +from ccpi.framework import AcquisitionData, ImageData, DataContainer +from ccpi.optimisation.ops import PowerMethodNonsquare +from ccpi.astra.processors import AstraForwardProjector, AstraBackProjector, \ + AstraForwardProjectorMC, AstraBackProjectorMC + +class AstraProjectorSimple(Operator): + """ASTRA projector modified to use DataSet and geometry.""" + def __init__(self, geomv, geomp, device): + super(AstraProjectorSimple, self).__init__() + + # Store volume and sinogram geometries. + self.sinogram_geometry = geomp + self.volume_geometry = geomv + + self.fp = AstraForwardProjector(volume_geometry=geomv, + sinogram_geometry=geomp, + proj_id=None, + device=device) + + self.bp = AstraBackProjector(volume_geometry=geomv, + sinogram_geometry=geomp, + proj_id=None, + device=device) + + # Initialise empty for singular value. + self.s1 = None + + def direct(self, IM): + self.fp.set_input(IM) + out = self.fp.get_output() + return out + + def adjoint(self, DATA): + self.bp.set_input(DATA) + out = self.bp.get_output() + return out + + #def delete(self): + # astra.data2d.delete(self.proj_id) + + def get_max_sing_val(self): + self.s1, sall, svec = PowerMethodNonsquare(self,10) + return self.s1 + + def size(self): + # Only implemented for 2D + return ( (self.sinogram_geometry.angles.size, \ + self.sinogram_geometry.pixel_num_h), \ + (self.volume_geometry.voxel_num_x, \ + self.volume_geometry.voxel_num_y) ) + + def create_image_data(self): + inputsize = self.size()[1] + return DataContainer(numpy.random.randn(inputsize[0], + inputsize[1])) + + +class AstraProjectorMC(Operator): + """ASTRA Multichannel projector""" + def __init__(self, geomv, geomp, device): + super(AstraProjectorMC, self).__init__() + + # Store volume and sinogram geometries. + self.sinogram_geometry = geomp + self.volume_geometry = geomv + + self.fp = AstraForwardProjectorMC(volume_geometry=geomv, + sinogram_geometry=geomp, + proj_id=None, + device=device) + + self.bp = AstraBackProjectorMC(volume_geometry=geomv, + sinogram_geometry=geomp, + proj_id=None, + device=device) + + # Initialise empty for singular value. + self.s1 = None + + def direct(self, IM): + self.fp.set_input(IM) + out = self.fp.get_output() + return out + + def adjoint(self, DATA): + self.bp.set_input(DATA) + out = self.bp.get_output() + return out + + #def delete(self): + # astra.data2d.delete(self.proj_id) + + def get_max_sing_val(self): + 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 + return ( (self.sinogram_geometry.angles.size, \ + self.sinogram_geometry.pixel_num_h), \ + (self.volume_geometry.voxel_num_x, \ + self.volume_geometry.voxel_num_y) ) + + def create_image_data(self): + inputsize = self.size()[1] + return DataContainer(numpy.random.randn(self.volume_geometry.channels, + inputsize[0], + inputsize[1])) + diff --git a/Wrappers/Python/ccpi/astra/processors.py b/Wrappers/Python/ccpi/astra/processors.py new file mode 100755 index 0000000..16c1f78 --- /dev/null +++ b/Wrappers/Python/ccpi/astra/processors.py @@ -0,0 +1,205 @@ +from ccpi.framework import DataSetProcessor, ImageData, AcquisitionData +from ccpi.astra.utils import convert_geometry_to_astra +import astra + + +class AstraForwardProjector(DataSetProcessor): + '''AstraForwardProjector + + Forward project ImageData to AcquisitionData using ASTRA proj_id. + + Input: ImageData + Parameter: proj_id + Output: AcquisitionData + ''' + + def __init__(self, + volume_geometry=None, + sinogram_geometry=None, + proj_id=None, + device='cpu'): + kwargs = { + 'volume_geometry' : volume_geometry, + 'sinogram_geometry' : sinogram_geometry, + 'proj_id' : proj_id, + 'device' : device + } + + #DataSetProcessor.__init__(self, **kwargs) + super(AstraForwardProjector, self).__init__(**kwargs) + + self.set_ImageGeometry(volume_geometry) + self.set_AcquisitionGeometry(sinogram_geometry) + + # Set up ASTRA Volume and projection geometry, not to be stored in self + vol_geom, proj_geom = convert_geometry_to_astra(self.volume_geometry, + self.sinogram_geometry) + + # ASTRA projector, to be stored + if device == 'cpu': + # Note that 'line' only one option + if self.sinogram_geometry.geom_type == 'parallel': + self.set_projector(astra.create_projector('line', proj_geom, vol_geom) ) + elif self.sinogram_geometry.geom_type == 'cone': + self.set_projector(astra.create_projector('line_fanflat', proj_geom, vol_geom) ) + else: + NotImplemented + elif device == 'gpu': + self.set_projector(astra.create_projector('cuda', proj_geom, vol_geom) ) + else: + NotImplemented + + def check_input(self, dataset): + if dataset.number_of_dimensions == 3 or\ + dataset.number_of_dimensions == 2: + return True + else: + raise ValueError("Expected input dimensions is 2 or 3, got {0}"\ + .format(dataset.number_of_dimensions)) + + def set_projector(self, proj_id): + self.proj_id = proj_id + + def set_ImageGeometry(self, volume_geometry): + self.volume_geometry = volume_geometry + + def set_AcquisitionGeometry(self, sinogram_geometry): + self.sinogram_geometry = sinogram_geometry + + def process(self): + IM = self.get_input() + DATA = AcquisitionData(geometry=self.sinogram_geometry) + #sinogram_id, DATA = astra.create_sino( IM.as_array(), + # self.proj_id) + sinogram_id, DATA.array = astra.create_sino(IM.as_array(), + self.proj_id) + astra.data2d.delete(sinogram_id) + #return AcquisitionData(array=DATA, geometry=self.sinogram_geometry) + return DATA + +class AstraBackProjector(DataSetProcessor): + '''AstraBackProjector + + Back project AcquisitionData to ImageData using ASTRA proj_id. + + Input: AcquisitionData + Parameter: proj_id + Output: ImageData + ''' + + def __init__(self, + volume_geometry=None, + sinogram_geometry=None, + proj_id=None, + device='cpu'): + kwargs = { + 'volume_geometry' : volume_geometry, + 'sinogram_geometry' : sinogram_geometry, + 'proj_id' : proj_id, + 'device' : device + } + + #DataSetProcessor.__init__(self, **kwargs) + super(AstraBackProjector, self).__init__(**kwargs) + + self.set_ImageGeometry(volume_geometry) + self.set_AcquisitionGeometry(sinogram_geometry) + + # Set up ASTRA Volume and projection geometry, not to be stored in self + vol_geom, proj_geom = convert_geometry_to_astra(self.volume_geometry, + self.sinogram_geometry) + + # ASTRA projector, to be stored + if device == 'cpu': + # Note that 'line' only one option + if self.sinogram_geometry.geom_type == 'parallel': + self.set_projector(astra.create_projector('line', proj_geom, vol_geom) ) + elif self.sinogram_geometry.geom_type == 'cone': + self.set_projector(astra.create_projector('line_fanflat', proj_geom, vol_geom) ) + else: + NotImplemented + elif device == 'gpu': + self.set_projector(astra.create_projector('cuda', proj_geom, vol_geom) ) + else: + NotImplemented + + def check_input(self, dataset): + if dataset.number_of_dimensions == 3 or dataset.number_of_dimensions == 2: + return True + else: + raise ValueError("Expected input dimensions is 2 or 3, got {0}"\ + .format(dataset.number_of_dimensions)) + + def set_projector(self, proj_id): + self.proj_id = proj_id + + def set_ImageGeometry(self, volume_geometry): + self.volume_geometry = volume_geometry + + def set_AcquisitionGeometry(self, sinogram_geometry): + self.sinogram_geometry = sinogram_geometry + + def process(self): + DATA = self.get_input() + IM = ImageData(geometry=self.volume_geometry) + rec_id, IM.array = astra.create_backprojection(DATA.as_array(), + self.proj_id) + astra.data2d.delete(rec_id) + return IM + +class AstraForwardProjectorMC(AstraForwardProjector): + '''AstraForwardProjector Multi channel + + Forward project ImageData to AcquisitionDataSet using ASTRA proj_id. + + Input: ImageDataSet + Parameter: proj_id + Output: AcquisitionData + ''' + def check_input(self, dataset): + if dataset.number_of_dimensions == 2 or \ + dataset.number_of_dimensions == 3 or \ + dataset.number_of_dimensions == 4: + return True + else: + raise ValueError("Expected input dimensions is 2 or 3, got {0}"\ + .format(dataset.number_of_dimensions)) + def process(self): + IM = self.get_input() + #create the output AcquisitionData + DATA = AcquisitionData(geometry=self.sinogram_geometry) + + for k in range(DATA.geometry.channels): + sinogram_id, DATA.as_array()[k] = astra.create_sino(IM.as_array()[k], + self.proj_id) + astra.data2d.delete(sinogram_id) + return DATA + +class AstraBackProjectorMC(AstraBackProjector): + '''AstraBackProjector Multi channel + + Back project AcquisitionData to ImageData using ASTRA proj_id. + + Input: AcquisitionData + Parameter: proj_id + Output: ImageData + ''' + def check_input(self, dataset): + if dataset.number_of_dimensions == 2 or \ + dataset.number_of_dimensions == 3 or \ + dataset.number_of_dimensions == 4: + return True + else: + raise ValueError("Expected input dimensions is 2 or 3, got {0}"\ + .format(dataset.number_of_dimensions)) + def process(self): + DATA = self.get_input() + + IM = ImageData(geometry=self.volume_geometry) + + for k in range(IM.geometry.channels): + rec_id, IM.as_array()[k] = astra.create_backprojection( + DATA.as_array()[k], + self.proj_id) + astra.data2d.delete(rec_id) + return IM
\ No newline at end of file diff --git a/Wrappers/Python/ccpi/astra/astra_utils.py b/Wrappers/Python/ccpi/astra/utils.py index 9f8fe46..9f8fe46 100644..100755 --- a/Wrappers/Python/ccpi/astra/astra_utils.py +++ b/Wrappers/Python/ccpi/astra/utils.py diff --git a/Wrappers/Python/conda-recipe/meta.yaml b/Wrappers/Python/conda-recipe/meta.yaml index 12e83e2..ca91ed0 100755 --- a/Wrappers/Python/conda-recipe/meta.yaml +++ b/Wrappers/Python/conda-recipe/meta.yaml @@ -19,7 +19,7 @@ requirements: - python - numpy - scipy - - ccpi-common + - ccpi-framework - astra-toolbox about: diff --git a/Wrappers/Python/wip/demo_sophiabeads.py b/Wrappers/Python/wip/demo_sophiabeads.py index e3c7f3a..33ce9ec 100755 --- a/Wrappers/Python/wip/demo_sophiabeads.py +++ b/Wrappers/Python/wip/demo_sophiabeads.py @@ -3,8 +3,8 @@ from ccpi.io.reader import XTEKReader import numpy as np import matplotlib.pyplot as plt from ccpi.framework import ImageGeometry, AcquisitionGeometry, AcquisitionData, ImageData -from ccpi.astra.astra_ops import AstraProjectorSimple -from ccpi.reconstruction.algs import CGLS +from ccpi.astra.ops import AstraProjectorSimple +from ccpi.optimisation.algs import CGLS # Set up reader object and read the data datareader = XTEKReader("C:/Users/mbbssjj2/Documents/SophiaBeads_256_averaged/SophiaBeads_256_averaged.xtekct") diff --git a/Wrappers/Python/wip/simple_demo_astra.py b/Wrappers/Python/wip/simple_demo_astra.py index 3365504..26f3ff4 100755 --- a/Wrappers/Python/wip/simple_demo_astra.py +++ b/Wrappers/Python/wip/simple_demo_astra.py @@ -2,9 +2,9 @@ #sys.path.append("..") from ccpi.framework import ImageData , ImageGeometry, AcquisitionGeometry -from ccpi.reconstruction.algs import FISTA, FBPD, CGLS -from ccpi.reconstruction.funcs import Norm2sq, Norm1 , TV2D -from ccpi.astra.astra_ops import AstraProjectorSimple +from ccpi.optimisation.algs import FISTA, FBPD, CGLS +from ccpi.optimisation.funcs import Norm2sq, Norm1 , TV2D +from ccpi.astra.ops import AstraProjectorSimple import numpy as np import matplotlib.pyplot as plt diff --git a/Wrappers/Python/wip/simple_mc_demo.py b/Wrappers/Python/wip/simple_mc_demo.py index f77a678..7369d8f 100755 --- a/Wrappers/Python/wip/simple_mc_demo.py +++ b/Wrappers/Python/wip/simple_mc_demo.py @@ -2,9 +2,9 @@ #sys.path.append("..") from ccpi.framework import ImageData, AcquisitionData, ImageGeometry, AcquisitionGeometry -from ccpi.reconstruction.algs import FISTA -from ccpi.reconstruction.funcs import Norm2sq, Norm1 -from ccpi.astra.astra_ops import AstraProjectorMC +from ccpi.optimisation.algs import FISTA +from ccpi.optimisation.funcs import Norm2sq, Norm1 +from ccpi.astra.ops import AstraProjectorMC import numpy import matplotlib.pyplot as plt |