diff options
| -rwxr-xr-x | Wrappers/Python/ccpi/plugins/operators/CCPiProjectorSimple.py (renamed from Wrappers/Python/ccpi/plugins/ops.py) | 31 | ||||
| -rw-r--r-- | Wrappers/Python/ccpi/plugins/operators/__init__.py | 1 | ||||
| -rw-r--r-- | Wrappers/Python/ccpi/plugins/processors/__init__.py | 4 | ||||
| -rwxr-xr-x | Wrappers/Python/ccpi/plugins/processors/processors.py (renamed from Wrappers/Python/ccpi/plugins/processors.py) | 102 | ||||
| -rw-r--r-- | Wrappers/Python/setup.py | 3 |
5 files changed, 81 insertions, 60 deletions
diff --git a/Wrappers/Python/ccpi/plugins/ops.py b/Wrappers/Python/ccpi/plugins/operators/CCPiProjectorSimple.py index f2e0dd3..eedcac3 100755 --- a/Wrappers/Python/ccpi/plugins/ops.py +++ b/Wrappers/Python/ccpi/plugins/operators/CCPiProjectorSimple.py @@ -18,7 +18,7 @@ # limitations under the License.
import numpy
-from ccpi.optimisation.ops import Operator, PowerMethodNonsquare
+from ccpi.optimisation.operators import Operator, LinearOperator
from ccpi.framework import ImageData, DataContainer , \
ImageGeometry, AcquisitionGeometry
from ccpi.plugins.processors import CCPiBackwardProjector, \
@@ -39,8 +39,7 @@ class CCPiProjectorSimple(Operator): raise TypeError('Can only handle parallel beam')
# set-up the geometries if compatible
- geoms = setupCCPiGeometries(geomv.voxel_num_x,geomv.voxel_num_y,
- geomv.voxel_num_z, geomp.angles, 0)
+ geoms = setupCCPiGeometries(geomv, geomp, 0)
vg = ImageGeometry(voxel_num_x=geoms['output_volume_x'],
@@ -74,7 +73,8 @@ class CCPiProjectorSimple(Operator): self.bp = CCPiBackwardProjector(image_geometry=vg,
acquisition_geometry=pg,
- output_axes_order=['horizontal_x','horizontal_y','vertical'])
+ output_axes_order=[ ImageGeometry.HORIZONTAL_X,
+ ImageGeometry.HORIZONTAL_Y, ImageGeometry.VERTICAL])
# Initialise empty for singular value.
self.s1 = None
@@ -85,17 +85,21 @@ class CCPiProjectorSimple(Operator): return True
def direct(self, image_data, out=None):
- self.fp.set_input(image_data)
+ self.fp.set_input(image_data.subset(dimensions=['horizontal_x','horizontal_y','vertical']))
if out is None:
out = self.fp.get_output()
+ #return out.subset(dimensions=[AcquisitionGeometry.ANGLE , AcquisitionGeometry.VERTICAL ,
+ # AcquisitionGeometry.HORIZONTAL])
return out
else:
out.fill(self.fp.get_output())
def adjoint(self, acquisition_data, out=None):
- self.bp.set_input(acquisition_data)
+ self.bp.set_input(acquisition_data.subset(dimensions=['angle','vertical','horizontal']))
if out is None:
out = self.bp.get_output()
+ #return out.subset(dimensions=[ImageGeometry.VERTICAL, ImageGeometry.HORIZONTAL_Y,
+ # ImageGeometry.HORIZONTAL_X])
return out
else:
out.fill(self.bp.get_output())
@@ -103,8 +107,13 @@ class CCPiProjectorSimple(Operator): #def delete(self):
# astra.data2d.delete(self.proj_id)
- def get_max_sing_val(self):
- a = PowerMethodNonsquare(self,10)
+ def norm(self):
+ x0 = self.domain_geometry().allocate(ImageGeometry.RANDOM,
+ #dimension_labels=['horizontal_x', 'horizontal_y','vertical']
+ )
+ print (x0.dimension_labels)
+
+ a = LinearOperator.PowerMethod(self,10, x0)
self.s1 = a[0]
return self.s1
@@ -133,8 +142,8 @@ class CCPiProjectorSimple(Operator): self.vg.center_x,
self.vg.center_y,
self.vg.center_z,
- self.vg.channels,
- ['horizontal_x','horizontal_y','vertical'] )
+ self.vg.channels
+ )
def range_geometry(self):
return AcquisitionGeometry(self.ag.geom_type,
@@ -147,5 +156,5 @@ class CCPiProjectorSimple(Operator): self.ag.dist_source_center,
self.ag.dist_center_detector,
self.ag.channels,
- ['angle','vertical','horizontal'])
+ )
diff --git a/Wrappers/Python/ccpi/plugins/operators/__init__.py b/Wrappers/Python/ccpi/plugins/operators/__init__.py new file mode 100644 index 0000000..f5ffd0c --- /dev/null +++ b/Wrappers/Python/ccpi/plugins/operators/__init__.py @@ -0,0 +1 @@ +from .CCPiProjectorSimple import CCPiProjectorSimple
\ No newline at end of file diff --git a/Wrappers/Python/ccpi/plugins/processors/__init__.py b/Wrappers/Python/ccpi/plugins/processors/__init__.py new file mode 100644 index 0000000..e688671 --- /dev/null +++ b/Wrappers/Python/ccpi/plugins/processors/__init__.py @@ -0,0 +1,4 @@ +from .processors import setupCCPiGeometries
+from .processors import CCPiForwardProjector
+from .processors import CCPiBackwardProjector
+from .processors import AcquisitionDataPadder
diff --git a/Wrappers/Python/ccpi/plugins/processors.py b/Wrappers/Python/ccpi/plugins/processors/processors.py index d31c03b..84de79f 100755 --- a/Wrappers/Python/ccpi/plugins/processors.py +++ b/Wrappers/Python/ccpi/plugins/processors/processors.py @@ -22,46 +22,50 @@ from ccpi.framework import DataProcessor, AcquisitionData,\ from ccpi.reconstruction.parallelbeam import alg as pbalg
import numpy
-def setupCCPiGeometries(voxel_num_x, voxel_num_y, voxel_num_z, angles, counter):
+def setupCCPiGeometries(ig, ag, counter):
+ Phantom_ccpi = ig.allocate(dimension_labels=[ImageGeometry.HORIZONTAL_X,
+ ImageGeometry.HORIZONTAL_Y,
+ ImageGeometry.VERTICAL])
- vg = ImageGeometry(voxel_num_x=voxel_num_x,voxel_num_y=voxel_num_y, voxel_num_z=voxel_num_z)
- Phantom_ccpi = ImageData(geometry=vg,
- dimension_labels=['horizontal_x','horizontal_y','vertical'])
- #.subset(['horizontal_x','horizontal_y','vertical'])
- # ask the ccpi code what dimensions it would like
+ voxel_per_pixel = 1
+ angles = ag.angles
+ geoms = pbalg.pb_setup_geometry_from_image(Phantom_ccpi.as_array(),
+ angles,
+ voxel_per_pixel )
- voxel_per_pixel = 1
- geoms = pbalg.pb_setup_geometry_from_image(Phantom_ccpi.as_array(),
- angles,
- voxel_per_pixel )
-
- pg = AcquisitionGeometry('parallel',
- '3D',
- angles,
- geoms['n_h'], 1.0,
- geoms['n_v'], 1.0 #2D in 3D is a slice 1 pixel thick
- )
-
- center_of_rotation = Phantom_ccpi.get_dimension_size('horizontal_x') / 2
- ad = AcquisitionData(geometry=pg,dimension_labels=['angle','vertical','horizontal'])
- geoms_i = pbalg.pb_setup_geometry_from_acquisition(ad.as_array(),
- angles,
- center_of_rotation,
- voxel_per_pixel )
-
- counter+=1
-
- if counter < 4:
- if (not ( geoms_i == geoms )):
- print ("not equal and {0}".format(counter))
- X = max(geoms['output_volume_x'], geoms_i['output_volume_x'])
- Y = max(geoms['output_volume_y'], geoms_i['output_volume_y'])
- Z = max(geoms['output_volume_z'], geoms_i['output_volume_z'])
- return setupCCPiGeometries(X,Y,Z,angles, counter)
+ pg = AcquisitionGeometry('parallel',
+ '3D',
+ angles,
+ geoms['n_h'], 1.0,
+ geoms['n_v'], 1.0 #2D in 3D is a slice 1 pixel thick
+ )
+
+ center_of_rotation = Phantom_ccpi.get_dimension_size('horizontal_x') / 2
+ #ad = AcquisitionData(geometry=pg,dimension_labels=['angle','vertical','horizontal'])
+ ad = pg.allocate(dimension_labels=[AcquisitionGeometry.ANGLE,
+ AcquisitionGeometry.VERTICAL,
+ AcquisitionGeometry.HORIZONTAL])
+ geoms_i = pbalg.pb_setup_geometry_from_acquisition(ad.as_array(),
+ angles,
+ center_of_rotation,
+ voxel_per_pixel )
+
+ counter+=1
+
+ if counter < 4:
+ print (geoms, geoms_i)
+ if (not ( geoms_i == geoms )):
+ print ("not equal and {} {} {}".format(counter, geoms['output_volume_z'], geoms_i['output_volume_z']))
+ X = max(geoms['output_volume_x'], geoms_i['output_volume_x'])
+ Y = max(geoms['output_volume_y'], geoms_i['output_volume_y'])
+ Z = max(geoms['output_volume_z'], geoms_i['output_volume_z'])
+ return setupCCPiGeometries(X,Y,Z,angles, counter)
+ else:
+ print ("happy now {} {} {}".format(counter, geoms['output_volume_z'], geoms_i['output_volume_z']))
+
+ return geoms
else:
- return geoms
- else:
- return geoms_i
+ return geoms_i
class CCPiForwardProjector(DataProcessor):
@@ -102,7 +106,7 @@ class CCPiForwardProjector(DataProcessor): raise ValueError("Expected input dimensions is 2 or 3, got {0}"\
.format(dataset.number_of_dimensions))
- def process(self):
+ def process(self, out=None):
volume = self.get_input()
volume_axes = volume.get_data_axes_order(new_order=self.default_image_axes_order)
@@ -115,12 +119,15 @@ class CCPiForwardProjector(DataProcessor): pixels = pbalg.pb_forward_project(volume.as_array(),
self.acquisition_geometry.angles,
pixel_per_voxel)
- out = AcquisitionData(geometry=self.acquisition_geometry,
- label_dimensions=self.default_acquisition_axes_order)
- out.fill(pixels)
+
+ out = self.acquisition_geometry.allocate(
+ dimension_labels=self.output_axes_order)
out_axes = out.get_data_axes_order(new_order=self.output_axes_order)
if not out_axes == [0,1,2]:
- out.array = numpy.transpose(out.array, out_axes)
+ print ("transpose")
+ pixels = numpy.transpose(pixels, out_axes)
+ out.fill(pixels)
+
return out
else:
raise ValueError('Cannot process cone beam')
@@ -146,7 +153,8 @@ class CCPiBackwardProjector(DataProcessor): output_axes_order=None):
if output_axes_order is None:
# default ccpi projector image storing order
- output_axes_order = ['horizontal_x','horizontal_y','vertical']
+ #output_axes_order = ['horizontal_x','horizontal_y','vertical']
+ output_axes_order = ['vertical', 'horizontal_y','horizontal_x',]
kwargs = {
'image_geometry' : image_geometry,
'acquisition_geometry' : acquisition_geometry,
@@ -165,7 +173,7 @@ class CCPiBackwardProjector(DataProcessor): raise ValueError("Expected input dimensions is 2 or 3, got {0}"\
.format(dataset.number_of_dimensions))
- def process(self):
+ def process(self, out=None):
projections = self.get_input()
projections_axes = projections.get_data_axes_order(new_order=self.default_acquisition_axes_order)
if not projections_axes == [0,1,2]:
@@ -185,9 +193,7 @@ class CCPiBackwardProjector(DataProcessor): center_of_rotation,
pixel_per_voxel
)
- out = ImageData(geometry=self.image_geometry,
- dimension_labels=self.default_image_axes_order)
-
+ out = self.image_geometry.allocate()
out_axes = out.get_data_axes_order(new_order=self.output_axes_order)
if not out_axes == [0,1,2]:
back = numpy.transpose(back, out_axes)
@@ -231,7 +237,7 @@ class AcquisitionDataPadder(DataProcessor): raise ValueError("Expected input dimensions is 2 or 3, got {0}"\
.format(dataset.number_of_dimensions))
- def process(self):
+ def process(self, out=None):
projections = self.get_input()
w = projections.get_dimension_size('horizontal')
delta = w - 2 * self.center_of_rotation
diff --git a/Wrappers/Python/setup.py b/Wrappers/Python/setup.py index f4c42e8..ad12d1c 100644 --- a/Wrappers/Python/setup.py +++ b/Wrappers/Python/setup.py @@ -32,7 +32,8 @@ if cil_version == '': setup( name="ccpi-plugins", version=cil_version, - packages=['ccpi' , 'ccpi.plugins'], + packages=['ccpi' , 'ccpi.plugins', 'ccpi.plugins.operators', + 'ccpi.plugins.processors'], # Project uses reStructuredText, so ensure that the docutils get # installed or upgraded on the target machine |