CL reconstruction artifact

[zeminLi66]

testExperiment.2025-12-12.11-50-40.mp4

experimentResult.2025-12-12.11-51-09.mp4

I am currently conducting a cone-beam CT (CBCT) 3D reconstruction experiment. I generated projection data from a simulated 3D model and then used an iterative algorithm for the reconstruction.
Methodology and Setup:
I performed the reconstruction using the SART (Simultaneous Algebraic Reconstruction Technique) algorithm. Crucially, the forward and backward projection methods I employed are sourced from the TIGRE toolbox:
Forward Projection: Siddon algorithm.
Backward Projection: Backprojection 2 algorithm.
Observed Problem:
The resulting reconstructed volume, visible in the output video, contains significant image artifacts, specifically numerous distinct circular artifacts localized in the lower region of the reconstructed volume.
Key Observations (Artifact Mitigation):
I discovered that these circular artifacts are reliably eliminated under two independent conditions:
Simultaneously reducing both the voxel size and the detector pixel size to half of their original dimensions.
Simultaneously doubling the distance from the radiation source to the object center (DSO) and the distance from the object center to the detector (DOD). (This action significantly reduces the effective cone angle.)
Hypothesis and Inquiry:
Based on the fact that increasing the geometric distances (which reduces the cone angle) or increasing the sampling rate (halving the pixel/voxel size) solves the issue, I suspect the artifacts are primarily caused by an excessive cone angle, leading to insufficient sampling or potential truncation/interpolation errors near the edges.
I would be grateful for your expert insight on the following points:
What is the precise fundamental reason for the appearance of these circular artifacts in this specific cone-beam geometry (especially given the use of TIGRE's Siddon/Backprojection 2)? Is my suspicion about the excessive cone angle or sampling density correct?
Aside from physically changing the geometry or substantially increasing the resolution, are there algorithmic solutions, projection modifications, or correction methods that can be implemented to solve or mitigate these artifacts within the SART framework?
Thank you very much for your time and assistance.

[zeminLi66]

https://github.com/user-attachments/assets/b70b2e7a-f53a-411f-bfb7-22e81f231a70
I conducted an additional experiment, as demonstrated in the video: I extracted a 2D image slice from the 3D model after each iteration. I used a total of 127 projection images and performed 127 iterations. However, these artifacts only began to appear after the 75th iteration.

[zeminLi66]

The projection method is shown in the figure

[AnderBiguri]

Interesting.

Can you share the geometry too? The standard pixel size to avoid aliasing artifacts should be detector*DSO/DSD. Is this the case?

[zeminLi66]

I used the same data set for both forward projection and reconstruction. The reconstruction result is perfect if DSO and DOD are not multiplied by 0.33.
According to the formula, I set the voxel size to 0.0167, but there are still artifacts，It shouldn't be that problem

[AnderBiguri]

Does this only happen with SART?
Are you using your own C++ implementation of SART? In any case, try replacing the matrices V and W by a scalar with their mean.
There has been edge cases where these matrices cause strange numerical errors. I wonder if that is what we are seeing. Mostly because your strange artifacts only appear in very specific areas, I wonder if numerical precission error is the cause.

[zeminLi66]

I only implemented SART, I haven't implemented anything else yet, so I'm not sure. I manually implemented TIGRE's SART algorithm in C++, which allowed me to achieve acceleration and quickly obtain the reconstruction results. You're right. If I change the size or reconstruct another object, as long as the object or detector has a larger area or is closer to the source, these artifacts will appear, and the position and shape of each occurrence are basically the same， The effect is essentially the same as the Python DEMO. However, when updating the 3D model using backprojection in the final iteration step, even without dividing the backprojected value by V(the voxel weight), artifacts still remain. The result is shown in the figure.

[AnderBiguri]

@zeminLi66 very interesting. Try to completely remove that V in all iterations. Its just a minor weighting that is slightly different per projection, but has very limited influence in the reconstruction. Its average value should be geo.dVoxel.

In essence, it has the total length of the ray over the image. It has happened in the past that when this intersection is minuscule for some pixels, then inv_V is very very large, and it causes numerical errors. I wonder if that is what you are seing

[zeminLi66]

I have been studying the TIGRE code. If I recall correctly, V seems to only take effect at the end when updating the 3D model—the back projection is multiplied by the inverse of V and then added to the result. I removed this V in the final function, and it had a noticeable impact on the imaging: the image became quite blurry, and that small ring of artifacts remained but also turned blury. I have also tried to save each forward projection during the iteration process, and I found that there is no problem with the forward projection. It is because the 3D model has artifacts first, and then the forward projection appears. So I suspect that it may be a problem with the back projection? Of course, I cannot completely rule out whether it is a problem with the V parameter.Unfortunately, I am at home now, but I can show you the results without V this and forward projection result Monday when I return to work. In fact, I have always been curious about the mathematical principle behind this V. I would really appreciate it if you could explain it to me.

[AnderBiguri]

Yes, just to be clear, V has to be replaced by the average of V, not just removed, to have a small effect. The same with W. They should have a minimal effect if replaced by their average value. Could be either V or W, but its worth removing them to see if its happening.

Mathematically V is just a weight to make all updates of each pixel have the same importance. Its a "update importance" normalization factor + general scaling factor (the average).

[zeminLi66]

I understand, I will try it out later. If it works, I will get back to you

[zeminLi66]

experimentResult.2025-12-15.11-11-56.mp4

I tried to force the values of V to 0.015 (geo. dVoxel), and the results are shown in the video. I also tried to change the value of W, but the result was still not ideal.

[zeminLi66]

https://github.com/user-attachments/assets/dc67cab4-3e8e-4bbf-9b7f-6be8cb6aa7d9
If I use the default backprojection (not backprojection2，W and V have not been changed to mean values), the result will be as shown in the video above.

[zeminLi66]

The issue has been resolved. The root cause was the back-projection. I applied some constraints to the weights, which provided a perfect fix. Haha, thank you for your previous response.

[zeminLi66]

I have identified the root cause of the issue: during the back-projection process in TIGRE, the calculation of lsq (prior to the weight calculation) does not account for Euler ZYZ rotations. Consequently, if a user specifies alpha, theta, and psi angles, the lsq value becomes inaccurate, leading to incorrect weight calculations.

[AnderBiguri]

@zeminLi66 ah I see. There is no clear solution, as the weigths are FDK weigths, which should only work for circular rotation, mathematically speaking. But this should be counterbalanced by the V matrix using the same weight, right?