Time-varying random phase maps by decimated subband synthesis; drag-to-zoom, scrollbar pan and zoom buttons for the waveforms
12 changed files+602−152
README.mdmodified+3−1View file
@@ -10,7 +10,9 @@ The source sound A is split into overlapping Gaussian frequency bands, the short
1010
1111 where z_k is the complex subband signal and phi_k is the phase map evaluated at the band center. Every subband envelope |z_k(t)| is unchanged, so A and B have the same spectrogram at the analysis resolution, but their waveforms differ. Because phi_k does not depend on time, the operation is a linear time-invariant filter with response G(f) = sum_k exp(i phi_k) H_k(f), and the app computes it that way, with one FFT pair, in a Web Worker.
1212
13-The controls are the window width sigma (0.3 ms to 50 ms), the phase map (smoothed random with adjustable smoothing and strength, a pure delay as a control, or a chirp whose group delay grows linearly with frequency), and the source. The app plots the map, the magnitude response and group delay of the equivalent filter, the waveforms of A and B with a zoom, and the spectrograms of A and B on a shared scale. A and B can be played singly or in sequence with the same gain.
13+The random map can also vary in time. It is then a random field phi(f, t) smoothed in frequency over a chosen number of bands and in time over a chosen time scale, and B is synthesized band by band as y(t) = Re sum_k exp(i phi_k(t)) z_k(t). Each band is narrow, so its subband signal is computed at a reduced sample rate from its slice of the spectrum, rotated, and transformed back, which keeps the cost small even with thousands of bands. A phase that changes in time is frequency modulation: band k shifts by phi_k'(t) / 2 pi, and for a field of standard deviation s and time scale T the typical shift is s / (2 pi T), which the app reports. The expectation is that slow variation (seconds) is inaudible and fast variation (tens of milliseconds) is heard as warble or roughness, with the transition near the pitch threshold of about 1 Hz for low tones.
14+
15+The controls are the window width sigma (0.3 ms to 50 ms), the phase map (smoothed random with adjustable smoothing in frequency and in time, strength, and seed; a pure delay as a control; or a chirp whose group delay grows linearly with frequency), and the source. The app plots the map over time and frequency, the magnitude response and group delay of the equivalent filter (at the middle of the clip when the map varies), the waveforms of A and B with a zoom, and the spectrograms of A and B on a shared scale. A and B can be played singly or in sequence with the same gain.
1416
1517 The design and the reasoning behind it are described in [docs/spectrogram-phase-proposal.md](../phase-perception/docs/spectrogram-phase-proposal.md) in the first app's repository.
1618
src/bank.tsmodified+130−0View file
@@ -166,3 +166,133 @@ export function spectrogram(x: Float32Array, fs: number, sigmaSec: number, fmax:
166166 }
167167 return { data, rows, cols, durationSec: n / fs, fmax: rows * ((nb + 1) / rows) * binHz }
168168 }
169+
170+/** S(f) = sum_k W(f - f_k) on the bins 0..L/2, the normalizer that makes the bank sum to one. */
171+export function bandSum(centers: Float64Array, sigmaSec: number, L: number, fs: number): Float64Array {
172+ const sf = sigmaF(sigmaSec)
173+ const binHz = fs / L
174+ const nb = L / 2 + 1
175+ const S = new Float64Array(nb)
176+ const reach = 5 * sf
177+ const inv2s2 = 1 / (2 * sf * sf)
178+ for (let k = 0; k < centers.length; k++) {
179+ const fc = centers[k]
180+ const i0 = Math.max(0, Math.ceil((fc - reach) / binHz))
181+ const i1 = Math.min(nb - 1, Math.floor((fc + reach) / binHz))
182+ for (let i = i0; i <= i1; i++) {
183+ const d = i * binHz - fc
184+ S[i] += Math.exp(-d * d * inv2s2)
185+ }
186+ }
187+ return S
188+}
189+
190+const scratch = new Map<number, { re: Float64Array; im: Float64Array }>()
191+const scratchFor = (m: number) => {
192+ let s = scratch.get(m)
193+ if (!s) {
194+ s = { re: new Float64Array(m), im: new Float64Array(m) }
195+ scratch.set(m, s)
196+ }
197+ s.re.fill(0)
198+ s.im.fill(0)
199+ return s
200+}
201+
202+/**
203+ * y = Re sum_k exp(i phi_k(t)) z_k(t) with a phase that varies in time, by
204+ * exact subband synthesis. Each band is narrow, so its subband signal is
205+ * computed at a reduced sample rate: the band's slice of the spectrum is
206+ * demodulated to baseband with a small inverse FFT, rotated sample by sample,
207+ * transformed back, and added into the output spectrum at its original bins.
208+ * The small transform is sized to hold the band plus the frequency shift the
209+ * rotation can impose (maxRateHz), so nothing aliases.
210+ *
211+ * phiAt(k, times) must return the phase of band k at the given times, in
212+ * seconds from the start of the padded buffer (the signal starts at padSec).
213+ */
214+export function applyTimeVarying(
215+ x: Float32Array,
216+ fs: number,
217+ sigmaSec: number,
218+ centers: Float64Array,
219+ padSec: number,
220+ phiAt: (k: number, times: Float64Array) => Float64Array,
221+ maxRateHz: number,
222+): Float32Array {
223+ const n = x.length
224+ const pad = Math.ceil(padSec * fs)
225+ const L = nextPow2(n + 2 * pad)
226+ const half = L / 2
227+ const binHz = fs / L
228+ const sf = sigmaF(sigmaSec)
229+ const reach = 5 * sf
230+ const inv2s2 = 1 / (2 * sf * sf)
231+ const S = bandSum(centers, sigmaSec, L, fs)
232+
233+ const re = new Float64Array(L)
234+ const im = new Float64Array(L)
235+ re.set(x, pad)
236+ fftPlan(L).transform(re, im)
237+ // Positive half of the spectrum. The rotated subbands are one-sided, and the real part of a
238+ // one-sided signal has half its spectrum on each side, so the positive half of X (not the
239+ // doubled analytic spectrum) is what the Hermitian mirror below needs.
240+ const xr = new Float64Array(half + 1)
241+ const xi = new Float64Array(half + 1)
242+ for (let i = 0; i <= half; i++) {
243+ xr[i] = re[i]
244+ xi[i] = im[i]
245+ }
246+
247+ const yr = new Float64Array(L)
248+ const yi = new Float64Array(L)
249+ const extra = Math.ceil((1.5 * maxRateHz) / binHz) + 16
250+ for (let k = 0; k < centers.length; k++) {
251+ const fc = centers[k]
252+ const i0 = Math.max(0, Math.ceil((fc - reach) / binHz))
253+ const i1 = Math.min(half, Math.floor((fc + reach) / binHz))
254+ if (i1 < i0) continue
255+ const M = Math.min(L, nextPow2(i1 - i0 + 1 + 2 * extra))
256+ const ic = Math.round(fc / binHz)
257+ const { re: sr, im: si } = scratchFor(M)
258+ for (let i = i0; i <= i1; i++) {
259+ const d = i * binHz - fc
260+ const w = S[i] > 0 ? Math.exp(-d * d * inv2s2) / S[i] : 0
261+ const j = (((i - ic) % M) + M) % M
262+ sr[j] = xr[i] * w
263+ si[j] = xi[i] * w
264+ }
265+ const plan = fftPlan(M)
266+ plan.transform(sr, si, true) // z_k demodulated, sampled every L/M samples
267+ const times = new Float64Array(M)
268+ for (let m = 0; m < M; m++) times[m] = (m * L) / (M * fs)
269+ const phi = phiAt(k, times)
270+ for (let m = 0; m < M; m++) {
271+ const c = Math.cos(phi[m])
272+ const s = Math.sin(phi[m])
273+ const a = sr[m]
274+ const b = si[m]
275+ sr[m] = a * c - b * s
276+ si[m] = a * s + b * c
277+ }
278+ plan.transform(sr, si)
279+ for (let j = 0; j < M; j++) {
280+ const off = j < M / 2 ? j : j - M
281+ const i = ic + off
282+ if (i < 0 || i > half) continue
283+ yr[i] += sr[j]
284+ yi[i] += si[j]
285+ }
286+ }
287+ // Hermitian mirror so the output is real
288+ yi[0] = 0
289+ yi[half] = 0
290+ for (let i = 1; i < half; i++) {
291+ yr[L - i] = yr[i]
292+ yi[L - i] = -yi[i]
293+ }
294+ fftPlan(L).transform(yr, yi, true)
295+ const y = new Float32Array(n)
296+ for (let i = 0; i < n; i++) y[i] = yr[pad + i]
297+ return y
298+}
src/components/AnalysisPanel.tsxmodified+57−21View file
@@ -1,5 +1,7 @@
1-import type { PhaseMapSpec } from '../phasemap'
1+import { UNIFORM_STD, type PhaseMapSpec } from '../phasemap'
2+import { phaseColor } from '../plot'
23 import type { RenderResult } from '../types'
4+import { ImagePlot } from './ImagePlot'
35 import { LinePlot } from './LinePlot'
46
57 type Props = {
@@ -13,9 +15,9 @@ type Props = {
1315
1416 const SIGMA_LOG_MIN = Math.log10(0.3)
1517 const SIGMA_LOG_MAX = Math.log10(50)
18+const TIME_LOG_MIN = 0 // 1 ms
19+const TIME_LOG_MAX = Math.log10(3000)
1620 const sig3 = (x: number) => Number(x.toPrecision(3))
17-
18-const MAP_COLOR = '#d9641e'
1921 const THRESHOLD_MS = 1
2022
2123 /** Window width, phase map choice and parameters, and the plots that describe the resulting filter. */
@@ -23,21 +25,25 @@ export function AnalysisPanel({ sigmaMs, setSigmaMs, map, setMap, fmax, result }
2325 const sigmaF = 1000 / (2 * Math.PI * sigmaMs)
2426 const setKind = (kind: PhaseMapSpec['kind']) => {
2527 if (kind === map.kind) return
26- if (kind === 'random') setMap({ kind, smoothBands: 8, seed: 1, strength: 1 })
28+ if (kind === 'random') setMap({ kind, smoothBands: 8, timeScaleMs: null, seed: 1, strength: 1 })
2729 else if (kind === 'delay') setMap({ kind, delayMs: sig3(sigmaMs / 4) })
2830 else setMap({ kind, delayAtNyquistMs: sig3(sigmaMs / 4) })
2931 }
32+ const timeVarying = map.kind === 'random' && map.timeScaleMs !== null
33+ // a Gaussian-smoothed field with std s and time scale T has a derivative of std s / T
34+ const fmHz = map.kind === 'random' && map.timeScaleMs !== null ? (map.strength * UNIFORM_STD) / (2 * Math.PI * (map.timeScaleMs / 1000)) : 0
3035
3136 return (
3237 <section className="card">
3338 <h2>Filter bank and phase map</h2>
3439 <p className="hint">
3540 The sound is split into overlapping Gaussian bands of width sigma in time (and 1 / (2 pi sigma) in frequency). Sound B is made by
36- rotating the phase of each band by phi(f) and adding the bands back together, which leaves every band's envelope unchanged. For a
37- fixed map this is an all-pass filter; its group delay, minus the slope of the phase, is what smears transients. Peaks of group delay
38- beyond about {THRESHOLD_MS} ms are reported to be audible on impulsive sounds. Because the map is sampled at the band centers, it
39- must change by less than about a radian from one band to the next; the magnitude plot shows when it does not. For a pure delay this
40- means delays shorter than sigma, so the achievable group delay grows with the window width.
41+ rotating the phase of each band by phi(f, t) and adding the bands back together, which leaves every band's envelope unchanged. A map
42+ that is fixed in time is an all-pass filter; its group delay, minus the slope of the phase in frequency, is what smears transients,
43+ and peaks beyond about {THRESHOLD_MS} ms are reported to be audible on impulsive sounds. A map that changes in time shifts the
44+ frequency of each band by the rate of change of its phase, which is audible on steady sounds once the shift approaches the pitch
45+ threshold. Because the map is sampled at the band centers, it must change by less than about a radian from one band to the next; the
46+ magnitude plot shows when it does not.
4147 </p>
4248 <div className="global-controls">
4349 <label>
@@ -81,7 +87,7 @@ export function AnalysisPanel({ sigmaMs, setSigmaMs, map, setMap, fmax, result }
8187 {map.kind === 'random' && (
8288 <>
8389 <label>
84- <span>Smoothing</span>
90+ <span>Smoothing in frequency</span>
8591 <input type="range" min={1} max={40} step={1} value={map.smoothBands} onChange={(e) => setMap({ ...map, smoothBands: Number(e.target.value) })} />
8692 <span className="slider-unit">
8793 {map.smoothBands} bands ({sig3(map.smoothBands * sigmaF)} Hz)
@@ -113,6 +119,41 @@ export function AnalysisPanel({ sigmaMs, setSigmaMs, map, setMap, fmax, result }
113119 </label>
114120 )}
115121 </div>
122+ {map.kind === 'random' && (
123+ <div className="global-controls">
124+ <label>
125+ <input type="checkbox" checked={map.timeScaleMs === null} onChange={(e) => setMap({ ...map, timeScaleMs: e.target.checked ? null : 200 })} />
126+ <span>Fixed in time</span>
127+ </label>
128+ {map.timeScaleMs !== null && (
129+ <label>
130+ <span>Smoothing in time</span>
131+ <input
132+ type="range"
133+ className="wide-range"
134+ min={TIME_LOG_MIN}
135+ max={TIME_LOG_MAX}
136+ step={0.01}
137+ value={Math.min(TIME_LOG_MAX, Math.max(TIME_LOG_MIN, Math.log10(map.timeScaleMs)))}
138+ onChange={(e) => setMap({ ...map, timeScaleMs: sig3(10 ** Number(e.target.value)) })}
139+ />
140+ <input
141+ type="number"
142+ className="slider-value"
143+ min={0.1}
144+ step="any"
145+ value={map.timeScaleMs}
146+ onChange={(e) => {
147+ const v = Number(e.target.value)
148+ if (Number.isFinite(v) && v > 0) setMap({ ...map, timeScaleMs: v })
149+ }}
150+ />
151+ <span className="slider-unit">ms</span>
152+ <span className="info">typical frequency shift {sig3(fmHz)} Hz; the window is {sigmaMs} ms</span>
153+ </label>
154+ )}
155+ </div>
156+ )}
116157
117158 {result && result.maxMagDeviationDb > 1 && (
118159 <div className="warning">
@@ -123,23 +164,18 @@ export function AnalysisPanel({ sigmaMs, setSigmaMs, map, setMap, fmax, result }
123164 {result && (
124165 <div className="response-plots">
125166 <div>
126- <div className="plot-label">Phase map phi(f), wrapped to (-180°, 180°]</div>
127- <LinePlot
128- series={[{ x: result.map.f, y: Array.from(result.map.y, (v) => (v * 180) / Math.PI), color: MAP_COLOR, dots: result.map.f.length <= 400, width: 1 }]}
129- xRange={[0, fmax]}
130- yRange={[-190, 190]}
131- xLabel="frequency (Hz)"
132- yLabel="degrees"
133- zeroLine
134- />
167+ <div className="plot-label">
168+ Phase map phi(f, t){result.timeVarying ? '' : ', fixed in time'} <span className="phase-legend" /> -180° to 180°
169+ </div>
170+ <ImagePlot image={result.phaseField} color={phaseColor} label="" labelColor="#000" />
135171 </div>
136172 <div>
137- <div className="plot-label">Magnitude response of the equivalent filter (should stay near 0 dB)</div>
173+ <div className="plot-label">Magnitude response of the equivalent filter{timeVarying ? ' at mid-clip' : ''} (should stay near 0 dB)</div>
138174 <LinePlot series={[{ x: result.magDb.f, y: result.magDb.y, color: '#1f5fa8' }]} xRange={[0, fmax]} yRange={[-30, 3]} xLabel="frequency (Hz)" yLabel="dB" zeroLine />
139175 </div>
140176 <div>
141177 <div className="plot-label">
142- Group delay of the equivalent filter (largest {result.maxGroupDelayMs.toFixed(2)} ms)
178+ Group delay{timeVarying ? ' at mid-clip' : ''} (largest {result.maxGroupDelayMs.toFixed(2)} ms)
143179 </div>
144180 <LinePlot
145181 series={[{ x: result.groupDelayMs.f, y: result.groupDelayMs.y, color: '#1f5fa8' }]}
src/components/SpectrogramPlot.tsx →src/components/ImagePlot.tsxrenamed+27−27View file
@@ -1,25 +1,26 @@
11 import { useCallback, useMemo } from 'react'
22 import { useCanvas } from '../hooks/useCanvas'
3-import { AXIS, FONT, FRAME, formatTick, magma, niceTicks } from '../plot'
4-import type { SpectrogramData } from '../types'
3+import { AXIS, FONT, FRAME, formatTick, niceTicks } from '../plot'
4+import type { ImageData2D } from '../types'
55
66 type Props = {
7- spec: SpectrogramData
8- /** dB below the shared maximum shown as black */
9- dbRange: number
7+ image: ImageData2D
8+ /** maps a data value to rgb */
9+ color: (v: number) => [number, number, number]
1010 label: string
11- color: string
11+ labelColor: string
1212 /** time range to outline, e.g. the waveform view */
1313 highlight?: [number, number]
14+ className?: string
1415 }
1516
1617 const M = { l: 48, r: 12, t: 8, b: 26 }
1718
18-/** Draws a dB spectrogram as an image with time and frequency axes. */
19-export function SpectrogramPlot({ spec, dbRange, label, color, highlight }: Props) {
20- // Rasterize once per spectrogram, then scale to the canvas on every draw.
21- const image = useMemo(() => {
22- const { rows, cols, db } = spec
19+/** Draws a time-frequency image (spectrogram, phase field) with axes. */
20+export function ImagePlot({ image, color, label, labelColor, highlight, className }: Props) {
21+ // Rasterize once per image, then scale to the canvas on every draw.
22+ const raster = useMemo(() => {
23+ const { rows, cols, data } = image
2324 const off = document.createElement('canvas')
2425 off.width = cols
2526 off.height = rows
@@ -30,8 +31,7 @@ export function SpectrogramPlot({ spec, dbRange, label, color, highlight }: Prop
3031 for (let r = 0; r < rows; r++) {
3132 const yImg = rows - 1 - r // row 0 (0 Hz) at the bottom
3233 for (let c = 0; c < cols; c++) {
33- const v = db[r * cols + c]
34- const [R, G, B] = magma(1 + v / dbRange)
34+ const [R, G, B] = color(data[r * cols + c])
3535 const o = 4 * (yImg * cols + c)
3636 d[o] = R
3737 d[o + 1] = G
@@ -41,23 +41,23 @@ export function SpectrogramPlot({ spec, dbRange, label, color, highlight }: Prop
4141 }
4242 ictx.putImageData(img, 0, 0)
4343 return off
44- }, [spec, dbRange])
44+ }, [image, color])
4545
4646 const draw = useCallback(
4747 (ctx: CanvasRenderingContext2D, w: number, h: number) => {
4848 const pw = w - M.l - M.r
4949 const ph = h - M.t - M.b
5050 ctx.font = FONT
51- if (image) {
51+ if (raster) {
5252 ctx.imageSmoothingEnabled = false
53- ctx.drawImage(image, M.l, M.t, pw, ph)
53+ ctx.drawImage(raster, M.l, M.t, pw, ph)
5454 }
5555 ctx.strokeStyle = FRAME
5656 ctx.lineWidth = 1
5757 ctx.strokeRect(M.l + 0.5, M.t + 0.5, pw, ph)
5858 if (highlight) {
59- const x0 = M.l + (Math.max(0, highlight[0]) / spec.durationSec) * pw
60- const x1 = M.l + (Math.min(spec.durationSec, highlight[1]) / spec.durationSec) * pw
59+ const x0 = M.l + (Math.max(0, highlight[0]) / image.durationSec) * pw
60+ const x1 = M.l + (Math.min(image.durationSec, highlight[1]) / image.durationSec) * pw
6161 ctx.strokeStyle = 'rgba(255, 255, 255, 0.85)'
6262 ctx.setLineDash([4, 3])
6363 ctx.strokeRect(Math.round(x0) + 0.5, M.t + 0.5, Math.max(2, Math.round(x1 - x0)), ph)
@@ -67,16 +67,16 @@ export function SpectrogramPlot({ spec, dbRange, label, color, highlight }: Prop
6767 ctx.fillStyle = AXIS
6868 ctx.textAlign = 'right'
6969 ctx.textBaseline = 'middle'
70- for (const v of niceTicks(0, spec.fmax, 4)) {
71- if (v > spec.fmax) continue
72- const y = M.t + ph - (v / spec.fmax) * ph
70+ for (const v of niceTicks(0, image.fmax, 4)) {
71+ if (v > image.fmax) continue
72+ const y = M.t + ph - (v / image.fmax) * ph
7373 ctx.fillText(formatTick(v), M.l - 5, y)
7474 }
7575 ctx.textAlign = 'center'
7676 ctx.textBaseline = 'top'
77- for (const v of niceTicks(0, spec.durationSec, 6)) {
78- if (v > spec.durationSec) continue
79- const x = M.l + (v / spec.durationSec) * pw
77+ for (const v of niceTicks(0, image.durationSec, 6)) {
78+ if (v > image.durationSec) continue
79+ const x = M.l + (v / image.durationSec) * pw
8080 ctx.fillText(formatTick(v), x, M.t + ph + 4)
8181 }
8282 ctx.textBaseline = 'bottom'
@@ -88,14 +88,14 @@ export function SpectrogramPlot({ spec, dbRange, label, color, highlight }: Prop
8888 ctx.fillText('Hz', 0, 0)
8989 ctx.restore()
9090
91- ctx.fillStyle = color
91+ ctx.fillStyle = labelColor
9292 ctx.font = 'bold 12px system-ui, sans-serif'
9393 ctx.textAlign = 'left'
9494 ctx.textBaseline = 'top'
9595 ctx.fillText(label, M.l + 6, M.t + 5)
9696 },
97- [image, spec, label, color, highlight],
97+ [raster, image, label, labelColor, highlight],
9898 )
9999 const ref = useCanvas(draw)
100- return <canvas ref={ref} className="plot spectrogram-plot" />
100+ return <canvas ref={ref} className={`plot ${className ?? 'image-plot'}`} />
101101 }
src/components/PlaybackPanel.tsxmodified+42−12View file
@@ -1,8 +1,8 @@
1-import { useCallback, useEffect, useState } from 'react'
2-import type { View } from '../hooks/useZoomPan'
3-import { formatSeconds } from '../plot'
1+import { useCallback, useEffect, useRef, useState } from 'react'
2+import { clampView, type View } from '../hooks/useZoomPan'
3+import { dbColor, formatSeconds } from '../plot'
44 import type { RenderResult, Source } from '../types'
5-import { SpectrogramPlot } from './SpectrogramPlot'
5+import { ImagePlot } from './ImagePlot'
66 import { WaveformPlot } from './WaveformPlot'
77
88 type Props = {
@@ -17,17 +17,40 @@ type Props = {
1717 const COLOR_A = '#1f5fa8'
1818 const COLOR_B = '#d9641e'
1919 const DB_RANGE = 70
20+const specColor = dbColor(DB_RANGE)
2021
2122 /** Play A and B, and compare their waveforms and spectrograms. */
2223 export function PlaybackPanel({ source, result, stale, playing, play, stop }: Props) {
2324 const fs = source.fs
2425 const duration = source.samples.length / fs
26+ const minSpan = 16 / fs
2527 const [view, setView] = useState<View>({ start: 0, span: duration })
2628 // A new clip resets the view; the same clip re-rendered with a new map keeps it.
2729 useEffect(() => {
2830 setView({ start: 0, span: source.samples.length / source.fs })
2931 }, [source])
3032 const onView = useCallback((v: View) => setView(v), [])
33+ const zoomBy = (factor: number) => {
34+ const c = view.start + view.span / 2
35+ const span = view.span * factor
36+ setView(clampView({ start: c - span / 2, span }, duration, minSpan))
37+ }
38+
39+ // A native horizontal scrollbar pans the view: its content is as many times wider than the plot as the clip is longer than the view.
40+ const scrollRef = useRef<HTMLDivElement>(null)
41+ useEffect(() => {
42+ const el = scrollRef.current
43+ if (!el) return
44+ const target = (view.start / duration) * el.scrollWidth
45+ if (Math.abs(el.scrollLeft - target) > 1) el.scrollLeft = target
46+ }, [view, duration])
47+ const onScroll = () => {
48+ const el = scrollRef.current
49+ if (!el || el.scrollWidth <= el.clientWidth) return
50+ const start = (el.scrollLeft / el.scrollWidth) * duration
51+ if (Math.abs(start - view.start) <= duration / el.scrollWidth) return // our own programmatic scroll
52+ setView(clampView({ start, span: view.span }, duration, minSpan))
53+ }
3154
3255 let yMax = 0
3356 for (let i = 0; i < source.samples.length; i++) yMax = Math.max(yMax, Math.abs(source.samples[i]))
@@ -43,7 +66,7 @@ export function PlaybackPanel({ source, result, stale, playing, play, stop }: Pr
4366 <p className="hint">
4467 A is the source; B has the phase map applied. Both are played with the same gain{stale ? ' (B is being re-rendered)' : ''}. The
4568 waveforms differ by {Math.round(100 * result.relDiff)}% of A's RMS, while the spectrograms below, computed at the analysis
46- resolution, should look the same.
69+ resolution, should look the same{result.timeVarying ? ' when the map changes slowly compared with the window' : ''}.
4770 </p>
4871 <div className="button-row">
4972 <button className="primary" onClick={() => play('A', ['A'])} disabled={busy || stale}>
@@ -65,22 +88,29 @@ export function PlaybackPanel({ source, result, stale, playing, play, stop }: Pr
6588 </div>
6689
6790 <div className="view-controls">
68- <span className="info">Scroll to zoom, drag to pan, double-click to reset. Both waveforms share the view.</span>
91+ <span className="info">Drag across a waveform to zoom into that range. The scroll wheel also zooms; the bar below the waveforms pans.</span>
6992 <span className="spacer" />
70- <span className="info">
71- {zoomed ? `showing ${formatSeconds(view.span)} from ${formatSeconds(view.start)}` : `full clip, ${formatSeconds(duration)}`}
72- </span>
93+ <span className="info">{zoomed ? `${formatSeconds(view.span)} from ${formatSeconds(view.start)}` : `full clip, ${formatSeconds(duration)}`}</span>
94+ <button className="small" onClick={() => zoomBy(0.5)} disabled={view.span <= minSpan}>
95+ Zoom in
96+ </button>
97+ <button className="small" onClick={() => zoomBy(2)} disabled={!zoomed}>
98+ Zoom out
99+ </button>
73100 <button className="small" onClick={() => setView({ start: 0, span: duration })} disabled={!zoomed}>
74- Reset zoom
101+ Full clip
75102 </button>
76103 </div>
77104 <WaveformPlot samples={source.samples} fs={fs} view={view} onView={onView} yMax={yMax} label="A" color={COLOR_A} />
78105 <WaveformPlot samples={result.y} fs={fs} view={view} onView={onView} yMax={yMax} label="B" color={COLOR_B} />
106+ <div className="scrollbar" ref={scrollRef} onScroll={onScroll} aria-label="pan">
107+ <div style={{ width: `${(100 * duration) / view.span}%` }} />
108+ </div>
79109
80110 <div className="plot-label">Spectrograms at the analysis resolution, {DB_RANGE} dB range, shared scale{zoomed ? '; the box marks the waveform view' : ''}</div>
81111 <div className="spectrograms">
82- <SpectrogramPlot spec={result.specA} dbRange={DB_RANGE} label="A" color="#fff" highlight={highlight} />
83- <SpectrogramPlot spec={result.specB} dbRange={DB_RANGE} label="B" color="#fff" highlight={highlight} />
112+ <ImagePlot image={result.specA} color={specColor} label="A" labelColor="#fff" highlight={highlight} />
113+ <ImagePlot image={result.specB} color={specColor} label="B" labelColor="#fff" highlight={highlight} />
84114 </div>
85115 </section>
86116 )
src/components/WaveformPlot.tsxmodified+10−4View file
@@ -1,4 +1,4 @@
1-import { useCallback } from 'react'
1+import { useCallback, useRef } from 'react'
22 import { useCanvas } from '../hooks/useCanvas'
33 import { useZoomPan, type View } from '../hooks/useZoomPan'
44 import { AXIS, FONT, FRAME, GRID, formatTick, niceTicks } from '../plot'
@@ -19,7 +19,7 @@ const MARGINS = { left: M.l, right: M.r }
1919 const MIN_SPAN_SAMPLES = 16
2020
2121 /**
22- * One waveform with wheel zoom, drag pan, and double-click reset. Draws a
22+ * One waveform with drag-to-zoom, wheel zoom, and double-click reset. Draws a
2323 * min/max envelope when there are many samples per pixel, a line otherwise.
2424 */
2525 export function WaveformPlot({ samples, fs, view, onView, yMax, label, color }: Props) {
@@ -125,6 +125,12 @@ export function WaveformPlot({ samples, fs, view, onView, yMax, label, color }:
125125 [samples, fs, startSec, spanSec, yMax, label, color],
126126 )
127127 const ref = useCanvas(draw)
128- useZoomPan(ref, view, duration, MIN_SPAN_SAMPLES / fs, MARGINS, onView)
129- return <canvas ref={ref} className="plot waveform-plot" />
128+ const overlay = useRef<HTMLDivElement>(null)
129+ useZoomPan(ref, overlay, view, duration, MIN_SPAN_SAMPLES / fs, MARGINS, onView)
130+ return (
131+ <div className="plot-wrap">
132+ <canvas ref={ref} className="plot waveform-plot" />
133+ <div ref={overlay} className="select-overlay" style={{ top: M.t, bottom: M.b, display: 'none' }} />
134+ </div>
135+ )
130136 }
src/hooks/useZoomPan.tsmodified+50−17View file
@@ -10,32 +10,59 @@ export function clampView(v: View, duration: number, minSpan: number): View {
1010 }
1111
1212 /**
13- * Mouse-wheel zoom about the cursor, drag to pan, double-click to reset, on a
14- * canvas whose plot area spans horizontally from `left` to `width - right`.
15- * The listeners are registered once and read the latest view through a ref,
16- * so a drag survives the re-renders it causes.
13+ * Zoom interactions on a canvas whose plot area spans horizontally from
14+ * `left` to `width - right`: drag across a range to zoom into it (the range is
15+ * shown in the overlay element while dragging), wheel to zoom about the cursor
16+ * or, with horizontal scrolling, to pan, and double-click to reset. The
17+ * listeners are registered once and read the latest view through a ref.
1718 */
18-export function useZoomPan(ref: RefObject<HTMLCanvasElement | null>, view: View, duration: number, minSpan: number, margins: { left: number; right: number }, onView: (v: View) => void): void {
19+export function useZoomPan(
20+ ref: RefObject<HTMLCanvasElement | null>,
21+ overlayRef: RefObject<HTMLDivElement | null>,
22+ view: View,
23+ duration: number,
24+ minSpan: number,
25+ margins: { left: number; right: number },
26+ onView: (v: View) => void,
27+): void {
1928 const viewRef = useRef(view)
2029 viewRef.current = view
2130 const onViewRef = useRef(onView)
2231 onViewRef.current = onView
23- const drag = useRef<{ x: number; start: number } | null>(null)
32+ const drag = useRef<{ x0: number; x1: number } | null>(null)
2433
2534 useEffect(() => {
2635 const el = ref.current
2736 if (!el) return
2837 const plotWidth = () => el.clientWidth - margins.left - margins.right
38+ const plotX = (clientX: number) => clientX - el.getBoundingClientRect().left - margins.left
2939 const timeAt = (clientX: number) => {
3040 const v = viewRef.current
31- const rect = el.getBoundingClientRect()
32- return v.start + ((clientX - rect.left - margins.left) / plotWidth()) * v.span
41+ return v.start + (plotX(clientX) / plotWidth()) * v.span
3342 }
3443 const emit = (v: View) => onViewRef.current(clampView(v, duration, minSpan))
44+ const showSelection = () => {
45+ const ov = overlayRef.current
46+ const d = drag.current
47+ if (!ov) return
48+ if (!d) {
49+ ov.style.display = 'none'
50+ return
51+ }
52+ const a = Math.max(0, Math.min(d.x0, d.x1))
53+ const b = Math.min(plotWidth(), Math.max(d.x0, d.x1))
54+ ov.style.display = 'block'
55+ ov.style.left = `${margins.left + a}px`
56+ ov.style.width = `${Math.max(1, b - a)}px`
57+ }
3558
3659 const onWheel = (e: WheelEvent) => {
3760 e.preventDefault()
3861 const v = viewRef.current
62+ if (Math.abs(e.deltaX) > Math.abs(e.deltaY)) {
63+ emit({ start: v.start + (e.deltaX / plotWidth()) * v.span, span: v.span })
64+ return
65+ }
3966 const factor = Math.exp(e.deltaY * 0.0015)
4067 const t = timeAt(e.clientX)
4168 const span = v.span * factor
@@ -43,26 +70,32 @@ export function useZoomPan(ref: RefObject<HTMLCanvasElement | null>, view: View,
4370 }
4471 const onDown = (e: PointerEvent) => {
4572 if (e.button !== 0) return
46- drag.current = { x: e.clientX, start: viewRef.current.start }
73+ const x = plotX(e.clientX)
74+ drag.current = { x0: x, x1: x }
4775 el.setPointerCapture(e.pointerId)
48- el.style.cursor = 'grabbing'
76+ showSelection()
4977 }
5078 const onMove = (e: PointerEvent) => {
5179 const d = drag.current
5280 if (!d) return
53- const v = viewRef.current
54- const dt = ((e.clientX - d.x) / plotWidth()) * v.span
55- emit({ start: d.start - dt, span: v.span })
81+ d.x1 = plotX(e.clientX)
82+ showSelection()
5683 }
5784 const onUp = (e: PointerEvent) => {
58- if (!drag.current) return
85+ const d = drag.current
86+ if (!d) return
5987 drag.current = null
88+ showSelection()
6089 if (el.hasPointerCapture(e.pointerId)) el.releasePointerCapture(e.pointerId)
61- el.style.cursor = 'grab'
90+ const a = Math.min(d.x0, d.x1)
91+ const b = Math.max(d.x0, d.x1)
92+ if (b - a < 4) return
93+ const v = viewRef.current
94+ const pw = plotWidth()
95+ emit({ start: v.start + (a / pw) * v.span, span: ((b - a) / pw) * v.span })
6296 }
6397 const onDouble = () => emit({ start: 0, span: duration })
6498
65- el.style.cursor = 'grab'
6699 el.style.touchAction = 'none'
67100 el.addEventListener('wheel', onWheel, { passive: false })
68101 el.addEventListener('pointerdown', onDown)
@@ -78,5 +111,5 @@ export function useZoomPan(ref: RefObject<HTMLCanvasElement | null>, view: View,
78111 el.removeEventListener('pointercancel', onUp)
79112 el.removeEventListener('dblclick', onDouble)
80113 }
81- }, [ref, duration, minSpan, margins])
114+ }, [ref, overlayRef, duration, minSpan, margins])
82115 }
src/index.cssmodified+46−1View file
@@ -79,10 +79,55 @@ body {
7979 margin-top: 6px;
8080 }
8181
82-.spectrogram-plot {
82+.image-plot {
8383 height: 240px;
8484 }
8585
86+.response-plots .image-plot {
87+ height: 170px;
88+}
89+
90+.plot-wrap {
91+ position: relative;
92+}
93+
94+.waveform-plot {
95+ cursor: crosshair;
96+}
97+
98+.select-overlay {
99+ position: absolute;
100+ background: rgba(31, 95, 168, 0.18);
101+ border-left: 1px solid var(--accent);
102+ border-right: 1px solid var(--accent);
103+ pointer-events: none;
104+}
105+
106+.scrollbar {
107+ overflow-x: auto;
108+ overflow-y: hidden;
109+ height: 18px;
110+ margin: 0 12px 4px 48px;
111+}
112+
113+.scrollbar > div {
114+ height: 1px;
115+}
116+
117+.phase-legend {
118+ display: inline-block;
119+ width: 90px;
120+ height: 9px;
121+ vertical-align: middle;
122+ margin: 0 4px 0 8px;
123+ border-radius: 2px;
124+ background: linear-gradient(90deg, hsl(0 60% 55%), hsl(60 60% 55%), hsl(120 60% 55%), hsl(180 60% 55%), hsl(240 60% 55%), hsl(300 60% 55%), hsl(360 60% 55%));
125+}
126+
127+.global-controls input[type='checkbox'] {
128+ accent-color: var(--accent);
129+}
130+
86131 .spectrograms {
87132 display: grid;
88133 grid-template-columns: repeat(auto-fit, minmax(420px, 1fr));
src/phasemap.tsmodified+128−34View file
@@ -1,8 +1,10 @@
11 /**
22 * Phase maps: the angle phi_k applied to band k, as a function of the band
3- * center frequency. The maps must vary smoothly from one band to the next,
4- * otherwise neighbouring bands partly cancel where they overlap and the
5- * magnitude response acquires notches.
3+ * center frequency and, for the random map, optionally of time. Maps must
4+ * vary smoothly from one band to the next, otherwise neighbouring bands
5+ * partly cancel where they overlap and the magnitude response acquires
6+ * notches. A map that varies in time is frequency modulation of each band:
7+ * the instantaneous frequency shifts by (d phi / dt) / 2 pi.
68 */
79
810 export type PhaseMapSpec =
@@ -10,6 +12,8 @@ export type PhaseMapSpec =
1012 kind: 'random'
1113 /** Gaussian smoothing width across bands */
1214 smoothBands: number
15+ /** Gaussian smoothing width in time; null means the map does not change in time */
16+ timeScaleMs: number | null
1317 seed: number
1418 /** multiplies the map; 1 gives the spread of uniform random phases */
1519 strength: number
@@ -25,7 +29,9 @@ export type PhaseMapSpec =
2529 delayAtNyquistMs: number
2630 }
2731
28-export const defaultPhaseMap: PhaseMapSpec = { kind: 'random', smoothBands: 8, seed: 1, strength: 1 }
32+export const defaultPhaseMap: PhaseMapSpec = { kind: 'random', smoothBands: 8, timeScaleMs: null, seed: 1, strength: 1 }
33+
34+export const isTimeVarying = (spec: PhaseMapSpec): boolean => spec.kind === 'random' && spec.timeScaleMs !== null
2935
3036 /** Small seeded generator (mulberry32) so a map can be reproduced from its seed. */
3137 export function mulberry32(seed: number): () => number {
@@ -39,47 +45,135 @@ export function mulberry32(seed: number): () => number {
3945 }
4046 }
4147
42-const UNIFORM_STD = Math.PI / Math.sqrt(3)
48+/** Standard deviation of a phase drawn uniformly from (-pi, pi]. */
49+export const UNIFORM_STD = Math.PI / Math.sqrt(3)
4350
44-function smoothedRandom(n: number, smoothBands: number, seed: number, strength: number): Float64Array {
45- const rnd = mulberry32(seed)
46- const raw = new Float64Array(n)
47- for (let i = 0; i < n; i++) raw[i] = (rnd() * 2 - 1) * Math.PI
48- if (smoothBands <= 0) return raw.map((v) => v * strength)
49- const radius = Math.ceil(4 * smoothBands)
51+/**
52+ * A random phase field on a coarse grid of (band, time) cells, smoothed with
53+ * Gaussians in both directions and interpolated bilinearly. Row r holds band
54+ * index r * dK and column c holds time c * dT.
55+ */
56+export type PhaseField = {
57+ values: Float32Array
58+ nK: number
59+ nT: number
60+ dK: number
61+ dT: number
62+ /** largest |d phi / dt| / 2 pi over the grid, in Hz: the largest frequency shift the field imposes */
63+ maxRateHz: number
64+}
65+
66+/** Gaussian smoothing along one axis of a row-major nK x nT grid, with the kernel renormalized at the edges. */
67+function smoothAxis(v: Float32Array, nK: number, nT: number, alongTime: boolean, sigmaCells: number): void {
68+ if (sigmaCells < 0.3) return
69+ const radius = Math.ceil(4 * sigmaCells)
5070 const kern = new Float64Array(2 * radius + 1)
51- for (let j = -radius; j <= radius; j++) kern[j + radius] = Math.exp(-0.5 * (j / smoothBands) ** 2)
52- const out = new Float64Array(n)
53- for (let i = 0; i < n; i++) {
54- let s = 0
55- let wsum = 0
56- for (let j = -radius; j <= radius; j++) {
57- const m = i + j
58- if (m < 0 || m >= n) continue
59- const w = kern[j + radius]
60- s += w * raw[m]
61- wsum += w
71+ for (let j = -radius; j <= radius; j++) kern[j + radius] = Math.exp(-0.5 * (j / sigmaCells) ** 2)
72+ const len = alongTime ? nT : nK
73+ const lines = alongTime ? nK : nT
74+ const stride = alongTime ? 1 : nT
75+ const lineStride = alongTime ? nT : 1
76+ const tmp = new Float64Array(len)
77+ for (let l = 0; l < lines; l++) {
78+ const base = l * lineStride
79+ for (let i = 0; i < len; i++) {
80+ let s = 0
81+ let wsum = 0
82+ const j0 = Math.max(-radius, -i)
83+ const j1 = Math.min(radius, len - 1 - i)
84+ for (let j = j0; j <= j1; j++) {
85+ const w = kern[j + radius]
86+ s += w * v[base + (i + j) * stride]
87+ wsum += w
88+ }
89+ tmp[i] = s / wsum
6290 }
63- out[i] = s / wsum
91+ for (let i = 0; i < len; i++) v[base + i * stride] = tmp[i]
6492 }
65- // Smoothing shrinks the spread; rescale so the map keeps the std of uniform phases.
93+}
94+
95+export function randomField(nBands: number, durationSec: number, smoothBands: number, timeScaleSec: number | null, seed: number, strength: number): PhaseField {
96+ const sb = Math.max(0, smoothBands)
97+ // grid spacing: about a third of each smoothing scale, bounded so the grid stays small
98+ const dK = Math.max(1, Math.floor(sb / 3), Math.ceil(nBands / 1024))
99+ const nK = Math.ceil(Math.max(0, nBands - 1) / dK) + 1
100+ let nT = 1
101+ let dT = durationSec
102+ if (timeScaleSec !== null) {
103+ dT = Math.max(timeScaleSec / 3, durationSec / 4096)
104+ nT = Math.ceil(durationSec / dT) + 1
105+ }
106+ const v = new Float32Array(nK * nT)
107+ const rnd = mulberry32(seed)
108+ for (let i = 0; i < v.length; i++) v[i] = (rnd() * 2 - 1) * Math.PI
109+ smoothAxis(v, nK, nT, false, sb / dK)
110+ if (timeScaleSec !== null) smoothAxis(v, nK, nT, true, timeScaleSec / dT)
111+ // rescale so the field keeps the spread of uniform phases, times the strength
66112 let mean = 0
67- for (let i = 0; i < n; i++) mean += out[i]
68- mean /= n
69- let v = 0
70- for (let i = 0; i < n; i++) v += (out[i] - mean) ** 2
71- const std = Math.sqrt(v / Math.max(1, n - 1))
113+ for (let i = 0; i < v.length; i++) mean += v[i]
114+ mean /= v.length
115+ let s2 = 0
116+ for (let i = 0; i < v.length; i++) s2 += (v[i] - mean) ** 2
117+ const std = Math.sqrt(s2 / Math.max(1, v.length - 1))
72118 const scale = std > 0 ? (strength * UNIFORM_STD) / std : 0
73- for (let i = 0; i < n; i++) out[i] = (out[i] - mean) * scale
119+ for (let i = 0; i < v.length; i++) v[i] = (v[i] - mean) * scale
120+ let maxRate = 0
121+ if (nT > 1) {
122+ for (let r = 0; r < nK; r++) {
123+ for (let c = 0; c + 1 < nT; c++) {
124+ const d = Math.abs(v[r * nT + c + 1] - v[r * nT + c])
125+ if (d > maxRate) maxRate = d
126+ }
127+ }
128+ maxRate /= 2 * Math.PI * dT
129+ }
130+ return { values: v, nK, nT, dK, dT, maxRateHz: maxRate }
131+}
132+
133+/** The field for band k as a function of the grid time index (interpolated between grid rows). */
134+export function fieldRow(field: PhaseField, k: number): Float64Array {
135+ const { values, nK, nT, dK } = field
136+ const r = Math.min(nK - 1, k / dK)
137+ const r0 = Math.floor(r)
138+ const r1 = Math.min(nK - 1, r0 + 1)
139+ const w = r - r0
140+ const out = new Float64Array(nT)
141+ for (let c = 0; c < nT; c++) out[c] = (1 - w) * values[r0 * nT + c] + w * values[r1 * nT + c]
74142 return out
75143 }
76144
77-/** phi_k for every band center. */
145+/** Linear interpolation of a grid row at arbitrary times (seconds). */
146+export function sampleRow(row: Float64Array, dT: number, times: ArrayLike<number>, out: Float64Array): Float64Array {
147+ const nT = row.length
148+ for (let i = 0; i < times.length; i++) {
149+ if (nT === 1) {
150+ out[i] = row[0]
151+ continue
152+ }
153+ const c = Math.min(nT - 1, Math.max(0, times[i] / dT))
154+ const c0 = Math.floor(c)
155+ const c1 = Math.min(nT - 1, c0 + 1)
156+ const w = c - c0
157+ out[i] = (1 - w) * row[c0] + w * row[c1]
158+ }
159+ return out
160+}
161+
162+export function fieldAt(field: PhaseField, k: number, tSec: number): number {
163+ const row = fieldRow(field, k)
164+ return sampleRow(row, field.dT, [tSec], new Float64Array(1))[0]
165+}
166+
167+/** phi_k for every band center, for maps that do not vary in time. */
78168 export function phaseMap(spec: PhaseMapSpec, centers: Float64Array, fs: number): Float64Array {
79169 const n = centers.length
80170 switch (spec.kind) {
81- case 'random':
82- return smoothedRandom(n, spec.smoothBands, spec.seed, spec.strength)
171+ case 'random': {
172+ const field = randomField(n, 1, spec.smoothBands, null, spec.seed, spec.strength)
173+ const out = new Float64Array(n)
174+ for (let k = 0; k < n; k++) out[k] = fieldAt(field, k, 0)
175+ return out
176+ }
83177 case 'delay': {
84178 const tau = spec.delayMs / 1000
85179 return centers.map((f) => -2 * Math.PI * f * tau)
@@ -94,7 +188,7 @@ export function phaseMap(spec: PhaseMapSpec, centers: Float64Array, fs: number):
94188 }
95189
96190 /** Largest band-to-band slope of the map, as a group delay in seconds. Used to size the FFT padding. */
97-export function maxGroupDelaySec(phi: Float64Array, df: number): number {
191+export function maxGroupDelaySec(phi: ArrayLike<number>, df: number): number {
98192 let m = 0
99193 for (let i = 1; i < phi.length; i++) {
100194 let d = Math.abs(phi[i] - phi[i - 1])
src/plot.tsmodified+29−0View file
@@ -48,3 +48,32 @@ export const FONT = '11px system-ui, sans-serif'
4848 export const AXIS = '#666'
4949 export const GRID = '#e6e8eb'
5050 export const FRAME = '#cfd3d8'
51+
52+function hslToRgb(h: number, s: number, l: number): [number, number, number] {
53+ const c = (1 - Math.abs(2 * l - 1)) * s
54+ const hp = (((h % 360) + 360) % 360) / 60
55+ const x = c * (1 - Math.abs((hp % 2) - 1))
56+ let r = 0
57+ let g = 0
58+ let b = 0
59+ if (hp < 1) [r, g, b] = [c, x, 0]
60+ else if (hp < 2) [r, g, b] = [x, c, 0]
61+ else if (hp < 3) [r, g, b] = [0, c, x]
62+ else if (hp < 4) [r, g, b] = [0, x, c]
63+ else if (hp < 5) [r, g, b] = [x, 0, c]
64+ else [r, g, b] = [c, 0, x]
65+ const m = l - c / 2
66+ return [255 * (r + m), 255 * (g + m), 255 * (b + m)]
67+}
68+
69+/** Cyclic colormap for an angle in radians: hue goes once around the circle. */
70+export function phaseColor(phi: number): [number, number, number] {
71+ const h = ((phi + Math.PI) / (2 * Math.PI)) * 360
72+ return hslToRgb(h, 0.6, 0.55)
73+}
74+
75+/** Spectrogram colormap: dB value (<= 0) to rgb, black at -range. */
76+export const dbColor =
77+ (range: number) =>
78+ (v: number): [number, number, number] =>
79+ magma(1 + v / range)
src/types.tsmodified+13−11View file
@@ -17,15 +17,15 @@ export type RenderRequest = {
1717 fmax: number
1818 maxCols: number
1919 maxRows: number
20- /** number of points in the response and map curves */
20+ /** number of points in the response curves */
2121 curvePoints: number
2222 }
2323
2424 export type Curve = { f: Float32Array; y: Float32Array }
2525
26-export type SpectrogramData = {
27- /** rows x cols in dB relative to a shared maximum; row 0 is 0 Hz */
28- db: Float32Array
26+/** A rows x cols image over time (columns) and frequency (rows, row 0 at 0 Hz). */
27+export type ImageData2D = {
28+ data: Float32Array
2929 rows: number
3030 cols: number
3131 durationSec: number
@@ -42,17 +42,19 @@ export type RenderResult = {
4242 relDiff: number
4343 bands: number
4444 sigmaF: number
45- /** the requested map at the band centers, wrapped to (-pi, pi] */
46- map: Curve
47- /** magnitude of the effective filter in dB */
45+ /** whether the map changed over the clip */
46+ timeVarying: boolean
47+ /** the phase map over time and frequency, in radians wrapped to (-pi, pi] */
48+ phaseField: ImageData2D
49+ /** magnitude in dB of the filter equivalent to the map (at the middle of the clip if it varies) */
4850 magDb: Curve
49- /** group delay of the effective filter in ms */
51+ /** group delay in ms of that filter */
5052 groupDelayMs: Curve
51- /** largest |group delay| in ms over the plotted range */
5253 maxGroupDelayMs: number
5354 /** largest |magnitude response| in dB over the plotted range; large values mean the map is too rough for the window */
5455 maxMagDeviationDb: number
55- specA: SpectrogramData
56- specB: SpectrogramData
56+ /** spectrograms in dB relative to a shared maximum */
57+ specA: ImageData2D
58+ specB: ImageData2D
5759 elapsedMs: number
5860 }
src/worker.tsmodified+67−24View file
@@ -1,11 +1,12 @@
11 /**
22 * Web Worker: renders the phase-rotated sound B from the source A, and the
3- * curves and spectrograms shown in the app. One message in, one message out.
3+ * curves and images shown in the app. One message in, one message out.
44 */
5-import { applyFilter, bandCenters, groupDelay, sigmaF, spectrogram, type Spectrogram } from './bank'
6-import { maxGroupDelaySec, phaseMap, wrapAngle } from './phasemap'
5+import { applyFilter, applyTimeVarying, bandCenters, effectiveFilter, groupDelay, sigmaF, spectrogram, type Spectrogram } from './bank'
6+import { nextPow2 } from './fft'
7+import { fieldAt, fieldRow, isTimeVarying, maxGroupDelaySec, phaseMap, randomField, sampleRow, wrapAngle } from './phasemap'
78 import { peak, rms } from './sources'
8-import type { Curve, RenderRequest, RenderResult, SpectrogramData } from './types'
9+import type { Curve, ImageData2D, RenderRequest, RenderResult } from './types'
910
1011 /** Sample an array defined on FFT bins onto `points` evenly spaced frequencies in [0, fmax]. */
1112 function sampleCurve(values: ArrayLike<number>, binHz: number, fmax: number, points: number): Curve {
@@ -20,21 +21,73 @@ function sampleCurve(values: ArrayLike<number>, binHz: number, fmax: number, poi
2021 return { f, y }
2122 }
2223
23-function toDb(s: Spectrogram, ref: number): SpectrogramData {
24+function toDb(s: Spectrogram, ref: number): ImageData2D {
2425 const db = new Float32Array(s.data.length)
2526 for (let i = 0; i < db.length; i++) db[i] = 20 * Math.log10(Math.max(s.data[i], 1e-12) / ref)
26- return { db, rows: s.rows, cols: s.cols, durationSec: s.durationSec, fmax: s.fmax }
27+ return { data: db, rows: s.rows, cols: s.cols, durationSec: s.durationSec, fmax: s.fmax }
2728 }
2829
30+const FIELD_ROWS = 160
31+const FIELD_COLS = 320
32+
2933 function render(req: RenderRequest): RenderResult {
3034 const t0 = performance.now()
3135 const { samples: x, fs, fmax } = req
36+ const n = x.length
3237 const sigma = req.sigmaMs / 1000
3338 const centers = bandCenters(fs, sigma)
3439 const df = sigmaF(sigma)
35- const phi = phaseMap(req.map, centers, fs)
36- const padSec = 8 * sigma + maxGroupDelaySec(phi, df) + 0.02
37- const { y, g } = applyFilter(x, fs, sigma, phi, centers, padSec)
40+ const kMax = Math.min(centers.length - 1, Math.floor(fmax / df))
41+ const timeVarying = isTimeVarying(req.map)
42+
43+ let y: Float32Array
44+ let phiRef: Float64Array // the map used for the response plots
45+ let field: ImageData2D
46+ let L: number
47+ if (timeVarying && req.map.kind === 'random') {
48+ // group delay of a smoothed random map stays below the window width, so the padding need not depend on it
49+ const padSec = 10 * sigma + 0.05
50+ const pad = Math.ceil(padSec * fs)
51+ L = nextPow2(n + 2 * pad)
52+ const rf = randomField(centers.length, L / fs, req.map.smoothBands, req.map.timeScaleMs! / 1000, req.map.seed, req.map.strength)
53+ const rows = new Map<number, Float64Array>()
54+ const rowFor = (k: number) => {
55+ let r = rows.get(k)
56+ if (!r) {
57+ r = fieldRow(rf, k)
58+ rows.set(k, r)
59+ }
60+ return r
61+ }
62+ y = applyTimeVarying(x, fs, sigma, centers, padSec, (k, times) => sampleRow(rowFor(k), rf.dT, times, new Float64Array(times.length)), rf.maxRateHz)
63+ const tMid = (pad + n / 2) / fs
64+ phiRef = new Float64Array(centers.length)
65+ for (let k = 0; k < centers.length; k++) phiRef[k] = fieldAt(rf, k, tMid)
66+ // the field over the clip, for display
67+ const fr = Math.min(FIELD_ROWS, kMax + 1)
68+ const data = new Float32Array(fr * FIELD_COLS)
69+ const times = new Float64Array(FIELD_COLS)
70+ for (let c = 0; c < FIELD_COLS; c++) times[c] = (pad + ((c + 0.5) * n) / FIELD_COLS) / fs
71+ const tmp = new Float64Array(FIELD_COLS)
72+ for (let r = 0; r < fr; r++) {
73+ const k = Math.round((r * kMax) / Math.max(1, fr - 1))
74+ sampleRow(rowFor(k), rf.dT, times, tmp)
75+ for (let c = 0; c < FIELD_COLS; c++) data[r * FIELD_COLS + c] = wrapAngle(tmp[c])
76+ }
77+ field = { data, rows: fr, cols: FIELD_COLS, durationSec: n / fs, fmax: kMax * df }
78+ } else {
79+ const phi = phaseMap(req.map, centers, fs)
80+ const padSec = 8 * sigma + maxGroupDelaySec(phi, df) + 0.02
81+ const out = applyFilter(x, fs, sigma, phi, centers, padSec)
82+ y = out.y
83+ phiRef = phi
84+ L = Math.round(fs / out.g.binHz)
85+ const fr = Math.min(FIELD_ROWS, kMax + 1)
86+ const data = new Float32Array(fr)
87+ for (let r = 0; r < fr; r++) data[r] = wrapAngle(phi[Math.round((r * kMax) / Math.max(1, fr - 1))])
88+ field = { data, rows: fr, cols: 1, durationSec: n / fs, fmax: kMax * df }
89+ }
90+ const g = effectiveFilter(centers, phiRef, sigma, L, fs)
3891
3992 // shared gain: target RMS from A, but never let either sound clip
4093 const pk = Math.max(peak(x), peak(y))
@@ -43,7 +96,7 @@ function render(req: RenderRequest): RenderResult {
4396
4497 let d = 0
4598 let e = 0
46- for (let i = 0; i < x.length; i++) {
99+ for (let i = 0; i < n; i++) {
47100 d += (y[i] - x[i]) ** 2
48101 e += x[i] * x[i]
49102 }
@@ -55,7 +108,6 @@ function render(req: RenderRequest): RenderResult {
55108 const tau = groupDelay(g)
56109 const tauMs = new Float64Array(nb)
57110 for (let i = 0; i < nb; i++) tauMs[i] = tau[i] * 1000
58- const gdCurve = sampleCurve(tauMs, g.binHz, fmax, req.curvePoints)
59111 let maxGd = 0
60112 let maxMagDev = 0
61113 const lastBin = Math.min(nb - 1, Math.round(fmax / g.binHz))
@@ -64,16 +116,6 @@ function render(req: RenderRequest): RenderResult {
64116 maxMagDev = Math.max(maxMagDev, Math.abs(magDb[i]))
65117 }
66118
67- // the map itself, wrapped, at the band centers up to fmax (subsampled if there are many bands)
68- const kMax = Math.min(centers.length - 1, Math.floor(fmax / df))
69- const step = Math.max(1, Math.ceil((kMax + 1) / req.curvePoints))
70- const mapF: number[] = []
71- const mapY: number[] = []
72- for (let k = 0; k <= kMax; k += step) {
73- mapF.push(centers[k])
74- mapY.push(wrapAngle(phi[k]))
75- }
76-
77119 const sA = spectrogram(x, fs, sigma, fmax, req.maxCols, req.maxRows)
78120 const sB = spectrogram(y, fs, sigma, fmax, req.maxCols, req.maxRows)
79121 let ref = 1e-12
@@ -87,9 +129,10 @@ function render(req: RenderRequest): RenderResult {
87129 relDiff,
88130 bands: centers.length,
89131 sigmaF: df,
90- map: { f: Float32Array.from(mapF), y: Float32Array.from(mapY) },
132+ timeVarying,
133+ phaseField: field,
91134 magDb: sampleCurve(magDb, g.binHz, fmax, req.curvePoints),
92- groupDelayMs: gdCurve,
135+ groupDelayMs: sampleCurve(tauMs, g.binHz, fmax, req.curvePoints),
93136 maxGroupDelayMs: maxGd,
94137 maxMagDeviationDb: maxMagDev,
95138 specA: toDb(sA, ref),
@@ -100,6 +143,6 @@ function render(req: RenderRequest): RenderResult {
100143
101144 self.onmessage = (e: MessageEvent<RenderRequest>) => {
102145 const res = render(e.data)
103- const transfer = [res.y.buffer, res.specA.db.buffer, res.specB.db.buffer]
146+ const transfer = [res.y.buffer, res.specA.data.buffer, res.specB.data.buffer, res.phaseField.data.buffer]
104147 ;(self as unknown as Worker).postMessage(res, transfer)
105148 }