---
jupytext:
  text_representation:
    extension: .md
    format_name: myst
    format_version: 0.13
    jupytext_version: 1.10.3
kernelspec:
  display_name: Python 3 (ipykernel)
  language: python
  name: python3
mystnb:
  execution_timeout: 180
---

(ch:test-tomography)=
# Practical: Electron Tomography

## Introduction

Unlike single-particle analysis, **electron tomography** (cryo-ET) reconstructs a 3D volume from images of a single specimen tilted to many angles by physically rotating the sample stage. We work through the problem in 2D: the unknown object is a 2D image, and we acquire 1D projections at various tilt angles, then reconstruct the 2D image from those projections.

**Axis convention:** The electron beam travels along $z$. The tilt axis is $y$ (pointing through the imaging plane). Tilting the sample about $y$ rotates the $xz$-plane, so each tilt angle $\alpha$ produces a 1D projection along $x$ of the rotated image.

A typical cryo-ET experiment collects tilt images from $-60°$ to $+60°$ in steps of $1°$–$3°$, limited by grid-bar occlusion and increased sample thickness at high tilt.

```{admonition} Interactive elements
:class: tip
Click **Live Code** in the top toolbar to activate the kernel, then expand each **Show code** toggle and click ▶ to run the cell.
```

### Helper functions

```{code-cell} ipython3
:tags: [hide-input]
import io, warnings
import numpy as np
import matplotlib
matplotlib.use('agg')
import matplotlib.pyplot as plt
from scipy.ndimage import rotate as nd_rotate, gaussian_filter
from ipywidgets import IntSlider, FloatSlider, VBox, HBox, Layout, Output
from IPython.display import display, Image, clear_output
warnings.filterwarnings('ignore')

_sl  = Layout(width="420px")
_sty = {"description_width": "160px"}

def fig2img(fig, dpi=100):
    buf = io.BytesIO()
    fig.savefig(buf, format='png', dpi=dpi, bbox_inches='tight')
    buf.seek(0); plt.close(fig)
    return Image(data=buf.read())

# ── test image ────────────────────────────────────────────────────────────────

def make_2D_test_image(size=100):
    """
    2D phantom representing a cross-section of a biological specimen:
      - ring  = microtubule cross-section
      - lines = membrane sheets
      - dot   = protein complex
    Pure-numpy implementation (no skimage needed).
    """
    y, x = np.mgrid[0:size, 0:size].astype(float)
    img  = np.zeros((size, size))

    # Ring (microtubule): hollow circle, top-left quadrant
    cx1, cy1 = size*0.33, size*0.33
    R1 = np.sqrt((x - cx1)**2 + (y - cy1)**2)
    img += 0.9 * ((R1 > size*0.12) & (R1 < size*0.18))

    # Horizontal membrane sheet
    cy2 = size*0.68
    img += 0.85 * (np.abs(y - cy2) < size*0.018)

    # Vertical membrane sheet
    cx3 = size*0.67
    img += 0.85 * (np.abs(x - cx3) < size*0.018) * ((y > size*0.4) & (y < size*0.9))

    # Dot (protein complex), top-right quadrant
    cx4, cy4 = size*0.72, size*0.35
    R4 = np.sqrt((x - cx4)**2 + (y - cy4)**2)
    img += 0.75 * (R4 < size*0.06)

    return np.clip(img, 0, 1)[::-1]

# ── Radon transform ───────────────────────────────────────────────────────────

def radon(im, angles):
    """
    Compute 1D projections of im at each tilt angle.
    Rotation is applied to the image; projection is the column sum.
    Returns sinogram of shape (len(angles), im.shape[1]).
    """
    sino = np.zeros((len(angles), im.shape[1]))
    for k, ang in enumerate(angles):
        rotated    = nd_rotate(im, ang, reshape=False, cval=0)
        sino[k, :] = rotated.sum(axis=0)
    return sino

# ── backprojection ────────────────────────────────────────────────────────────

def repeat_line(line, N):
    """Repeat 1D line along a new axis N times → shape (N, len(line))."""
    return np.tile(line[np.newaxis, :], (N, 1))

def backproject(sinogram, angles):
    """
    Simple (unfiltered) backprojection.
    Each 1D projection is repeated into a 2D slab, rotated, and summed.
    """
    n = sinogram.shape[1]
    recon = np.zeros((n, n))
    for proj, ang in zip(sinogram, angles):
        slab    = repeat_line(proj / n, n)
        recon  += nd_rotate(slab, -ang, reshape=False, cval=0)
    return recon / max(len(angles), 1)

# ── filtered backprojection ───────────────────────────────────────────────────

def ramp_filter(sinogram):
    """Apply ramp filter |ν| to each projection in the sinogram."""
    n       = sinogram.shape[1]
    freqs   = np.fft.rfftfreq(n)
    ramp    = np.abs(freqs)
    sino_f  = np.zeros_like(sinogram)
    for i, proj in enumerate(sinogram):
        fft_proj    = np.fft.rfft(proj)
        sino_f[i]   = np.fft.irfft(fft_proj * ramp, n)
    return sino_f

def filtered_backproject(sinogram, angles):
    """Filtered backprojection with ramp filter."""
    sino_f = ramp_filter(sinogram)
    return backproject(sino_f, angles)

# ── Fourier power spectrum ────────────────────────────────────────────────────

def power_spectrum(im):
    F   = np.fft.fft2(im)
    F   = np.fft.fftshift(F)
    return np.log(np.abs(F)**2 + 1)

print("Helper functions loaded.")
```

---

## The 2D test image

We work with a 2D phantom that mimics a biological cross-section: a hollow **ring** (microtubule), two **line segments** (membrane sheets), and a small **dot** (protein complex).

```{code-cell} ipython3
:tags: [hide-input]
im = make_2D_test_image(150)

fig, axes = plt.subplots(1, 2, figsize=(10, 4.5))
axes[0].imshow(im, cmap='gray', origin='lower')
axes[0].set_xlabel('x [pixels]'); axes[0].set_ylabel('z [pixels]')
axes[0].set_title('2D test image (150×150)')
axes[1].imshow(power_spectrum(im), cmap='inferno', origin='lower')
axes[1].set_title('Power spectrum (log scale)')
for ax in axes: ax.axis('on')
plt.tight_layout()
display(fig2img(fig))
```

---

## Tilting the sample and acquiring projections

To acquire a tilt series, the specimen is physically rotated about the $y$-axis. At each tilt angle $\alpha$, the transmitted electrons form a 1D projection (the column sum of the rotated image). The interactive below shows the specimen at a chosen tilt angle together with the resulting 1D projection profile.

```{code-cell} ipython3
:tags: [hide-input]
tilt_sl = IntSlider(value=0, min=-70, max=70, step=1,
                     description="Tilt angle (°)", style=_sty, layout=_sl)
out_tilt = Output()

def _update_tilt(_=None):
    ang     = tilt_sl.value
    rotated = nd_rotate(im, ang, reshape=False, cval=0)
    proj    = rotated.sum(axis=0) / im.shape[0]

    fig, axes = plt.subplots(1, 3, figsize=(13, 4))
    axes[0].imshow(im,      cmap='gray', origin='lower'); axes[0].axis('off')
    axes[0].set_title("Original specimen")
    axes[1].imshow(rotated, cmap='gray', origin='lower'); axes[1].axis('off')
    axes[1].set_title(f"Specimen at α={ang}°")
    axes[2].plot(proj, np.arange(len(proj)), color='steelblue', linewidth=1.4)
    axes[2].set_xlim(left=0)
    axes[2].set_xlabel("Projection intensity"); axes[2].set_ylabel("x [pixel]")
    axes[2].set_title("1D projection (column sum)")
    axes[2].grid(alpha=0.3)
    plt.tight_layout()
    with out_tilt:
        clear_output(wait=True); display(fig2img(fig))

tilt_sl.observe(_update_tilt, names='value')
display(VBox([tilt_sl, out_tilt]))
_update_tilt()
```

---

## The sinogram

Collecting projections over the full tilt range and stacking them row by row produces the **sinogram** — so named because a point object traces a sinusoidal curve across it as the tilt angle varies.

```{code-cell} ipython3
:tags: [hide-input]
angles_full = np.arange(-60, 61, 3)
sino_full   = radon(im, angles_full)

fig, axes = plt.subplots(1, 2, figsize=(12, 4))
axes[0].imshow(sino_full, cmap='gray', origin='lower', aspect='auto',
               extent=[0, im.shape[1], angles_full[0], angles_full[-1]])
axes[0].set_xlabel('x [pixels]'); axes[0].set_ylabel('Tilt angle (°)')
axes[0].set_title(f'Sinogram — {len(angles_full)} tilts, ±60°')

# Show sinusoidal trajectory of a single bright pixel
cy_test, cx_test = 100, 50
axes[1].plot(angles_full,
             [nd_rotate(im, a, reshape=False)[cy_test, cx_test] for a in angles_full],
             color='steelblue', lw=1.4)
axes[1].set_xlabel('Tilt angle (°)'); axes[1].set_ylabel('Projection value')
axes[1].set_title(f'Projection value of pixel ({cx_test},{cy_test}) vs tilt')
axes[1].grid(alpha=0.3)
plt.tight_layout()
display(fig2img(fig))
```

---

## Task 1 – Backprojection

The simplest reconstruction method is **backprojection**: for each projection at angle $\alpha$, repeat the 1D line into a 2D slab, rotate it back by $-\alpha$, and sum over all angles.

$$
A_\text{BP}(x,z) = \sum_\alpha p_\alpha\!\left(x\cos\alpha + z\sin\alpha\right)
$$

This is correct in principle but overweights low-spatial-frequency components (every point in the image receives contributions from every projection), producing a blurry, star-shaped reconstruction artifact.

```{admonition} Task 1
:class: note
1. Take a 1D projection `projection`.
2. Repeat it $N$ times along a new axis using `repeat_line` to get a 2D slab of shape $(N, N)$.
3. Rotate the slab back by $-\alpha$ using `nd_rotate`.
4. Sum over all angles. Divide by the number of angles.
```

```{code-cell} ipython3
:tags: [hide-input]
max_tilt_sl = IntSlider(value=60, min=10, max=70, step=5,
                         description="Max tilt (°)",   style=_sty, layout=_sl)
tinc_sl     = IntSlider(value=3,  min=1,  max=15, step=1,
                         description="Tilt increment (°)", style=_sty, layout=_sl)
out_bp = Output()

def _update_bp(_=None):
    angs   = np.arange(-max_tilt_sl.value, max_tilt_sl.value+1, tinc_sl.value)
    sino   = radon(im, angs)
    recon  = backproject(sino, angs)
    ps_rec = power_spectrum(recon)

    fig, axes = plt.subplots(1, 4, figsize=(16, 4))
    axes[0].imshow(im,    cmap='gray', origin='lower'); axes[0].axis('off')
    axes[0].set_title("True image")
    axes[1].imshow(sino,  cmap='gray', origin='lower', aspect='auto')
    axes[1].set_title(f"Sinogram ({len(angs)} tilts)")
    axes[1].set_xlabel("x"); axes[1].set_ylabel("Tilt index")
    axes[2].imshow(recon, cmap='gray', origin='lower'); axes[2].axis('off')
    axes[2].set_title(f"Backprojection\n(±{max_tilt_sl.value}°, Δ={tinc_sl.value}°)")
    axes[3].imshow(ps_rec, cmap='inferno', origin='lower'); axes[3].axis('off')
    axes[3].set_title("Fourier space of reconstruction")
    plt.tight_layout()
    with out_bp:
        clear_output(wait=True); display(fig2img(fig))

for w in [max_tilt_sl, tinc_sl]:
    w.observe(_update_bp, names='value')
display(VBox([max_tilt_sl, tinc_sl, out_bp]))
_update_bp()
```

---

## Task 2 – Filtered backprojection

The blurring artifact of plain backprojection is caused by the over-sampling of low frequencies: every projection adds a constant along its direction, contributing to the low-frequency region of Fourier space more than to high frequencies. The fix is to pre-filter each projection with a **ramp filter** $|\nu|$ before backprojecting. This exactly counteracts the $1/|\nu|$ weighting introduced by backprojection.

$$
A_\text{FBP}(x,z) = \sum_\alpha \bigl[p_\alpha * h\bigr]\!\left(x\cos\alpha + z\sin\alpha\right), \qquad \hat h(\nu) = |\nu|
$$

```{admonition} Task 2
:class: note
1. Take a 1D projection.
2. Compute its Fourier transform using `np.fft.rfft`.
3. Multiply by the ramp filter $|\nu|$ (use `np.fft.rfftfreq` for the frequency axis).
4. Inverse-Fourier-transform back.
5. Use the filtered projections for backprojection.
```

```{code-cell} ipython3
:tags: [hide-input]
mt_sl2   = IntSlider(value=60, min=10, max=70, step=5,
                      description="Max tilt (°)",       style=_sty, layout=_sl)
tinc_sl2 = IntSlider(value=3,  min=1,  max=15, step=1,
                      description="Tilt increment (°)", style=_sty, layout=_sl)
out_fbp = Output()

def _update_fbp(_=None):
    angs  = np.arange(-mt_sl2.value, mt_sl2.value+1, tinc_sl2.value)
    sino  = radon(im, angs)
    bp    = backproject(sino, angs)
    fbp   = filtered_backproject(sino, angs)
    ps_bp = power_spectrum(bp); ps_fbp = power_spectrum(fbp)

    fig, axes = plt.subplots(2, 3, figsize=(14, 8))
    axes[0,0].imshow(im,  cmap='gray', origin='lower'); axes[0,0].axis('off')
    axes[0,0].set_title("True image")
    axes[0,1].imshow(bp,  cmap='gray', origin='lower'); axes[0,1].axis('off')
    axes[0,1].set_title(f"Backprojection (BP)\n(±{mt_sl2.value}°, Δ={tinc_sl2.value}°)")
    axes[0,2].imshow(fbp, cmap='gray', origin='lower'); axes[0,2].axis('off')
    axes[0,2].set_title("Filtered BP (FBP)")
    axes[1,0].axis('off')
    axes[1,1].imshow(ps_bp,  cmap='inferno', origin='lower'); axes[1,1].axis('off')
    axes[1,1].set_title("BP Fourier space")
    axes[1,2].imshow(ps_fbp, cmap='inferno', origin='lower'); axes[1,2].axis('off')
    axes[1,2].set_title("FBP Fourier space")
    plt.suptitle("Comparison: BP vs filtered BP", fontsize=10)
    plt.tight_layout()
    with out_fbp:
        clear_output(wait=True); display(fig2img(fig))

for w in [mt_sl2, tinc_sl2]:
    w.observe(_update_fbp, names='value')
display(VBox([mt_sl2, tinc_sl2, out_fbp]))
_update_fbp()
```

---

## The Fourier Slice Theorem

The **projection-slice theorem** (also called the Fourier slice theorem) states that the 1D Fourier transform of a projection at angle $\alpha$ equals a central slice through the 2D Fourier transform of the object at the same angle:

$$
\mathcal{F}_1\{p_\alpha\}(\nu) = \mathcal{F}_2\{A\}(\nu\cos\alpha,\, \nu\sin\alpha)
$$

Each projection therefore fills one line (passing through the origin) in 2D Fourier space. To reconstruct the full 3D (or 2D) object, we need enough projections to cover Fourier space densely.

The interactive below shows, for a chosen tilt angle:
1. The rotated specimen and its 1D projection
2. The 2D Fourier transform of the specimen with the corresponding central slice highlighted
3. A verification that the 1D FT of the projection matches the highlighted central slice

```{code-cell} ipython3
:tags: [hide-input]
from scipy.ndimage import map_coordinates

fst_sl  = IntSlider(value=0, min=-70, max=70, step=1,
                     description="Tilt angle (°)", style=_sty, layout=_sl)
out_fst = Output()

_n_fst  = im.shape[0]
_F2_obj = np.fft.fftshift(np.fft.fft2(im))

def _extract_central_slice(F2, angle_deg):
    """Extract central slice from 2D FT at given angle via bilinear interpolation."""
    n    = F2.shape[0]
    ang  = np.deg2rad(angle_deg)
    t    = np.linspace(-n//2, n//2, n)
    xs   = t * np.cos(ang) + n//2
    ys   = t * np.sin(ang) + n//2
    re   = map_coordinates(F2.real, [ys, xs], order=1, mode='nearest')
    im_p = map_coordinates(F2.imag, [ys, xs], order=1, mode='nearest')
    return re + 1j*im_p

def _update_fst(_=None):
    ang  = fst_sl.value
    rot  = nd_rotate(im, ang, reshape=False, cval=0)
    proj = rot.sum(axis=0) / _n_fst

    # 1D FT of projection
    fft1d      = np.fft.fftshift(np.fft.fft(proj))
    fft1d_amp  = np.log(np.abs(fft1d) + 1)

    # Central slice from 2D FT
    cslice     = _extract_central_slice(_F2_obj, ang)
    cslice_amp = np.log(np.abs(cslice) + 1)

    # Draw slice line on 2D FT image
    ps2 = np.log(np.abs(_F2_obj) + 1)
    n   = _n_fst
    t   = np.linspace(-n//2, n//2, n)
    ang_r = np.deg2rad(ang)
    xs_line = t * np.cos(ang_r) + n//2
    ys_line = t * np.sin(ang_r) + n//2

    fig, axes = plt.subplots(1, 4, figsize=(16, 4))
    axes[0].imshow(rot, cmap='gray', origin='lower'); axes[0].axis('off')
    axes[0].set_title(f"Specimen at α={ang}°")

    axes[1].plot(proj, np.arange(len(proj)), color='steelblue', lw=1.4)
    axes[1].set_xlabel("Projection"); axes[1].set_ylabel("x")
    axes[1].set_title("1D projection"); axes[1].grid(alpha=0.3)

    axes[2].imshow(ps2, cmap='inferno', origin='lower')
    axes[2].plot(xs_line, ys_line, 'w-', linewidth=1.5, alpha=0.8)
    axes[2].set_title(f"2D FT + central slice\nat {ang}°"); axes[2].axis('off')

    axes[3].plot(fft1d_amp,  label='FT of projection', color='steelblue', lw=1.4)
    axes[3].plot(cslice_amp, label='Central slice of 2D FT', color='red', lw=1.4, ls='--')
    axes[3].set_title("Fourier slice theorem\nverification")
    axes[3].set_xlabel("Frequency index"); axes[3].legend(fontsize=8); axes[3].grid(alpha=0.3)

    plt.tight_layout()
    with out_fst:
        clear_output(wait=True); display(fig2img(fig))

fst_sl.observe(_update_fst, names='value')
display(VBox([fst_sl, out_fst]))
_update_fst()
```

---

## Missing wedge

In practice the tilt range is limited to roughly $\pm 60°$ because:
- **Grid bars** obstruct the beam at high tilt
- **Increased effective sample thickness** at high tilt degrades image quality

This means the central region of Fourier space is never filled — the **missing wedge**. Features oriented along the tilt axis (perpendicular to the tilt direction) are most affected: they appear elongated in the reconstruction.

The interactive below shows the Fourier coverage at different tilt ranges. The missing wedge is the dark triangular gap.

```{code-cell} ipython3
:tags: [hide-input]
mt_mw  = IntSlider(value=60, min=10, max=90, step=5,
                    description="Max tilt (°)",       style=_sty, layout=_sl)
di_mw  = IntSlider(value=3,  min=1,  max=10, step=1,
                    description="Tilt increment (°)", style=_sty, layout=_sl)
out_mw = Output()

def _update_mw(_=None):
    n    = im.shape[0]
    angs = np.arange(-mt_mw.value, mt_mw.value+1, di_mw.value)
    sino = radon(im, angs)
    fbp  = filtered_backproject(sino, angs)

    # Coverage map: mark filled lines in Fourier space
    coverage = np.zeros((n, n))
    t = np.linspace(-n//2, n//2, n).astype(int)
    for ang in angs:
        ang_r = np.deg2rad(ang)
        xs = np.clip((t * np.cos(ang_r) + n//2).astype(int), 0, n-1)
        ys = np.clip((t * np.sin(ang_r) + n//2).astype(int), 0, n-1)
        coverage[ys, xs] = 1
    coverage_s = gaussian_filter(coverage, 1.5)

    fig, axes = plt.subplots(1, 4, figsize=(16, 4))
    axes[0].imshow(im,         cmap='gray',    origin='lower'); axes[0].axis('off')
    axes[0].set_title("True image")
    axes[1].imshow(fbp,        cmap='gray',    origin='lower'); axes[1].axis('off')
    axes[1].set_title(f"FBP (±{mt_mw.value}°, Δ={di_mw.value}°)")
    axes[2].imshow(coverage_s, cmap='Blues',   origin='lower'); axes[2].axis('off')
    axes[2].set_title(f"Fourier coverage\n({len(angs)} tilts)")

    # Draw missing-wedge wedge lines on FBP
    mid = n//2
    for sign in [+1, -1]:
        ang_r = np.deg2rad(sign * mt_mw.value + 90)
        t_line = np.linspace(-n//2, n//2, n)
        axes[3].plot(t_line*np.cos(ang_r)+mid, t_line*np.sin(ang_r)+mid,
                     'r--', lw=1.2)
    axes[3].imshow(power_spectrum(fbp), cmap='inferno', origin='lower')
    axes[3].set_title("FBP Fourier space\n(red = missing wedge boundary)")
    axes[3].axis('off')
    plt.suptitle(f"Missing wedge: tilt range ±{mt_mw.value}°", fontsize=10)
    plt.tight_layout()
    with out_mw:
        clear_output(wait=True); display(fig2img(fig))

for w in [mt_mw, di_mw]:
    w.observe(_update_mw, names='value')
display(VBox([mt_mw, di_mw, out_mw]))
_update_mw()
```

---

## Task 3 – The Crowther criterion

To fully reconstruct an object of diameter $D$ at resolution $r$, how many tilt angles $m$ do we need? The answer follows from the Fourier slice theorem and the sampling theorem.

Each projection fills one central line in 2D Fourier space. At resolution $r$, we need the outermost Fourier shell (radius $1/r$) to be densely sampled. Two adjacent central lines at the edge of the Fourier disk must be separated by at most $1/D$ (the Nyquist frequency across the object). The angular separation between adjacent projections at radius $1/r$ is therefore:

$$
\Delta\alpha \leq \frac{1/D}{1/r} = \frac{r}{D}
$$

Since the total range is $\pi$ radians (projections at $\alpha$ and $\alpha+\pi$ are identical), the number of tilts needed is:

$$
\boxed{m \geq \frac{\pi}{\Delta\alpha} = \frac{\pi D}{r}}
$$

This is the **Crowther criterion**. It tells us that to achieve resolution $r$ from an object of diameter $D$, we need at least $\pi D/r$ projections — independent of whether we are doing tomography or SPA (in SPA, particles in random orientations naturally fill the full angular range).

```{code-cell} ipython3
:tags: [hide-input]
# Interactive: show how many tilts are needed for different D and r
_sl_cr = Layout(width="440px")
D_sl = IntSlider(value=10, min=2, max=50, step=1,
                  description="Diameter D (nm)", style=_sty, layout=_sl_cr)
r_sl = FloatSlider(value=1.0, min=0.5, max=5.0, step=0.5,
                    description="Resolution r (nm)", style=_sty, layout=_sl_cr)
out_cr = Output()

def _update_cr(_=None):
    D     = D_sl.value; r = r_sl.value
    m_min = int(np.ceil(np.pi * D / r))
    D_arr = np.linspace(2, 50, 200)
    r_arr = np.array([0.5, 1.0, 2.0, 3.0, 5.0])

    fig, axes = plt.subplots(1, 2, figsize=(12, 4.5))
    for ri in r_arr:
        axes[0].plot(D_arr, np.pi * D_arr / ri, lw=1.6, label=f"r={ri} nm")
    axes[0].axvline(D, color='k', ls='--', lw=1.2)
    axes[0].axhline(m_min, color='red', ls=':', lw=1.2,
                    label=f"m={m_min} (D={D} nm, r={r} nm)")
    axes[0].set_xlabel("Object diameter D (nm)"); axes[0].set_ylabel("Min. tilts m")
    axes[0].set_title("Crowther criterion: m ≥ πD/r")
    axes[0].legend(fontsize=8); axes[0].grid(alpha=0.3)

    # Show what m tilts actually cover
    angs_cr = np.linspace(-90, 90, m_min, endpoint=False)
    n_cr    = 80
    cov_cr  = np.zeros((n_cr, n_cr))
    t_idx   = np.arange(n_cr) - n_cr//2
    for ang in angs_cr:
        ang_r = np.deg2rad(ang)
        xs = np.clip((t_idx*np.cos(ang_r)+n_cr//2).astype(int), 0, n_cr-1)
        ys = np.clip((t_idx*np.sin(ang_r)+n_cr//2).astype(int), 0, n_cr-1)
        cov_cr[ys, xs] = 1
    axes[1].imshow(gaussian_filter(cov_cr, 1), cmap='Blues', origin='lower')
    axes[1].set_title(f"Fourier coverage with m={m_min} tilts\n"
                      f"(D={D} nm, r={r} nm → full coverage at radius D/r={D//r:.0f})")
    axes[1].axis('off')
    plt.tight_layout()
    with out_cr:
        clear_output(wait=True); display(fig2img(fig))

for w in [D_sl, r_sl]:
    w.observe(_update_cr, names='value')
display(VBox([D_sl, r_sl, out_cr]))
_update_cr()
```

---

## Summary

| Concept | Key idea |
|---------|----------|
| Tilt series | Physical rotation of specimen; each tilt gives one 1D projection |
| Sinogram | Stack of 1D projections; point objects trace sinusoids |
| Radon transform | Mathematical projection operator $p_\alpha(t) = \int A(t\cos\alpha - s\sin\alpha, t\sin\alpha + s\cos\alpha)\,ds$ |
| Backprojection | Smear each projection back along its direction and sum |
| Ramp filter | $\hat h(\nu) = |\nu|$; corrects low-frequency overweighting of backprojection |
| Fourier slice theorem | $\mathcal{F}_1\{p_\alpha\}(\nu) = \mathcal{F}_2\{A\}(\nu\cos\alpha, \nu\sin\alpha)$ |
| Missing wedge | Triangular gap in Fourier coverage from limited tilt range |
| Crowther criterion | $m \geq \pi D / r$ tilts for resolution $r$ from object of diameter $D$ |