| SOCR ≫ | BPAD1 Website ≫ | BPAD GitHub ≫ |
Biomedical Physics with Applications to Disease (BPAD)
Every modality so far has been non-ionizing. This one is not. X-ray imaging buys penetration, speed, and geometric fidelity by delivering photons energetic enough to strip electrons from atoms. Medical tomography is organized around getting the most diagnostic information from the fewest such photons. Dose is not a footnote to this chapter, it is one of its governing variables.
X-ray radiography and computed tomography rest on a single physical principle: photons are attenuated as they pass through tissue, and different tissues attenuate differently. Plain radiography records one two-dimensional projection; CT collects many projections from many angles and computationally reconstructs a three-dimensional map of attenuation.
The conceptual sequence that organizes the chapter is
\[\underbrace{\text{produce photons}}_{\text{tube, spectrum}} \rightarrow \underbrace{\text{shape the beam}}_{\text{kVp, mAs, filtration}} \rightarrow \underbrace{\text{interact in tissue}}_{\text{photoelectric, Compton}} \rightarrow \underbrace{\text{measure projections}}_{\text{line integrals}}\]
\[\rightarrow \underbrace{\text{reconstruct}}_{\text{Radon, FBP, iterative}} \rightarrow \underbrace{\text{calibrate and display}}_{\text{HU, windowing}} \rightarrow \underbrace{\text{justify the dose}}_{\text{ALARA}}\]
Table 1 maps every earlier result this chapter depends on to the section that uses it. The chapter serves three uses. As a course text, each section develops the physics from first principles and anchors it in a clinical decision. As a reference, sections are largely self-contained. As classroom material, every section ends with a slide-ready Section summary and a Checkpoint question, and three computational laboratories build the signature skills: forward projection and the sinogram, filtered back projection, and artifact simulation.
Chapter-level clinical puzzle. A patient with acute abdominal pain needs a diagnosis in minutes. A radiologist must find a 200 µm microcalcification in dense breast tissue. An interventionalist must see a 2 mm cerebral aneurysm through the skull. A trauma team needs to know whether a metal implant is loose. And in every one of these the physician must be able to answer: was this dose justified, and could the same answer have been obtained with less? Four imaging tasks, one physics, and a fifth question that no other chapter in this book has to ask.
Photon story for the chapter. Follow one photon. It is born when an electron accelerated through 120 kV slams into a tungsten anode and is deflected by a nucleus, radiating away part of its kinetic energy. It leaves the tube, passes through an aluminium filter that removes its lower-energy siblings, and enters the patient. Inside, it faces two fates: it can be absorbed photoelectrically, its whole energy transferred to a bound electron, or it can Compton scatter, losing part of its energy and changing direction. If it survives unscattered to the detector it carries information about the total attenuation along one straight line. Millions of such lines, from hundreds of angles, are what a Fourier transform turns back into a picture of the inside of a person.
Key idea. X-ray and CT images are not photographs. They are maps of the linear attenuation coefficient \(\mu(x,y,z)\), which itself is a function of tissue density, effective atomic number, and photon energy. Everything clinical follows from that dependence, and everything about safety follows from the ionizing energy required to measure it.
| # | Objective | Section |
|---|---|---|
| 1 | Explain bremsstrahlung and characteristic X-ray production, and predict how the spectrum changes with kVp, mAs, and filtration. | 5.1 |
| 2 | Compute radiation yield and explain the anode heat problem quantitatively. | 5.1 |
| 3 | Distinguish photoelectric absorption, Compton scattering, coherent scattering, and pair production, and state the energy regime of each. | 5.2 |
| 4 | Use the \(Z^3/E^3\) photoelectric scaling and the near-\(Z\)-independence of Compton to predict how contrast changes with kVp. | 5.2 |
| 5 | Explain K-edges quantitatively and justify the choice of iodine and barium as contrast agents. | 5.2 |
| 6 | Distinguish linear from mass attenuation coefficients, apply the mixture rule, and compute half-value layers. | 5.3 |
| 7 | Explain beam hardening from a polychromatic spectrum and predict cupping and streak artifacts. | 5.3 |
| 8 | Quantify the effect of kVp, mAs, filtration, collimation, and grids on contrast, noise, and dose. | 5.4 |
| 9 | Define absorbed, equivalent, and effective dose, interpret CTDI\(_{\text{vol}}\) and DLP, and compare examination doses with natural background. | 5.5 |
| 10 | Distinguish stochastic from deterministic radiation effects and apply justification and optimization. | 5.5 |
| 11 | Explain why a projection is a line integral and how projections are organized into a sinogram. | 5.7 |
| 12 | State the Radon transform and derive the central slice theorem. | 5.7, 5.8 |
| 13 | Explain why unfiltered back projection blurs, why the ramp filter corrects it, and why practical kernels are band-limited. | 5.9 |
| 14 | Apply the angular sampling requirement \(N_{\text{views}} \gtrsim (\pi/2)N_{\text{det}}\) and recognize view-aliasing streaks. | 5.9 |
| 15 | Derive Hounsfield units from \(\mu\) and reproduce the standard clinical HU table. | 5.10 |
| 16 | Apply window level and width, and explain why the same data supports several displays. | 5.10 |
| 17 | Identify CT artifacts and match each to its physical cause and correction. | 5.11 |
| 18 | Explain dual-energy CT material decomposition and connect it to multi-wavelength unmixing elsewhere in BPAD. | 5.12 |
| 19 | Relate resolution, contrast, noise, and dose through the detectability condition. | 5.13 |
| 20 | Implement forward projection, filtered back projection, and artifact simulation in R. |
5.15 |
| Symbol | Meaning | Units |
|---|---|---|
| \(E,\ h\nu\) | photon energy | keV |
| \(V\), kVp | tube voltage; peak kilovoltage | kV |
| mAs | tube current \(\times\) exposure time | mA s |
| \(Z,\ Z_{\text{eff}}\) | atomic number; effective atomic number of a mixture | — |
| \(\mu\) | linear attenuation coefficient | cm\(^{-1}\) |
| \(\mu/\rho\) | mass attenuation coefficient | cm\(^2\) g\(^{-1}\) |
| \(\mu_{\text{en}}/\rho\) | mass energy-absorption coefficient | cm\(^2\) g\(^{-1}\) |
| \(\rho\) | physical density | g cm\(^{-3}\) |
| \(\tau,\ \sigma_{\text{C}},\ \sigma_{\text{R}}\) | photoelectric, Compton, coherent cross-sections | cm\(^2\) g\(^{-1}\) |
| \(x_{1/2}\) | half-value layer | cm (or mm Al) |
| \(\Phi,\ \Psi\) | photon fluence; energy fluence | cm\(^{-2}\); keV cm\(^{-2}\) |
| \(D\) | absorbed dose | Gy \(=\) J kg\(^{-1}\) |
| \(H,\ E_{\text{eff}}\) | equivalent dose; effective dose | Sv |
| \(w_R,\ w_T\) | radiation and tissue weighting factors | — |
| CTDI\(_{\text{vol}}\), DLP | volume CT dose index; dose-length product | mGy; mGy cm |
| \(P(\theta,s)\) | projection (line integral) at angle \(\theta\), detector position \(s\) | cm\(^{-1}\) cm |
| \(R[f]\) | Radon transform of \(f\) | — |
| \(k\) | spatial frequency | cm\(^{-1}\) |
| \(|k|\) | ramp filter | cm\(^{-1}\) |
| \(N_{\text{det}},\ N_{\text{views}}\) | detector channels; projection angles | — |
| HU | Hounsfield unit | — |
| \(L,\ W\) | window level; window width | HU |
| \(w\) | voxel width | cm |
| \(\delta\mu\) | attenuation difference between lesion and background | cm\(^{-1}\) |
| DQE | detective quantum efficiency | — |
| Earlier result | Chapter | Role here | Section |
|---|---|---|---|
| Exponential decay / first-order ODE | 1, 2, 3 | Beer–Lambert attenuation, half-value layer | 5.3 |
| Photon energy \(E = hc/\lambda\) | 2 | why X-rays ionize and optical photons do not | 5.1 |
| Beer–Lambert law and \(\mu_a\) | 2 | identical mathematics, different interaction physics | 5.3 |
| Multi-wavelength linear unmixing | 2, 3 | dual-energy CT material decomposition | 5.12 |
| Fourier transform, 1D and 2D | 1 | central slice theorem, ramp filter, FBP | 5.8, 5.9 |
| Convolution and windowing | 1, 4 | ramp kernel; band-limited reconstruction kernels | 5.9 |
| Sampling and the Nyquist criterion | 1, 3, 4 | angular and detector sampling; view aliasing | 5.9 |
| Poisson counting statistics, \(\sqrt{N}\) | 1, 6 | quantum noise, dose–noise trade-off | 5.4, 5.13 |
| Point spread function | 1, 3, 4, 7 | focal-spot and detector blur | 5.4 |
| Linear systems \(Af = p\), least squares | 1, 8 | iterative reconstruction | 5.9 |
| Rose criterion / detectability | 1 | lesion visibility and required fluence | 5.13 |
| Back projection in photoacoustics | 3 | same inverse-problem structure, acoustic instead of X-ray | 5.9 |
| Safety indices (MI/TI, SAR) | 3, 4 | the ionizing analogue: absorbed and effective dose | 5.5 |
An X-ray tube accelerates electrons from a heated cathode through a potential difference \(V\) onto a high-atomic-number anode, almost always tungsten. Each electron arrives with kinetic energy \(eV\), and the Duane–Hunt limit sets the maximum photon energy that can be produced:
\[\begin{equation} E_{\max} = eV, \qquad \lambda_{\min} = \frac{hc}{eV} . \tag{1} \end{equation}\]
A tube operated at 120 kVp therefore emits nothing above 120 keV. Note carefully that “kVp” is a peak voltage: the spectrum is continuous below it, and the mean photon energy of a filtered clinical beam is only about one third to one half of the peak.
Bremsstrahlung (“braking radiation”) is produced when an electron is deflected by the Coulomb field of a nucleus and radiates away part of its kinetic energy. Because the deflection can be arbitrarily gentle or violent, the emitted photon can carry any energy from nearly zero up to \(eV\), producing a continuous spectrum. To a good approximation the unfiltered intensity per unit energy follows Kramers’ law,
\[\begin{equation} I(E) \propto Z\,(E_{\max} - E), \qquad 0 \le E \le E_{\max}, \tag{2} \end{equation}\]
a straight line falling to zero at the Duane–Hunt limit. The observed spectrum differs because low-energy photons are strongly removed by the anode itself, the tube window, and added filtration.
If an incident electron ejects an inner-shell electron, an outer-shell electron falls into the vacancy and emits a photon of energy equal to the level difference:
\[\begin{equation} E_{\text{char}} = E_{\text{outer}} - E_{\text{inner}} . \tag{3} \end{equation}\]
Because atomic levels are discrete and element-specific, these appear as sharp lines. For tungsten the K\(\alpha\) lines lie at 57.98 and 59.32 keV and K\(\beta\) near 67.2 keV, and they appear only when the tube voltage exceeds the tungsten K-shell binding energy of 69.5 keV. A tube run at 60 kVp produces bremsstrahlung but no tungsten K lines at all.
The fraction of incident electron energy converted into X-rays is the radiation yield, approximately
\[\begin{equation} Y \approx 9\times10^{-10}\, Z\, V \qquad (V \text{ in volts}), \tag{4} \end{equation}\]
which for tungsten (\(Z = 74\)) at 120 kV gives \(Y \approx 0.8\%\). More than 99% of the electron energy becomes heat in the anode. This single number explains rotating anodes, molybdenum-graphite anode discs, oil cooling, focal-spot sizing, and the exposure limits that force a radiographer to wait between exposures.
## Aluminium mass attenuation, cm^2/g, interpolated on a log-log grid
al_E <- c(10, 15, 20, 30, 40, 50, 60, 80, 100, 150)
al_murho <- c(26.23, 7.955, 3.441, 1.128, 0.5685, 0.3681, 0.2778, 0.2018, 0.1704, 0.1378)
mu_al <- function(E) exp(approx(log(al_E), log(al_murho), log(pmax(E, 10)), rule = 2)$y)*2.70
kramers <- function(E, kVp) pmax(kVp - E, 0) # Eq. (kramers), Z folded into scale
filt <- function(E, t_mm) exp(-mu_al(E)*t_mm/10) # t in mm Al
Eg <- seq(10, 140, by = 0.5)
raw120 <- kramers(Eg, 120)
fil120 <- raw120*filt(Eg, 2.5)
sc <- max(raw120)
p_filt <- ggplot(rbind(
data.frame(E = Eg, I = raw120/sc, q = "unfiltered (Kramers)"),
data.frame(E = Eg, I = fil120/sc, q = "after 2.5 mm Al")),
aes(E, I, colour = q)) +
geom_line(linewidth = 0.95) +
scale_colour_manual(values = bpad_pal[c(7,1)]) +
labs(x = "Photon energy (keV)", y = "Relative photon number",
title = "Filtration removes the low-energy tail",
subtitle = "120 kVp")
kv <- c(70, 100, 140)
spec_df <- do.call(rbind, lapply(kv, function(k) {
s <- kramers(Eg, k)*filt(Eg, 2.5)
data.frame(E = Eg, I = s/sc, kVp = factor(paste0(k, " kVp"),
levels = paste0(kv, " kVp")))
}))
p_kvp <- ggplot(spec_df, aes(E, I, colour = kVp)) +
geom_line(linewidth = 0.95) +
geom_segment(data = data.frame(E = c(59.3, 67.2), y0 = 0,
y1 = c(0.34, 0.20)),
aes(x = E, xend = E, y = y0, yend = y1),
inherit.aes = FALSE, colour = bpad_pal[2], linewidth = 1.1) +
annotate("text", x = 63, y = 0.40, label = "W K lines\n(only above 69.5 kVp)",
size = 2.7, colour = bpad_pal[2]) +
scale_colour_manual(values = bpad_pal[c(3,1,4)]) +
labs(x = "Photon energy (keV)", y = "Relative photon number",
title = "Raising kVp raises endpoint and mean energy")
bpad_grid(p_filt, p_kvp, ncol = 2)Figure 1: Model X-ray tube spectra. Left: the unfiltered Kramers continuum, Eq. (kramers), together with the effect of 2.5 mm of aluminium filtration; filtration removes the low-energy photons that would otherwise deposit skin dose without reaching the detector, and in doing so raises the mean beam energy. Right: filtered spectra at three tube voltages, with the tungsten characteristic K lines that appear only once the tube voltage exceeds the 69.5 keV K-edge. Raising kVp raises both the endpoint and the mean energy; raising mAs would scale these curves vertically without changing their shape at all.
mean_E <- function(kVp, tmm) {
s <- kramers(Eg, kVp)*filt(Eg, tmm); sum(s*Eg)/sum(s)
}
knitr::kable(data.frame(
kVp = c(70, 100, 120, 140),
mean_E_unfiltered = round(sapply(c(70,100,120,140), mean_E, tmm = 0), 1),
mean_E_2p5mm_Al = round(sapply(c(70,100,120,140), mean_E, tmm = 2.5), 1),
ratio_mean_to_peak = round(sapply(c(70,100,120,140), mean_E, tmm = 2.5)/
c(70,100,120,140), 2),
yield_percent = round(100*9e-10*74*c(70,100,120,140)*1e3, 2)),
col.names = c("kVp","Mean E, unfiltered (keV)","Mean E, 2.5 mm Al (keV)",
"Mean / peak","Radiation yield (%)"),
caption = "Mean photon energy is roughly half the peak kilovoltage, and radiation yield is under one percent: the anode is overwhelmingly a heater that happens to emit some X-rays.")| kVp | Mean E, unfiltered (keV) | Mean E, 2.5 mm Al (keV) | Mean / peak | Radiation yield (%) |
|---|---|---|---|---|
| 70 | 29.8 | 40.9 | 0.58 | 0.47 |
| 100 | 39.8 | 52.3 | 0.52 | 0.67 |
| 120 | 46.5 | 59.6 | 0.50 | 0.80 |
| 140 | 53.2 | 66.8 | 0.48 | 0.93 |
Worked Example 5.1 (the anode heat problem). A chest radiograph uses 120 kVp and 5 mAs. The electrical energy delivered is
\[W = V \times I \times t = (120\times10^{3}\ \text{V})(5\times10^{-3}\ \text{A s}) = 600\ \text{J}.\]
With a yield of \(Y \approx 0.8\%\), only about 4.8 J leaves as X-rays and roughly 595 J is deposited as heat in a focal spot a millimetre or so across. That is a power density comparable to a welding arc, which is why the anode rotates at up to 10,000 rpm to spread the load over a track rather than a point. A CT scan delivering hundreds of mAs per rotation makes the heat problem, not the physics, the binding engineering constraint.
Section 5.1 summary.
Checkpoint 5.1. A tube is operated at 60 kVp. Will the tungsten K\(\alpha\) lines appear in the spectrum? Explain using Eq. (3), and state what changes if the anode is molybdenum (\(E_K = 20.0\) keV) instead — the choice made in mammography.
This section is the physical foundation of everything that follows. Contrast in X-ray imaging exists because different tissues attenuate differently, and they attenuate differently because of two competing interaction mechanisms with completely different dependences on atomic number and photon energy. Without this, statements like “higher kVp reduces contrast” are rules to memorize; with it, they are consequences.
| Interaction | Energy regime | \(Z\) dependence | Effect on the image |
|---|---|---|---|
| Coherent (Rayleigh) scattering | low (\(<30\) keV) | \(\sim Z^2/E^2\) | small-angle scatter; a few percent of events |
| Photoelectric absorption | low to mid | \(\propto Z^3/E^3\) | the source of contrast; photon vanishes |
| Compton scattering | mid to high | \(\propto\) electron density, nearly \(Z\)-independent | the source of scatter degradation |
| Pair production | \(>1.022\) MeV | \(\propto Z^2\) | irrelevant at diagnostic energies |
Only the middle two matter in the 20–140 keV diagnostic range.
A photon is completely absorbed by a bound electron, which is ejected with kinetic energy \(E - E_{\text{binding}}\). The interaction requires a bound electron, so it is strongly favoured when the photon energy is just above a shell binding energy and when the atom holds electrons tightly. The mass attenuation cross-section scales approximately as
\[\begin{equation} \frac{\tau}{\rho} \;\propto\; \frac{Z^{3}}{E^{3}} . \tag{5} \end{equation}\]
Both powers matter enormously. The \(Z^3\) makes bone (calcium, \(Z = 20\)) and iodine (\(Z = 53\)) stand out dramatically from soft tissue (\(Z_{\text{eff}} \approx 7.4\)). The \(E^{-3}\) makes the effect collapse as the beam hardens: doubling the photon energy cuts photoelectric absorption roughly eightfold.
A photon scatters off an essentially free outer electron, transferring part of its energy and changing direction by \(\theta\). Conservation of energy and momentum give the Compton relation
\[\begin{equation} E' = \frac{E}{1 + \dfrac{E}{m_ec^{2}}(1 - \cos\theta)}, \qquad m_ec^{2} = 511\ \text{keV} . \tag{6} \end{equation}\]
Because the electron is treated as free, the cross-section depends on electron density per gram, not on atomic number. Electrons per gram is \(N_A Z/A \approx N_A/2\) for essentially every element except hydrogen, so Compton attenuation per gram is nearly the same for fat, muscle, and bone. It therefore carries almost no tissue-discriminating information, and the scattered photon that reaches the detector from the wrong direction actively destroys information (Section 5.4.4).
The two mechanisms cross over at an energy that depends on \(Z\): around 25 keV in soft tissue, around 40 keV in bone, and much higher in iodine. Below the crossover, contrast is excellent and dose is high; above it, penetration is excellent and contrast is poor.
Ef <- exp(seq(log(20), log(150), length.out = 200))
interp <- function(col) exp(approx(log(atten$E_keV), log(atten[[col]]), log(Ef))$y)
df_mu <- rbind(
data.frame(E = Ef, mu = interp("water"), tissue = "Soft tissue (water)"),
data.frame(E = Ef, mu = interp("adipose"), tissue = "Adipose"),
data.frame(E = Ef, mu = interp("bone"), tissue = "Cortical bone"))
p_mu <- ggplot(df_mu, aes(E, mu, colour = tissue)) +
geom_line(linewidth = 0.95) +
scale_x_log10() + scale_y_log10() +
scale_colour_manual(values = bpad_pal[c(5,1,2)]) +
labs(x = "Photon energy (keV, log scale)",
y = expression(mu/rho*" (cm"^2*"/g, log scale)"),
title = "Mass attenuation coefficients")
ratio <- data.frame(E = Ef, r = interp("bone")/interp("water"))
p_rat <- ggplot(ratio, aes(E, r)) +
geom_hline(yintercept = 1, linetype = "dotted", colour = "grey50") +
geom_line(linewidth = 1, colour = bpad_pal[2]) +
annotate("rect", xmin = 24, xmax = 32, ymin = 0.9, ymax = 5.2,
fill = bpad_pal[3], alpha = 0.15) +
annotate("text", x = 28, y = 4.6, label = "mammography", size = 2.8,
colour = bpad_pal[3]) +
annotate("rect", xmin = 50, xmax = 80, ymin = 0.9, ymax = 5.2,
fill = bpad_pal[1], alpha = 0.12) +
annotate("text", x = 64, y = 4.0, label = "CT / general\nradiography",
size = 2.8, colour = bpad_pal[1]) +
scale_x_log10() +
labs(x = "Photon energy (keV, log scale)",
y = "Bone / soft-tissue attenuation ratio",
title = "Contrast collapses as the beam hardens")
bpad_grid(p_mu, p_rat, ncol = 2)Figure 2: Why kVp controls contrast. Left: measured mass attenuation coefficients for soft tissue and cortical bone; both fall with energy, but bone falls far faster because its photoelectric component carries the Z-cubed advantage that disappears above about 60 keV. Right: the resulting bone-to-soft-tissue attenuation ratio, which collapses from nearly 5 at 20 keV to essentially 1 at 150 keV. This single curve explains why mammography uses 26 to 30 kVp, why chest radiography uses 120 kVp to deliberately suppress rib contrast, and why CT must trade contrast for penetration.
knitr::kable(data.frame(
E_keV = atten$E_keV,
water = atten$water, bone = atten$bone,
bone_over_water = round(atten$bone/atten$water, 2),
adipose_over_water = round(atten$adipose/atten$water, 3)),
col.names = c("E (keV)","Water mu/rho","Bone mu/rho","Bone / water",
"Adipose / water"),
caption = "Contrast between tissues falls monotonically with photon energy. At 150 keV bone and soft tissue attenuate almost identically per gram, and the only remaining contrast comes from density.")| E (keV) | Water mu/rho | Bone mu/rho | Bone / water | Adipose / water |
|---|---|---|---|---|
| 20 | 0.8096 | 4.0010 | 4.94 | 0.778 |
| 30 | 0.3756 | 1.3310 | 3.54 | 0.844 |
| 40 | 0.2683 | 0.6655 | 2.48 | 0.881 |
| 50 | 0.2269 | 0.4242 | 1.87 | 0.910 |
| 60 | 0.2059 | 0.3148 | 1.53 | 0.925 |
| 80 | 0.1837 | 0.2229 | 1.21 | 0.936 |
| 100 | 0.1707 | 0.1855 | 1.09 | 0.938 |
| 150 | 0.1505 | 0.1480 | 0.98 | 0.944 |
cat(sprintf("Bone-to-soft-tissue ratio: %.2f at 30 keV, %.2f at 100 keV -- a %.0f%% loss.\n",
atten$bone[2]/atten$water[2], atten$bone[7]/atten$water[7],
100*(1 - (atten$bone[7]/atten$water[7] - 1)/(atten$bone[2]/atten$water[2] - 1))))## Bone-to-soft-tissue ratio: 3.54 at 30 keV, 1.09 at 100 keV -- a 97% loss.
cat(sprintf("Adipose-to-water differs by only %.1f%% at 30 keV and %.1f%% at 100 keV:\n",
100*(1 - atten$adipose[2]/atten$water[2]),
100*(1 - atten$adipose[7]/atten$water[7])))## Adipose-to-water differs by only 15.6% at 30 keV and 6.2% at 100 keV:
## this is the entire difficulty of mammography, and why it operates near 28 kVp.
From the equation to the protocol. Equation (5) is a design rule. Mammography must separate adipose from glandular tissue, whose mass attenuation coefficients differ by only about 16% at 30 keV and by 6% at 100 keV. There is no receiver gain that recovers a difference the physics did not create, so mammography operates at 26–30 kVp with a molybdenum or rhodium anode, accepting high skin dose and poor penetration in exchange for the \(Z^3/E^3\) contrast advantage. Chest radiography faces the opposite problem — ribs obscuring lung — and deliberately uses 120 kVp to suppress bone contrast. Same equation, opposite protocol.
The photoelectric cross-section rises steeply just above the binding energy of a shell, because a photon that could not eject a K-shell electron suddenly can. This produces a K-edge: a discontinuous jump in attenuation.
The clinical consequence is decisive. Iodine’s K-edge lies at 33.2 keV and barium’s at 37.4 keV, both comfortably inside the diagnostic spectrum, so both become extraordinarily attenuating exactly where clinical beams have plenty of photons. Calcium’s K-edge, at 4.0 keV, is far below any usable photon energy, so bone’s contrast comes from the smooth \(Z^3/E^3\) tail rather than from an edge.
p_ke <- ggplot(iodine, aes(E_keV, murho)) +
geom_line(linewidth = 0.95, colour = bpad_pal[2]) +
geom_point(size = 1.5, colour = bpad_pal[2]) +
geom_line(data = data.frame(E_keV = Ef, murho = interp("water")),
aes(E_keV, murho), colour = bpad_pal[1], linewidth = 0.95) +
geom_vline(xintercept = k_edges["Iodine"], linetype = "dotted", colour = "grey40") +
annotate("text", x = 36, y = 45, hjust = 0, size = 2.9, colour = "grey30",
label = "K-edge\n33.2 keV") +
annotate("text", x = 100, y = 0.35, label = "soft tissue", size = 2.9,
colour = bpad_pal[1]) +
annotate("text", x = 100, y = 4.0, label = "iodine", size = 2.9,
colour = bpad_pal[2]) +
scale_x_log10() + scale_y_log10() +
labs(x = "Photon energy (keV, log scale)",
y = expression(mu/rho*" (cm"^2*"/g, log)"),
title = "Iodine versus soft tissue")
Ei <- iodine$E_keV
rat_i <- iodine$murho/exp(approx(log(atten$E_keV), log(atten$water), log(Ei), rule = 2)$y)
p_ri <- ggplot(data.frame(E = Ei, r = rat_i), aes(E, r)) +
geom_line(linewidth = 0.95, colour = bpad_pal[4]) +
geom_point(size = 1.6, colour = bpad_pal[4]) +
geom_vline(xintercept = k_edges["Iodine"], linetype = "dotted", colour = "grey40") +
scale_x_log10() + scale_y_log10() +
labs(x = "Photon energy (keV, log scale)",
y = "Iodine / soft-tissue ratio (log)",
title = "Contrast-agent advantage per gram")
bpad_grid(p_ke, p_ri, ncol = 2)Figure 3: The iodine K-edge and why contrast agents work. Left: the mass attenuation coefficient of iodine jumps by a factor of about 5.5 across its K-edge at 33.2 keV, and remains far above soft tissue throughout the diagnostic range. Right: the iodine-to-soft-tissue attenuation ratio per gram, which exceeds 50 across most of the diagnostic band. This is why a few milligrams of iodine per millilitre transforms an invisible vessel into a bright one, and why dual-energy CT can identify iodine specifically by imaging on either side of the edge.
knitr::kable(data.frame(
element = names(k_edges), K_edge_keV = as.numeric(k_edges),
in_diagnostic_beam = ifelse(k_edges > 15 & k_edges < 120, "yes", "no"),
use = c("bone mineral (edge unusable)", "vascular / CT contrast",
"gastrointestinal contrast", "MRI contrast; CT K-edge imaging",
"anode material (sets K lines)")),
col.names = c("Element","K-edge (keV)","Inside the diagnostic beam?","Role"),
caption = "K-edge energies. Iodine and barium are chosen not merely because Z is high but because their K-edges land inside the usable photon spectrum.")| Element | K-edge (keV) | Inside the diagnostic beam? | Role | |
|---|---|---|---|---|
| Calcium | Calcium | 4.04 | no | bone mineral (edge unusable) |
| Iodine | Iodine | 33.17 | yes | vascular / CT contrast |
| Barium | Barium | 37.44 | yes | gastrointestinal contrast |
| Gadolinium | Gadolinium | 50.24 | yes | MRI contrast; CT K-edge imaging |
| Tungsten | Tungsten | 69.53 | yes | anode material (sets K lines) |
cat(sprintf("Iodine K-edge jump: mu/rho goes from %.2f to %.2f cm^2/g, a factor of %.1f\n",
6.10, 33.40, 33.40/6.10))## Iodine K-edge jump: mu/rho goes from 6.10 to 33.40 cm^2/g, a factor of 5.5
cat(sprintf("across roughly 0.4 keV. At 50 keV iodine attenuates %.0f times more per gram\n",
12.20/0.2269))## across roughly 0.4 keV. At 50 keV iodine attenuates 54 times more per gram
## than soft tissue, which is why milligram quantities produce vivid contrast.
Two mechanisms, two consequences. Photoelectric absorption is the friend: it depends steeply on \(Z\) and \(E\), it creates essentially all diagnostic contrast, and the photon is removed cleanly from the beam. Compton scattering is the enemy: it is nearly \(Z\)-independent so it carries almost no tissue information, and the scattered photon may still reach the detector, adding a background that dilutes contrast. Raising kVp shifts the balance from friend to enemy — which is exactly why penetration and contrast trade against each other.
Section 5.2 summary.
Checkpoint 5.2. A radiologist wants to see rib fractures and switches from 120 kVp to 70 kVp. Using Eq. (5) and the bone-to-soft-tissue ratio figure, predict what happens to bone contrast, to patient dose, and to the required mAs, and state which of the three the radiologist is most likely to have forgotten.
For a monoenergetic beam through a uniform absorber, the fractional loss per unit thickness is constant, giving the first-order ODE \(dI/dx = -\mu I\) and hence
\[\begin{equation} I(x) = I_0\,e^{-\mu x}, \qquad \mu = -\frac{1}{x}\ln\!\frac{I}{I_0}, \qquad x_{1/2} = \frac{\ln 2}{\mu} . \tag{7} \end{equation}\]
Equation (7) is the same equation as optical absorption in Chapter 2 and acoustic attenuation in Chapter 3. Only the interaction physics differs.
The linear coefficient \(\mu\) (cm\(^{-1}\)) depends on density and therefore on physical state; the mass coefficient \(\mu/\rho\) (cm\(^2\) g\(^{-1}\)) does not. Water and steam have identical \(\mu/\rho\) but wildly different \(\mu\). Tabulated data are always mass coefficients, and the conversion is simply \(\mu = (\mu/\rho)\rho\).
For a mixture or compound, mass coefficients combine by mass fraction:
\[\begin{equation} \left(\frac{\mu}{\rho}\right)_{\text{mix}} = \sum_i w_i \left(\frac{\mu}{\rho}\right)_i , \qquad \sum_i w_i = 1 . \tag{8} \end{equation}\]
Equation (8) is what lets us compute the attenuation of blood carrying a known iodine concentration, or of bone as a calcium-hydroxyapatite–collagen–water mixture.
A common error. Mass fractions, not volume or molar fractions, enter Eq. (8). Mixing them up produces errors of tens of percent in contrast-agent calculations, and always in the direction of overestimating the agent’s effect.
For two adjacent paths of the same thickness \(x\) through materials with coefficients \(\mu_1\) and \(\mu_2\), the transmitted intensities differ, and the subject contrast is
\[\begin{equation} C = \frac{I_1 - I_2}{I_1} = 1 - e^{-(\mu_2-\mu_1)x} \approx (\mu_2 - \mu_1)\,x \quad\text{for small } (\Delta\mu)x . \tag{9} \end{equation}\]
Contrast therefore grows with the difference in attenuation coefficients and with thickness — but \(\Delta\mu\) itself collapses with energy (Section 5.2.4), which is the whole trade.
Real beams are polychromatic. Because low-energy photons are removed preferentially (Eq. (5)), the surviving beam becomes progressively more penetrating with depth: its mean energy rises and its effective attenuation coefficient falls. The attenuation is therefore not a single exponential, and \(-\ln(I/I_0)\) is not proportional to thickness — which breaks the line-integral assumption on which CT reconstruction rests.
## Water mass attenuation, interpolated from the chapter table
mu_water_E <- function(E) exp(approx(log(atten$E_keV), log(atten$water),
log(pmin(pmax(E, 20), 150)), rule = 2)$y)*1.0
Eb <- seq(20, 120, by = 0.5)
spec <- pmax(120 - Eb, 0)*exp(-mu_al(Eb)*0.25) # 120 kVp, 2.5 mm Al
spec <- spec/sum(spec)
xs <- seq(0, 30, by = 0.25)
poly_I <- sapply(xs, function(x) sum(spec*Eb*exp(-mu_water_E(Eb)*x)))
poly_I <- poly_I/poly_I[1]
mono_I <- exp(-mu_water_E(60)*xs)
Ebar <- sapply(xs, function(x) {
w <- spec*Eb*exp(-mu_water_E(Eb)*x); sum(w*Eb)/sum(w) })
mu_eff <- c(NA, -diff(log(poly_I))/diff(xs))
p_bh1 <- ggplot(rbind(
data.frame(x = xs, I = mono_I, q = "monoenergetic 60 keV"),
data.frame(x = xs, I = poly_I, q = "polychromatic 120 kVp")),
aes(x, I, colour = q)) +
geom_line(linewidth = 0.95) +
scale_y_log10() + scale_colour_manual(values = bpad_pal[c(7,2)]) +
labs(x = "Water thickness (cm)", y = "Transmitted energy fluence (log)",
title = "Polychromatic beams bend")
p_bh2 <- ggplot(data.frame(x = xs, E = Ebar), aes(x, E)) +
geom_line(linewidth = 0.95, colour = bpad_pal[1]) +
labs(x = "Water thickness (cm)", y = "Mean photon energy (keV)",
title = "The beam hardens with depth")
p_bh3 <- ggplot(data.frame(x = xs, m = mu_eff), aes(x, m)) +
geom_line(linewidth = 0.95, colour = bpad_pal[3], na.rm = TRUE) +
labs(x = "Water thickness (cm)",
y = expression("effective "*mu*" (cm"^-1*")"),
title = "Effective attenuation falls")
bpad_grid(p_bh1, p_bh2, p_bh3, ncol = 3)Figure 4: Beam hardening computed from a real spectrum rather than assumed. Left: transmitted energy fluence through water for a monoenergetic 60 keV beam and for a filtered 120 kVp polychromatic beam; the polychromatic curve bends upward on a logarithmic axis, meaning it attenuates less and less per centimetre. Centre: the mean photon energy of the surviving beam, which hardens from 68 keV at the surface to 84 keV at 30 cm. Right: the resulting effective attenuation coefficient, which falls by about 15 percent across that path. Because CT assumes a fixed coefficient along each ray, this deviation is exactly what produces cupping and streak artifacts.
cat(sprintf("mean beam energy: %.1f keV at the surface -> %.1f keV at 30 cm\n",
Ebar[1], Ebar[length(Ebar)]))## mean beam energy: 68.3 keV at the surface -> 84.0 keV at 30 cm
cat(sprintf("effective mu: %.4f /cm near the surface -> %.4f /cm at 30 cm (a %.0f%% fall)\n",
mu_eff[2], mu_eff[length(mu_eff)],
100*(1 - mu_eff[length(mu_eff)]/mu_eff[2])))## effective mu: 0.2170 /cm near the surface -> 0.1847 /cm at 30 cm (a 15% fall)
## A CT reconstruction that assumes a single mu per ray must therefore mis-assign
## attenuation, producing the cupping artifact demonstrated in Section 5.11.
Section 5.3 summary.
Checkpoint 5.3. A tissue has \(\mu/\rho = 0.19\) cm\(^2\)g\(^{-1}\) at 70 keV and density 1.06 g cm\(^{-3}\). Compute \(\mu\), the half-value layer, and the fraction transmitted through 20 cm. Then explain why the same \(\mu/\rho\) in lung (\(\rho = 0.26\)) gives a completely different image.
| Control | Physical effect | Contrast | Noise | Dose |
|---|---|---|---|---|
| kVp \(\uparrow\) | raises endpoint and mean energy | \(\downarrow\) (photoelectric collapses) | \(\downarrow\) | \(\uparrow\) (roughly as kVp\(^2\)–kVp\(^3\)) |
| mAs \(\uparrow\) | more photons, same spectrum | unchanged | \(\downarrow\) as \(1/\sqrt{\text{mAs}}\) | \(\uparrow\) linearly |
| Filtration \(\uparrow\) | removes low-energy photons | slight \(\downarrow\) | slight \(\uparrow\) | \(\downarrow\) (skin dose) |
| Collimation \(\uparrow\) | smaller irradiated field | \(\uparrow\) (less scatter) | — | \(\downarrow\) |
The asymmetry is the point: mAs is the clean knob. It changes noise and dose in a known way and leaves contrast alone. kVp changes everything at once.
X-ray detection is photon counting, so the number of photons in a pixel is Poisson distributed (Chapter 1) with \(\sigma_N = \sqrt{N}\). The signal-to-noise ratio is therefore
\[\begin{equation} \mathrm{SNR} = \frac{N}{\sqrt{N}} = \sqrt{N} \;\propto\; \sqrt{\text{mAs}} , \tag{10} \end{equation}\]
Equation (10) is the single most consequential relation in radiographic dose management: halving the noise costs four times the dose.
set.seed(7)
n <- 90
gx <- matrix(seq(-1, 1, length.out = n), n, n)
gy <- t(gx)
obj <- 1 - 0.35*(gx^2 + gy^2 < 0.45^2) - 0.12*((gx-0.4)^2 + (gy+0.35)^2 < 0.13^2)
levels_mAs <- c(1, 4, 16, 64)
N0 <- 25
op <- par(mfrow = c(2, 2), mar = c(0.5, 0.5, 2.2, 0.5))
for (m in levels_mAs) {
cnt <- matrix(rpois(n*n, N0*m*obj), n, n)
image(cnt, col = gpal, axes = FALSE, asp = 1,
main = sprintf("%dx mAs", m)); box(col = "grey70")
}Figure 5: Quantum noise and the cost of reducing it. Left: simulated Poisson images of the same object at four exposure levels, illustrating that image quality improves only as the square root of dose. Right: measured noise versus relative mAs, following the inverse-square-root law of Eq. (quantum-noise); the dotted lines show that halving the noise requires quadrupling the exposure and therefore the dose.
par(op)
mgrid <- 2^seq(0, 7, by = 0.25)
noise <- sapply(mgrid, function(m) {
cnt <- matrix(rpois(n*n, N0*m*obj), n, n)/(N0*m)
sd(cnt[gx^2 + gy^2 > 0.7^2])
})
p_nd <- ggplot(data.frame(m = mgrid, s = noise/noise[1]), aes(m, s)) +
geom_line(linewidth = 0.95, colour = bpad_pal[1]) +
geom_point(size = 1.3, colour = bpad_pal[1]) +
geom_segment(aes(x = 1, xend = 4, y = 0.5, yend = 0.5), linetype = "dotted",
colour = bpad_pal[2]) +
geom_segment(aes(x = 4, xend = 4, y = 0.15, yend = 0.5), linetype = "dotted",
colour = bpad_pal[2]) +
annotate("text", x = 5, y = 0.62, hjust = 0, size = 3, colour = bpad_pal[2],
label = "half the noise\ncosts 4x the dose") +
scale_x_log10() + scale_y_log10() +
labs(x = "Relative mAs (log scale)", y = "Relative noise (log scale)",
title = "Noise falls only as one over the square root of dose")
print(p_nd)Figure 6: Quantum noise and the cost of reducing it. Left: simulated Poisson images of the same object at four exposure levels, illustrating that image quality improves only as the square root of dose. Right: measured noise versus relative mAs, following the inverse-square-root law of Eq. (quantum-noise); the dotted lines show that halving the noise requires quadrupling the exposure and therefore the dose.
cat(sprintf("fitted slope of log(noise) on log(mAs) = %.3f (theory -0.500)\n",
coef(lm(log(noise) ~ log(mgrid)))[2]))## fitted slope of log(noise) on log(mAs) = -0.500 (theory -0.500)
Filtration (typically 2.5 mm aluminium equivalent, or K-edge filters such as rhodium in mammography) removes photons too soft to traverse the patient. Those photons would deposit skin dose and contribute nothing to the image. Section 5.1’s figure showed the effect directly: filtration raises the mean beam energy at modest cost in contrast.
Collimation restricts the field to the anatomy of interest. It reduces integral dose and, because scatter is produced in the irradiated volume, reduces scatter roughly in proportion to the irradiated area.
Compton-scattered photons arrive at the detector from the wrong direction and add a spatially smooth background. Their effect is measured by the scatter-to-primary ratio (SPR), which can exceed 4:1 in an adult abdomen. Scatter degrades contrast by the dilution factor
\[\begin{equation} C_{\text{observed}} = \frac{C_{\text{primary}}}{1 + \mathrm{SPR}} . \tag{11} \end{equation}\]
An antiscatter grid of lead strips absorbs obliquely travelling photons. Its performance is described by the grid ratio \(r = h/D\) (strip height over interspace width) and the Bucky factor, the multiplicative increase in exposure required to restore the detector signal:
\[\begin{equation} B = \frac{\text{exposure with grid}}{\text{exposure without grid}} . \tag{12} \end{equation}\]
SPR <- seq(0, 6, by = 0.05)
p_sc <- ggplot(data.frame(s = SPR, c = 1/(1 + SPR)), aes(s, c)) +
geom_line(linewidth = 0.95, colour = bpad_pal[1]) +
geom_vline(xintercept = c(0.5, 4), linetype = "dotted", colour = "grey45") +
annotate("text", x = 0.6, y = 0.85, hjust = 0, size = 2.9, colour = "grey30",
label = "extremity") +
annotate("text", x = 4.1, y = 0.55, hjust = 0, size = 2.9, colour = "grey30",
label = "adult abdomen,\nno grid") +
labs(x = "Scatter-to-primary ratio", y = "Fraction of primary contrast retained",
title = "Scatter dilutes contrast")
grids <- data.frame(
ratio = c(0, 5, 8, 10, 12, 16),
SPR_out = c(4.0, 1.2, 0.7, 0.5, 0.4, 0.25),
bucky = c(1.0, 2.0, 3.0, 3.5, 4.0, 5.5))
grids$contrast_retained <- 1/(1 + grids$SPR_out)
p_gr <- ggplot(grids, aes(bucky, contrast_retained)) +
geom_path(linewidth = 0.9, colour = bpad_pal[2]) +
geom_point(size = 2.4, colour = bpad_pal[2]) +
geom_text(aes(label = ifelse(ratio == 0, "no grid", paste0("r = ", ratio))),
vjust = -0.9, size = 2.8) +
coord_cartesian(ylim = c(0, 1.05)) +
labs(x = "Bucky factor (relative dose)", y = "Fraction of contrast retained",
title = "Grids buy contrast with dose")
bpad_grid(p_sc, p_gr, ncol = 2)Figure 7: The scatter dilemma quantified. Left: observed contrast as a fraction of primary contrast, Eq. (scatter-contrast); at a scatter-to-primary ratio of 4, typical of an adult abdomen without a grid, four fifths of the available contrast is destroyed. Right: the grid trade-off; higher grid ratios remove more scatter and recover more contrast, but the Bucky factor means the patient dose must rise to compensate. Grids are therefore used for thick body parts and omitted for extremities, where scatter is low and the dose penalty is not worth paying.
knitr::kable(grids[, c("ratio","SPR_out","bucky","contrast_retained")],
col.names = c("Grid ratio","Residual SPR","Bucky factor","Contrast retained"),
digits = 2,
caption = "Representative grid performance. Moving from no grid to r = 10 recovers most of the contrast at roughly 3.5 times the dose, which is worthwhile for an abdomen and wasteful for a wrist.")| Grid ratio | Residual SPR | Bucky factor | Contrast retained |
|---|---|---|---|
| 0 | 4.00 | 1.0 | 0.20 |
| 5 | 1.20 | 2.0 | 0.45 |
| 8 | 0.70 | 3.0 | 0.59 |
| 10 | 0.50 | 3.5 | 0.67 |
| 12 | 0.40 | 4.0 | 0.71 |
| 16 | 0.25 | 5.5 | 0.80 |
Modern systems use flat-panel digital detectors: either indirect (a scintillator such as CsI converts X-rays to light, read out by a photodiode array) or direct (a photoconductor such as amorphous selenium produces charge directly). The figure of merit that combines absorption efficiency, conversion gain, and added electronic noise is the detective quantum efficiency
\[\begin{equation} \mathrm{DQE}(f) = \frac{\mathrm{SNR}_{\text{out}}^{2}(f)}{\mathrm{SNR}_{\text{in}}^{2}(f)} \le 1 , \tag{13} \end{equation}\]
a spatial-frequency-dependent measure of how much of the information carried by the incident photons survives detection. A detector with DQE 0.6 needs roughly 1.7 times the dose of an ideal one to achieve the same image quality — which is why detector upgrades are genuinely dose-reduction interventions.
Photon-counting detectors, now entering clinical CT, register individual photons and their energies rather than integrating total energy. They eliminate electronic noise, remove the energy weighting that suppresses low-energy (high-contrast) photons, and provide spectral information for free (Section 5.12).
Section 5.4 summary.
Checkpoint 5.4. A technologist doubles mAs to reduce noise in a noisy radiograph. By what factor does the noise fall, and by what factor does the dose rise? Would raising kVp instead have been better or worse, and on which axis?
Every previous imaging chapter in this book described a non-ionizing modality. This one does not, and the difference is not rhetorical: X-ray photons carry enough energy to eject electrons from atoms, break chemical bonds, and damage DNA. A biomedical physics treatment of CT that omits dosimetry is incomplete in the same way that a treatment of ultrasound without the mechanical index would be.
\[\begin{equation} D = \frac{dE}{dm}\ [\mathrm{Gy} = \mathrm{J\,kg^{-1}}], \qquad H_T = w_R\, D_T\ [\mathrm{Sv}], \qquad E_{\text{eff}} = \sum_T w_T\, H_T\ [\mathrm{Sv}] . \tag{14} \end{equation}\]
Absorbed dose \(D\) is pure physics: energy deposited per unit mass. Equivalent dose \(H\) weights it by radiation type through \(w_R\), which is 1 for X-rays and photons but 20 for alpha particles (relevant in Chapter 6). Effective dose \(E_{\text{eff}}\) weights each organ by its radiosensitivity \(w_T\) and sums, producing a single whole-body-equivalent number that permits comparison across very different examinations. For a monoenergetic beam,
\[\begin{equation} D = \Psi\left(\frac{\mu_{\text{en}}}{\rho}\right) = (h\nu)\,\Phi\left(\frac{\mu_{\text{en}}}{\rho}\right), \tag{15} \end{equation}\]
where \(\mu_{\text{en}}/\rho\) is the mass energy-absorption coefficient — distinct from \(\mu/\rho\), because scattered energy that leaves the volume is not deposited in it.
CT reports two numbers on every study.
CTDI\(_{\text{vol}}\) (mGy) is the average dose within the scanned slice, measured in a standard 16 cm (head) or 32 cm (body) phantom. It describes dose intensity and does not depend on scan length.
DLP \(=\) CTDI\(_{\text{vol}} \times\) scan length (mGy cm) describes the total energy imparted. Effective dose is estimated from it by a body-region conversion factor \(k\):
\[\begin{equation} E_{\text{eff}} \approx k \times \mathrm{DLP} . \tag{16} \end{equation}\]
CTDI\(_{\text{vol}}\) is not patient dose. It is a phantom-based output index of the scanner, not the dose to the individual in front of you. A small child scanned at an adult CTDI\(_{\text{vol}}\) receives a substantially higher organ dose than the number suggests, because the same radiation is deposited in less mass. Size-specific dose estimates (SSDE) exist precisely to correct this, and paediatric protocols must be adjusted rather than inherited.
exams <- data.frame(
exam = c("Chest radiograph (PA)", "Extremity radiograph", "Mammogram (2 view)",
"Abdominal radiograph", "Head CT", "Chest CT", "Chest CT (low dose)",
"Abdomen/pelvis CT", "CT angiography", "Multiphase abdomen CT"),
mSv = c(0.02, 0.001, 0.4, 0.7, 2.0, 7.0, 1.5, 10.0, 12.0, 25.0))
bkgd_yr <- 3.0
exams$days <- exams$mSv/(bkgd_yr/365)
exams$exam <- factor(exams$exam, levels = exams$exam[order(exams$mSv)])
p_d1 <- ggplot(exams, aes(exam, mSv)) +
geom_col(fill = bpad_pal[1]) +
geom_hline(yintercept = bkgd_yr, linetype = "dashed", colour = bpad_pal[2]) +
annotate("text", x = 2, y = bkgd_yr*1.6, size = 2.9, colour = bpad_pal[2],
label = "annual natural background") +
scale_y_log10() + coord_flip() +
labs(x = NULL, y = "Effective dose (mSv, log scale)",
title = "Typical effective doses")
p_d2 <- ggplot(exams, aes(exam, days)) +
geom_col(fill = bpad_pal[3]) +
geom_text(aes(label = ifelse(days >= 1, sprintf("%.0f d", days),
sprintf("%.1f d", days))),
hjust = -0.15, size = 2.7) +
scale_y_log10(limits = c(0.05, 6000)) + coord_flip() +
labs(x = NULL, y = "Equivalent days of natural background (log scale)",
title = "The same doses, in a unit patients understand")
bpad_grid(p_d1, p_d2, ncol = 2)Figure 8: Putting radiation dose in perspective. Left: typical effective doses for common examinations on a logarithmic scale, spanning three orders of magnitude from a chest radiograph to a multiphase abdominal CT; the dashed line marks the annual natural background dose. Right: the same doses expressed as the equivalent number of days of natural background radiation, which is the framing patients find most interpretable. A chest radiograph is a few days; an abdominal CT is a few years.
knitr::kable(data.frame(
examination = as.character(exams$exam),
effective_mSv = exams$mSv,
background_days = round(exams$days, 1),
chest_xray_equivalents = round(exams$mSv/0.02, 0)),
col.names = c("Examination","Effective dose (mSv)","Background days",
"Chest radiographs"),
caption = "Representative effective doses. Values vary widely with scanner, protocol, and patient size; they are indicative rather than prescriptive.")| Examination | Effective dose (mSv) | Background days | Chest radiographs |
|---|---|---|---|
| Chest radiograph (PA) | 0.020 | 2.4 | 1 |
| Extremity radiograph | 0.001 | 0.1 | 0 |
| Mammogram (2 view) | 0.400 | 48.7 | 20 |
| Abdominal radiograph | 0.700 | 85.2 | 35 |
| Head CT | 2.000 | 243.3 | 100 |
| Chest CT | 7.000 | 851.7 | 350 |
| Chest CT (low dose) | 1.500 | 182.5 | 75 |
| Abdomen/pelvis CT | 10.000 | 1216.7 | 500 |
| CT angiography | 12.000 | 1460.0 | 600 |
| Multiphase abdomen CT | 25.000 | 3041.7 | 1250 |
cat(sprintf("An abdomen/pelvis CT delivers about %.0f times the dose of a chest radiograph\n",
10/0.02))## An abdomen/pelvis CT delivers about 500 times the dose of a chest radiograph
## and roughly 3.3 years of natural background radiation.
Radiation effects divide cleanly, and confusing them causes both unnecessary alarm and unnecessary complacency.
| Stochastic | Deterministic (tissue reactions) | |
|---|---|---|
| Mechanism | DNA mutation leading to cancer | cell killing |
| Threshold | assumed none (linear no-threshold model) | yes, typically \(\gtrsim\) 0.5–2 Gy |
| Dose affects | probability of the effect | severity of the effect |
| Diagnostic relevance | the relevant concern for CT and radiography | relevant only in prolonged fluoroscopy and interventional work |
| Example | radiation-induced malignancy | skin erythema, epilation, cataract |
Diagnostic imaging operates in the stochastic regime. The linear no-threshold model used for radiation protection assumes risk is proportional to dose with no safe threshold; it is a conservative regulatory convention rather than a firmly established biological fact at low doses, and it should be described to students as such.
Radiation protection in medicine rests on two principles, in this order.
Justification. The examination must do more good than harm. The most effective dose reduction available is the scan that is not performed, or is replaced by ultrasound (Chapter 3) or MRI (Chapter 4).
Optimization (ALARA). Given that the study is justified, the dose must be as low as reasonably achievable consistent with answering the clinical question. In practice this means automatic exposure control, iterative or deep-learning reconstruction that tolerates lower photon counts, size-adapted paediatric protocols, restricting scan length, and avoiding unnecessary multiphase acquisitions.
Three chapters, three safety physics. Chapter 3 limits cavitation and heating through MI and TI, both governed by acoustic pressure and intensity. Chapter 4 limits RF heating through SAR \(\propto B_0^2\alpha^2\), and manages a permanent projectile hazard. Chapter 5 limits ionization through absorbed and effective dose. The first two have no cumulative lifetime risk model; this one does, which is why justification — asking whether to image at all — appears here and nowhere else.
Section 5.5 summary.
Checkpoint 5.5. A CT abdomen reports CTDI\(_{\text{vol}} = 12\) mGy over a 40 cm scan length, with \(k = 0.015\) mSv (mGy cm)\(^{-1}\). Compute the DLP and the estimated effective dose, express it in chest-radiograph equivalents, and state one protocol change that would halve it and what it would cost diagnostically.
A radiograph collapses a three-dimensional patient onto a two-dimensional detector: every structure along a ray contributes to the same pixel, so a lesion behind the heart, a rib, and the spine all superimpose. Two consequences follow. Depth information is lost entirely, and low-contrast structures are buried under the integrated attenuation of everything else along the ray.
CT solves both by measuring the same line integrals from many angles and reconstructing the attenuation map. The gain is dramatic: radiographic contrast resolution is limited to differences of a few percent, whereas CT routinely distinguishes tissues differing by less than 0.5% in \(\mu\) — about 5 HU.
Analogy. A radiograph is one shadow. CT uses many shadows from many directions to reconstruct the object that cast them. The mathematics of “reconstructing an object from its shadows” is the Radon transform, and its inverse is what a CT scanner computes several hundred times per second.
For a ray through a non-uniform object, the Beer–Lambert law generalizes to
\[\begin{equation} I = I_0\exp\!\left(-\int_{\text{ray}} \mu(x,y)\,ds\right), \tag{17} \end{equation}\]
and taking the negative logarithm linearizes the measurement:
\[\begin{equation} P(\theta,s) = -\ln\!\frac{I}{I_0} = \int_{\text{ray}} \mu(x,y)\,ds . \tag{18} \end{equation}\]
Equation (18) is the central step of CT: it converts a multiplicative, exponential measurement into an additive line integral, which is what makes linear reconstruction possible at all. It is also exactly the step that beam hardening (Section 5.3.4) invalidates, because with a polychromatic beam \(-\ln(I/I_0)\) is not a line integral of any single \(\mu\).
Formally, for a 2D function \(f(x,y)\),
\[\begin{equation} R[f](\theta,s) = \iint f(x,y)\,\delta(x\cos\theta + y\sin\theta - s)\,dx\,dy , \tag{19} \end{equation}\]
where the delta function selects the line \(x\cos\theta + y\sin\theta = s\). CT projection data are the Radon transform of \(\mu(x,y)\), and reconstruction is the inverse Radon transform.
Stacking projections over all angles gives the sinogram \(P(\theta, s)\): detector position horizontally, angle vertically. A point object at radius \(r\) and angle \(\phi\) projects to \(s = r\cos(\theta - \phi)\) — a sinusoid, which is where the name comes from.
## Analytic Radon transform of a Gaussian blob at (x0, y0) with width b:
## a Gaussian in s centred at x0 cos(theta) + y0 sin(theta), amplitude b*sqrt(pi).
blob_sino <- function(s, th, x0, y0, b, amp = 1)
amp*b*sqrt(pi)*exp(-((s - (x0*cos(th) + y0*sin(th)))^2)/b^2)
blob_obj <- function(x, y, x0, y0, b, amp = 1)
amp*exp(-((x - x0)^2 + (y - y0)^2)/b^2)
xs <- seq(-3, 3, length.out = 220)
ths <- seq(0, pi, length.out = 220)
Gs <- expand.grid(s = xs, th = ths)
Go <- expand.grid(x = xs, y = xs)
## --- single blob
Gs$v1 <- blob_sino(Gs$s, Gs$th, 1.4, 0.0, 0.35)
Go$v1 <- blob_obj (Go$x, Go$y, 1.4, 0.0, 0.35)
## --- three blobs at different radii
Gs$v2 <- blob_sino(Gs$s, Gs$th, 1.8, 0.0, 0.30, 1.0) +
blob_sino(Gs$s, Gs$th, -0.6, 1.0, 0.30, 0.8) +
blob_sino(Gs$s, Gs$th, 0.0, 0.0, 0.30, 0.6)
Go$v2 <- blob_obj (Go$x, Go$y, 1.8, 0.0, 0.30, 1.0) +
blob_obj (Go$x, Go$y, -0.6, 1.0, 0.30, 0.8) +
blob_obj (Go$x, Go$y, 0.0, 0.0, 0.30, 0.6)
mk_obj <- function(col, ttl)
ggplot(Go, aes(x, y, fill = .data[[col]])) + geom_raster() +
coord_equal() + scale_fill_gradient(low = "black", high = "white", guide = "none") +
labs(title = ttl, x = "x", y = "y")
mk_sin <- function(col, ttl)
ggplot(Gs, aes(s, th, fill = .data[[col]])) + geom_raster() +
scale_fill_gradient(low = "black", high = "white", guide = "none") +
scale_y_continuous(breaks = c(0, pi/2, pi), labels = c("0", "pi/2", "pi")) +
labs(title = ttl, x = "detector position s", y = "projection angle theta")
bpad_grid(mk_obj("v1", "Object: one off-centre blob"),
mk_sin("v1", "Sinogram: a single sinusoid"),
mk_obj("v2", "Object: three blobs"),
mk_sin("v2", "Sinogram: three superposed sinusoids"),
ncol = 2)Figure 9: Objects and their sinograms, computed rather than asserted. Top row: a single off-centre Gaussian blob, and its analytic sinogram, which is a pure sinusoid whose amplitude equals the object’s distance from the rotation centre and whose phase encodes its angular position. Bottom row: three blobs at different radii, and the superposition of three sinusoids that results. Reading a sinogram is a learnable skill: the amplitude of each trace tells you how far from the centre a structure lies, and a structure at the exact centre produces a straight vertical line.
## Sinogram reading rules, verifiable from the figure:
## - amplitude of a trace = distance of the structure from the rotation centre
## - phase of a trace = angular position of the structure
## - a straight vertical line = a structure exactly at the centre
## - the blob at radius 1.8 traces s = 1.8 cos(theta), peaking at s = 1.8
Why \(\theta\) spans only \([0,\pi)\). A parallel-beam projection at \(\theta\) and at \(\theta + \pi\) measure the same set of lines traversed in opposite directions and therefore carry identical information: \(P(\theta + \pi, s) = P(\theta, -s)\). Displaying several multiples of \(2\pi\) of angle merely repeats the data and obscures this symmetry. Fan-beam and helical geometries break the symmetry and do require a full rotation plus the fan angle.
The central slice theorem (or Fourier slice theorem) explains why projections suffice. Let \(P_\theta(s)\) be the projection at angle \(\theta\). Then
\[\begin{equation} \mathcal{F}_1\{P_\theta\}(k) = \mathcal{F}_2\{f\}(k\cos\theta,\; k\sin\theta) . \tag{20} \end{equation}\]
In words: the one-dimensional Fourier transform of a projection equals a radial line, through the origin at angle \(\theta\), of the two-dimensional Fourier transform of the object. Each projection therefore hands us one line of the object’s 2D spectrum, and enough lines determine the whole spectrum, whose inverse transform is the image.
The proof is three lines. Writing the projection as a line integral and transforming,
\[\mathcal{F}_1\{P_\theta\}(k) = \int P_\theta(s)e^{-i2\pi ks}ds = \iint f(x,y)\,e^{-i2\pi k(x\cos\theta + y\sin\theta)}\,dx\,dy,\]
which is precisely \(\mathcal{F}_2\{f\}\) evaluated at \((u,v) = (k\cos\theta, k\sin\theta)\). The delta function of Eq. (19) has done all the work by converting the one-dimensional exponential into the two-dimensional one.
Np <- 128
axp <- (seq_len(Np) - (Np + 1)/2)/(Np/2)
gxp <- matrix(axp, Np, Np); gyp <- matrix(axp, Np, Np, byrow = TRUE)
ph <- matrix(0, Np, Np)
ph[gxp^2/0.8^2 + gyp^2/0.95^2 <= 1] <- 0.3
ph[(gxp - 0.25)^2 + (gyp + 0.15)^2 <= 0.18^2] <- 0.7
ph[(gxp + 0.3)^2 + (gyp - 0.25)^2 <= 0.12^2] <- 0.55
fftshift2 <- function(M) {
n <- nrow(M); m <- ncol(M); h <- n %/% 2; k <- m %/% 2
M[c((h+1):n, 1:h), c((k+1):m, 1:k)]
}
F2 <- fftshift2(fft(ph))
F2m <- log(1 + Mod(F2))
dfF <- expand.grid(u = seq_len(Np), v = seq_len(Np)); dfF$m <- as.vector(F2m)
cenp <- Np/2 + 1
ang <- c(20, 60, 110)*pi/180
lines_df <- do.call(rbind, lapply(ang, function(a)
data.frame(u = cenp + c(-1,1)*40*cos(a), v = cenp + c(-1,1)*40*sin(a),
g = sprintf("%.0f deg", a*180/pi))))
p_f2 <- ggplot(dfF, aes(u, v, fill = m)) + geom_raster() +
scale_fill_viridis_c(guide = "none") + coord_equal() +
geom_path(data = lines_df, aes(u, v, group = g), inherit.aes = FALSE,
colour = "white", linewidth = 0.6) +
labs(title = "2D Fourier magnitude with radial slices", x = NULL, y = NULL)
## Numerical check at one angle
theta0 <- 60*pi/180
nDetC <- 181; SmaxC <- sqrt(2)
Xv <- as.vector(gxp); Yv <- as.vector(gyp); Vv <- as.vector(ph)
sC <- Xv*cos(theta0) + Yv*sin(theta0)
bC <- pmin(pmax(round((sC + SmaxC)/(2*SmaxC)*(nDetC - 1)) + 1, 1), nDetC)
aggC <- tapply(Vv, bC, sum)
projC <- numeric(nDetC); projC[as.integer(names(aggC))] <- aggC
F1 <- Mod(fftshift2(matrix(fft(projC), nDetC, 1))[, 1])
kk <- (seq_len(nDetC) - (nDetC %/% 2 + 1))
## corresponding radial slice of F2
rad <- seq(-30, 30)
uu <- cenp + rad*cos(theta0); vv <- cenp + rad*sin(theta0)
slice <- Mod(F2[cbind(pmin(pmax(round(uu),1),Np), pmin(pmax(round(vv),1),Np))])
cmp <- rbind(
data.frame(k = kk[abs(kk) <= 30]*(nDetC/Np)^0, y = F1[abs(kk) <= 30]/max(F1),
src = "1D FT of the projection"),
data.frame(k = rad, y = slice/max(slice), src = "radial slice of the 2D FT"))
p_cs <- ggplot(cmp, aes(k, y, colour = src)) +
geom_line(linewidth = 0.85) +
scale_y_log10() + scale_colour_manual(values = bpad_pal[c(1,2)]) +
labs(x = "spatial frequency index", y = "normalized magnitude (log)",
title = "Theorem verified at 60 degrees")
bpad_grid(p_f2, p_cs, ncol = 2)Figure 10: The central slice theorem verified numerically. Left: the 2D Fourier magnitude of a phantom, with three radial lines marked at the angles of three projections. Right: for one angle, the 1D Fourier transform of the projection is overlaid on the corresponding radial slice extracted from the 2D transform; the two agree to within interpolation error, confirming Eq. (central-slice). The consequence is that projections fill the frequency plane radially, so samples are dense near the origin and sparse at high frequency, which is exactly the imbalance the ramp filter corrects.
Sections 5.7–5.8 summary.
Checkpoint 5.6. In a sinogram, one trace is a flat horizontal line at \(s = 0\) for all angles, and another is a sinusoid of amplitude 8 cm. What can you say about the position of each structure, and which would be more affected by a small patient rotation?
Back projection smears each projection back along the direction it came from:
\[\begin{equation} f_b(x,y) = \int_0^{\pi} P(\theta,\; x\cos\theta + y\sin\theta)\,d\theta . \tag{21} \end{equation}\]
This is not the inverse Radon transform. The central slice theorem shows why. Projections fill the frequency plane along radial lines, so the density of samples at radius \(k\) falls as \(1/k\) — the plane is oversampled near the origin and undersampled at high frequency. Back projecting without correction therefore weights low frequencies too heavily, producing a reconstruction that is the true image convolved with a \(1/r\) blur.
Writing the inverse 2D Fourier transform in polar coordinates introduces the Jacobian \(|k|\,dk\,d\theta\), and that Jacobian is the correction:
\[\begin{equation} f(x,y) = \int_0^{\pi}\left[\int_{-\infty}^{\infty} \tilde P_\theta(k)\,|k|\,e^{i2\pi ks}\,dk\right]_{s = x\cos\theta + y\sin\theta} d\theta . \tag{22} \end{equation}\]
The factor \(|k|\) in Eq. (22) is the ramp filter. Its origin is purely geometric — it compensates for the radial sampling density — but its consequence is practical: it amplifies high frequencies, and therefore noise, without bound. Real systems use band-limited kernels with a cutoff \(k_{\max}\):
\[\begin{equation} H(k) = |k|\,W(k), \qquad W_{\text{Ram-Lak}} = \mathbb{1}(|k|\le k_{\max}), \quad W_{\text{Shepp-Logan}} = \mathrm{sinc}\!\left(\frac{k}{2k_{\max}}\right), \quad W_{\text{Hann}} = \tfrac12\!\left(1 + \cos\frac{\pi k}{k_{\max}}\right). \tag{23} \end{equation}\]
Kernel choice is the single most visible reconstruction decision a radiographer makes: sharp kernels for lung and bone, smooth kernels for brain and liver, and the trade is always resolution against noise.
The detector sets the highest measurable spatial frequency; the number of views must be sufficient that the azimuthal sampling at that frequency is no coarser than the radial sampling. This gives the standard requirement
\[\begin{equation} N_{\text{views}} \;\gtrsim\; \frac{\pi}{2}\,N_{\text{det}} . \tag{24} \end{equation}\]
Violating Eq. (24) produces view-aliasing streaks radiating from high-contrast structures — the CT analogue of the Nyquist violations met in Chapters 3 and 4.
## --- phantom
N <- 96
ax <- (seq_len(N) - (N + 1)/2)/(N/2)
gx <- matrix(ax, N, N); gy <- matrix(ax, N, N, byrow = TRUE)
disc <- function(M, x0, y0, r, val) { M[(gx - x0)^2 + (gy - y0)^2 <= r^2] <- val; M }
phantom <- matrix(0, N, N)
phantom <- disc(phantom, 0.00, 0.00, 0.80, 0.30)
phantom <- disc(phantom, 0.00, 0.00, 0.74, 0.15)
phantom <- disc(phantom, -0.30, 0.10, 0.14, 0.50)
phantom <- disc(phantom, 0.32, -0.12, 0.16, 0.60)
phantom <- disc(phantom, 0.00, -0.40, 0.07, 0.80)
nDet <- 141; Smax <- sqrt(2)
det <- seq(-Smax, Smax, length.out = nDet)
Xv <- as.vector(gx); Yv <- as.vector(gy); Vv <- as.vector(phantom)
forward <- function(thetas) {
sino <- matrix(0, nDet, length(thetas))
for (k in seq_along(thetas)) {
s <- Xv*cos(thetas[k]) + Yv*sin(thetas[k])
b <- pmin(pmax(round((s + Smax)/(2*Smax)*(nDet - 1)) + 1L, 1L), nDet)
agg <- tapply(Vv, b, sum)
col <- numeric(nDet); col[as.integer(names(agg))] <- agg
sino[, k] <- col
}
sino
}
ramp_filter <- function(sino) {
f <- 0:(nDet - 1); f[f > nDet/2] <- f[f > nDet/2] - nDet
rp <- abs(f)/max(abs(f))
apply(sino, 2, function(p) Re(fft(fft(p)*rp, inverse = TRUE))/nDet)
}
backproject <- function(proj, thetas) {
recon <- matrix(0, N, N)
for (k in seq_along(thetas)) {
ct <- cos(thetas[k]); st <- sin(thetas[k])
for (j in seq_len(N)) {
s <- ax*ct + ax[j]*st
idx <- (s + Smax)/(2*Smax)*(nDet - 1) + 1
i0 <- pmin(pmax(floor(idx), 1L), nDet - 1L); fr <- idx - i0
recon[, j] <- recon[, j] + proj[i0, k]*(1 - fr) + proj[i0 + 1L, k]*fr
}
}
recon*pi/length(thetas)
}
thetas <- seq(0, pi, length.out = 180)
sino <- forward(thetas)
rec_bp <- backproject(sino, thetas)
rec_fbp <- backproject(ramp_filter(sino), thetas)
th16 <- seq(0, pi, length.out = 16)
rec_16 <- backproject(ramp_filter(forward(th16)), th16)
view_counts <- c(4, 8, 16, 24, 32, 48, 64, 96, 128, 180, 256)
corr <- sapply(view_counts, function(nv) {
th <- seq(0, pi, length.out = nv)
cor(as.vector(phantom), as.vector(backproject(ramp_filter(forward(th)), th)))
})
need <- pi/2*nDet
op <- par(mfrow = c(2, 3), mar = c(2.2, 2.2, 2.4, 0.8))
image(ax, ax, phantom, col = gpal, axes = FALSE, asp = 1, main = "Phantom (truth)")
image(seq_along(thetas), det, t(sino), col = gpal, axes = FALSE,
main = "Sinogram", xlab = "", ylab = "")
image(ax, ax, rec_bp, col = gpal, axes = FALSE, asp = 1,
main = "Unfiltered back projection")
image(ax, ax, rec_fbp, col = gpal, axes = FALSE, asp = 1,
main = "Filtered back projection (180 views)")
image(ax, ax, rec_16, col = gpal, axes = FALSE, asp = 1,
main = "FBP, only 16 views: streaks")
par(mar = c(4.2, 4.2, 2.4, 0.8))
plot(view_counts, corr, type = "b", pch = 16, log = "x", col = bpad_pal[1],
xlab = "number of views", ylab = "correlation with truth",
main = "Angular sampling")
abline(v = need, lty = 2, col = bpad_pal[2])
text(need, min(corr) + 0.05, sprintf("(pi/2)N_det = %.0f", need),
pos = 2, cex = 0.75, col = bpad_pal[2])Figure 11: Reconstruction from projections. Top row: numerical phantom, its sinogram, and unfiltered back projection, which recovers the object’s location but is severely blurred by the one-over-r point spread implied by radial sampling. Bottom row: ramp-filtered back projection, which restores sharp edges; reconstruction from only 16 views, showing the characteristic view-aliasing streaks radiating from high-contrast structures; and the correlation with truth as a function of view count, which saturates near the value predicted by Eq. (view-sampling), marked by the dashed line.
par(op)
cat(sprintf("correlation with truth: unfiltered BP = %.3f, filtered BP = %.3f\n",
cor(as.vector(phantom), as.vector(rec_bp)),
cor(as.vector(phantom), as.vector(rec_fbp))))## correlation with truth: unfiltered BP = 0.771, filtered BP = 0.965
cat(sprintf("with %d detector channels, Eq. (view-sampling) requires about %.0f views\n",
nDet, need))## with 141 detector channels, Eq. (view-sampling) requires about 221 views
knitr::kable(data.frame(views = view_counts, correlation = round(corr, 4)),
caption = "Reconstruction fidelity versus the number of projection angles. Improvement is steep below the sampling requirement and negligible above it: acquiring more views past that point costs dose and time for nothing.")| views | correlation |
|---|---|
| 4 | 0.2624 |
| 8 | 0.4375 |
| 16 | 0.6678 |
| 24 | 0.7927 |
| 32 | 0.8546 |
| 48 | 0.9112 |
| 64 | 0.9357 |
| 96 | 0.9450 |
| 128 | 0.9623 |
| 180 | 0.9654 |
| 256 | 0.9752 |
Iterative reconstruction writes the problem as a linear system \(Af = p\), where \(A\) encodes ray paths, and solves it by repeatedly comparing predicted with measured projections:
\[\begin{equation} f^{(n+1)} = f^{(n)} + \lambda A^{\mathsf T}\!\left(p - Af^{(n)}\right), \tag{25} \end{equation}\]
which is exactly the least-squares gradient step of Chapter 1. Its advantage is that the statistical model (Poisson noise) and physical constraints (non-negativity, beam hardening, scatter) can be built in, which is why iterative methods enabled the substantial dose reductions of the last fifteen years.
Deep-learning reconstruction learns either the mapping from projections to images or a denoising step applied to a conventional reconstruction. It is now clinically deployed, and Chapter 8’s cautions apply in full: a network can produce an image that is plausible, smooth, and confidently wrong, because a learned prior can hallucinate structure not supported by the measured data. Validation must therefore test whether lesions that were present survive and whether lesions that were absent fail to appear.
| Method | Idea | Strength | Limitation |
|---|---|---|---|
| Back projection | smear projections back | intuitive | \(1/r\) blur; not the inverse transform |
| Filtered back projection | apply \(|k|\) first | fast, analytic, predictable | noise amplification; needs full sampling |
| Iterative | solve \(Af = p\) with a noise model | lower dose, handles artifacts | compute cost; can look “plastic” |
| Deep learning | learn the mapping or the denoiser | strongest dose reduction | hallucination risk; needs task-based validation |
Section 5.9 summary.
Checkpoint 5.7. A scanner with 900 detector channels acquires 400 views per rotation. Is it adequately sampled? If not, what artifact would you expect, where would it be worst, and what are the two ways to fix it?
Reconstruction yields \(\mu(x,y)\) in cm\(^{-1}\), whose numerical value depends on the beam spectrum and therefore on the scanner. Normalizing to water removes most of that dependence:
\[\begin{equation} \mathrm{HU} = 1000\,\frac{\mu_{\text{tissue}} - \mu_{\text{water}}}{\mu_{\text{water}}} . \tag{26} \end{equation}\]
By construction Eq. (26) puts water at 0 HU and air at \(-1000\) HU. The scale is deliberately fine: one HU is a 0.1% change in \(\mu\), so clinically meaningful soft-tissue differences of 10–40 HU correspond to attenuation differences of 1–4%.
mu_water70 <- 0.1937*1.00
tissue_hu <- data.frame(
tissue = c("Air","Lung","Fat","Water","White matter","Grey matter",
"Muscle","Blood","Liver","Trabecular bone","Cortical bone"),
murho = c(0.1866, 0.1937, 0.1815, 0.1937, 0.1940, 0.1943,
0.1946, 0.1944, 0.1948, 0.2160, 0.2861),
rho = c(0.00120, 0.26, 0.92, 1.00, 1.035, 1.038,
1.05, 1.06, 1.06, 1.35, 1.92))
tissue_hu$mu <- tissue_hu$murho*tissue_hu$rho
tissue_hu$HU <- 1000*(tissue_hu$mu - mu_water70)/mu_water70
tissue_hu$tissue <- factor(tissue_hu$tissue, levels = tissue_hu$tissue)
tissue_hu$idx <- seq_len(nrow(tissue_hu))
ggplot(tissue_hu, aes(idx, HU)) +
annotate("rect", xmin = 0.4, xmax = 11.6, ymin = 20, ymax = 70,
fill = bpad_pal[3], alpha = 0.25) +
geom_col(fill = bpad_pal[1], width = 0.7) +
geom_hline(yintercept = 0, colour = "grey40") +
geom_text(aes(label = round(HU)),
vjust = ifelse(tissue_hu$HU >= 0, -0.4, 1.3), size = 2.8) +
annotate("text", x = 6, y = 330, size = 2.9, colour = bpad_pal[3],
label = "the entire soft-tissue diagnostic range: 20 to 70 HU") +
scale_x_continuous(breaks = tissue_hu$idx, labels = tissue_hu$tissue) +
coord_cartesian(ylim = c(-1150, 2100)) +
theme(axis.text.x = element_text(angle = 40, hjust = 1)) +
labs(x = NULL, y = "Hounsfield units",
title = "HU computed from mu/rho and density at 70 keV")Figure 12: Hounsfield units derived from first principles. Bars show HU computed from tabulated mass attenuation coefficients and physical densities at a 70 keV effective energy, using Eq. (hounsfield); the values reproduce the standard clinical table without being taken from it. Note the enormous dynamic range, from -1000 for air to nearly 1900 for cortical bone, against which the entire diagnostic soft-tissue range from about 20 to 70 HU is a narrow band, shaded here. That compression is precisely why windowing is not optional.
knitr::kable(data.frame(
tissue = as.character(tissue_hu$tissue),
murho = tissue_hu$murho, density = tissue_hu$rho,
mu = round(tissue_hu$mu, 4), HU_computed = round(tissue_hu$HU),
HU_clinical = c(-1000, -700, -100, 0, 30, 40, 45, 55, 60, 300, 1800)),
col.names = c("Tissue","mu/rho (cm2/g)","Density (g/cm3)","mu (1/cm)",
"HU computed","HU clinical reference"),
caption = "Computed HU values against standard clinical references. The agreement across four orders of magnitude in mu confirms that HU is a physically derived quantity, not an empirical convention.")| Tissue | mu/rho (cm2/g) | Density (g/cm3) | mu (1/cm) | HU computed | HU clinical reference |
|---|---|---|---|---|---|
| Air | 0.1866 | 0.0012 | 0.0002 | -999 | -1000 |
| Lung | 0.1937 | 0.2600 | 0.0504 | -740 | -700 |
| Fat | 0.1815 | 0.9200 | 0.1670 | -138 | -100 |
| Water | 0.1937 | 1.0000 | 0.1937 | 0 | 0 |
| White matter | 0.1940 | 1.0350 | 0.2008 | 37 | 30 |
| Grey matter | 0.1943 | 1.0380 | 0.2017 | 41 | 40 |
| Muscle | 0.1946 | 1.0500 | 0.2043 | 55 | 45 |
| Blood | 0.1944 | 1.0600 | 0.2061 | 64 | 55 |
| Liver | 0.1948 | 1.0600 | 0.2065 | 66 | 60 |
| Trabecular bone | 0.2160 | 1.3500 | 0.2916 | 505 | 300 |
| Cortical bone | 0.2861 | 1.9200 | 0.5493 | 1836 | 1800 |
## A 1 HU change corresponds to a 0.10% change in mu.
cat(sprintf("Grey and white matter differ by only %.0f HU, i.e. %.1f%% in attenuation:\n",
tissue_hu$HU[6] - tissue_hu$HU[5],
(tissue_hu$HU[6] - tissue_hu$HU[5])/10))## Grey and white matter differ by only 5 HU, i.e. 0.5% in attenuation:
## detecting that difference is what forced CT to reach 0.5% contrast resolution.
A display offers roughly 256 grey levels and the eye distinguishes fewer; CT data span about 4000 HU. Windowing maps a selected HU range onto the full grey scale:
\[\begin{equation} g(\mathrm{HU}) = \mathrm{clip}\!\left(\frac{\mathrm{HU} - (L - W/2)}{W},\,0,\,1\right), \tag{27} \end{equation}\]
with level \(L\) the centre and width \(W\) the range. Narrow windows maximize contrast within a narrow HU band and saturate everything outside it.
Nw <- 220
axw <- seq(-1, 1, length.out = Nw)
gxw <- matrix(axw, Nw, Nw); gyw <- matrix(axw, Nw, Nw, byrow = TRUE)
r2 <- gxw^2 + gyw^2
hu <- matrix(-1000, Nw, Nw)
hu[r2 <= 0.85^2] <- 40
hu[r2 <= 0.85^2 & r2 >= 0.80^2] <- 700
hu[(gxw + 0.30)^2 + (gyw - 0.20)^2 <= 0.12^2] <- -100
hu[(gxw - 0.35)^2 + (gyw + 0.10)^2 <= 0.10^2] <- 60
hu[(gxw - 0.05)^2 + (gyw + 0.35)^2 <= 0.06^2] <- 20
apply_window <- function(hu, L, W) pmin(pmax((hu - (L - W/2))/W, 0), 1)
wins <- list(c(40, 400, "Soft tissue\nL=40 W=400"),
c(-600, 1500, "Lung\nL=-600 W=1500"),
c(400, 1800, "Bone\nL=400 W=1800"),
c(35, 40, "Narrow stroke\nL=35 W=40"))
op <- par(mfrow = c(1, 4), mar = c(0.6, 0.6, 3.0, 0.6))
for (w in wins) {
image(axw, axw, apply_window(hu, as.numeric(w[1]), as.numeric(w[2])),
col = gpal, axes = FALSE, asp = 1, main = w[3], cex.main = 0.95)
box(col = "grey70")
}Figure 13: One data set, four displays. A synthetic HU slice containing air, fat, soft tissue, blood, a bony rim, and a subtle 20 HU lesion is shown in four standard windows. The bony rim dominates the bone window and the lesion is invisible there; the lesion appears only in the soft-tissue and narrow stroke windows. Nothing about the underlying data changes between panels: windowing is a display transformation, and choosing the wrong one can hide the finding entirely.
par(op)
lesion <- 20; bg <- 40
knitr::kable(data.frame(
window = c("Soft tissue","Lung","Bone","Narrow stroke"),
L = c(40, -600, 400, 35), W = c(400, 1500, 1800, 40),
displayed_contrast = sprintf("%.1f%%",
100*abs(apply_window(lesion, c(40,-600,400,35), c(400,1500,1800,40)) -
apply_window(bg, c(40,-600,400,35), c(400,1500,1800,40))))),
col.names = c("Window","Level (HU)","Width (HU)","Displayed lesion contrast"),
caption = "The same 20 HU lesion against 40 HU background, displayed four ways. The narrow stroke window renders it at 50 percent grey-scale contrast; the bone window renders it at about 1 percent, which is invisible.")| Window | Level (HU) | Width (HU) | Displayed lesion contrast |
|---|---|---|---|
| Soft tissue | 40 | 400 | 5.0% |
| Lung | -600 | 1500 | 1.3% |
| Bone | 400 | 1800 | 1.1% |
| Narrow stroke | 35 | 40 | 50.0% |
From the equation to the reading room. Equation (27) shows that displayed contrast is \(\Delta\mathrm{HU}/W\). A 20 HU lesion is rendered at 5% of the grey scale in a 400 HU soft-tissue window and at 1.1% in an 1800 HU bone window — below visual threshold. Narrow “stroke windows” (\(W \approx 30\)–40 HU) exist precisely to make early ischaemic grey-white differentiation visible, and a radiologist who reviews a head CT only in a standard window can miss an early infarct that the data fully contain.
Section 5.10 summary.
Checkpoint 5.8. A 15 HU lesion sits in 45 HU liver. Compute the displayed contrast in a \(W = 400\) window and in a \(W = 150\) window. What is the cost of the narrower window, and which structures would you lose?
Every CT artifact is a modelling assumption failing, and Table 10 organizes them accordingly.
| Artifact | Assumption violated | Mechanism | Correction |
|---|---|---|---|
| Beam hardening (cupping) | monoenergetic beam | effective \(\mu\) falls with depth | polynomial correction, filtration, dual-energy |
| Beam hardening (streaks) | monoenergetic beam | dense structures harden the beam unequally | iterative correction, higher kVp |
| Metal artifact | finite dynamic range, monoenergetic beam | photon starvation plus extreme hardening | MAR algorithms, dual-energy, higher kVp |
| View aliasing (streaks) | adequate angular sampling | too few views, Eq. (24) | more views, interpolation |
| Ring artifact | uniform detector response | a mis-calibrated detector channel | detector calibration, ring filters |
| Motion | object static during acquisition | inconsistent projections | faster rotation, gating, motion correction |
| Partial volume | uniform voxel content | voxel spans two tissues | thinner slices, isotropic acquisition |
| Scatter | photons travel straight | scattered photons mis-assigned | collimation, scatter correction |
| Photon starvation | adequate counts on every ray | too few photons through thick paths | tube-current modulation, iterative recon |
## Uniform water cylinder
Nc <- 96
axc <- (seq_len(Nc) - (Nc + 1)/2)/(Nc/2)
gxc <- matrix(axc, Nc, Nc); gyc <- matrix(axc, Nc, Nc, byrow = TRUE)
cyl <- matrix(0, Nc, Nc); cyl[gxc^2 + gyc^2 <= 0.8^2] <- 1
cyl2 <- cyl
cyl2[(gxc - 0.45)^2 + (gyc)^2 <= 0.12^2] <- 4 # dense insert
cyl2[(gxc + 0.45)^2 + (gyc)^2 <= 0.12^2] <- 4
nD <- 141; Sm <- sqrt(2)
Xc <- as.vector(gxc); Yc <- as.vector(gyc)
thc <- seq(0, pi, length.out = 160)
## Polychromatic forward projection: attenuate the SPECTRUM, then take -ln(I/I0),
## which is exactly what a scanner does and exactly what breaks the line-integral model.
Es <- seq(30, 120, by = 3)
wts <- pmax(120 - Es, 0)*exp(-mu_al(Es)*0.25); wts <- wts/sum(wts)
muE <- mu_water_E(Es) # water, cm^-1 per keV
poly_forward <- function(V) {
sino <- matrix(0, nD, length(thc))
for (k in seq_along(thc)) {
s <- Xc*cos(thc[k]) + Yc*sin(thc[k])
b <- pmin(pmax(round((s + Sm)/(2*Sm)*(nD - 1)) + 1L, 1L), nD)
agg <- tapply(V, b, sum)
L <- numeric(nD); L[as.integer(names(agg))] <- agg # path length in "units"
Itot <- sapply(L, function(len) sum(wts*Es*exp(-muE*len*0.35)))
sino[, k] <- -log(Itot/sum(wts*Es))
}
sino
}
bp2 <- function(proj) {
f <- 0:(nD - 1); f[f > nD/2] <- f[f > nD/2] - nD
rp <- abs(f)/max(abs(f))
pf <- apply(proj, 2, function(p) Re(fft(fft(p)*rp, inverse = TRUE))/nD)
R <- matrix(0, Nc, Nc)
for (k in seq_along(thc)) {
ct <- cos(thc[k]); st <- sin(thc[k])
for (j in seq_len(Nc)) {
s <- axc*ct + axc[j]*st
idx <- (s + Sm)/(2*Sm)*(nD - 1) + 1
i0 <- pmin(pmax(floor(idx), 1L), nD - 1L); fr <- idx - i0
R[, j] <- R[, j] + pf[i0, k]*(1 - fr) + pf[i0 + 1L, k]*fr
}
}
R*pi/length(thc)
}
rec_cup <- bp2(poly_forward(as.vector(cyl)))
rec_strk <- bp2(poly_forward(as.vector(cyl2)))
op <- par(mfrow = c(3, 2), mar = c(2.2, 2.4, 2.4, 0.8))
image(axc, axc, cyl, col = gpal, axes = FALSE, asp = 1, main = "Uniform cylinder (truth)")
image(axc, axc, cyl2, col = gpal, axes = FALSE, asp = 1, main = "Cylinder + dense inserts")
image(axc, axc, rec_cup, col = gpal, axes = FALSE, asp = 1,
main = "Polychromatic recon: cupping")
image(axc, axc, rec_strk, col = gpal, axes = FALSE, asp = 1,
main = "Polychromatic recon: streaks")
par(mar = c(4.2, 4.2, 2.4, 0.8))
mid <- Nc/2
prof <- rec_cup[, mid]; prof <- prof/max(prof)
plot(axc, prof, type = "l", lwd = 2, col = bpad_pal[2],
xlab = "position", ylab = "relative reconstructed value",
main = "Profile: the cup")
abline(h = max(prof), lty = 3, col = "grey50")
prof2 <- rec_strk[, mid]; prof2 <- prof2/max(prof2)
plot(axc, prof2, type = "l", lwd = 2, col = bpad_pal[1],
xlab = "position", ylab = "relative reconstructed value",
main = "Profile through the inserts")Figure 14: Two artifacts produced by breaking a modelling assumption. Left column: a uniform water cylinder, and the cupping artifact that appears when the projections are generated with a polychromatic beam but reconstructed as though the beam were monoenergetic; the profile through the centre shows the characteristic depression, which can reach tens of Hounsfield units and mimic pathology. Right column: the same phantom with two dense inserts, and the dark streaks between them produced by the same mechanism operating unequally along different rays.
par(op)
inner <- prof[abs(axc) < 0.1]; outer <- prof[abs(axc) > 0.55 & abs(axc) < 0.75]
cat(sprintf("cupping depth: centre %.3f vs periphery %.3f, a %.1f%% depression\n",
mean(inner), mean(outer), 100*(1 - mean(inner)/mean(outer))))## cupping depth: centre 0.751 vs periphery 0.825, a 8.9% depression
## A uniform object has been reconstructed as non-uniform purely because the
## polychromatic measurement is not the line integral the algorithm assumed.
Because photoelectric and Compton attenuation depend differently on energy (Section 5.2), attenuation measured at two energies contains enough information to separate them, and hence to identify materials rather than merely measure their attenuation. Writing \(\mu\) as a two-component decomposition,
\[\begin{equation} \mu(E) \;\approx\; a_{\text{PE}}\,f_{\text{PE}}(E) \;+\; a_{\text{C}}\,f_{\text{C}}(E), \tag{28} \end{equation}\]
two measurements at distinct effective energies give two equations in the two unknowns \(a_{\text{PE}}, a_{\text{C}}\) — a \(2\times2\) linear system solved exactly as in Chapter 1.
One inverse problem, four chapters. This is the fourth appearance of the same structure. Chapter 2 unmixed oxy- and deoxyhaemoglobin from two optical wavelengths. Chapter 3 generalized it to multispectral photoacoustics. Chapter 4 solved a two-unknown system for tensor components. Here two X-ray energies separate photoelectric from Compton attenuation. In every case the conditioning of the system depends on choosing measurement channels where the basis functions differ most — which for dual-energy CT means separating the two effective spectra as far as possible, typically 80 and 140 kVp.
Clinically this yields virtual monoenergetic images (reconstructed as though the beam were monoenergetic, which suppresses beam hardening and metal artifact), iodine maps and virtual non-contrast images (separating iodine from soft tissue), material-specific imaging such as uric acid for gout and calcium for bone-marrow oedema, and improved low-energy contrast without the dose penalty of actually acquiring at low kVp.
Four quantities are coupled, and no protocol improves all of them at once.
If an object of diameter \(L\) is imaged with voxel width \(w\), it spans \(n = L/w\) voxels and its circular cross-section contains about
\[\begin{equation} \frac{A_{\text{object}}}{A_{\text{voxel}}} = \frac{\pi L^{2}/4}{w^{2}} = \frac{\pi n^{2}}{4} \tag{29} \end{equation}\]
voxels. Halving \(w\) quadruples the in-plane voxel count and, in 3D, multiplies it by eight — while each voxel now intercepts a correspondingly smaller share of the photons.
Detecting a lesion requires its contrast to exceed the noise by some factor \(k\) (the Rose criterion uses \(k \approx 3\)–5). Carrying the Poisson statistics through the reconstruction gives the scaling
\[\begin{equation} \Phi \;\gtrsim\; \frac{k^{2}}{w^{4}\,(\delta\mu)^{2}} , \tag{30} \end{equation}\]
where \(\Phi\) is the required photon fluence, \(w\) the voxel width, and \(\delta\mu\) the attenuation difference between lesion and background.
Where the fourth power comes from, and why it is brutal. The signal from a lesion of size \(w\) is proportional to the extra attenuation it produces, \(\delta\mu \cdot w\), measured over an area \(w^2\) that collects \(\Phi w^2\) photons. The noise on that measurement is \(\sqrt{\Phi w^{2}}\). Requiring signal-to-noise \(\ge k\) gives \(\delta\mu\, w \sqrt{\Phi w^{2}} \gtrsim k\), i.e. \(\Phi \gtrsim k^2/(w^4 (\delta\mu)^2)\). The consequence is unforgiving: halving the voxel width costs sixteen times the fluence, and halving the contrast costs four times. Dimensionally the expression is consistent, since \([k^2/(w^4\delta\mu^2)] = 1/L^{2}\), the units of fluence.
wg <- 10^seq(log10(0.02), log10(0.5), length.out = 200) # cm
dmu_levels <- c(0.002, 0.005, 0.02) # 1/cm
kk <- 4
det_df <- do.call(rbind, lapply(dmu_levels, function(d)
data.frame(w = wg, phi = kk^2/(wg^4*d^2),
dmu = factor(sprintf("%.3f /cm", d),
levels = sprintf("%.3f /cm", dmu_levels)))))
p_dt <- ggplot(det_df, aes(w*10, phi, colour = dmu)) +
geom_line(linewidth = 0.95) +
scale_x_log10() + scale_y_log10() +
scale_colour_manual(values = bpad_pal[c(3,1,2)]) +
labs(x = "Voxel width (mm, log scale)",
y = expression("Required fluence (arbitrary, log scale)"),
title = "Fluence needed to detect a lesion")
grid_d <- expand.grid(w = 10^seq(log10(0.03), log10(0.4), length.out = 120),
d = 10^seq(log10(0.001), log10(0.05), length.out = 120))
grid_d$logphi <- log10(kk^2/(grid_d$w^4*grid_d$d^2))
p_iso <- ggplot(grid_d, aes(w*10, d, fill = logphi)) +
geom_raster() +
geom_contour(data = grid_d, aes(w*10, d, z = logphi),
breaks = seq(2, 12, by = 2), colour = "white",
linewidth = 0.45, inherit.aes = FALSE) +
scale_x_log10() + scale_y_log10() + scale_fill_viridis_c() +
theme(legend.position = "right", legend.title = element_text(size = 8)) +
labs(x = "Voxel width (mm, log)", y = expression(delta*mu*" (1/cm, log)"),
fill = "log10 fluence", title = "Iso-dose contours")
bpad_grid(p_dt, p_iso, ncol = 2)Figure 15: The cost of resolution and contrast. Left: required photon fluence as a function of voxel width at three lesion contrasts, from Eq. (detectability), on logarithmic axes; the slope of minus four means that every factor of two in resolution costs a factor of sixteen in dose. Right: the same relation as an iso-fluence map over the resolution-contrast plane, with contours of constant required dose. Any diagonal movement toward the lower-left corner, which is where small low-contrast lesions live, is extraordinarily expensive.
knitr::kable(data.frame(
change = c("halve voxel width", "quarter voxel width",
"halve lesion contrast", "halve both"),
fluence_factor = c(2^4, 4^4, 2^2, 2^4*2^2)),
col.names = c("Change","Fluence (and dose) multiplier"),
caption = "Consequences of Eq. (detectability). Detecting a lesion half the size and half the contrast requires 64 times the dose, which is why sub-millimetre low-contrast CT is fundamentally hard rather than merely unrefined.")| Change | Fluence (and dose) multiplier |
|---|---|
| halve voxel width | 16 |
| quarter voxel width | 256 |
| halve lesion contrast | 4 |
| halve both | 64 |
From the clinical question to the protocol.
| Question | Choice | Physical reason | Cost |
|---|---|---|---|
| Is there a fracture? | radiograph, moderate kVp | bone-soft-tissue contrast is intrinsically large | 2D superposition |
| Is there a microcalcification? | mammography, 26–30 kVp, Mo/Rh | maximize \(Z^3/E^3\) contrast between adipose and glandular | high skin dose, thin part only |
| Where exactly is the lesion? | CT | depth localization via many projections | ~500 chest radiographs of dose |
| Is a vessel patent? | CT angiography with iodine | K-edge at 33.2 keV sits in the beam | contrast reaction, renal load |
| Is the beam hardening or is it pathology? | dual-energy CT | photoelectric/Compton separation | dose or hardware |
| Is there haemorrhage in the brain? | non-contrast CT, narrow window | 60–80 HU blood against 40 HU brain | needs \(W \approx 40\) to be visible |
| Is the implant loose? | CT with MAR, higher kVp | reduce photon starvation and hardening | residual artifact |
| Could this be answered without radiation? | ultrasound (Ch. 3) or MRI (Ch. 4) | justification precedes optimization | time, availability, contrast |
Goal. Compute transmitted intensities, half-value layers, and subject contrast directly from tabulated coefficients, and show numerically that contrast and penetration trade against each other.
E_lab <- exp(seq(log(20), log(150), length.out = 300))
mu_w <- exp(approx(log(atten$E_keV), log(atten$water), log(E_lab), rule = 2)$y)*1.00
mu_b <- exp(approx(log(atten$E_keV), log(atten$bone), log(E_lab), rule = 2)$y)*1.92
x_body <- 20; x_bone <- 2
trans <- exp(-mu_w*x_body)
contr <- 1 - exp(-(mu_b - mu_w)*x_bone)
p_l1 <- ggplot() +
geom_line(aes(E_lab, trans/max(trans), colour = "penetration (20 cm tissue)"),
linewidth = 0.95) +
geom_line(aes(E_lab, contr, colour = "bone contrast (2 cm)"), linewidth = 0.95) +
scale_x_log10() + scale_colour_manual(values = bpad_pal[c(1,2)]) +
labs(x = "Photon energy (keV, log scale)", y = "Normalized value",
title = "Neither curve alone has an optimum")
## Dose-normalized figure of merit.
## CNR ~ contrast * sqrt(N_detected); N_detected ~ N0 * transmission.
## Entrance dose ~ N0 * E * (mu_en/rho), which tracks mu/rho at these energies.
## Hence FOM = contrast * sqrt(transmission / (E * mu_w)).
fom_of <- function(xb) {
tr <- exp(-mu_w*xb)
cc <- 1 - exp(-(mu_b - mu_w)*x_bone)
f <- cc*sqrt(tr/(E_lab*mu_w))
f/max(f)
}
fom_df <- do.call(rbind, lapply(c(5, 20, 30), function(xb)
data.frame(E = E_lab, f = fom_of(xb),
body = factor(sprintf("%d cm thick", xb),
levels = sprintf("%d cm thick", c(5, 20, 30))))))
opt <- do.call(rbind, lapply(c(5, 20, 30), function(xb)
data.frame(E = E_lab[which.max(fom_of(xb))], f = 1,
body = factor(sprintf("%d cm thick", xb),
levels = sprintf("%d cm thick", c(5, 20, 30))))))
p_l2 <- ggplot(fom_df, aes(E, f, colour = body)) +
geom_line(linewidth = 0.95) +
geom_point(data = opt, size = 2.4) +
geom_text(data = opt, aes(label = sprintf("%.0f keV", E)),
vjust = -0.9, size = 2.8, show.legend = FALSE) +
scale_x_log10() + scale_colour_manual(values = bpad_pal[c(3,1,2)]) +
coord_cartesian(ylim = c(0, 1.15)) +
labs(x = "Photon energy (keV, log scale)",
y = "CNR per unit entrance dose (normalized)",
title = "The optimum, and how it moves with body size")
bpad_grid(p_l1, p_l2, ncol = 2)Figure 16: Lab 1. Left: penetration through 20 cm of soft tissue rises monotonically with photon energy while bone contrast falls monotonically, so neither alone identifies an optimum. Right: a dose-normalized figure of merit, contrast-to-noise per unit entrance dose, which does peak; for a 20 cm body it peaks near 50 keV, and the optimum shifts downward for thinner parts and upward for thicker ones. This is the physical origin of the diagnostic energy window and of the rule that thicker patients need higher kVp.
Esel <- c(20, 30, 50, 70, 100, 150)
knitr::kable(data.frame(
E_keV = Esel,
mu_tissue = round(approx(E_lab, mu_w, Esel, rule = 2)$y, 4),
HVL_cm = round(log(2)/approx(E_lab, mu_w, Esel, rule = 2)$y, 2),
transmitted_20cm = signif(exp(-approx(E_lab, mu_w, Esel, rule = 2)$y*20), 3),
bone_contrast_2cm = round(approx(E_lab, contr, Esel, rule = 2)$y, 3)),
col.names = c("E (keV)","mu tissue (1/cm)","HVL (cm)","Transmitted through 20 cm",
"Bone contrast over 2 cm"),
caption = "At 20 keV essentially nothing penetrates an abdomen; at 150 keV plenty penetrates but bone is nearly invisible. The usable window is the compromise between them.")| E (keV) | mu tissue (1/cm) | HVL (cm) | Transmitted through 20 cm | Bone contrast over 2 cm |
|---|---|---|---|---|
| 20 | 0.8096 | 0.86 | 1.00e-07 | 1.000 |
| 30 | 0.3759 | 1.84 | 5.44e-04 | 0.987 |
| 50 | 0.2269 | 3.05 | 1.07e-02 | 0.691 |
| 70 | 0.1937 | 3.58 | 2.08e-02 | 0.461 |
| 100 | 0.1707 | 4.06 | 3.29e-02 | 0.310 |
| 150 | 0.1505 | 4.61 | 4.93e-02 | 0.235 |
for (xb in c(5, 20, 30))
cat(sprintf("body thickness %2d cm -> optimum near %3.0f keV\n",
xb, E_lab[which.max(fom_of(xb))]))## body thickness 5 cm -> optimum near 40 keV
## body thickness 20 cm -> optimum near 50 keV
## body thickness 30 cm -> optimum near 59 keV
cat(sprintf("at 20 keV only %.2e of the beam survives 20 cm, so an unusable dose would\n",
exp(-approx(E_lab, mu_w, 20, rule = 2)$y*20)))## at 20 keV only 9.29e-08 of the beam survives 20 cm, so an unusable dose would
## be required to form an image: low-kVp imaging is restricted to thin body parts.
Try it yourself. Change x_body from 20 cm to 5 cm (an extremity) and re-run. The
optimum shifts downward in energy, which is exactly why extremity and mammographic
protocols use far lower kVp than abdominal ones. Then set x_bone to 0.2 cm to model a
microcalcification and observe how much harder the task becomes.
The full demonstration appears in Section 5.9. Use that code as a starting point for the following investigations.
Try it yourself.
rpois on the transmitted
intensity, then take the negative logarithm) and observe how the ramp filter amplifies it.Goal. Quantify how window settings determine whether a lesion is visible, and connect this to the detectability condition of Section 5.13.
Wg <- seq(20, 2000, by = 5)
lesion_hu <- 20; bg_hu <- 45
disp_contrast <- abs(lesion_hu - bg_hu)/Wg
p_w1 <- ggplot(data.frame(W = Wg, c = 100*disp_contrast), aes(W, c)) +
annotate("rect", xmin = 20, xmax = 2000, ymin = 0.5, ymax = 5,
fill = bpad_pal[2], alpha = 0.15) +
geom_line(linewidth = 0.95, colour = bpad_pal[1]) +
geom_vline(xintercept = c(40, 400, 1800), linetype = "dotted", colour = "grey45") +
annotate("text", x = 45, y = 48, hjust = 0, size = 2.7, colour = "grey30",
label = "stroke") +
annotate("text", x = 420, y = 20, hjust = 0, size = 2.7, colour = "grey30",
label = "soft tissue") +
annotate("text", x = 1200, y = 12, hjust = 0, size = 2.7, colour = "grey30",
label = "bone") +
scale_x_log10() + scale_y_log10() +
labs(x = "Window width (HU, log scale)",
y = "Displayed contrast (% of grey scale, log)",
title = "Narrow windows amplify contrast")
noise_hu <- c(3, 10, 25)
cnr <- do.call(rbind, lapply(noise_hu, function(s)
data.frame(W = Wg, cnr = abs(lesion_hu - bg_hu)/s,
sigma = factor(sprintf("%d HU noise", s),
levels = sprintf("%d HU noise", noise_hu)))))
p_w2 <- ggplot(cnr, aes(W, cnr, colour = sigma)) +
geom_line(linewidth = 0.95) +
geom_hline(yintercept = 4, linetype = "dashed", colour = "grey40") +
annotate("text", x = 1500, y = 4.8, size = 2.8, colour = "grey30",
label = "Rose criterion k = 4") +
scale_x_log10() + scale_y_log10() +
scale_colour_manual(values = bpad_pal[c(3,1,2)]) +
labs(x = "Window width (HU, log scale)", y = "Contrast-to-noise ratio (log)",
title = "Windowing cannot create information")
bpad_grid(p_w1, p_w2, ncol = 2)Figure 17: Lab 3. Left: displayed grey-scale contrast of a 20 HU lesion against a 45 HU background as a function of window width; the relation is exactly the inverse proportionality of Eq. (windowing), and the shaded band marks the region below roughly 5 percent where the difference becomes hard to perceive. Right: the effect of adding realistic image noise; with 10 HU of noise, narrowing the window amplifies the noise just as fast as the signal, so contrast-to-noise ratio is flat and no window setting can rescue a lesion that the dose did not resolve.
## A 20 HU lesion in 45 HU background is displayed at:
for (W in c(40, 150, 400, 1800))
cat(sprintf(" W = %4d HU -> %.1f%% of the grey scale%s\n", W,
100*abs(lesion_hu - bg_hu)/W,
ifelse(100*abs(lesion_hu-bg_hu)/W < 5, " (below visual threshold)", "")))## W = 40 HU -> 62.5% of the grey scale
## W = 150 HU -> 16.7% of the grey scale
## W = 400 HU -> 6.2% of the grey scale
## W = 1800 HU -> 1.4% of the grey scale (below visual threshold)
##
## But the contrast-to-noise ratio does not depend on W at all: windowing rescales
## signal and noise identically. It reveals information already present; it cannot
## supply information the photons did not carry.
X-ray imaging and CT rest on a single measurement — the attenuation of a photon beam along a line — and on a single inversion — recovering a two- or three-dimensional map from its line integrals. Everything else in this chapter elaborates one of those two ideas or manages their cost.
Three themes recur.
Contrast is interaction physics, not image processing. The photoelectric cross-section scales as \(Z^3/E^3\) and the Compton cross-section barely depends on \(Z\) at all. That single asymmetry explains why kVp controls contrast, why bone-to-soft-tissue contrast collapses from 3.5 to 1.1 between 30 and 100 keV, why mammography operates at 28 kVp and chest radiography at 120, why iodine and barium were chosen from the whole periodic table, and why dual-energy CT can identify materials rather than merely measure them. A student who holds Eq. (5) in mind can derive the protocol; one who does not must memorize it.
Reconstruction is a Fourier statement with a sampling requirement. The central slice theorem says each projection is one radial line of the object’s spectrum. Radial sampling oversamples low frequencies, which is why plain back projection blurs and why the ramp filter — the polar Jacobian — is the exact correction. The same geometry sets the angular sampling requirement \(N_{\text{views}} \gtrsim (\pi/2)N_{\text{det}}\), whose violation produces streaks in exactly the way that Nyquist violations produce aliasing in Chapters 3 and 4. And the whole edifice rests on \(-\ln(I/I_0)\) being a line integral, which polychromatic beam hardening quietly breaks — the origin of cupping and metal artifact.
Dose is a governing variable, not an afterthought. This is the only modality in this book whose photons ionize. Detectability scales as \(\Phi \gtrsim k^2/(w^4\delta\mu^2)\), so halving the voxel width costs sixteen times the dose and halving the contrast costs four. Noise falls only as \(1/\sqrt{\text{mAs}}\). An abdominal CT delivers roughly three years of natural background radiation. These are not reasons to avoid CT — the diagnostic benefit usually dominates by a wide margin — but they are reasons why justification, the question of whether to image at all, belongs in a physics chapter rather than a policy appendix.
Several directions are moving quickly. Photon-counting detectors eliminate electronic noise, remove the energy weighting that suppresses the most contrast-rich photons, and deliver spectral information intrinsically. Deep-learning reconstruction has enabled substantial dose reduction and brings the hallucination and generalization risks that Chapter 8 treats at length. Dual-energy and spectral CT convert an attenuation image into a material map. And automatic exposure control with size-specific dosimetry is steadily closing the gap between the dose a protocol delivers and the dose the clinical question actually requires.
| Term | Meaning |
|---|---|
| ALARA | As Low As Reasonably Achievable; the optimization principle for dose |
| Beam hardening | Rise in mean beam energy with depth as soft photons are removed |
| Bremsstrahlung | Continuous X-ray spectrum from decelerating electrons |
| Bucky factor | Exposure increase required when a grid is used |
| Characteristic X-ray | Discrete photon from an inner-shell atomic transition |
| Compton scattering | Photon scatter off a quasi-free electron; nearly \(Z\)-independent |
| CTDI\(_{\text{vol}}\) | Volume CT dose index; phantom-based scanner output, not patient dose |
| Central slice theorem | 1D transform of a projection equals a radial slice of the 2D spectrum |
| Cupping | Beam-hardening artifact in which a uniform object appears darker centrally |
| DLP | Dose-length product; CTDI\(_{\text{vol}}\) times scan length |
| DQE | Detective quantum efficiency; fraction of incident-photon information retained |
| Dual-energy CT | Material decomposition from two effective spectra |
| Duane–Hunt limit | Maximum photon energy \(E_{\max} = eV\) |
| Effective dose | Tissue-weighted whole-body-equivalent dose, in sieverts |
| Filtered back projection | Ramp-filtered analytic reconstruction |
| Half-value layer | Thickness halving beam intensity, \(\ln 2/\mu\) |
| Hounsfield unit | Attenuation normalized to water; 1 HU is 0.1% in \(\mu\) |
| K-edge | Discontinuous attenuation rise at a shell binding energy |
| Line integral | \(-\ln(I/I_0)\); the quantity CT reconstructs from |
| Linear attenuation coefficient | \(\mu\) (cm\(^{-1}\)); density-dependent |
| Mass attenuation coefficient | \(\mu/\rho\) (cm\(^2\) g\(^{-1}\)); density-independent |
| Photoelectric absorption | Complete photon absorption; \(\propto Z^3/E^3\) |
| Radon transform | Map from an object to its complete set of line integrals |
| Ramp filter | \(|k|\) weighting; the polar Jacobian of the inverse transform |
| Rose criterion | Detectability threshold, signal-to-noise of roughly 4–5 |
| Scatter-to-primary ratio | Ratio of scattered to primary photons at the detector |
| Sinogram | Projection data indexed by angle and detector position |
| Stochastic effect | Probabilistic radiation effect (cancer); no assumed threshold |
| View aliasing | Streak artifact from too few projection angles |
| Window level / width | Display centre and range mapped onto the grey scale |
A tube operates at 100 kVp and 200 mAs. (a) What is the maximum photon energy? (b) Compute the electrical energy delivered. (c) Using \(Y \approx 9\times10^{-10}ZV\) with \(Z = 74\), find the X-ray energy produced and the heat deposited. (d) Comment on the anode design consequence.
Soft tissue has \(\mu/\rho = 0.2059\) cm\(^2\)g\(^{-1}\) at 60 keV, density 1.06 g cm\(^{-3}\). (a) Find \(\mu\). (b) Find the half-value layer. (c) Find the fraction transmitted through 25 cm. (d) How many HVLs is 25 cm?
Using the chapter’s table, compute the bone-to-soft-tissue attenuation ratio at 30 keV and at 100 keV. (a) By what factor does the excess attenuation (\(\mu_b/\mu_w - 1\)) fall? (b) Explain the clinical consequence for rib visualization on a chest radiograph.
Assume \(\tau/\rho \propto Z^3/E^3\). (a) By what factor does photoelectric absorption change when photon energy doubles? (b) By what factor does it differ between iodine (\(Z = 53\)) and soft tissue (\(Z_{\text{eff}} = 7.4\)) at the same energy? (c) Comment on why milligram quantities of iodine are visible.
A solution is 5% iodine and 95% water by mass. At 50 keV, \(\mu/\rho\) is 12.2 cm\(^2\)g\(^{-1}\) for iodine and 0.227 cm\(^2\)g\(^{-1}\) for water. (a) Compute the mixture’s \(\mu/\rho\). (b) With density 1.05 g cm\(^{-3}\), compute \(\mu\). (c) By what factor is it more attenuating than blood (\(\mu \approx 0.24\) cm\(^{-1}\))?
A radiograph has unacceptable noise. (a) By what factor must mAs increase to halve the noise? (b) What happens to dose? (c) If instead kVp is raised, what happens to noise, contrast, and dose, and why is this not equivalent?
An abdominal radiograph has SPR \(= 4\). (a) What fraction of the primary contrast survives? (b) A grid reduces SPR to 0.5 with a Bucky factor of 3.5. What fraction survives now? (c) By what factor did contrast improve, and at what dose cost? (d) Would you use the grid for a wrist (SPR \(= 0.3\))?
A CT abdomen reports CTDI\(_{\text{vol}} = 14\) mGy over a 35 cm scan. Using \(k = 0.015\) mSv (mGy cm)\(^{-1}\): (a) compute DLP; (b) estimate effective dose; (c) express it in chest-radiograph equivalents (0.02 mSv) and in days of natural background (3 mSv/yr).
A ray passes through 5 cm of fat (\(\mu = 0.18\) cm\(^{-1}\)), 12 cm of soft tissue (\(\mu = 0.21\)), and 3 cm of bone (\(\mu = 0.48\)). (a) Compute \(-\ln(I/I_0)\). (b) What single uniform \(\mu\) over the 20 cm path would give the same reading? (c) Explain why one projection cannot separate the three tissues.
In a sinogram, a trace follows \(s = 6\cos(\theta - 40^\circ)\) cm. (a) Where is the structure? (b) What trace would a structure at the isocentre produce? (c) What would two structures at the same radius but opposite sides look like?
A scanner has 800 detector channels. (a) How many views are required by Eq. (24)? (b) If only 300 are acquired, what artifact appears and where is it worst? (c) Name two remedies and their costs.
At 70 keV, water has \(\mu = 0.1937\) cm\(^{-1}\). (a) Compute HU for a tissue with \(\mu = 0.2064\). (b) Compute \(\mu\) for a lesion at \(-60\) HU and identify it. (c) What percentage change in \(\mu\) does 1 HU represent?
A 25 HU lesion sits in 50 HU liver, with 8 HU image noise. (a) Compute displayed contrast at \(W = 400\) and \(W = 100\). (b) Compute the contrast-to-noise ratio in each. (c) Explain why (a) changes and (b) does not, and state what would actually improve detectability.
Using \(\Phi \gtrsim k^2/(w^4\delta\mu^2)\): (a) by what factor must fluence rise to detect a lesion of half the diameter at the same contrast? (b) Half the contrast at the same size? (c) Both? (d) Comment on the feasibility of routine sub-millimetre low-contrast CT.
A CT of a uniform water phantom shows the centre 40 HU lower than the periphery. (a) Name the artifact and its mechanism. (b) Why does raising kVp reduce it? (c) Why does dual-energy CT reduce it more fundamentally? (d) Why might a radiologist mistake it for pathology?
A 25-year-old presents with right lower quadrant pain. CT would deliver about 10 mSv. (a) Estimate the background-equivalent exposure. (b) Name two alternative modalities and the physical property each measures. (c) State the principle that governs choosing among them, and explain why it is not the same as ALARA.
| Problem | Computed answer |
|---|---|
| 1(a) | 100 keV |
| 1(b) | 20000 J |
| 1(c) | X-rays 133.2 J, heat 19867 J (yield 0.67%) |
| 2(a) | 0.2183 /cm |
| 2(b) | 3.18 cm |
| 2(c) | 0.0043 |
| 2(d) | 7.9 HVLs |
| 3(a) | excess falls from 2.54 to 0.09, a factor of 29 |
| 4(a) | falls by 2^3 = 8 |
| 4(b) | (53/7.4)^3 = 367 |
| 5(a) | 0.826 cm^2/g |
| 5(b) | 0.867 /cm |
| 5(c) | 3.6 x |
| 6(a) | mAs x4; dose x4 |
| 7(a) | 20% |
| 7(b) | 67% |
| 7(c) | contrast x3.3 at 3.5x the dose |
| 8(a) | 490 mGy cm |
| 8(b) | 7.35 mSv |
| 8(c) | 368 chest radiographs; 894 days of background |
| 9(a) | 4.86 |
| 9(b) | 0.2430 /cm |
| 10(a) | radius 6 cm, azimuth 40 degrees |
| 11(a) | 1257 views |
| 12(a) | 66 HU |
| 12(b) | mu = 0.1821 /cm -> fat |
| 12(c) | 0.1% |
| 13(a) | 6.25% at W=400; 25.0% at W=100 |
| 13(b) | CNR = 3.12 in both |
| 14(a) | 16x |
| 14(b) | 4x |
| 14(c) | 64x |
| 16(a) | 3.3 years of background |
From the table, bone/water is 3.54 at 30 keV and 1.09 at 100 keV. The ratio falls by a factor of 3.3, but the diagnostically relevant quantity is the excess attenuation \(\mu_b/\mu_w - 1\), which falls from 2.54 to 0.09 — a factor of 29. That is why a 120 kVp chest radiograph deliberately renders ribs faintly: the technique is chosen to suppress bone so that lung parenchyma behind it becomes visible. The same physics, run backwards, is why a rib series uses lower kVp.