changed default values

2 files changed, 50 insertions, 31 deletions

diff --git a/Wrappers/Python/ccpi/reconstruction/geoms.py b/Wrappers/Python/ccpi/reconstruction/geoms.py
index 9f81fa0..a887c0c 100644
--- a/Wrappers/Python/ccpi/reconstruction/geoms.py
+++ b/Wrappers/Python/ccpi/reconstruction/geoms.py
@@ -2,12 +2,12 @@
 class VolumeGeometry:
     
     def __init__(self, 
-                 voxel_num_x=None, 
-                 voxel_num_y=None, 
-                 voxel_num_z=None, 
-                 voxel_size_x=None, 
-                 voxel_size_y=None, 
-                 voxel_size_z=None, 
+                 voxel_num_x=0, 
+                 voxel_num_y=0, 
+                 voxel_num_z=0, 
+                 voxel_size_x=1, 
+                 voxel_size_y=1, 
+                 voxel_size_z=1, 
                  center_x=0, 
                  center_y=0, 
                  center_z=0, 
@@ -37,10 +37,16 @@ class VolumeGeometry:
         return self.center_y + 0.5*self.voxel_num_y*self.voxel_size_y
         
     def getMinZ(self):
-        return self.center_z - 0.5*self.voxel_num_z*self.voxel_size_z
+        if not voxel_num_z == 0:
+            return self.center_z - 0.5*self.voxel_num_z*self.voxel_size_z
+        else:
+            return 0
         
     def getMaxZ(self):
-        return self.center_z + 0.5*self.voxel_num_z*self.voxel_size_z
+        if not voxel_num_z == 0:
+            return self.center_z + 0.5*self.voxel_num_z*self.voxel_size_z
+        else:
+            return 0
         
     
 class SinogramGeometry:
@@ -49,9 +55,9 @@ class SinogramGeometry:
                  geom_type, 
                  dimension, 
                  angles, 
-                 pixel_num_h=None, 
+                 pixel_num_h=0, 
                  pixel_size_h=1, 
-                 pixel_num_v=None, 
+                 pixel_num_v=0, 
                  pixel_size_v=1, 
                  dist_source_center=None, 
                  dist_center_detector=None, 
diff --git a/Wrappers/Python/test/simple_mc_demo.py b/Wrappers/Python/wip/simple_mc_demo.py
index a6959f9..54c9134 100644..100755
--- a/Wrappers/Python/test/simple_mc_demo.py
+++ b/Wrappers/Python/wip/simple_mc_demo.py
@@ -11,44 +11,57 @@ from ccpi.reconstruction.ops import *
 from ccpi.reconstruction.astra_ops import *
 from ccpi.reconstruction.geoms import *
 
-import numpy as np
+import numpy
 import matplotlib.pyplot as plt
 
-test_case = 2   # 1=parallel2D, 2=cone2D
+test_case = 1   # 1=parallel2D, 2=cone2D
 
 # Set up phantom
 N = 128
-
+M = 100
 numchannels = 3
 
-x = np.zeros((N,N,1,numchannels))
+vg = VolumeGeometry(voxel_num_x=N,voxel_num_y=N,channels=numchannels)
+Phantom = VolumeData(geometry=vg)
+
+x = Phantom.as_array()
+x[0 , round(N/4):round(3*N/4) , round(N/4):round(3*N/4)  ] = 1.0
+x[0 , round(N/8):round(7*N/8) , round(3*N/8):round(5*N/8)] = 2.0
+
+x[1 , round(N/4):round(3*N/4) , round(N/4):round(3*N/4)  ] = 0.7
+x[1 , round(N/8):round(7*N/8) , round(3*N/8):round(5*N/8)] = 1.2
+
+x[2 , round(N/4):round(3*N/4) , round(N/4):round(3*N/4)  ] = 1.5
+x[2 , round(N/8):round(7*N/8) , round(3*N/8):round(5*N/8)] = 2.2
+
+#x = numpy.zeros((N,N,1,numchannels))
 
-x[round(N/4):round(3*N/4),round(N/4):round(3*N/4),:,0] = 1.0
-x[round(N/8):round(7*N/8),round(3*N/8):round(5*N/8),:,0] = 2.0
+#x[round(N/4):round(3*N/4),round(N/4):round(3*N/4),:,0] = 1.0
+#x[round(N/8):round(7*N/8),round(3*N/8):round(5*N/8),:,0] = 2.0
 
-x[round(N/4):round(3*N/4),round(N/4):round(3*N/4),:,1] = 0.7
-x[round(N/8):round(7*N/8),round(3*N/8):round(5*N/8),:,1] = 1.2
+#x[round(N/4):round(3*N/4),round(N/4):round(3*N/4),:,1] = 0.7
+#x[round(N/8):round(7*N/8),round(3*N/8):round(5*N/8),:,1] = 1.2
 
-x[round(N/4):round(3*N/4),round(N/4):round(3*N/4),:,2] = 1.5
-x[round(N/8):round(7*N/8),round(3*N/8):round(5*N/8),:,2] = 2.2
+#x[round(N/4):round(3*N/4),round(N/4):round(3*N/4),:,2] = 1.5
+#x[round(N/8):round(7*N/8),round(3*N/8):round(5*N/8),:,2] = 2.2
 
 f, axarr = plt.subplots(1,numchannels)
 for k in numpy.arange(3):
-    axarr[k].imshow(x[:,:,0,k],vmin=0,vmax=2.5)
+    axarr[k].imshow(x[k],vmin=0,vmax=2.5)
 plt.show()
 
-vg = VolumeGeometry(N,N,None, 1,1,None,channels=numchannels)
+#vg = VolumeGeometry(N,N,None, 1,1,None,channels=numchannels)
 
 
-Phantom = VolumeData(x,geometry=vg)
+#Phantom = VolumeData(x,geometry=vg)
 
 # Set up measurement geometry
 angles_num = 20; # angles number
 
 if test_case==1:
-    angles = np.linspace(0,np.pi,angles_num,endpoint=False)
+    angles = numpy.linspace(0,numpy.pi,angles_num,endpoint=False)
 elif test_case==2:
-    angles = np.linspace(0,2*np.pi,angles_num,endpoint=False)
+    angles = numpy.linspace(0,2*numpy.pi,angles_num,endpoint=False)
 else:
     NotImplemented
 
@@ -76,7 +89,7 @@ elif test_case==2:
                           channels=numchannels)
 
 # ASTRA operator using volume and sinogram geometries
-Aop = AstraProjectorMC(vg, pg, 'gpu')
+Aop = AstraProjectorMC(vg, pg, 'cpu')
 
 
 # Try forward and backprojection
@@ -84,21 +97,21 @@ b = Aop.direct(Phantom)
 
 fb, axarrb = plt.subplots(1,numchannels)
 for k in numpy.arange(3):
-    axarrb[k].imshow(b.as_array()[:,:,0,k],vmin=0,vmax=250)
+    axarrb[k].imshow(b.as_array()[k],vmin=0,vmax=250)
 plt.show()
 
 out2 = Aop.adjoint(b)
 
 fo, axarro = plt.subplots(1,numchannels)
 for k in range(3):
-    axarro[k].imshow(out2.as_array()[:,:,0,k],vmin=0,vmax=3500)
+    axarro[k].imshow(out2.as_array()[k],vmin=0,vmax=3500)
 plt.show()
 
 # Create least squares object instance with projector and data.
 f = Norm2sq(Aop,b,c=0.5)
 
 # Initial guess
-x_init = VolumeData(np.zeros(x.shape),geometry=vg)
+x_init = VolumeData(numpy.zeros(x.shape),geometry=vg)
 
 # FISTA options
 opt = {'tol': 1e-4, 'iter': 200}
@@ -109,7 +122,7 @@ x_fista0, it0, timing0, criter0 = FISTA(x_init, f, None, opt)
 
 ff0, axarrf0 = plt.subplots(1,numchannels)
 for k in numpy.arange(3):
-    axarrf0[k].imshow(x_fista0.as_array()[:,:,0,k],vmin=0,vmax=2.5)
+    axarrf0[k].imshow(x_fista0.as_array()[k],vmin=0,vmax=2.5)
 plt.show()
 
 plt.semilogy(criter0)
@@ -126,7 +139,7 @@ x_fista1, it1, timing1, criter1 = FISTA(x_init, f, g0, opt)
 
 ff1, axarrf1 = plt.subplots(1,numchannels)
 for k in numpy.arange(3):
-    axarrf1[k].imshow(x_fista1.as_array()[:,:,0,k],vmin=0,vmax=2.5)
+    axarrf1[k].imshow(x_fista1.as_array()[k],vmin=0,vmax=2.5)
 plt.show()
 
 plt.semilogy(criter1)