Initial commit

This commit is contained in:
Peter Stockings
2025-05-17 09:21:16 +10:00
commit 86f7a69b37
7 changed files with 765 additions and 0 deletions

0
.static Normal file
View File

108
CTAccordion .js Normal file
View File

@@ -0,0 +1,108 @@
export const steps = [
{
title: "1. Upload / Input Image",
content: `You upload or use a default grayscale image. This is treated like a 2D slice of a physical object.
The input image is a 2D function:
$ f(x, y) $
This represents how much X-rays are absorbed at each point.`,
},
{
title: "2. Radon Transform (Generating the Sinogram)",
content: `We rotate a virtual X-ray beam around the image and compute line integrals at each angle — simulating how X-rays pass through.
\\[
R(\\theta, s) = \\int_{-\\infty}^{\\infty} f(s \\cos\\theta - t \\sin\\theta,\\ s \\sin\\theta + t \\cos\\theta)\\ dt
\\]
Where: <br>
- $\\theta$ = projection angle <br>
- $s$ = offset from center
![Radon transform](https://upload.wikimedia.org/wikipedia/commons/9/93/Radon_transform_sinogram.gif)`,
},
{
title: "3. Sinogram",
content: `You now have a 2D image where:
- X-axis = detector position
- Y-axis = angle
Each row is a projection.
Math:
$\text{Sinogram} = \{ R(\theta_1, s), R(\theta_2, s), \dots, R(\theta_n, s) \}$
![Sinogram](https://www.researchgate.net/profile/Samuel-Asante-2/publication/299856137/figure/fig3/AS:348226420002817@1460035056245/A-Shepp-Logan-Phantom-and-reconstructed-Image-Sinogram-a-Original-image-b-radon.png)`,
},
{
title: "4. Optional: Apply Ramp Filter (FBP)",
content: `If enabled, we sharpen each projection before back-projecting by amplifying high-frequency content.
Math:
\\[
P(\\omega) = \\mathcal{F}[R(\\theta, s)] \\\\
P'(\\omega) = |\\omega| \\cdot P(\\omega) \\\\
\\text{Filtered Projection} = \\mathcal{F}^{-1}[P'(\\omega)]
\\]
![Ramp filter graph](https://www.researchgate.net/publication/346858231/figure/fig4/AS:967035467091970@1607570630558/Applying-Ramp-filter-to-a-sinogram-preserves-high-frequency-features-The-filter-is.png)`,
},
{
title: "5. Back Projection",
content: `Each (filtered) projection is \"smeared\" back into the image space along its angle.
\\[
f'(x, y) = \\int_0^{\\pi} R'(\\theta,\\ x \\cos\\theta + y \\sin\\theta)\\ d\\theta
\\]`,
},
{
title: "6. Final Reconstruction",
content: `After all angles are added up, you get a reconstructed image resembling the original.
\[
f'(x, y) \approx f(x, y)
\]`,
},
];
export const StepAccordion = {
view(vnode) {
const index = vnode.attrs.index;
const step = steps[index];
const expanded = vnode.state.expanded ?? false;
return m("div", { class: "border-b border-gray-300 py-1 my-2" }, [
m(
"button",
{
class:
"w-full text-left font-bold text-lg text-gray-800 focus:outline-none",
onclick: () => {
vnode.state.expanded = !vnode.state.expanded;
},
},
step.title
),
vnode.state.expanded &&
m(
"div",
{
class: "mt-2 text-gray-700 whitespace-pre-wrap text-sm",
onupdate: () => {
if (window.MathJax) window.MathJax.typesetPromise();
},
},
m.trust(markdownToHTML(step.content))
),
]);
},
};
// You must provide a markdownToHTML() function or use a library like marked.js
function markdownToHTML(text) {
return text
.replace(/\n/g, "<br>")
.replace(
/!\[(.*?)\]\((.*?)\)/g,
'<img alt="$1" src="$2" class="my-2 rounded shadow max-w-full">'
);
}

326
UploadImageComponent.js Normal file
View File

@@ -0,0 +1,326 @@
import {
generateSinogram,
reconstructImageFromSinogram,
convertToGrayscale,
} from "./sinogram.js";
import { StepAccordion } from "./CTAccordion .js";
export const UploadImageComponent = {
hasLoadedInitialImage: false,
angleCount: 180,
imageUrl:
"https://upload.wikimedia.org/wikipedia/commons/e/e5/Shepp_logan.png",
sinogramUrl: null,
reconstructedUrl: null,
defaultImageUrl:
"https://upload.wikimedia.org/wikipedia/commons/e/e5/Shepp_logan.png",
reconstructionFrames: [],
currentFrameIndex: 0,
renderMode: "grayscale", // or "heatmap"
useFBP: true,
drawAngleOverlay(theta) {
const canvas = this.overlayCanvas;
if (!canvas || !this.imageElement) return;
const ctx = canvas.getContext("2d");
const w = canvas.width;
const h = canvas.height;
const cx = w / 2;
const cy = h / 2;
const len = Math.max(w, h);
ctx.clearRect(0, 0, w, h);
const dx = len * Math.cos(theta);
const dy = len * Math.sin(theta);
ctx.strokeStyle = "rgba(255,0,0,0.8)";
ctx.lineWidth = 2;
ctx.beginPath();
ctx.moveTo(cx - dx, cy - dy);
ctx.lineTo(cx + dx, cy + dy);
ctx.stroke();
},
isOverlayReady() {
return (
this.overlayCanvas &&
this.imageElement &&
this.imageElement.complete &&
this.imageElement.naturalWidth > 0
);
},
async loadAndProcess(url, isUploaded = false) {
this.imageUrl = url;
this.sinogramUrl = "loading";
this.reconstructedUrl = null;
m.redraw();
this.imageUrl = url;
let finalUrl = url;
if (isUploaded) {
finalUrl = await convertToGrayscale(url);
this.imageUrl = finalUrl;
}
this.sinogramUrl = await generateSinogram(
finalUrl,
this.angleCount,
this.drawAngleOverlay.bind(this)
);
m.redraw();
this.reconstructionFrames = [];
this.currentFrameIndex = 0;
this.reconstructedUrl = await reconstructImageFromSinogram(
this.sinogramUrl,
undefined,
(angle, frameUrl) => {
this.reconstructionFrames.push(frameUrl);
this.currentFrameIndex = this.reconstructionFrames.length - 1;
this.reconstructedUrl = frameUrl;
m.redraw();
},
this.renderMode,
this.useFBP
);
},
oninit() {
this.loadAndProcessDebounced = debounce((url) => {
this.loadAndProcess(url);
}, 300);
},
view() {
return m(
"div",
{
class: "flex flex-col items-center min-h-screen bg-gray-100 py-10 px-4",
},
[
// Header
m("header", { class: "mb-10 text-center" }, [
m(
"h1",
{ class: "text-4xl font-bold text-gray-800 mb-2" },
"Sinogram Generator"
),
m(
"p",
{ class: "text-gray-600 text-lg" },
"Upload a grayscale image to simulate CT scan projections"
),
]),
m(StepAccordion, { index: 0 }),
// Upload Box
m(
"div",
{
class:
"w-full max-w-lg border-4 border-dashed border-gray-400 bg-white rounded-xl p-6 text-center hover:bg-gray-50 cursor-pointer transition",
ondragover: (e) => e.preventDefault(),
ondrop: (e) => {
e.preventDefault();
const file = e.dataTransfer.files[0];
if (file && file.type.startsWith("image/")) {
const url = URL.createObjectURL(file);
this.loadAndProcess(url, true);
}
},
onclick: () => document.getElementById("fileInput").click(),
},
[
m(
"p",
{ class: "text-gray-500" },
"Click or drag a grayscale image here"
),
m("input", {
id: "fileInput",
type: "file",
class: "hidden",
accept: "image/*",
onchange: (e) => {
const file = e.target.files[0];
if (file && file.type.startsWith("image/")) {
const url = URL.createObjectURL(file);
this.loadAndProcess(url);
}
},
}),
]
),
m(StepAccordion, { index: 1 }),
// Image Preview
m("div", { class: "relative mt-6 w-full max-w-md" }, [
m("img", {
src: this.imageUrl,
class: "rounded shadow max-w-full h-auto mx-auto",
onload: (e) => {
this.imageElement = e.target;
// Only start once both image and canvas are ready
if (this.isOverlayReady() && !this.hasLoadedInitialImage) {
this.hasLoadedInitialImage = true;
this.loadAndProcess(this.imageUrl);
}
},
}),
m("canvas", {
width: this.imageElement?.width || 0,
height: this.imageElement?.height || 0,
style: "position:absolute; top:0; left:0; pointer-events:none;",
oncreate: ({ dom }) => {
this.overlayCanvas = dom;
// Trigger load if image was already ready
if (this.isOverlayReady() && !this.hasLoadedInitialImage) {
this.hasLoadedInitialImage = true;
this.loadAndProcess(this.imageUrl);
}
},
}),
]),
// Angle Slider
m("div", { class: "mt-6 w-full max-w-md" }, [
m(
"label",
{ class: "block text-sm font-medium text-gray-700 mb-1" },
`Number of angles: ${this.angleCount}`
),
m("input", {
type: "range",
min: 5,
max: 360,
value: this.angleCount,
step: 1,
class: "w-full",
oninput: (e) => {
this.angleCount = parseInt(e.target.value, 10);
this.loadAndProcessDebounced(this.imageUrl); // reprocess with new angle count
},
}),
]),
m(StepAccordion, { index: 2 }),
// Sinogram
this.sinogramUrl &&
m("div", { class: "mt-10 w-full max-w-md text-center" }, [
m(
"h2",
{ class: "text-xl font-semibold text-gray-700 mb-4" },
"Generated Sinogram"
),
this.sinogramUrl === "loading"
? m("p", "Processing...")
: m("img", {
src: this.sinogramUrl,
alt: "Sinogram",
class: "rounded shadow max-w-full h-auto mx-auto",
}),
]),
m(StepAccordion, { index: 3 }),
// Reconstructed
this.reconstructedUrl &&
m("div", { class: "mt-10 w-full max-w-md text-center" }, [
m(
"h2",
{ class: "text-xl font-semibold text-gray-700 mb-4" },
"Reconstructed Image (Back Projection)"
),
m("img", {
src: this.reconstructionFrames[this.currentFrameIndex],
alt: "Reconstructed",
class: "rounded shadow max-w-full h-auto mx-auto",
}),
m("div", { class: "mt-6 w-full max-w-md text-center" }, [
m(
"label",
{ class: "text-sm text-gray-600 mr-2" },
"Render style:"
),
m(
"select",
{
value: this.renderMode,
onchange: (e) => {
this.renderMode = e.target.value;
this.loadAndProcess(this.imageUrl); // re-render using selected mode
},
},
[
m("option", { value: "heatmap" }, "Heatmap"),
m("option", { value: "grayscale" }, "Grayscale"),
]
),
]),
m("div", { class: "mt-4 w-full max-w-md text-left" }, [
m("label", [
m("input", {
type: "checkbox",
checked: this.useFBP,
onchange: (e) => {
this.useFBP = e.target.checked;
this.loadAndProcess(this.imageUrl); // regenerate with or without FBP
},
}),
m(
"span",
{ class: "ml-2 text-gray-700" },
"Use Filtered Back Projection (Ramp)"
),
]),
]),
this.reconstructionFrames.length > 1 &&
m("div", { class: "mt-4" }, [
m("input", {
type: "range",
min: 0,
max: this.reconstructionFrames.length - 1,
value: this.currentFrameIndex,
step: 1,
oninput: (e) => {
this.currentFrameIndex = +e.target.value;
this.reconstructedUrl =
this.reconstructionFrames[this.currentFrameIndex];
},
}),
m(
"p",
{ class: "text-sm text-gray-500 mt-1" },
`Angle ${this.currentFrameIndex + 1} / ${
this.reconstructionFrames.length
}`
),
]),
m(StepAccordion, { index: 4 }),
]),
]
);
},
};
function debounce(fn, delay) {
let timeout;
return (...args) => {
clearTimeout(timeout);
timeout = setTimeout(() => fn(...args), delay);
};
}

3
app.js Normal file
View File

@@ -0,0 +1,3 @@
import { UploadImageComponent } from "./UploadImageComponent.js";
m.mount(document.body, UploadImageComponent);

92
fbp.js Normal file
View File

@@ -0,0 +1,92 @@
// Helper: Next power of 2
function nextPow2(n) {
return 1 << (32 - Math.clz32(n - 1));
}
// Basic Cooley-Tukey FFT (real/imag arrays)
function fft1D(re, im) {
const N = re.length;
if (N <= 1) return;
// Bit-reversal permutation
const rev = new Uint32Array(N);
let logN = Math.log2(N);
for (let i = 0; i < N; i++) {
let x = i,
y = 0;
for (let j = 0; j < logN; j++) {
y <<= 1;
y |= x & 1;
x >>= 1;
}
rev[i] = y;
}
for (let i = 0; i < N; i++) {
if (i < rev[i]) {
[re[i], re[rev[i]]] = [re[rev[i]], re[i]];
[im[i], im[rev[i]]] = [im[rev[i]], im[i]];
}
}
// FFT
for (let s = 1; s <= logN; s++) {
const m = 1 << s;
const m2 = m >> 1;
const wAngle = (-2 * Math.PI) / m;
const cosW = Math.cos(wAngle);
const sinW = Math.sin(wAngle);
for (let k = 0; k < N; k += m) {
let wr = 1,
wi = 0;
for (let j = 0; j < m2; j++) {
const tRe = wr * re[k + j + m2] - wi * im[k + j + m2];
const tIm = wr * im[k + j + m2] + wi * re[k + j + m2];
const uRe = re[k + j];
const uIm = im[k + j];
re[k + j] = uRe + tRe;
im[k + j] = uIm + tIm;
re[k + j + m2] = uRe - tRe;
im[k + j + m2] = uIm - tIm;
const tempWr = wr;
wr = wr * cosW - wi * sinW;
wi = tempWr * sinW + wi * cosW;
}
}
}
}
// Inverse FFT
function ifft1D(re, im) {
// Conjugate
for (let i = 0; i < re.length; i++) im[i] = -im[i];
fft1D(re, im);
// Normalize and re-conjugate
const N = re.length;
for (let i = 0; i < N; i++) {
re[i] /= N;
im[i] = -im[i] / N;
}
}
// Apply ramp filter in frequency domain
export function applyRampFilter(projection) {
const N = nextPow2(projection.length);
const re = new Float32Array(N);
const im = new Float32Array(N);
for (let i = 0; i < projection.length; i++) {
re[i] = projection[i];
}
fft1D(re, im);
for (let i = 0; i < N / 2; i++) {
const freq = i / N;
re[i] *= freq;
im[i] *= freq;
re[N - i - 1] *= freq;
im[N - i - 1] *= freq;
}
ifft1D(re, im);
return Array.from(re.slice(0, projection.length));
}

27
index.html Normal file
View File

@@ -0,0 +1,27 @@
<!DOCTYPE html>
<html lang="en">
<head>
<meta charset="UTF-8" />
<meta name="viewport" content="width=device-width, initial-scale=1.0" />
<title>Sinogram App</title>
<script src="https://cdn.tailwindcss.com"></script>
<script src="https://unpkg.com/mithril/mithril.js"></script>
<script>
window.MathJax = {
tex: {
inlineMath: [['$', '$'], ['\\(', '\\)']],
displayMath: [['\\[', '\\]']],
},
svg: { fontCache: 'global' }
};
</script>
<script src="https://cdn.jsdelivr.net/npm/mathjax@3/es5/tex-svg.js"></script>
</head>
<body class="bg-gray-100">
<script type="module" src="./app.js"></script>
</body>
</html>

209
sinogram.js Normal file
View File

@@ -0,0 +1,209 @@
import { applyRampFilter } from "./fbp.js";
export async function generateSinogram(
imageUrl,
angles = 180,
drawAngleCallback = null
) {
const image = await loadImage(imageUrl);
const size = Math.max(image.width, image.height);
const projections = [];
const canvas = Object.assign(document.createElement("canvas"), {
width: size,
height: size,
});
const ctx = canvas.getContext("2d");
for (let angle = 0; angle < angles; angle++) {
const theta = (angle * Math.PI) / angles;
// 🔁 Call visual overlay for this angle
if (drawAngleCallback) drawAngleCallback(theta);
// (Optional: add delay for animation)
await new Promise((r) => setTimeout(r, 0.01));
// Clear canvas
ctx.clearRect(0, 0, size, size);
// Transform and draw rotated image
ctx.save();
ctx.translate(size / 2, size / 2);
ctx.rotate(theta);
ctx.drawImage(image, -image.width / 2, -image.height / 2);
ctx.restore();
// Read pixel data
const { data } = ctx.getImageData(0, 0, size, size);
// Sum brightness vertically (simulate X-ray projection)
const projection = [];
for (let x = 0; x < size; x++) {
let sum = 0;
for (let y = 0; y < size; y++) {
const i = (y * size + x) * 4;
const gray = data[i]; // red channel (since grayscale)
sum += gray;
}
projection.push(sum / size); // normalize
}
projections.push(projection);
}
// Create sinogram canvas
const sinogramCanvas = Object.assign(document.createElement("canvas"), {
width: size,
height: angles,
});
const sinCtx = sinogramCanvas.getContext("2d");
const imgData = sinCtx.createImageData(size, angles);
for (let y = 0; y < angles; y++) {
for (let x = 0; x < size; x++) {
const val = projections[y][x];
const i = (y * size + x) * 4;
imgData.data[i + 0] = val;
imgData.data[i + 1] = val;
imgData.data[i + 2] = val;
imgData.data[i + 3] = 255;
}
}
sinCtx.putImageData(imgData, 0, 0);
return sinogramCanvas.toDataURL();
}
export async function reconstructImageFromSinogram(
sinogramUrl,
size = 256,
onFrame = null,
renderMode = "heatmap",
useFBP = true
) {
const sinogramImage = await loadImage(sinogramUrl);
const canvas = Object.assign(document.createElement("canvas"), {
width: sinogramImage.width,
height: sinogramImage.height,
});
const ctx = canvas.getContext("2d");
ctx.drawImage(sinogramImage, 0, 0);
const sinogramData = ctx.getImageData(
0,
0,
sinogramImage.width,
sinogramImage.height
).data;
size = sinogramImage.width; // match size to sinogram resolution
const outputCanvas = Object.assign(document.createElement("canvas"), {
width: size,
height: size,
});
const outputCtx = outputCanvas.getContext("2d");
const accum = new Float32Array(size * size);
const center = size / 2;
const angles = sinogramImage.height;
const width = sinogramImage.width;
for (let angle = 0; angle < angles; angle++) {
const theta = (angle * Math.PI) / angles;
let projection = [];
for (let x = 0; x < width; x++) {
const i = (angle * width + x) * 4;
projection.push(sinogramData[i]);
}
if (useFBP) {
projection = applyRampFilter(projection);
}
for (let y = 0; y < size; y++) {
for (let x = 0; x < size; x++) {
const x0 = x - center;
const y0 = center - y; // flip y
const s = Math.round(
x0 * Math.cos(theta) + y0 * Math.sin(theta) + width / 2
);
if (s >= 0 && s < width) {
accum[y * size + x] += projection[s];
}
}
}
if (onFrame) {
// normalize and draw current frame
let maxVal = 0;
for (let i = 0; i < accum.length; i++) {
if (accum[i] > maxVal) maxVal = accum[i];
}
const imageData = outputCtx.createImageData(size, size);
for (let i = 0; i < accum.length; i++) {
let val = accum[i] / maxVal;
val = Math.min(1, Math.max(0, val));
let r, g, b;
if (renderMode === "grayscale") {
const gray = Math.round(val * 255);
r = g = b = gray;
} else {
[r, g, b] = getHeatmapColor(val);
}
imageData.data[i * 4 + 0] = r;
imageData.data[i * 4 + 1] = g;
imageData.data[i * 4 + 2] = b;
imageData.data[i * 4 + 3] = 255;
}
outputCtx.putImageData(imageData, 0, 0);
await new Promise((r) => setTimeout(r, 1));
onFrame(angle, outputCanvas.toDataURL());
}
}
return outputCanvas.toDataURL();
}
// Heatmap mapping: blue → green → yellow → red
function getHeatmapColor(value) {
const r = Math.min(255, Math.max(0, 255 * Math.min(1, 4 * (value - 0.75))));
const g = Math.min(255, Math.max(0, 255 * (4 * Math.abs(value - 0.5) - 1)));
const b = Math.min(255, Math.max(0, 255 * (1 - 4 * value)));
return [r, g, b];
}
function loadImage(src) {
return new Promise((resolve) => {
const img = new Image();
img.crossOrigin = "anonymous";
img.onload = () => resolve(img);
img.src = src;
});
}
export async function convertToGrayscale(imageUrl) {
const image = await loadImage(imageUrl);
const canvas = Object.assign(document.createElement("canvas"), {
width: image.width,
height: image.height,
});
const ctx = canvas.getContext("2d");
// Draw original image
ctx.drawImage(image, 0, 0);
// Get pixel data
const imageData = ctx.getImageData(0, 0, canvas.width, canvas.height);
const data = imageData.data;
// Convert to grayscale: set R, G, B to luminance
for (let i = 0; i < data.length; i += 4) {
const r = data[i];
const g = data[i + 1];
const b = data[i + 2];
const luminance = 0.299 * r + 0.587 * g + 0.114 * b;
data[i] = data[i + 1] = data[i + 2] = luminance;
}
ctx.putImageData(imageData, 0, 0);
return canvas.toDataURL();
}