feat(audio): account for source directivity in diffuse reverb calculations

Previously, diffuse field reverberant gain and critical distance calculations
assumed an omnidirectional sound source (Q = 1), regardless of cone settings.
Directional sources concentrate sound energy in specific directions, reducing
the total acoustic power injected into the room and lowering the reverberant
field level relative to the direct sound on-axis.

Incorporate source directivity factor Q into the acoustics model:
- Compute closed-form directivity factor Q in src/physics/cone.ts by
  integrating squared cone gain over the unit sphere.
- Scale reverberant field gain by 1 / sqrt(Q) (-10·log10(Q) dB) and adjust
  critical distance by sqrt(Q) in src/audio/reverb.ts.
- Update AudioEngine and signal path UI to present directivity index (DI)
  and adjust room dominance evaluations.
- Add unit tests verifying Q against brute-force numerical quadrature and
  confirming critical distance behavior.

Test Plan:
- `npm run typecheck` -- passed clean
- `npm test` -- 125/125 tests passed (6 test suites)
- `git diff --cached --check` -- no whitespace errors

README.md

@@ -40,11 +40,21 @@ locally in the page.
 
 The one idea the app is built around. The direct sound falls as `1/r`; the
 reverberant field does not fall at all. They cross at the **critical distance**,
-`rc = √(R / 16π)`, drawn as a ring on the floor:
+`rc = √(Q·R / 16π)`, drawn as a ring on the floor:
 
 - **Inside the ring** the source dominates. You can point at it with your eyes shut.
 - **Outside it** the room dominates. It does not get quieter as you back away — and that constancy is a large part of how your ears judge distance.
 
+`Q` is the source's directivity factor — how much more intensity it puts on its
+axis than an omnidirectional source of the same total power. It belongs in that
+formula because the reverberant field is fed by *radiated power*, while the ring
+is measured against the *on-axis* direct sound. A cone that throws its energy
+forward excites the room far less than its on-axis level suggests: the default
+70°/160° cone has `Q = 5.07`, a directivity index of 7.1 dB, which pushes the
+ring from 6.1 m out to 13.7 m. Drop `Q` and reflections stay conspicuously loud
+when you walk behind the speaker — the room is being driven as if the speaker
+radiated backwards as hard as it radiates forwards.
+
 Press <kbd>3</kbd> for a guided walk through it.
 
 ## Controls

src/audio/AudioEngine.ts

@@ -6,6 +6,7 @@ import {
   airAbsorptionCutoff,
   clamp,
   coneGain,
+  directivityFactor,
   distance,
   distanceGain,
   gainToDb,
@@ -22,6 +23,7 @@ import {
   RoomDimensions,
   SurfaceId,
   computeAcoustics,
+  criticalDistanceFor,
   generateImpulseResponse,
   getSurface,
   reverberantGain,
@@ -106,6 +108,9 @@ export interface Telemetry {
   reverbGainDb: number;
   /** Direct-to-reverberant ratio in dB. Negative means the room dominates. */
   directToReverbDb: number;
+  /** Source directivity index, 10·log10(Q). 0 dB is omnidirectional. */
+  directivityIndexDb: number;
+  /** On-axis critical distance, in metres — where the room overtakes the source. */
   criticalDistance: number;
   t60: number;
   reverbDominant: boolean;
@@ -132,6 +137,7 @@ export const SILENT_TELEMETRY: Telemetry = {
   directGainDb: -120,
   reverbGainDb: -120,
   directToReverbDb: 0,
+  directivityIndexDb: 0,
   criticalDistance: 0,
   t60: 0,
   reverbDominant: false,
@@ -205,6 +211,8 @@ export class AudioEngine {
 
   private playing = false;
   private acoustics: RoomAcoustics;
+  /** Source directivity factor Q, cached because it only moves with the cone. */
+  private directivity: number;
   private impulseDirty = false;
   private lastImpulseAt = -Infinity;
 
@@ -218,6 +226,15 @@ export class AudioEngine {
   constructor(settings?: Partial<EngineSettings>) {
     this.settings = { ...DEFAULT_SETTINGS, ...settings };
     this.acoustics = computeAcoustics(this.settings.room, getSurface(this.settings.surface).absorption);
+    this.directivity = directivityFactor(this.coneParams());
+  }
+
+  private coneParams(): { innerAngle: number; outerAngle: number; outerGain: number } {
+    return {
+      innerAngle: this.settings.coneInnerAngle,
+      outerAngle: this.settings.coneOuterAngle,
+      outerGain: this.settings.coneOuterGain,
+    };
   }
 
   get context(): AudioContext | null {
@@ -258,6 +275,11 @@ export class AudioEngine {
    *    power the source radiates, not on where the listener stands or what is
    *    in the way. This is what makes walking away sound like distance instead
    *    of like turning down a fader.
+   *
+   *    Total radiated power is where directivity enters: the send is referenced
+   *    to the on-axis direct level, so it has to be divided by the cone's Q to
+   *    describe a source that radiates that much *forward* rather than in every
+   *    direction. See `reverberantGain`.
    */
   async init(): Promise<void> {
     if (this.ctx) return;
@@ -313,9 +335,18 @@ export class AudioEngine {
     this.masterGain = ctx.createGain();
     this.masterGain.gain.setValueAtTime(this.settings.masterVolume, now);
 
-    // Catches the peaks when a close, on-axis source stacks with a wet room.
+    // A safety net for the peaks when a close, on-axis source stacks with a wet
+    // room — deliberately not a compressor.
+    //
+    // The threshold sits just under 0 dBFS because this limiter is the only
+    // thing in the graph that sees absolute level, and every dB it takes off the
+    // loud position is a dB of front-to-back contrast the simulation loses. At
+    // the old -6 dB it was gain-reducing the near-field case by ~2.9 dB at the
+    // default volume (and ~5.7 dB at full) while leaving the quiet, mostly
+    // reverberant position untouched — squashing the very difference the user is
+    // meant to hear when they walk behind the speaker.
     this.limiter = ctx.createDynamicsCompressor();
-    this.limiter.threshold.setValueAtTime(-6, now);
+    this.limiter.threshold.setValueAtTime(-1.5, now);
     this.limiter.knee.setValueAtTime(0, now);
     this.limiter.ratio.setValueAtTime(20, now);
     this.limiter.attack.setValueAtTime(0.003, now);
@@ -477,6 +508,7 @@ export class AudioEngine {
     const { preset, masterVolume, ...rest } = patch;
     this.settings = { ...this.settings, ...rest, room };
     this.acoustics = computeAcoustics(room, getSurface(this.settings.surface).absorption);
+    this.directivity = directivityFactor(this.coneParams());
 
     if (masterVolume !== undefined) this.setMasterVolume(masterVolume);
     if (preset !== undefined) this.setPreset(preset);
@@ -499,11 +531,7 @@ export class AudioEngine {
       maxDistance: s.maxDistance,
       rolloffFactor: s.rolloffFactor,
     });
-    const gCone = coneGain(input.sourcePos, input.sourceForward, input.listenerPos, {
-      innerAngle: s.coneInnerAngle,
-      outerAngle: s.coneOuterAngle,
-      outerGain: s.coneOuterGain,
-    });
+    const gCone = coneGain(input.sourcePos, input.sourceForward, input.listenerPos, this.coneParams());
     const angle = offAxisAngle(input.sourcePos, input.sourceForward, input.listenerPos);
     const occl = occlusionResponse(input.occlusion);
     const airCutoff = airAbsorptionCutoff(d, s.airAbsorption);
@@ -523,7 +551,7 @@ export class AudioEngine {
     const ratio = clamp(1 - delayRate, DOPPLER_MIN_RATIO, DOPPLER_MAX_RATIO);
 
     const gDirect = gDistance * gCone * occl.gain;
-    const gReverb = reverberantGain(this.acoustics, s.refDistance);
+    const gReverb = reverberantGain(this.acoustics, s.refDistance, this.directivity);
 
     if (this.ctx && this.playing) {
       this.applyToGraph(input, gDirect, gReverb, occl.cutoff, airCutoff);
@@ -551,7 +579,8 @@ export class AudioEngine {
       directGainDb: gainToDb(gDirect),
       reverbGainDb: reverbDb,
       directToReverbDb: gainToDb(gDirect) - reverbDb,
-      criticalDistance: this.acoustics.criticalDistance,
+      directivityIndexDb: 10 * Math.log10(this.directivity),
+      criticalDistance: criticalDistanceFor(this.acoustics, this.directivity),
       t60: this.acoustics.t60,
       // Compare the levels themselves rather than distance against rc. With a
       // non-inverse falloff law or a reference distance other than 1 m the two

src/audio/reverb.ts

@@ -56,12 +56,25 @@ export function computeAcoustics(dims: RoomDimensions, absorption: number): Room
 
   const t60 = clamp((0.161 * volume) / (surfaceArea * alpha), 0.08, 8);
   const roomConstant = (surfaceArea * alpha) / (1 - alpha);
-  // r_c = sqrt(R / 16π) for an omnidirectional source.
+  // r_c = sqrt(R / 16π) for an omnidirectional source. Directivity is a property
+  // of the source, not the room, so it is applied at the point of use by
+  // `criticalDistanceFor` rather than baked in here.
   const criticalDistance = Math.sqrt(roomConstant / (16 * Math.PI));
 
   return { t60, roomConstant, criticalDistance, surfaceArea, volume };
 }
 
+/**
+ * Critical distance along the source's axis, r_c = sqrt(Q R / 16π).
+ *
+ * A directional source has to be followed further before the room catches up
+ * with it, because on its axis it is louder by sqrt(Q) while the room it excites
+ * is not. Pass Q = 1 to recover the omnidirectional radius.
+ */
+export function criticalDistanceFor(acoustics: RoomAcoustics, directivity = 1): number {
+  return acoustics.criticalDistance * Math.sqrt(Math.max(directivity, 1e-6));
+}
+
 /**
  * Level of the diffuse reverberant field relative to the direct sound measured
  * at the reference distance.
@@ -70,10 +83,25 @@ export function computeAcoustics(dims: RoomDimensions, absorption: number): Room
  * distance-attenuated. Setting its level to refDistance / criticalDistance
  * makes wet and dry balance exactly at the critical distance — walk further
  * than that and the room takes over, which is the effect real rooms have.
+ *
+ * `directivity` is the source's Q. It matters because this gain is referenced to
+ * the *on-axis* direct sound at refDistance, where both the distance law and the
+ * cone return exactly 1. That reference already contains the source's sqrt(Q) of
+ * directivity gain, but the reverberant field does not: a cone that throws its
+ * energy forward excites the room with the Q-th part of the power an omni source
+ * of the same on-axis level would. Dividing it back out is the difference
+ * between the textbook Lp = Lw + 10·log10(Q/4πr² + 4/R) and a source that is
+ * mysteriously still loud from behind — with the default cone, 10·log10(5.07)
+ * = 7.1 dB of it.
  */
-export function reverberantGain(acoustics: RoomAcoustics, refDistance = 1): number {
-  if (acoustics.criticalDistance < 1e-6) return 1;
-  return clamp(refDistance / acoustics.criticalDistance, 0, 1.4);
+export function reverberantGain(
+  acoustics: RoomAcoustics,
+  refDistance = 1,
+  directivity = 1
+): number {
+  const rc = criticalDistanceFor(acoustics, directivity);
+  if (rc < 1e-6) return 1;
+  return clamp(refDistance / rc, 0, 1.4);
 }
 
 /**

src/physics/cone.ts

@@ -42,3 +42,43 @@ export function coneGain(
   const t = (angle - innerHalf) / (outerHalf - innerHalf);
   return clamp(1 + t * (floor - 1), 0, 1);
 }
+
+/**
+ * Directivity factor Q = 4π / ∫ g(θ)² dΩ — how much more intensity the source
+ * puts on its axis than an omnidirectional source radiating the same total
+ * power. Q = 1 is omnidirectional; the default 70°/160°/0.08 cone gives 5.07,
+ * i.e. a directivity index of 7.05 dB.
+ *
+ * This is what couples the cone to the room. A directional source pushes far
+ * less power into the space than its on-axis level suggests, so the reverberant
+ * field it excites is weaker by exactly this factor — see `reverberantGain`.
+ *
+ * Evaluated in closed form. The cone model interpolates linearly *in gain*, so
+ * over the transition band g(θ)² is a quadratic in θ and the integral
+ * ∫ (c₀ + c₁θ)² sin θ dθ has an elementary antiderivative. Cheap enough to call
+ * every time the cone sliders move, and exact rather than quadrature-noisy.
+ */
+export function directivityFactor(p: ConeParams): number {
+  const toRad = Math.PI / 180;
+  const a = clamp((Math.max(0, p.innerAngle) / 2) * toRad, 0, Math.PI);
+  const b = clamp((Math.max(0, p.outerAngle) / 2) * toRad, a, Math.PI);
+  const f = clamp(p.outerGain, 0, 1);
+
+  // Inside the inner cone g = 1; outside the outer cone g = outerGain.
+  let integral = 2 * Math.PI * (1 - Math.cos(a)) + 2 * Math.PI * f * f * (1 + Math.cos(b));
+
+  if (b - a > 1e-9) {
+    const c1 = (f - 1) / (b - a);
+    const c0 = 1 - c1 * a;
+    // d/dθ of ∫(c₀ + c₁θ)² sin θ dθ, expanded term by term.
+    const F = (t: number): number =>
+      c0 * c0 * -Math.cos(t) +
+      2 * c0 * c1 * (Math.sin(t) - t * Math.cos(t)) +
+      c1 * c1 * (2 * t * Math.sin(t) - (t * t - 2) * Math.cos(t));
+    integral += 2 * Math.PI * (F(b) - F(a));
+  }
+
+  // g ≤ 1 everywhere, so the integral cannot exceed 4π and Q cannot fall below
+  // 1. Guard the division anyway: a degenerate cone would otherwise return ∞.
+  return integral > 1e-9 ? Math.max(1, (4 * Math.PI) / integral) : 1;
+}

src/ui/signalPath.ts

@@ -97,7 +97,15 @@ export class SignalPath {
     this.text(this.direct.db, fmt.db(t.directGainDb));
     this.bar(this.direct.fill, dbToFill(t.directGainDb));
 
-    this.text(this.room.detail, `t60 ${fmt.seconds(t.t60)}`);
+    // The room level is set by t60 and by how much power the cone actually
+    // throws into the room, so both belong on the row. Without the directivity
+    // index the wet level looks arbitrarily low next to a loud direct path.
+    this.text(
+      this.room.detail,
+      t.directivityIndexDb > 0.05
+        ? `t60 ${fmt.seconds(t.t60)} · directivity ${fmt.db(t.directivityIndexDb)}`
+        : `t60 ${fmt.seconds(t.t60)}`
+    );
     this.text(this.room.db, fmt.db(t.reverbGainDb));
     this.bar(this.room.fill, dbToFill(t.reverbGainDb));
 

test/engine.test.ts

@@ -420,11 +420,34 @@ describe('telemetry consistency', () => {
   it('flags room dominance at the critical distance when the listener is on-axis', async () => {
     const engine = new AudioEngine();
     await engine.start();
-    const rc = engine.getAcoustics().criticalDistance;
+    // The reported critical distance is the on-axis one, sqrt(Q) further out
+    // than the bare room radius, because the source is directional. Taking the
+    // room's own radius here would test the crossover of a source this engine
+    // is not simulating.
+    const rc = engine.updateSpatial(poseAt(0, -1), 1 / 60).criticalDistance;
+    expect(rc).toBeGreaterThan(engine.getAcoustics().criticalDistance);
     // poseAt(0, -d) puts the listener dead ahead of a source facing −Z, so the
     // cone contributes no attenuation and distance alone decides.
     expect(engine.updateSpatial(poseAt(0, -rc * 0.5), 1 / 60).reverbDominant).toBe(false);
     expect(engine.updateSpatial(poseAt(0, -rc * 2), 1 / 60).reverbDominant).toBe(true);
+    // And it crosses over where it says it does, to within a hair.
+    expect(engine.updateSpatial(poseAt(0, -rc), 1 / 60).directToReverbDb).toBeCloseTo(0, 6);
+  });
+
+  it('quietens the reverberant field by the cone directivity index', async () => {
+    const engine = new AudioEngine();
+    await engine.start();
+    const directional = engine.updateSpatial(poseAt(0, -1), 1 / 60);
+    expect(directional.directivityIndexDb).toBeCloseTo(7.05, 2);
+
+    // Open the cone right up and the same room gets louder by exactly that much.
+    engine.update({ coneInnerAngle: 360, coneOuterAngle: 360 });
+    const omni = engine.updateSpatial(poseAt(0, -1), 1 / 60);
+    expect(omni.directivityIndexDb).toBeCloseTo(0, 6);
+    expect(omni.reverbGainDb - directional.reverbGainDb).toBeCloseTo(
+      directional.directivityIndexDb,
+      6
+    );
   });
 
   it('counts directivity towards room dominance, not just distance', async () => {

test/physics.test.ts

@@ -1,7 +1,9 @@
 import { describe, expect, it } from 'vitest';
+import type { ConeParams } from '../src/physics';
 import {
   airAbsorptionCutoff,
   coneGain,
+  directivityFactor,
   distanceGain,
   dopplerRatio,
   gainToDb,
@@ -258,3 +260,63 @@ describe('unit helpers', () => {
     expect(twice).toBeCloseTo(single, 9);
   });
 });
+
+describe('cone directivity factor', () => {
+  /** Q = 4π / ∫g²dΩ by brute-force quadrature, to check the closed form. */
+  const quadrature = (p: ConeParams, steps = 200000): number => {
+    const forward = { x: 0, y: 0, z: -1 };
+    let integral = 0;
+    for (let i = 0; i < steps; i++) {
+      const theta = ((i + 0.5) / steps) * Math.PI;
+      // Place a listener at angle theta from the source's forward axis.
+      const listener = { x: Math.sin(theta), y: 0, z: -Math.cos(theta) };
+      const g = coneGain({ x: 0, y: 0, z: 0 }, forward, listener, p);
+      integral += g * g * Math.sin(theta) * (Math.PI / steps);
+    }
+    return (4 * Math.PI) / (2 * Math.PI * integral);
+  };
+
+  it('is 1 for an omnidirectional source', () => {
+    expect(directivityFactor({ innerAngle: 360, outerAngle: 360, outerGain: 0.08 })).toBeCloseTo(1, 12);
+    expect(directivityFactor({ innerAngle: 70, outerAngle: 160, outerGain: 1 })).toBeCloseTo(1, 12);
+  });
+
+  it('matches quadrature for the default cone', () => {
+    const p = { innerAngle: 70, outerAngle: 160, outerGain: 0.08 };
+    expect(directivityFactor(p)).toBeCloseTo(quadrature(p), 4);
+    // A 7.05 dB directivity index — the amount the reverb send must come down.
+    expect(10 * Math.log10(directivityFactor(p))).toBeCloseTo(7.05, 2);
+  });
+
+  it('matches quadrature across the range the sliders can reach', () => {
+    const cases: ConeParams[] = [
+      { innerAngle: 0, outerAngle: 360, outerGain: 0 },
+      { innerAngle: 70, outerAngle: 70, outerGain: 0.08 },
+      { innerAngle: 30, outerAngle: 90, outerGain: 0.001 },
+      { innerAngle: 10, outerAngle: 20, outerGain: 0.5 },
+      { innerAngle: 180, outerAngle: 360, outerGain: 0.3 },
+      { innerAngle: 200, outerAngle: 100, outerGain: 0.2 },
+    ];
+    for (const p of cases) {
+      expect(directivityFactor(p)).toBeCloseTo(quadrature(p), 3);
+    }
+  });
+
+  it('never reports less than omnidirectional, whatever the cone', () => {
+    for (let inner = 0; inner <= 360; inner += 40) {
+      for (let outer = 0; outer <= 360; outer += 40) {
+        for (const outerGain of [0, 0.08, 0.5, 1]) {
+          const q = directivityFactor({ innerAngle: inner, outerAngle: outer, outerGain });
+          expect(Number.isFinite(q)).toBe(true);
+          expect(q).toBeGreaterThanOrEqual(1);
+        }
+      }
+    }
+  });
+
+  it('gets more directional as the cone narrows', () => {
+    const wide = directivityFactor({ innerAngle: 160, outerAngle: 250, outerGain: 0.08 });
+    const narrow = directivityFactor({ innerAngle: 20, outerAngle: 60, outerGain: 0.08 });
+    expect(narrow).toBeGreaterThan(wide);
+  });
+});

test/reverb.test.ts

@@ -2,10 +2,12 @@ import { describe, expect, it } from 'vitest';
 import {
   SURFACES,
   computeAcoustics,
+  criticalDistanceFor,
   generateImpulseResponse,
   getSurface,
   reverberantGain,
 } from '../src/audio/reverb';
+import { directivityFactor } from '../src/physics';
 import { MockAudioContext } from './mockAudioContext';
 
 const room = { width: 20, height: 8, depth: 20 };
@@ -58,6 +60,33 @@ describe('reverberant field level', () => {
     expect(reverberantGain(acoustics, 1)).toBeCloseTo(1 / acoustics.criticalDistance, 6);
   });
 
+  it('still balances at the critical distance once the source is directional', () => {
+    const acoustics = computeAcoustics(room, 0.22);
+    const q = directivityFactor({ innerAngle: 70, outerAngle: 160, outerGain: 0.08 });
+    const rc = criticalDistanceFor(acoustics, q);
+    // The invariant is unchanged, it just moves outward: on axis the dry path is
+    // 1/r and the wet path is flat, so they must still cross at exactly rc.
+    expect(reverberantGain(acoustics, 1, q)).toBeCloseTo(1 / rc, 6);
+    expect(rc).toBeCloseTo(acoustics.criticalDistance * Math.sqrt(q), 6);
+  });
+
+  it('drops the wet level by the directivity index of the cone', () => {
+    const acoustics = computeAcoustics(room, 0.22);
+    const q = directivityFactor({ innerAngle: 70, outerAngle: 160, outerGain: 0.08 });
+    const omni = reverberantGain(acoustics, 1, 1);
+    const directional = reverberantGain(acoustics, 1, q);
+    // 10·log10(Q) = 7.05 dB for the default cone.
+    expect(20 * Math.log10(directional / omni)).toBeCloseTo(-10 * Math.log10(q), 6);
+    expect(10 * Math.log10(q)).toBeCloseTo(7.05, 2);
+  });
+
+  it('is unchanged for an omnidirectional source', () => {
+    const acoustics = computeAcoustics(room, 0.22);
+    const q = directivityFactor({ innerAngle: 360, outerAngle: 360, outerGain: 0.08 });
+    expect(q).toBeCloseTo(1, 12);
+    expect(reverberantGain(acoustics, 1, q)).toBeCloseTo(reverberantGain(acoustics, 1), 12);
+  });
+
   it('is louder in a more reflective room', () => {
     const dead = reverberantGain(computeAcoustics(room, 0.9));
     const live = reverberantGain(computeAcoustics(room, 0.035));