Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Half-Fan Weighting Artifacts #92

Open
lc82111 opened this issue Aug 9, 2024 · 12 comments
Open

Half-Fan Weighting Artifacts #92

lc82111 opened this issue Aug 9, 2024 · 12 comments

Comments

@lc82111
Copy link

lc82111 commented Aug 9, 2024

I observed a ring artifact in the overlapping region when using half-fan weighting. Are there any solutions?

image

image

cfg = get_config(sdd=990.805, sod=537.5, nAngles=720,
                 nBins=2, nPad=0, detPhyOffsetX=-138.23, detPhyOffsetY=-96,
                 volZDim=128, volXDim=512, volYDim=512, volXSpacing=1, volYSpacing=1, volZSpacing=1.28)
projs1, angles1 = readDcmProj(cfg, './20240717_504/AVG_504_120_30_air_cw.png', './20240717_504/504_120_30_cw/')
angles = np.linspace(-180, 180, cfg.nAngles).astype(np.float32)

leapct = tomographicModels()
leapct.set_conebeam(cfg.nAngles, cfg.detRows, cfg.detCols, cfg.detSpacing, cfg.detSpacing, cfg.detRowCenter, cfg.detColCenter, angles, cfg.SOD, cfg.SDD)
leapct.set_offsetScan(True)
leapct.set_truncatedScan(True)
# leapct.set_default_volume()
leapct.set_volume(cfg.volXDim, cfg.volYDim, cfg.volZDim, cfg.volXSpacing, cfg.volZSpacing, cfg.volOffsetX, cfg.volOffsetY, cfg.volOffsetZ)
leapct.print_parameters()

vol = leapct.allocateVolume(0)
leapct.FBP(projs1,vol)

print(vol.max(), vol.min(), vol.shape) # 0.030418077 -0.008564696 (128, 512, 512)
print(vol[40, 120:400, 120:400].max(), vol[40, 120:400, 120:400].min())  # 0.019195706 -0.00047323745
print(projs1[4].max(), projs1[4].min(), projs1[4].mean())  # 3.6975892 0.0 1.7670326
vol[vol<0] = 0


fig, axe = plt.subplots(2, 4, figsize=(20, 10))
# vmax = 0.018
vmax = 0.019
vmin = 0.015
_ = axe[0,0].imshow(vol[10, 100:412, 100:412], cmap='gray', vmin=vmin, vmax=vmax) #; fig.colorbar(_, ax=axe[0,0])
_ = axe[0,1].imshow(vol[20, 100:412, 100:412], cmap='gray', vmin=vmin, vmax=vmax)# ; fig.colorbar(_, ax=axe[0,1])
_ = axe[0,2].imshow(vol[40, 100:412, 100:412], cmap='gray', vmin=vmin, vmax=vmax)# ; fig.colorbar(_, ax=axe[0,2])
_ = axe[0,3].imshow(vol[50, 100:412, 100:412], cmap='gray', vmin=vmin, vmax=vmax)# ; fig.colorbar(_, ax=axe[0,3])
_ = axe[1,0].imshow(vol[60, 100:412, 100:412],  cmap='gray', vmax=vmax)# ; fig.colorbar(_, ax=axe[1,3])
_ = axe[1,1].imshow(vol[70, 100:412, 100:412],  cmap='gray', vmax=vmax)# ; fig.colorbar(_, ax=axe[1,3])
_ = axe[1,2].imshow(vol[80, 100:412, 100:412],  cmap='gray', vmax=vmax)# ; fig.colorbar(_, ax=axe[1,3])
_ = axe[1,3].imshow(vol[110, 100:412, 100:412],  cmap='gray', vmax=vmax)# ; fig.colorbar(_, ax=axe[1,3])

======== CT Cone-Beam Geometry ========
number of angles: 720
number of detector elements (rows, cols): 704 x 704
angular range: 360.500702 degrees
detector pixel size: 0.615057 mm x 0.615057 mm
center detector pixel: 507.583130, 576.243469
sod = 537.500000 mm
sdd = 990.804993 mm

======== CT Volume ========
number of voxels (x, y, z): 512 x 512 x 128
voxel size: 1.000000 mm x 1.000000 mm x 1.280000 mm
FOV: [-256.000000, 256.000000] x [-256.000000, 256.000000] x [-81.919999, 81.919999]

SDD: 990.805
SOD: 537.5
detColCenter: 576.2434642032332
detCols: 704
detPhyOffsetX: -138.23
detPhyOffsetY: -96
detPhySize: 433
detPixOffsetX: -224.74346420323323
detPixOffsetY: -156.08314087759814
detRowCenter: 507.5831408775981
detRows: 704
detSpacing: 0.6150568181818182
device: cuda:0
nAngles: 720
nBins: 2
nPad: 0
volOffsetX: 0
volOffsetY: 0
volOffsetZ: 0
volPhyXDim: 512
volPhyYDim: 512
volPhyZDim: 163.84
volXDim: 512
volXSpacing: 1
volYDim: 512
volYSpacing: 1
volZDim: 128
volZSpacing: 1.28

@kylechampley
Copy link
Collaborator

That's weird. I simulated data with your geometry and I did not get an artifact. If you are able to send your data to me, I can take a look it at.

@hws203
Copy link

hws203 commented Aug 9, 2024

@lc82111 I guess that your thinking ring artifact is a real ring of your target sample. As you can see the projection image, there is a circle tap on top of cylinder. Those border area of this tap may be a ring.

@lc82111
Copy link
Author

lc82111 commented Aug 10, 2024

@kylechampley Thanks for your kind help. The raw projections can be found at this link:

https://drive.google.com/file/d/1xb7X9mQo2kcIhWD_UGi40JXHCmd-6zQq/view?usp=sharing

@lc82111
Copy link
Author

lc82111 commented Aug 10, 2024

@lc82111 I guess that your thinking ring artifact is a real ring of your target sample. As you can see the projection image, there is a circle tap on top of cylinder. Those border area of this tap may be a ring.

@hws203 Thanks for the suggestion. However, I believe the ring artifact I'm observing is not a physical characteristic of the sample. The phantom is flawless, and the ring and central dark region are absent in the full-scan data.

@kylechampley
Copy link
Collaborator

@lc82111, the data in the link has no air scan projection. Could you send this and any other calibration files you have?

@kylechampley
Copy link
Collaborator

Well, I tried to reconstruct without the air scan data and I think I got something good. It is different from what you showed, but there is still that inner circle. I do think this inner circle is real. Here is my result.

image

And you can see that the inner cylinder is indeed in the raw data, so it is indeed real. See this projection:

image

@kylechampley
Copy link
Collaborator

Here is the script I used. I turned your dicom files into tiff files before running the reconstruction.

import numpy as np
from leap_preprocessing_algorithms import *
from leapctype import tomographicModels
leapct = tomographicModels()

angles = np.linspace(-180, 180, 720).astype(np.float32)
leapct.set_conebeam(720, 704*2, 704*2, 0.615057*0.5, 0.615057*0.5, 507.583130*2.0, 576.243469*2.0, angles, 537.5, 990.804993)
leapct.set_offsetScan(True)
leapct.set_truncatedScan(True)
leapct.set_default_volume(4.0) # just show a low resolution example
leapct.print_parameters()

g = leapct.load_projections(r'D:\tomography\issue92\raw.tif')
ROI = [30, 80, 950, 1000]
makeAttenuationRadiographs(leapct, g, None, None, ROI)

vol = leapct.allocateVolume(0)
leapct.FBP(g,vol)

leapct.display(vol)

@lc82111
Copy link
Author

lc82111 commented Aug 10, 2024

AVG_504_120_30_air_cw

Above is the averaged air scan projection in png format.

I will try your script and report back. Thanks for your help.

@lc82111
Copy link
Author

lc82111 commented Aug 10, 2024

image

image

To avoid ambiguity, I have added a yellow circle to highlight the artifact. Apologies for any confusion.

@lc82111
Copy link
Author

lc82111 commented Aug 10, 2024

This ring artifact appears to be more prominent in RTK reconstruction as below:

image

Furthermore, I believe the ring artifact is likely attributed to the RTK DDF (Displaced Detector Filter) weighting, as evidenced by the abrupt changes in the weighted image:

image

Therefore, I utilized LEAP CT to validate my hypothesis.

@kylechampley
Copy link
Collaborator

I believe this is a cone-beam artifact caused by the void cylinder at the bottom of that projection. Notice that the diameter of your dark circle is about the same size as the air cylinder at the bottom of the projection.

I think you also have some beam hardening artifacts here. If you know something about the spectra, you could use LEAP's BHC algorithms to help mitigate this.

@hws203
Copy link

hws203 commented Aug 12, 2024

I think so too.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants