Open GreameLee opened 1 month ago
Which line is the bug in?
img_fbp = tigre.algorithms.fbp(limited_projection,geo,limited_angle)
@GreameLee I am quite sure the error message is much longer, and that it has more information :) Can you paste the entire error trace?
when I run
limited_angle = np.linspace(0, np.pi/2, 180) limited_projection = tigre.Ax(img, geo, limited_angle) print(limited_projection.shape) img_tv = ossart_tv(limited_projection,geo,limited_angle,5) img_tv.shape img_fbp = tigre.algorithms.fbp(limited_projection,geo,limited_angle) img_fbp.shape
ValueError Traceback (most recent call last) Cell In[40], line 6 4 img_tv = ossart_tv(limited_projection,geo,limited_angle,5) 5 img_tv.shape ----> 6 img_fbp = tigre.algorithms.fbp(limited_projection,geo,limited_angle) 7 img_fbp.shape
File c:\Users\Haodong_Li\AppData\Local\anaconda3\envs\pytorch\lib\sitepackages\tigre\algorithms\single_pass_algorithms.py:195, in fbp(proj, geo, angles, *kwargs) 192 gpuids = kwargs["gpuids"] if "gpuids" in kwargs else None 193 proj_filt = filtering(copy.deepcopy(proj), geox, 194 angles, parker=False, verbose=verbose) --> 195 return Atb(proj_filt, geo, angles, gpuids=gpuids) geo.DSO / geo.DSD
ValueError: operands could not be broadcast together with shapes (1,256,256) (180,) '''
I can not reproduce the error.
Can you show a minimal example using the TIGRE head data that also trhows an error, so I can copy paste it?
import tigre
import numpy as np
from tigre.utilities import sample_loader
from tigre.utilities import CTnoise
import tigre.algorithms as algs
from scipy.io import loadmat
import matplotlib.pyplot as plt
from tigre.algorithms.iterative_recon_alg import IterativeReconAlg
from tigre.algorithms.iterative_recon_alg import decorator
from tigre.utilities.im_3d_denoise import im3ddenoise
class OSSART_TV(IterativeReconAlg):
__doc__ = (
"OSSART_TV solves Cone Beam CT image reconstruction using Oriented Subsets\n"
"Simultaneous Algebraic Reconxtruction Technique with TV regularization algorithm\n"
"OSSART_TV(PROJ,GEO,ALPHA,NITER,BLOCKSIZE=20,TVLAMBDA=50,TVITER=50) \n"
"solves the reconstruction problem using the projection data PROJ taken\n"
"over ALPHA angles, corresponding to the geometry descrived in GEO,\n"
"using NITER iterations.\n"
) + IterativeReconAlg.__doc__
def __init__(self, proj, geo, angles, niter, **kwargs):
self.blocksize = 20 if 'blocksize' not in kwargs else kwargs['blocksize']
self.tvlambda = 50 if 'tvlambda' not in kwargs else kwargs['tvlambda']
self.tviter = 50 if 'tviter' not in kwargs else kwargs['tviter']
# these two settings work well for nVoxel=[254,254,254]
IterativeReconAlg.__init__(self, proj, geo, angles, niter, **kwargs)
ossart_tv = decorator(OSSART_TV, name="ossart_tv")
geo = tigre.geometry()
# VARIABLE DESCRIPTION UNITS
# -------------------------------------------------------------------------------------
# Distances
geo.DSD = 1536 # Distance Source Detector (mm)
geo.DSO = 1000 # Distance Source Origin (mm)
# Image parameters
geo.nVoxel = np.array([1, 256, 256]) # number of voxels (vx)
geo.sVoxel = np.array([1, 256, 256]) # total size of the image (mm)
geo.dVoxel = geo.sVoxel / geo.nVoxel # size of each voxel (mm)
print(geo.dVoxel)
# Detector parameters
geo.nDetector = np.array([1, 512]) # number of pixels (px)
geo.dDetector = np.array([geo.dVoxel[0], 0.8]) # size of each pixel (mm)
geo.sDetector = geo.nDetector * geo.dDetector # total size of the detector (mm)
# Offsets
geo.offOrigin = np.array([0, 0, 0]) # Offset of image from origin (mm)
geo.offDetector = np.array([0, 0]) # Offset of Detector (mm)
# MAKE SURE THAT THE DETECTOR PIXELS SIZE IN V IS THE SAME AS THE IMAGE!
geo.mode = "parallel"
#%% Define angles of projection and load phatom image
angles = np.linspace(0, 2 * np.pi, 180)
head = sample_loader.load_head_phantom(geo.nVoxel)
limited_angle = np.linspace(0, np.pi/2, 180)
limited_projection = tigre.Ax(img, geo, limited_angle)
print(limited_projection.shape)
img_tv = ossart_tv(limited_projection,geo,limited_angle,5)
print(img_tv.shape)
img_fbp = tigre.algorithms.fbp(limited_projection,geo,limited_angle)
print(img_fbp.shape)
I see.
try modifying the line that errors to:
return Atb(proj_filt, geo, angles, gpuids=gpuids) * geo.DSO[0] / geo.DSD[0]
This may solve the issue.
But after I modified the code as your comment, when I run single fbp algorithm I can see this bug:
'''
Traceback (most recent call last):
File "tigre_ct.py", line 72, in
I guess there is the problem that you need delete gpuids after TIGRE using each CT reconstruction algorithm
Hi @GreameLee, I am currently traveling but I'll fix this ASAP. I think it's just a mix of vectors/ints there, should be an easy fix, but can't do it until I have full access to a computer with a gpu.
Hi @GreameLee, I am currently traveling but I'll fix this ASAP. I think it's just a mix of vectors/ints there, should be an easy fix, but can't do it until I have full access to a computer with a gpu.
I just try to fix this by myself, it worked if I changed the
return Atb(proj_filt, geo, angles, gpuids=gpuids) * geo.DSO[0] / geo.DSD[0]
to
if isinstance(geo.DSO, int) == False:
a = Atb(proj_filt, geo, angles, gpuids=gpuids)* geo.DSO[0] / geo.DSD[0]
else:
a = Atb(proj_filt, geo, angles, gpuids=gpuids)* geo.DSO / geo.DSD
return a
when I run:
it get the bug, when I run:
it goes well