Over the winter break, I was reviewing Fourier Transform and its various applications. Fundamentally, Fourier Transform is about a shift in perspective - viewing the same data with a different lens makes some operations easier or faster. This piqued my curiosity about what other fundamental data transformations that I should be aware of. Like everyone in 2025, I jumped to my chatbot of choice, Gemini, to learn more.
There are a many! Gemini shared some common ones such as Principal Component Analysis (PCA), Log Transform, some I learnt briefly in school such as Laplace Transforms and Z-Transforms, but also some I know nothing about such as Radon Transform and Wavelet Transform.
This blog post is about my experience learning about Radon Transform and its application to computed tomography (CT) scan.
tl;dr I share the steps I took to learn about Radon Transform and CT scan to build a toy X ray in Threejs. I also share some reflections on process of learning and writing with LLM. You can find the app here: https://github.com/moomou/ct-voxel/
Computed tomography (CT) scan is the application of shining X-rays through an object, recording the X-rays original and final energy level, rotated over the object 180 degrees to build a 3D representation of the invisible innards. Why X-ray? In short, X-ray have higher energy and smaller wavelength than normal, visible lights that allow X-ray to penetrat materials that normal lights cannot (ref).
Imagine shining an ideal, infinitesimally, all-penetrating ray through human body at different angles, as shown in the image below.
As the ray penetrates different part of the objects (imagine X-ray hitting human skin, soft tissues, bones, etc), the intensity of the ray changes as some energy gets absorbed. If we know the exact intensity of the ray at each point of interest, we can then use that information to build a density image along the ray path. However, we can only measure ray’s source intensity (I0) and final intensity (I1) at the detector.
Beer–Lambert law is the empirical relationship that describes the attentuation of radiation through a medium that relates these 2 quantities via
where d is the distance traveled and μ is a constant describing the optical density of medium. For a human body, the density μ is not going to be a constant. Different body parts have different density. So the actual equation is more complicated and requires the use of an integral to sum up the different density at each (x, y) coordinate of the ray path
Move I0 to the left side, take the natural log of both side, and handling the negative sign, we get
If we can somehow solve for μ, we can then obtain the density at each X and Y coordinate. As it turns out, the right hand side is the Radon transform of a line (ie the ray projection) in the 2D plane. Briefly, radon tranform is an integral transform that sums up the values along a line - ie it maps a line to a point in 2D. Because the cartesian equation of a line y = mx + b cannot deal with vertical lines, Radon transform describes a line using the angle from the origin (θ) and perpendicular distance from origin to the line (ρ). As always, the actual math is more complicated so interested readers should consult their favorite Chatbot or Wikipedia for more information.
The key takeaway here is armed with ray intensities I0 and I1, applying Beer-Lambert law results in Radon transforms of the density function traced by X rays in coordinates (θ, ρ); radon transform of such lines is also know as the sinograms. This means we can apply the inverse Radon transform to solve μ.
The classic method for reversing the Radon Transform is Filtered Back Projection (FBP). It relies on the Central Slice Theorem, which proves that the 1D Fourier Transform of an X-ray projection is identical to a radial slice of the object’s 2D Fourier Transform. Because these radial slices are sampled in polar coordinates, they overlap densely at the origin (low frequencies). To correct for this $1/r$ oversampling (where r is the distance away from the origin), a Ramp Filter is applied in the frequency domain. This ensures that when the data is transformed back to the spaital domain (aka as an image grid)—a process called Back-Projection—the final reconstruction is sharp.
Note - It took me a while to grok the mapping between the 1D FFT and 2D FFT until I reviewed 2D FFT. Here is one brief lecture note I found helpful.
With the basics theory finally out of the way, let’s build a toy X ray simulator! I have always wanted to dabble with voxels and this gave me a perfect opportunity and I chose to go with Threejs. As the focus is in implementing X simulation, I vibe-coded a basic project template with Threejs and HTML canvas rendering functions.
From there, I built out a simple 3D voxel world with size 128 and represent the world with 2 arrays
const SIZE = 128;
const HALF_SIZE = SIZE / 2;
// stores object ID
const data = new Uint8Array(SIZE * SIZE * SIZE);
// stores the density of the material
const densitys = new Float32Array(SIZE * SIZE * SIZE);
I then model a simple human skull using spheres and cubes
(function setupScene() {
// 1. The "Skull" (Hollow Sphere)
// We add a large sphere, then an "air" sphere inside it
addSphere([64, 64, 64], 60, 1.0);
addSphere([64, 64, 64], 56, 0.0); // Hollow out the inside
// 2. Internal "Organs" with different densities
addSphere([45, 64, 64], 15, 0.3); // Left lobe
addSphere([83, 64, 64], 15, 0.3); // Right lobe
// 3. A "Spine" (Central Cube)
addCube([64, 64, 64], 12, 1.2);
// 4. A Small "Foreign Object" (High contrast test)
addSphere([64, 90, 40], 5, 1.5);
// 5. STREAK ARTIFACT TEST (Metal Insert)
// A very high density object. In real CT, these cause "streak artifacts."
addCube([100, 80, 64], 6, 5.0);
// 6. The "Scanning Table" (Bottom Plane)
for (let z = 0; z < SIZE; z++) {
for (let x = 0; x < SIZE; x++) {
const idx = get1DIndex(x, 5, z);
data[idx] = 1;
densitys[idx] = 0.8;
}
}
})()
Now the more interesting part. Here is the function to shoot X ray through the voxel world, one Y slice at a time through the XZ plane by summing up the density and applying the Beer-Lambert’s law.
function take2DXray(angle: number) {
// angle to radian
const rad = (Math.PI * angle) / 180;
const detector = new Float32Array(SIZE * SIZE);
for (let v = 0; v < SIZE; v++) { // v for y
for (let u = 0; u < SIZE; u++) { // u for x
let density = 0;
// note the range here is -64 to 64 because we want the X ray to rotate
// around the center of the voxel world but our coordinate system is from 0 to 128
for (let d = -64; d < 64; d++) { // d for z;
const x = Math.floor(
HALF_SIZE + (u - HALF_SIZE) * Math.cos(rad) - d * Math.sin(rad)
);
const z = Math.floor(
HALF_SIZE + (u - HALF_SIZE) * Math.sin(rad) + d * Math.cos(rad)
);
const y = v;
if (x >= 0 && x < SIZE && y >= 0 && y < SIZE && z >= 0 && z < SIZE) {
const voxel = getVoxel(x, y, z);
density += voxel.density;
}
}
const I = I0 * Math.exp(-density * 0.05); // 0.05 is the "attenuation coefficient"
detector[u + v * SIZE] = I;
}
}
return detector;
}
const xray2DByDeg: Array<Float32Array> = [];
// note we only need 180 degrees because 360 deg = 0 deg, 270 deg = 90 deg, etc.
for (let deg = 0; deg < 180; deg++) {
xray2DByDeg.push(take2DXray(deg));
}
xray2DByDeg contains the full projection data that allows us to build the 3D volume. A sinogram is a slice in this dataset at specific Y.
function getSinogramData(yValue: number) {
const sinogram = new Float32Array(180 * SIZE);
for (let deg = 0; deg < 180; deg++) {
const xray2d = xray2DByDeg[deg];
// now grab the density at Y
for (let u = 0; u < SIZE; u++) {
const observed = xray2d[u + yValue * SIZE];
const calculated = -Math.log(Math.max(observed, 0.001) / I0);
sinogram[u + deg * SIZE] = calculated;
}
}
return sinogram;
}
Given a sinogram, we can then construct the CT scan at a specific Y by doing the back projection
function backProjectWithRampFilter(y: number) {
const reconSlice = new Float32Array(SIZE * SIZE);
const sinogram = getSinogramData(y);
const filteredSinogram = new Float32Array(sinogram.length);
// apply ramp filter
// which is a V shaped filter - ie high freq gets thru
for (let deg = 0; deg < 180; deg++) {
const start = deg * SIZE;
const row = sinogram.slice(start, start + SIZE);
filteredSinogram.set(applyFFTRampFilter(row), start);
}
for (let deg = 0; deg < 180; deg++) {
const rad = (Math.PI * deg) / 180;
for (let z = 0; z < SIZE; z++) {
for (let x = 0; x < SIZE; x++) {
const xC = x - HALF_SIZE;
const zC = z - HALF_SIZE;
// inverse radon transform
const u = Math.floor(
HALF_SIZE + xC * Math.cos(rad) + zC * Math.sin(rad)
);
if (u >= 0 && u < SIZE) {
reconSlice[x + z * SIZE] += filteredSinogram[u + deg * SIZE];
}
}
}
}
return reconSlice;
}
This backprojection only builds the XZ slice (viewing the skull from top down).
To build the YZ slice, we can build the full 3D volume and then extract the YZ plane.
const full3DRecon = new Float32Array(SIZE * SIZE * SIZE);
for (let y = 0; y < SIZE; y++) {
const xzSlice = backProjectWithRampFilter(y);
for (let i = 0; i < SIZE * SIZE; i++) {
full3DRecon[y * (SIZE * SIZE) + i] = xzSlice[i];
}
}
function yzFrom3DRecon(xValue: number) {
const yzSlice = new Float32Array(SIZE * SIZE);
for (let z = 0; z < SIZE; z++) {
for (let y = 0; y < SIZE; y++) {
const volIdx = xValue + y * SIZE * SIZE + z * SIZE;
yzSlice[y + z * SIZE] = full3DRecon[volIdx];
}
}
return yzSlice;
}
Finally, putting everything together, here is the toy CT scan app.
You can find the full code here https://github.com/moomou/ct-voxel/
I am amazed at what I was able to build over a weekend with an interactive teacher like Gemini, having no prior knowledge about CT scan and Threejs. Almost felt like that I had superpower! However, such superpower is not without its downsides. During the course of my learning, I constantly asked myself what’s the point of building the Threejs app? Was I actually learning? Or was I simply reaping dopamine rewards without increasing my skills? Gemini generated a lot of code and definitely could have built the whole thing in one minute.
I was only able to alleviate some of that anxiety by deliberately avoiding some of Gemini’s eager help. Forcing myself to write the 2D X ray projection function allowed me to build a mental model for sinograms and ray projections. Only then was I able to onvince myself that I was learning. I am not envious of today’s students - it’s all too easy to cheat and avoid the hardwork.