updating version number

demo_leapctype/d24_circle_plus_line.py

@@ -42,8 +42,10 @@
 ct_circular.set_conebeam(numAngles, numRows, numCols, pixelSize, pixelSize, 0.5*(numRows-1), 0.5*(numCols-1)+10, ct_circular.setAngleArray(numAngles, 360.0), sod, sdd)
 
 # Now we define the line cone-beam data which requires that we use the modular-beam geometry
-numShifts = numRows//16
-shiftSpacing = pixelSize*16
+L = 16
+#L = 32
+numShifts = numRows//L
+shiftSpacing = pixelSize*L
 sourcePositions = np.ascontiguousarray(np.zeros((numShifts,3)).astype(np.float32), dtype=np.float32)
 moduleCenters = np.ascontiguousarray(np.zeros((numShifts,3)).astype(np.float32), dtype=np.float32)
 colVectors = np.ascontiguousarray(np.zeros((numShifts,3)).astype(np.float32), dtype=np.float32)
@@ -83,8 +85,11 @@
 
 
 # Specify simplified FORBILD head phantom
-ct_circular.set_FORBILD(f_true,True)
-
+#ct_circular.set_FORBILD(f_true,True)
+for i in range(11):
+    ct_circular.addObject(f_true, 4, [0.0, 0.0, 20.0*(i-5)], [100.0, 100.0, 5.0], 0.02)
+#ct_circular.display(f_true)
+#quit()
 
 # "Simulate" projection data
 ct_circular.project(g_circular,f_true)
@@ -95,11 +100,12 @@
 
 # Reconstruct the axial cone-beam data with FBP (FDK)
 ct_circular.FBP(g_circular, f)
+#ct_circular.display(f)
 
 # Now refine the reconstruction with the "line" CT data
 # This reconstruction should have no cone-beam artifacts
 # Try commenting out this line to see cone-beam artfiact in the reconstruction
-ct_line.LS(g_line, f, 10, True)
+ct_line.LS(g_line, f, 10)#, 'SQS', True)
 
 # Display the result
 ct_circular.display(f)

demo_leapctype/d99_speedTest.py

@@ -25,14 +25,14 @@
 # Set the scanner geometry
 #leapct.set_parallelbeam(numAngles, M, N, pixelSize, pixelSize, 0.5*(M-1), 0.5*(N-1), leapct.setAngleArray(numAngles, 360.0))
 #leapct.set_fanbeam(numAngles, M, N, pixelSize, pixelSize, 0.5*(M-1), 0.5*(N-1), leapct.setAngleArray(numAngles, 360.0), 1100, 1400)
-leapct.set_conebeam(numAngles, M, N, pixelSize, pixelSize, 0.5*(M-1), 0.5*(N-1), leapct.setAngleArray(numAngles, 360.0), 1100, 1400)
+leapct.set_conebeam(numAngles, M, N, pixelSize, pixelSize, 0.5*(M-1), 0.5*(N-1), leapct.setAngleArray(numAngles, 360.0), 1100, 1400, tiltAngle=0.2)
 
 # Set the volume parameters
 leapct.set_default_volume()
 leapct.print_parameters()
 
 # Set the backprojector model, 'SF' (the default setting), is more accurate, but 'VD' is faster
-#leapct.set_projector('VD')
+leapct.set_projector('VD')
 #leapct.convert_to_modularbeam()
 
 # Allocate space for the projections and the volume
@@ -55,3 +55,8 @@
 startTime = time.time()
 leapct.FBP(g,f)
 print('FBP time: ' + str(time.time()-startTime))
+
+# Filter Projections
+startTime = time.time()
+leapct.filterProjections(g)
+print('filter projections time: ' + str(time.time()-startTime))

setup.py

@@ -40,7 +40,7 @@
 
 setup(
     name='leapct',
-    version='1.24', 
+    version='1.25', 
     author='Kyle Champley, Hyojin Kim', 
     author_email='champley@gmail.com, hkim@llnl.gov', 
     description='LivermorE AI Projector for Computed Tomography (LEAPCT)', 

setup_AMD.py

@@ -113,7 +113,7 @@
 
 setup(
     name='leapct',
-    version='1.24', 
+    version='1.25', 
     author='Kyle Champley, Hyojin Kim', 
     author_email='champley@gmail.com, hkim@llnl.gov', 
     description='LivermorE AI Projector for Computed Tomography (LEAPCT)', 

setup_ctype.py

@@ -40,7 +40,7 @@
 
 setup(
     name='leapct',
-    version='1.24', 
+    version='1.25', 
     author='Kyle Champley, Hyojin Kim', 
     author_email='champley@gmail.com, hkim@llnl.gov', 
     description='LivermorE AI Projector for Computed Tomography (LEAPCT)', 

src/leapctype.py

@@ -857,7 +857,7 @@ def find_centerCol(self, g, iRow=-1, searchBounds=None):
            
         respectively.  For rays near the mid-plane, one can also use these cost functions for cone-parallel and cone-beam coordinates as well.
 
-        Note that it only works for parallel-, fan-, and cone-beam CT geometry types (i.e., everything but modular-beam)
+        Note that this only works for parallel-, fan-, and cone-beam CT geometry types (i.e., everything but modular-beam)
         and one may not get an accurate estimate if the projections are truncated on the right and/or left sides.
         If you have any bad edge detectors, these must be cropped out before running this algorithm.
         If this function does not return a good estimate, try changing the iRow parameter value or try using the
@@ -868,6 +868,7 @@ def find_centerCol(self, g, iRow=-1, searchBounds=None):
         Args:
             g (C contiguous float32 numpy array or torch tensor): projection data
             iRow (int): The detector row index to be used for the estimation, if no value is given, uses the row closest to the centerRow parameter
+            searchBounds (2-element array): optional argument to specify the interval for which to perform the search
             
         Returns:
             the error metric value
@@ -901,7 +902,7 @@ def find_tau(self, g, iRow=-1, searchBounds=None):
            
         For rays near the mid-plane, one can also use these cost functions for cone-beam coordinates as well.
 
-        Note that it only works for fan- and cone-beam CT geometry types
+        Note that this only works for fan- and cone-beam CT geometry types
         and one may not get an accurate estimate if the projections are truncated on the right and/or left sides.
         If you have any bad edge detectors, these must be cropped out before running this algorithm.
         If this function does not return a good estimate, try changing the iRow parameter value or try using the
@@ -912,6 +913,7 @@ def find_tau(self, g, iRow=-1, searchBounds=None):
         Args:
             g (C contiguous float32 numpy array or torch tensor): projection data
             iRow (int): The detector row index to be used for the estimation, if no value is given, uses the row closest to the centerRow parameter
+            searchBounds (2-element array): optional argument to specify the interval for which to perform the search
             
         Returns:
             the error metric value

src/parameters.cpp

@@ -942,7 +942,7 @@ void parameters::printAll()
 		if (tau != 0.0)
 			printf("tau = %f mm\n", tau);
 		if (geometry == CONE)
-			printf("gamma = %f degrees\n", tiltAngle);
+			printf("tiltAngle = %f degrees\n", tiltAngle);
 	}
 	else if (geometry == MODULAR)
 	{

src/tomographic_models.cpp

@@ -1512,6 +1512,11 @@ bool tomographicModels::doFBP(float* g, float* f, bool data_on_cpu)
 			params.offsetScan = false;
 	}
 	//*/
+	if (params.helicalPitch != 0.0 && params.tiltAngle != 0.0)
+	{
+		printf("Error: current implementation of helical FBP cannot handle detector tilt!\n");
+		return false;
+	}
 
 	//printf("v range: %f to %f\n", params.v(0), params.v(params.numRows-1));
 	if (data_on_cpu == true && FBP_multiGPU(g, f) == true)
@@ -2410,8 +2415,8 @@ bool tomographicModels::beam_hardening_heel_effect(float* g, float* anode_normal
 				float v = params.v(j,i) / sdd_cur;
 				float u = params.u(k,i) / sdd_cur;
 				float lineLength = sqrt(1.0 + u * u + v * v);
-				//float takeoffAngle = 90.0 - acos((u * anode_normal[0] + 1.0 * anode_normal[1] + v * anode_normal[2]) / lineLength) * 180.0 / PI;
-				float takeoffAngle = acos((u * anode_normal[0] + 1.0 * anode_normal[1] + v * anode_normal[2]) / lineLength) * 180.0 / PI;
+				float takeoffAngle = 90.0 - acos((u * anode_normal[0] + 1.0 * anode_normal[1] + v * anode_normal[2]) / lineLength) * 180.0 / PI;
+				//float takeoffAngle = acos((u * anode_normal[0] + 1.0 * anode_normal[1] + v * anode_normal[2]) / lineLength) * 180.0 / PI;
 				//if (i == 0 && j == params.numRows/2)
 				//	printf("takeoffAngle = %f (u=%f)\n", takeoffAngle, u);
 

src/tomographic_models.h

@@ -13,7 +13,7 @@
 #pragma once
 #endif
 
-#define LEAP_VERSION "1.24"
+#define LEAP_VERSION "1.25"
 
 /*
 #include <iostream>
@@ -35,7 +35,7 @@
 /**
  *  tomographicModels class
  * This is the main interface for LEAP.  All calls to forward project, backproject, FBP, filtering, noise filters come through this class.
- * The main job of this class is to set/get parameters, do error checks, and dispath jobs.  It contains almost no algorithm logic.
+ * The main job of this class is to set/get parameters, do error checks, and dispatch jobs.  It contains almost no algorithm logic.
  * In addition to the jobs listed above, this class is also responsible for divide jobs across multiple GPUs or dividing up GPU jobs so that
  * they fit into the available GPU memory.  Functions called from this class are either CPU based or single GPU based.
  */

unitTests/unit_tests.py

@@ -37,7 +37,7 @@
 #voxelScales = [2.0]
 #projection_methods = ['SF']
 #projection_methods = ['VD']
-#geometries = []
+geometries = []
 
 for ii in range(len(geometries)):
     igeom = geometries[ii]