1. 用强磁场（1.5T / 3T）把你体内水分子里的氢质子"排队"
2. 每个质子像陀螺一样绕磁场方向旋转（进动），转速 = 拉莫尔频率
   —— 转速由它所在位置的磁场强度决定
3. 打一发射频脉冲，把这群陀螺"推歪"
4. 陀螺们一边转一边发出电磁信号，线圈把它们的信号"加起来"收进来

import torch

# 三个"质子"，各自以不同频率旋转（拉莫尔进动），线圈只收到它们的合唱
fs, T = 256, 1.0                                   # 采样率、总时长
t = torch.arange(0, fs * T) / fs                   # 时间轴
freqs, amps = [3.0, 6.0, 9.0], [1.0, 0.6, 0.3]     # 频率 & 各自的"质子数量"

signal = sum(a * torch.cos(2 * torch.pi * f * t) for f, a in zip(freqs, amps))
print("线圈实际测到的只是一条乱七八糟的曲线:", signal[:5].round(decimals=3))

spec = torch.fft.rfft(signal).abs() / (len(t) / 2)  # 傅里叶：把合唱拆回独唱
peaks = torch.nonzero(spec > 0.05).flatten()
for k in peaks:
    print(f"  找回频率 {k.item():>2d} Hz, 强度 {spec[k].item():.2f}")
print("原始强度:", amps)

线圈实际测到的只是一条乱七八糟的曲线: tensor([1.9000, 1.8840, 1.8350, 1.7550, 1.6460])
  找回频率  3 Hz, 强度 1.00
  找回频率  6 Hz, 强度 0.60
  找回频率  9 Hz, 强度 0.30
原始强度: [1.0, 0.6, 0.3]

位置 x 处的磁场   = B0 + G·x
→ 该处质子的转速 = 正比于 (B0 + G·x)
→ 扫描一段时间后，位置 x 处的质子比中心多转了  k·x  的角度
                                              ↑
                                    k 由梯度强度和时间决定，机器说了算

import torch

N = 64
x = torch.arange(N).float()
body = torch.zeros(N); body[20:30] = 1.0; body[40:45] = 0.6   # 一维"人体"：两块组织

# 采一条 k 空间：第 k 步梯度让位置 x 的质子多转 2*pi*k*x/N 的相位，线圈把全身加起来
def acquire(body, k):
    phase = -2 * torch.pi * k * x / N
    return (body * torch.exp(1j * phase)).sum()      # 一个复数 = 一个采样点

ks = torch.arange(N)
kspace = torch.stack([acquire(body, k) for k in ks])  # 扫描仪真正存下来的东西
print("k 空间前 4 个采样点:", [f"{c:.2f}" for c in kspace[:4].tolist()])

recon = torch.fft.ifft(kspace).real                   # 逆傅里叶 = 还原图像
print("还原误差:", (recon - body).abs().max().item())
print("还原出的身体:", recon.round(decimals=2)[18:32])

k 空间前 4 个采样点: ['13.00+0.00j', '-8.77-3.98j', '-0.27+5.78j', '6.73-4.91j']
还原误差: 3.346924586367095e-06
还原出的身体: tensor([0., 0., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., -0., 0.])

如果热是一种"装在铁里的流体"：
  铁块里的热质是有限的 → 磨到一定程度就该"磨干"了
但实际上：
  只要摩擦不停，热就不停 → 热不可能是守恒的流体
  → 热必然与「运动」有关

N = 1,000,000 个采样点
  朴素算法: 10^12 次运算
  FFT:      2 × 10^7 次运算
  ——— 快了 5 万倍。算不动的事，变成了实时能干的事。

核爆:     一瞬间的各向同性挤压 → 短、脆、高频丰富、几乎没有长周期面波
天然地震: 断层错动（双力偶）    → 有拖很久的低频面波列

import torch
torch.manual_seed(0)

fs, T = 40.0, 120.0                                   # 地震台：40 Hz 采样，记录 2 分钟
t = torch.arange(int(fs * T)) / fs

def burst(t0, f, tau, amp):                           # 一段从 t0 开始、衰减时长 tau 的振动
    env = torch.where(t >= t0, torch.exp(-(t - t0) / tau), torch.zeros_like(t))
    return amp * env * torch.sin(2 * torch.pi * f * (t - t0))

# 核爆：一次性挤压，短、脆、高频丰富，几乎没有长周期面波
boom = burst(20, 4.0, 1.5, 1.0) + burst(20, 7.0, 0.8, 0.6)
# 天然地震：断层错动，有 P 波，更有拖很久的低频面波
quake = burst(20, 3.0, 1.2, 0.5) + burst(24, 0.08, 25.0, 1.0) + burst(26, 0.05, 30.0, 0.8)

def band_ratio(x):                                    # 高频能量 / 低频能量（mb:Ms 判据的玩具版）
    spec = torch.fft.rfft(x).abs() ** 2
    f = torch.fft.rfftfreq(len(x), 1 / fs)
    hi = spec[(f > 1.0) & (f < 10.0)].sum()
    lo = spec[(f > 0.02) & (f < 0.5)].sum()
    return (hi / lo).item()

for name, wave in [("核爆", boom), ("天然地震", quake)]:
    noise = 0.15 * torch.randn(len(t))
    x = wave + noise
    print(f"{name}: 时域峰值 {x.abs().max():.2f}  |  高频/低频能量比 = {band_ratio(x):9.2f}")

核爆: 时域峰值 1.46  |  高频/低频能量比 =     22.05
天然地震: 时域峰值 1.47  |  高频/低频能量比 =      0.15

import torch
torch.manual_seed(0)

N, f0 = 16384, 137                                    # 16384 个采样点，地球对面传来的一丝震动
n = torch.arange(N).float()
signal = 0.05 * torch.sin(2 * torch.pi * f0 * n / N)  # 振幅 0.05：比噪声小 20 倍

def peak_snr(x):                                      # 频域里，目标频率的峰 比 周围噪底 高多少倍
    spec = torch.fft.rfft(x).abs()
    peak = spec[f0]
    floor = torch.cat([spec[:f0 - 3], spec[f0 + 4:]]).median()
    return (peak / floor).item()

x1 = signal + torch.randn(N)                          # 单个台站：信号被噪声彻底淹没
print(f"时域信噪比        : {signal.std() / 1.0:.4f}   （示波器上什么都看不到）")
print(f"单台站频域峰信噪比: {peak_snr(x1):8.1f} 倍")

for m in [1, 4, 16, 64, 256]:                         # m 个台站：同一个信号 + 各自独立的噪声
    stack = sum(signal + torch.randn(N) for _ in range(m)) / m
    print(f"  {m:3d} 个台站叠加 -> 频域峰信噪比 {peak_snr(stack):8.1f} 倍")

时域信噪比        : 0.0354   （示波器上什么都看不到）
单台站频域峰信噪比:      4.7 倍
    1 个台站叠加 -> 频域峰信噪比      5.6 倍
    4 个台站叠加 -> 频域峰信噪比      7.2 倍
   16 个台站叠加 -> 频域峰信噪比     15.5 倍
   64 个台站叠加 -> 频域峰信噪比     31.2 倍
  256 个台站叠加 -> 频域峰信噪比     61.3 倍

import torch

def spin(z, angle):                      # 把复数 z 转 angle 弧度：z * e^{i*angle}
    return z * torch.exp(1j * torch.tensor(angle))

def mri_phase(rho, x, k):                # MRI：位置 x 的质子，在第 k 步梯度下转 k*x 度
    return spin(rho, k * x)

def rope(vec, m, theta=1.0):             # RoPE：第 m 个位置的词向量，转 m*theta 度
    return spin(vec, m * theta)

a = torch.tensor(1.0 + 0j)               # 一个"质子" / 一个词向量
for (m, n) in [(0, 3), (2, 5), (10, 13)]:
    mri  = (mri_phase(a, float(m), 1.0).conj() * mri_phase(a, float(n), 1.0)).real
    rp   = (rope(a, m).conj() * rope(a, n)).real
    print(f"位置对 ({m:2d},{n:2d})  距离={n-m}  MRI相干度={mri.item():.4f}  RoPE点积={rp.item():.4f}")

位置对 ( 0, 3)  距离=3  MRI相干度=-0.9900  RoPE点积=-0.9900
位置对 ( 2, 5)  距离=3  MRI相干度=-0.9900  RoPE点积=-0.9900
位置对 (10,13)  距离=3  MRI相干度=-0.9900  RoPE点积=-0.9900

import torch

# 造一张 64x64 的"断层图"：一个大圆（器官）+ 一个小圆（病灶）+ 一条细线（血管）
N = 64
yy, xx = torch.meshgrid(torch.arange(N).float(), torch.arange(N).float(), indexing='ij')
img = ((xx - 32) ** 2 + (yy - 32) ** 2 < 20 ** 2).float() * 0.5
img += ((xx - 40) ** 2 + (yy - 26) ** 2 < 5 ** 2).float() * 0.5
img[30:34, 12:52] += 0.4

K = torch.fft.fftshift(torch.fft.fft2(img))          # 全部 k 空间
r = ((xx - N // 2) ** 2 + (yy - N // 2) ** 2).sqrt()

for radius in [2, 4, 8, 16, 32]:                      # 只保留半径 radius 内的低频
    kept = K * (r <= radius)
    rec = torch.fft.ifft2(torch.fft.ifftshift(kept)).real
    energy = (kept.abs() ** 2).sum() / (K.abs() ** 2).sum()
    err = (rec - img).abs().mean()
    print(f"只采半径 {radius:2d} 的中心低频 -> 占总能量 {energy:.4%}, 图像误差 {err:.4f}")

只采半径  2 的中心低频 -> 占总能量 83.4581%, 图像误差 0.0940
只采半径  4 的中心低频 -> 占总能量 90.0067%, 图像误差 0.0707
只采半径  8 的中心低频 -> 占总能量 95.2790%, 图像误差 0.0449
只采半径 16 的中心低频 -> 占总能量 97.9025%, 图像误差 0.0258
只采半径 32 的中心低频 -> 占总能量 99.7135%, 图像误差 0.0104

MRI:   k 空间中心的低频 → 器官的大致轮廓
       k 空间边缘的高频 → 病灶边界、血管细节

RoPE:  转得慢的维度对   → 记住"隔了几十上百个词"的远距离关系
       转得快的维度对   → 分辨"相邻词"的细微先后

同一个设计哲学：一根针不够，得配一整套时针、分针、秒针。

import torch, math
import torch.nn.functional as F

N = 127; C = N // 2                                  # 用奇数边长，让"中心"正好落在一个像素上
yy, xx = torch.meshgrid(torch.arange(N).float(), torch.arange(N).float(), indexing='ij')
img = torch.exp(-((xx - C) ** 2 / 450 + (yy - C) ** 2 / 72))   # 一团沿 x 拉长的组织

def project(img, deg):                               # 站在 deg 角度看过去，把整块图压成一条线
    a = math.radians(deg)
    R = torch.tensor([[math.cos(a), -math.sin(a), 0.],
                      [math.sin(a),  math.cos(a), 0.]])[None]
    grid = F.affine_grid(R, (1, 1, N, N), align_corners=True)
    return F.grid_sample(img[None, None], grid, align_corners=True)[0, 0].sum(dim=0)

K2 = torch.fft.fftshift(torch.fft.fft2(torch.fft.ifftshift(img)))   # 完整的二维频域
t = torch.arange(N).float() - C

for deg in [0, 30, 45, 90]:
    p = project(img, deg)
    P = torch.fft.fftshift(torch.fft.fft(torch.fft.ifftshift(p)))   # 投影的一维傅里叶
    a = math.radians(deg)
    iy = (t * math.sin(a) + C).round().long()        # 频域中过原点、方向为 deg 的那条线
    ix = (t * math.cos(a) + C).round().long()
    line = K2[iy, ix]
    corr = (P.conj() * line).sum().abs() / (P.norm() * line.norm())
    print(f"{deg:3d}° 投影的一维FFT  ==  二维频域里同角度的那条线 ?  相似度 = {corr:.4f}")

  0° 投影的一维FFT  ==  二维频域里同角度的那条线 ?  相似度 = 1.0000
 30° 投影的一维FFT  ==  二维频域里同角度的那条线 ?  相似度 = 0.9947
 45° 投影的一维FFT  ==  二维频域里同角度的那条线 ?  相似度 = 0.9822
 90° 投影的一维FFT  ==  二维频域里同角度的那条线 ?  相似度 = 1.0000

import torch, math

N = 128; C = N // 2
yy, xx = torch.meshgrid(torch.arange(N).float(), torch.arange(N).float(), indexing='ij')
img = ((xx - C) ** 2 + (yy - C) ** 2 < 40 ** 2).float() * 0.5      # 器官
img += ((xx - 80) ** 2 + (yy - 52) ** 2 < 10 ** 2).float() * 0.5   # 病灶
img[60:68, 24:104] += 0.4                                          # 血管

K = torch.fft.fftshift(torch.fft.fft2(img))          # 上帝视角：完整的世界
t = torch.arange(-C, C).float()

def mask_of(n_angles):                                # n 个角度，各自在频域里划一条过原点的线
    m = torch.zeros(N, N, dtype=torch.bool)
    for a in [math.pi * i / n_angles for i in range(n_angles)]:
        iy = (t * math.sin(a) + C).round().long().clamp(0, N - 1)
        ix = (t * math.cos(a) + C).round().long().clamp(0, N - 1)
        m[iy, ix] = True
    return m

for n in [1, 2, 4, 8, 16, 32, 64, 128]:
    m = mask_of(n)
    rec = torch.fft.ifft2(torch.fft.ifftshift(K * m)).real
    cover = m.float().mean()
    err = (rec - img).abs().mean() / img.abs().mean()
    print(f"{n:3d} 个视角 -> 频域覆盖 {cover:6.1%}, 还原误差 {err:6.1%}")

  1 个视角 -> 频域覆盖   0.8%, 还原误差  94.3%
  2 个视角 -> 频域覆盖   1.6%, 还原误差  80.3%
  4 个视角 -> 频域覆盖   2.7%, 还原误差  47.4%
  8 个视角 -> 频域覆盖   5.7%, 还原误差  36.3%
 16 个视角 -> 频域覆盖  11.4%, 还原误差  24.0%
 32 个视角 -> 频域覆盖  21.7%, 还原误差  16.6%
 64 个视角 -> 频域覆盖  40.0%, 还原误差  10.6%
128 个视角 -> 频域覆盖  65.0%, 还原误差   6.2%

一条射线的走时 = 沿途所有介质速度的"积分"
              ↑
        这就是一次投影！（第 3 节那个"把整块图压扁成一条线"）

全球几千个台站 × 每年几万次地震 = 几百万条不同角度的射线
              ↓
        角度够多 → 反演出地球内部的三维速度结构

文本   = 从"语言"这个角度对世界做的一次投影
图像   = 从"视觉"这个角度对世界做的一次投影
音频   = 从"听觉"这个角度对世界做的一次投影
代码   = 从"可执行逻辑"这个角度对世界做的一次投影
视频   = 从"时间连续性"这个角度对世界做的一次投影
动作/触觉 = 从"物理交互"这个角度对世界做的一次投影

import torch
torch.manual_seed(0)

n, D, d = 200, 32, 4                      # 200 个"事物"，真实世界 32 维，每个模态只能看到 4 维
Z = torch.randn(n, D)                     # 柏拉图理想型：世界本身
def kernel(X):                            # 用"两两相似度"来刻画一个表征空间的形状
    Xn = X / X.norm(dim=1, keepdim=True)
    return Xn @ Xn.T

def align(A, B):                          # 两个表征空间的形状有多像
    a, b = A.flatten(), B.flatten()
    a, b = a - a.mean(), b - b.mean()
    return (a @ b / (a.norm() * b.norm())).item()

K_world = kernel(Z)
mods = [Z @ torch.randn(D, d) for _ in range(16)]      # 16 个模态：各自是世界的一次投影
print(f"单个模态 vs 世界: {align(kernel(mods[0]), K_world):.4f}")
for m in [1, 2, 4, 8, 16]:
    fused = sum(kernel(x) for x in mods[:m]) / m       # 把 m 个模态的"形状"叠起来
    print(f"融合 {m:2d} 个模态 -> 与真实世界的表征对齐度 {align(fused, K_world):.4f}")

单个模态 vs 世界: 0.3126
融合  1 个模态 -> 与真实世界的表征对齐度 0.3126
融合  2 个模态 -> 与真实世界的表征对齐度 0.4336
融合  4 个模态 -> 与真实世界的表征对齐度 0.5653
融合  8 个模态 -> 与真实世界的表征对齐度 0.6972
融合 16 个模态 -> 与真实世界的表征对齐度 0.8067

# 多头注意力：768 维 → 12 个头 × 64 维
k = k.view(B, T, self.n_head, C // self.n_head).transpose(1, 2)
# 每个头只"看"768 维空间中的一个 64 维切片 —— 一个头 = 一个角度
y = y.transpose(1, 2).contiguous().view(B, T, C)   # 12 个角度的观察拼回完整重建

python scripts/mri_to_multimodal_visualization.py