Add Tier-1 mastering metrics: Crest Factor, PSR, True Peak; LRA on LUFS

Three new metrics plug into the existing Metric registry. All inherit
the async compute path, so first-use is non-blocking and subsequent
switches hit the per-(file, metric) cache.

- CrestFactorMetric: 20*log10(peak/RMS) per 1 s window, 0.25 s hop.
  Uses cumsum-of-squares for O(N) RMS. Reference lines at 12 dB
  (roomy) and 6 dB (heavily limited).
- PSRMetric: sample-peak (dBFS) minus short-term LUFS over the same
  3 s window as LUFSMetric, so the two plots line up. Reference lines
  at 10 LU (good punch) and 4 LU (squashed).
- TruePeakMetric: ITU-R BS.1770 oversampled true peak via
  scipy.signal.resample_poly at 4x over 250 ms / 100 ms windows.
  data["integrated_tp_db"] carries the track-wide max.
- LUFSMetric now also computes Meter.loudness_range and surfaces it
  in the plot legend alongside the integrated value.
- Shared _to_dbfs helper floors 20*log10 at _EPS so silent windows
  don't explode.

Co-Authored-By: Claude Opus 4.7 <noreply@anthropic.com>

metrics.py

@@ -21,11 +21,21 @@ import matplotlib.colors as mcolors
 import matplotlib.cm as cm
 from matplotlib.figure import Figure
 import pyloudnorm as pyln
+from scipy import signal as scipy_signal
 
 from font_manager import safe_title
 from master_core import AudioFile
 
 
+# Small constant to keep 20*log10(...) from blowing up on perfect silence.
+_EPS = 1e-12
+
+
+def _to_dbfs(linear: np.ndarray | float) -> np.ndarray | float:
+  """Convert a linear magnitude to dBFS, floored at _EPS."""
+  return 20.0 * np.log10(np.maximum(linear, _EPS))
+
+
 class Metric(ABC):
   """A pluggable analysis metric."""
 
@@ -132,7 +142,7 @@ class WaveformMetric(Metric):
 
 
 class LUFSMetric(Metric):
-  """ITU-R BS.1770 loudness: short-term (3 s) time series + integrated value.
+  """ITU-R BS.1770 loudness: short-term (3 s) time series + integrated + LRA.
 
   Powered by pyloudnorm. The time series slides `meter.integrated_loudness`
   across the track because pyloudnorm doesn't expose a per-block series.
@@ -152,7 +162,6 @@ class LUFSMetric(Metric):
     sr = audio_file.sr
     meter = pyln.Meter(sr)
 
-    # pyloudnorm warns on clipping and on too-short audio; we handle both.
     with warnings.catch_warnings():
       warnings.simplefilter("ignore")
       integrated = self._safe_integrated(meter, y)
@@ -161,9 +170,9 @@ class LUFSMetric(Metric):
       hop_n = int(self.HOP_S * sr)
 
       if len(y) < window_n:
-        # Track shorter than 3 s — just one data point at the centre.
         times = np.array([len(y) / (2.0 * sr)])
         lufs = np.array([integrated if np.isfinite(integrated) else self.SILENCE_FLOOR])
+        lra = float("nan")
       else:
         n_windows = 1 + (len(y) - window_n) // hop_n
         lufs = np.empty(n_windows)
@@ -171,6 +180,10 @@ class LUFSMetric(Metric):
           start = i * hop_n
           lufs[i] = self._safe_integrated(meter, y[start:start + window_n])
         times = (np.arange(n_windows) * hop_n + window_n / 2.0) / sr
+        try:
+          lra = float(meter.loudness_range(y))
+        except (ValueError, FloatingPointError):
+          lra = float("nan")
 
     lufs = np.where(np.isfinite(lufs), lufs, self.SILENCE_FLOOR)
     lufs = np.clip(lufs, self.SILENCE_FLOOR, 0.0)
@@ -179,6 +192,7 @@ class LUFSMetric(Metric):
       "times": times,
       "lufs": lufs,
       "integrated": float(integrated),
+      "lra": lra,
     }
 
   @staticmethod
@@ -192,6 +206,7 @@ class LUFSMetric(Metric):
     times = data["times"]
     lufs = data["lufs"]
     integrated = data["integrated"]
+    lra = data.get("lra", float("nan"))
 
     fig = Figure(figsize=figsize, facecolor="white")
     ax = fig.add_subplot(111)
@@ -203,6 +218,10 @@ class LUFSMetric(Metric):
         label=f"Integrated: {integrated:.1f} LUFS",
       )
 
+    if np.isfinite(lra):
+      # Invisible plot entry to surface LRA in the legend without adding a line.
+      ax.plot([], [], " ", label=f"LRA: {lra:.1f} LU")
+
     # Streaming target reference (Spotify normalises to -14 LUFS).
     ax.axhline(-14.0, color="gray", linestyle=":", linewidth=0.8, alpha=0.6)
     ax.text(
@@ -221,11 +240,222 @@ class LUFSMetric(Metric):
     return fig
 
 
+class CrestFactorMetric(Metric):
+  """Crest factor = 20*log10(peak / RMS) per sliding window, in dB."""
+
+  id = "crest_factor"
+  display_name = "Crest Factor"
+
+  WINDOW_S = 1.0
+  HOP_S = 0.25
+
+  def compute(self, audio_file: AudioFile):
+    y = audio_file.y_mono.astype(np.float64, copy=False)
+    sr = audio_file.sr
+    window_n = int(self.WINDOW_S * sr)
+    hop_n = int(self.HOP_S * sr)
+
+    if len(y) < window_n:
+      times = np.array([len(y) / (2.0 * sr)])
+      peak = float(np.max(np.abs(y))) if len(y) else 0.0
+      rms = float(np.sqrt(np.mean(y * y))) if len(y) else 0.0
+      crest = 20.0 * np.log10(max(peak, _EPS) / max(rms, _EPS))
+      return {"times": times, "crest_db": np.array([crest])}
+
+    # RMS via cumulative-sum-of-squares (O(N)); peaks via sliding window view.
+    y2 = y * y
+    cumsum = np.concatenate(([0.0], np.cumsum(y2)))
+    n_windows = 1 + (len(y) - window_n) // hop_n
+    starts = np.arange(n_windows) * hop_n
+    ends = starts + window_n
+    mean_sq = (cumsum[ends] - cumsum[starts]) / window_n
+    rms = np.sqrt(np.maximum(mean_sq, _EPS))
+
+    abs_y = np.abs(y)
+    peaks = np.empty(n_windows)
+    for i in range(n_windows):
+      peaks[i] = np.max(abs_y[starts[i]:ends[i]])
+
+    crest_db = 20.0 * np.log10(np.maximum(peaks, _EPS) / rms)
+    times = (starts + window_n / 2.0) / sr
+    return {"times": times, "crest_db": crest_db}
+
+  def render(self, data, file_path, figsize=(10, 4)) -> Figure:
+    times = data["times"]
+    crest_db = data["crest_db"]
+
+    fig = Figure(figsize=figsize, facecolor="white")
+    ax = fig.add_subplot(111)
+    ax.plot(times, crest_db, color="#e09f3e", linewidth=1.4, label=f"Crest factor (1 s)")
+
+    # Rules of thumb: ~12 dB = roomy, ~6 dB = heavily limited.
+    ax.axhline(12.0, color="gray", linestyle=":", linewidth=0.8, alpha=0.6)
+    ax.text(times[-1], 12.0, "  12 dB", va="center", ha="left", fontsize=8, alpha=0.6)
+    ax.axhline(6.0, color="gray", linestyle=":", linewidth=0.8, alpha=0.6)
+    ax.text(times[-1], 6.0, "  6 dB (squashed)", va="center", ha="left", fontsize=8, alpha=0.6)
+
+    ax.set_ylim(0.0, 25.0)
+    ax.set_xlim(times[0], times[-1])
+    ax.set_ylabel("Crest factor (dB)")
+    ax.set_xlabel("Time (seconds)")
+    ax.set_title(safe_title(os.path.basename(file_path)))
+    ax.grid(True, alpha=0.3)
+    ax.legend(loc="lower right", fontsize=8)
+    fig.tight_layout()
+    return fig
+
+
+class PSRMetric(Metric):
+  """Peak-to-Short-term LUFS Ratio (sample-peak variant), in LU.
+
+  PSR = sample_peak_dBFS - short_term_LUFS  over the same 3 s windows used by
+  LUFSMetric. High PSR = punchy transients; low PSR = heavily limited.
+  """
+
+  id = "psr"
+  display_name = "PSR"
+
+  WINDOW_S = 3.0
+  HOP_S = 0.5
+  SILENCE_FLOOR = -70.0
+
+  def compute(self, audio_file: AudioFile):
+    y = audio_file.y_mono.astype(np.float64, copy=False)
+    sr = audio_file.sr
+    meter = pyln.Meter(sr)
+
+    window_n = int(self.WINDOW_S * sr)
+    hop_n = int(self.HOP_S * sr)
+
+    with warnings.catch_warnings():
+      warnings.simplefilter("ignore")
+      if len(y) < window_n:
+        times = np.array([len(y) / (2.0 * sr)])
+        peak_db = _to_dbfs(np.max(np.abs(y))) if len(y) else self.SILENCE_FLOOR
+        lufs = LUFSMetric._safe_integrated(meter, y)
+        psr = peak_db - lufs if np.isfinite(lufs) else 0.0
+        return {"times": times, "psr": np.array([psr])}
+
+      n_windows = 1 + (len(y) - window_n) // hop_n
+      abs_y = np.abs(y)
+      lufs_series = np.empty(n_windows)
+      peaks_db = np.empty(n_windows)
+      for i in range(n_windows):
+        start = i * hop_n
+        end = start + window_n
+        peaks_db[i] = _to_dbfs(np.max(abs_y[start:end]))
+        lufs_series[i] = LUFSMetric._safe_integrated(meter, y[start:end])
+      times = (np.arange(n_windows) * hop_n + window_n / 2.0) / sr
+
+    # PSR is meaningless where the loudness reading is below the absolute gate.
+    valid = np.isfinite(lufs_series) & (lufs_series > self.SILENCE_FLOOR)
+    psr = np.where(valid, peaks_db - lufs_series, np.nan)
+    return {"times": times, "psr": psr}
+
+  def render(self, data, file_path, figsize=(10, 4)) -> Figure:
+    times = data["times"]
+    psr = data["psr"]
+
+    fig = Figure(figsize=figsize, facecolor="white")
+    ax = fig.add_subplot(111)
+    ax.plot(times, psr, color="#7251b5", linewidth=1.4, label="PSR (3 s)")
+
+    # Ian Shepherd's rough thresholds.
+    ax.axhline(10.0, color="gray", linestyle=":", linewidth=0.8, alpha=0.6)
+    ax.text(times[-1], 10.0, "  10 LU (good punch)", va="center", ha="left", fontsize=8, alpha=0.6)
+    ax.axhline(4.0, color="gray", linestyle=":", linewidth=0.8, alpha=0.6)
+    ax.text(times[-1], 4.0, "  4 LU (squashed)", va="center", ha="left", fontsize=8, alpha=0.6)
+
+    ax.set_ylim(0.0, 25.0)
+    ax.set_xlim(times[0], times[-1])
+    ax.set_ylabel("PSR (LU)")
+    ax.set_xlabel("Time (seconds)")
+    ax.set_title(safe_title(os.path.basename(file_path)))
+    ax.grid(True, alpha=0.3)
+    ax.legend(loc="lower right", fontsize=8)
+    fig.tight_layout()
+    return fig
+
+
+class TruePeakMetric(Metric):
+  """ITU-R BS.1770 true peak via 4x polyphase oversampling, in dBTP.
+
+  Per-window true peak with a moderate hop so it renders quickly. Windows are
+  oversampled independently — slight edge under-detection at window boundaries
+  is masked by the 60% overlap.
+  """
+
+  id = "true_peak"
+  display_name = "True Peak"
+
+  WINDOW_S = 0.25
+  HOP_S = 0.1
+  OVERSAMPLE = 4
+
+  def compute(self, audio_file: AudioFile):
+    y = audio_file.y_mono.astype(np.float32, copy=False)
+    sr = audio_file.sr
+    window_n = int(self.WINDOW_S * sr)
+    hop_n = int(self.HOP_S * sr)
+
+    if len(y) < window_n:
+      y_up = scipy_signal.resample_poly(y, self.OVERSAMPLE, 1) if len(y) else np.zeros(1, dtype=np.float32)
+      peak_db = _to_dbfs(np.max(np.abs(y_up))) if len(y_up) else -70.0
+      return {
+        "times": np.array([len(y) / (2.0 * sr)]),
+        "tp_db": np.array([peak_db]),
+        "integrated_tp_db": float(peak_db),
+      }
+
+    n_windows = 1 + (len(y) - window_n) // hop_n
+    tp_db = np.empty(n_windows)
+    for i in range(n_windows):
+      start = i * hop_n
+      w = y[start:start + window_n]
+      w_up = scipy_signal.resample_poly(w, self.OVERSAMPLE, 1)
+      tp_db[i] = _to_dbfs(np.max(np.abs(w_up)))
+    times = (np.arange(n_windows) * hop_n + window_n / 2.0) / sr
+
+    integrated_tp_db = float(np.max(tp_db))
+    return {"times": times, "tp_db": tp_db, "integrated_tp_db": integrated_tp_db}
+
+  def render(self, data, file_path, figsize=(10, 4)) -> Figure:
+    times = data["times"]
+    tp_db = data["tp_db"]
+    integrated = data.get("integrated_tp_db", float("nan"))
+
+    fig = Figure(figsize=figsize, facecolor="white")
+    ax = fig.add_subplot(111)
+    ax.plot(times, tp_db, color="#c1121f", linewidth=1.0, label="True Peak (250 ms)")
+
+    # 0 dBTP = sample-level clip; -1 dBTP a common mastering ceiling.
+    ax.axhline(0.0, color="black", linestyle="--", linewidth=1.0, alpha=0.8)
+    ax.text(times[-1], 0.0, "  0 dBTP (clip)", va="center", ha="left", fontsize=8, alpha=0.7)
+    ax.axhline(-1.0, color="gray", linestyle=":", linewidth=0.8, alpha=0.6)
+    ax.text(times[-1], -1.0, "  -1 dBTP (typical ceiling)", va="center", ha="left", fontsize=8, alpha=0.6)
+
+    if np.isfinite(integrated):
+      ax.plot([], [], " ", label=f"Max: {integrated:.2f} dBTP")
+
+    ax.set_ylim(-30.0, 6.0)
+    ax.set_xlim(times[0], times[-1])
+    ax.set_ylabel("dBTP")
+    ax.set_xlabel("Time (seconds)")
+    ax.set_title(safe_title(os.path.basename(file_path)))
+    ax.grid(True, alpha=0.3)
+    ax.legend(loc="lower right", fontsize=8)
+    fig.tight_layout()
+    return fig
+
+
 METRICS: dict[str, Metric] = {
   m.id: m for m in (
     RMSPowerMetric(),
     WaveformMetric(),
     LUFSMetric(),
+    CrestFactorMetric(),
+    PSRMetric(),
+    TruePeakMetric(),
   )
 }
 DEFAULT_METRIC_ID = "rms_power"