Why Fourier Transform can be used in Tomography Reconstruction?

The Fourier transform is founded on Euler's formula:
$$e^{i\theta} = \cos(\theta) + i\sin(\theta) \tag{1}$$
This equation is, without a doubt, one of my favorites. It elegantly illustrates how the world, while often perceived as abstract or intangible, operates according to well-defined principles. Euler's formula shows that the real component can be extracted from the cosine function, representing the "true" part, while the sine function provides the imaginary component, which I interpret as noise. What fascinates me the most is that 'noise' shares the same fundamental nature as the true part, reflecting the inherent harmony between reality and its complex representations. When considering the relationship to angles, the transformation exhibits a form of periodic behavior. This periodicity indicates a repeating signal, which implies the presence of an underlying rule or pattern that can be leveraged repeatedly. By understanding and utilizing this repetition, we can potentially uncover consistent properties within the data that are valuable for analysis and signal processing.

In tomography reconstruction, the goal is to obtain cross-sectional images from projection data. Physically, this process is governed by the Beer-Lambert Law, which describes the relationship between cross-sections and projections:

$$ I = I_0 \cdot e^{- \mu x} \tag{2} $$

This equation forms the foundational principle of the experimental design. The system can be conceptualized as a sample processing framework, where the system itself can be represented as the material's absorption characteristics. In this case,The input is the X-ray signal, while the output is the projection image.Alternatively, the system can be viewed as the X-ray itself, with the 3D sample as the input and the resulting 2D projection as the output. A consideration worth exploring is whether this framework is time-invariant, as the X-ray intensity or absorption properties may vary over space and time.

In parallel beam tomography, each X-ray beam passes through the object and maps to a specific point on the projection image. The absorption along the beam’s path contributes to the overall reduction in X-ray intensity, influenced not only by the material’s properties but also by the total path length. This establishes an intrinsic connection between the cross-sectional image and the projection data, which can be uncovered using the Fourier transform. The Fourier transform, leveraging the exponential form expressed by Euler’s formula, multiplies the signal by a complex exponential function. When the cross-sectional signal is represented as 𝑔(𝑥, 𝑦), the continuous form of the 2D Fourier transform is provided below.
$$
G(f_X, f_Y) = \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} g(x, y) e^{-i 2 \pi (f_X x + f_Y y)} ,dx ,dy \tag{3}
$$
This function converts a 2D spatial signal into a 2D frequency signal, revealing inherent patterns and structures within the frequency domain. An intriguing question arises: Can a 2D spatial signal be directly transformed into a 1D frequency signal to enable potential compression or improve image quality? I believe this is achievable, supported by the projection-slice theorem, which offers a compelling reason for employing the Fourier transform in tomography.

Each point in the frequency domain is influenced by all corresponding points in the spatial domain, but the reverse is not always true. This raises an important question: Can an inverse transformation from the frequency domain fully reconstruct the original spatial signal, or does it mainly capture only the essential or high-quality components? By focusing on these significant elements, it may be possible to enhance the image or isolate specific features of interest.

• Projection-Slice theorem(from projected imaging to cross-section)

By rotating the setup, we change the direction of the X-ray entering the sample. Now, assuming that the X-ray enters along the y'-axis, we can extend the Beer-Lambert Law from Equation 2 to a more specific case. This case considers a continuous line from the cross-section 𝑔(𝑥,𝑦), which represents the desired outcome:
$$
I_{\theta}(x') =I_{0} e^{- \int_{-\infty}^{\infty} g(x' \cos \theta - y' \sin \theta, x' \sin \theta + y' \cos \theta)},dy' \tag{4}
$$
From the detection process, we obtain the intensity of the X-ray after passing through the sample. However, our primary focus is on the sample's absorption. To address this, we solve the function as follows:
$$
\ln \left( \frac{I_{\theta}(x')}{I_{0}} \right) \tag{5}
$$
To determine the sample's absorption along this line, we define it as the projection of the X-ray attenuation coefficient:
$$
p_{\theta}(x') = \int_{-\infty}^{\infty} g(x' \cos \theta - y' \sin \theta, x' \sin \theta + y' \cos \theta) , dy' \tag{6}
$$
This line represents an integral along the y'-axis from the 2D cross-section 𝑔(𝑥,𝑦),effectively compressing the information of a line into a single point. For a fixed angle, we can imagine slicing multiple parallel lines along the y'-axis from the sample (similar to cutting a pizza), with each line producing one data point on our projection detector, represented by 𝑥'. As we rotate, we investigate the continuity of the sample. At a fixed angle, the signal reduces to a 1D signal, which we then use fourier transform to the 1D siginal:
$$
P_{\theta}(f) = \int_{-\infty}^{\infty} p_{\theta}(x') e^{-i 2 \pi x' f} , dx' \tag{7}
$$
By introducing Equation 6 into Equation 7, we convert the 1D information back into 2D information, effectively uncompressing the data. This leads to the following equation:
$$
= \int_{-\infty}^{\infty} \left( \int_{-\infty}^{\infty} g(x' \cos \theta - y' \sin \theta, x' \sin \theta + y' \cos \theta) , dy' \right) e^{-i 2 \pi f x'} , dx'\tag{8}
$$
By transforming the plane from (x',y') to (x,y), we obtain the following:
$$
P_{\theta}(f) = \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} g(x, y) e^{-i 2 \pi (x \cos \theta + y \sin \theta) f} , dx , dy\tag{9}
$$
Equation 9 represents the Fourier transform of the 2D cross-section signal, but with the frequency expressed in polar coordinates,which is shown in equation 10. At a fixed angle, this frequency reduces to a 1D signal. This demonstrates that a 2D spatial signal can be converted into a 1D frequency signal, supporting the idea that we can transform information intensity across different domains by changing dimensions.
$$
P_{\theta}(f) = G(s \cos \theta, s \sin \theta)\tag{10}
$$

• Filtered backprojection(from cross-section to projected imaging)

For an unknown 2D cross-section g(x,y), we can apply the inverse 2D Fourier transform. This involves transforming the 2D cross-sectional signal to mathematically define it as follows:
$$
g(x, y) = \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} G(f_X, f_Y) e^{i 2 \pi (f_X x + f_Y y)} , df_X , df_Y \tag{11}
$$
Convert the frequency domain from rectangular coordinates to a polar coordinate system using the following method:
$$
f_X = s \cos \theta \
f_Y = s \sin \theta \
df_X , df_Y = |s| , ds , d\theta\tag{12}
$$
By substituting Equation 12 into Equation 11, we obtain the following result:
$$
g(x, y) = \int_0^{2\pi} d\theta \int_{-\infty}^{\infty} ds , |s| , G(s \cos \theta, s \sin \theta) e^{i 2 \pi s (x \cos \theta + y \sin \theta)}\tag{13}
$$
From the Projection-Slice Theorem，Equation 10, we derive the following result:
$$
g(x, y) = \int_0^{2\pi} d\theta \int_{-\infty}^{\infty} ds , |s| , P_{\theta}(s) e^{i 2 \pi s (x \cos \theta + y \sin \theta)}\tag{14}
$$
Thus, we obtain the corresponding relationship between the 2D cross-sectional signal and the 1D frequency signal, as shown in Equation 15. For each point on the 2D cross-section, the entire frequency domain must be considered in a polar coordinate system, with an appropriate scale factor. This scaling factor is essential when transforming from the frequency domain to the time or spatial domain, as it ensures consistency during the generalization process—particularly when converting from a polar frequency domain to a rectangular spatial domain.
$$
g(x, y) = \int_0^{2\pi} d\theta , f_{\theta}(x \cos \theta + y \sin \theta)\tag{15}
$$

The key question is: how many projections are needed? According to the Projection-Slice Theorem, if the input signal is circularly symmetric, only one projection is necessary. For a separable function, two projections suffice. However, for other types of input signals, a set of projections at various angles is required.

In the Projection-Slice approach, obtaining the 2D cross-sectional Fourier transform in polar coordinates is straightforward by applying the Fourier transform to the projections. However, converting from polar to rectangular coordinates can be computationally intensive due to the inversion performed on a rectangular grid. Additionally, this method requires waiting until all data is collected and transformed before performing the final inverse FFT.

On the other hand, the filtered back-projection method allows for real-time processing. For each angle, the back-projection can be performed, accumulated, and then followed by the next projection, all while convolving with the transfer function.

Please  forgive my ignorance and feel free to give suggenstions or any other ideas. Thank you.