1/**
2 * The solved field as the views need it, and how a complex value becomes a
3 * color: the pressure at phase t, Re(u exp(-i omega t)) = Re u cos wt + Im u sin wt,
4 * on a symmetric diverging scale, or the level 20 log10(|u| / ref) in dB.
5 */
6import type { Slice } from '../problem.ts';
8export type Mode = 'pressure' | 'level';
10export interface SliceField {
11 slice: Slice;
12 /** interleaved complex, NaN outside the room */
13 u: Float64Array;
14}
16export interface Field {
17 k: number;
18 beta: number;
19 x0: [number, number, number];
20 a: number;
21 p: number;
22 ngrid: number;
23 grid: Float64Array;
24 slices: SliceField[];
25 /** wall field at the nodes, interleaved complex */
26 uWall: Float64Array;
27 /** color limit for the pressure: 99th percentile of |u| more than half a wavelength from the source */
28 umax: number;
29 /** reference for the level: the median of |u| on the same points */
30 ref: number;
31}
33export const DB_MIN = -30;
34export const DB_MAX = 15;
36/** Color limits from the slices, as plot.py computes them. */
37export function colorLimits(slices: SliceField[], x0: ArrayLike<number>, k: number): { umax: number; ref: number } {
38 const vals: number[] = [];
39 const rmin = Math.PI / k;
40 for (const { slice, u } of slices) {
41 const n = slice.mask.length;
42 for (let q = 0; q < n; q++) {
43 if (!slice.mask[q]) continue;
44 const r = Math.hypot(slice.pts[q] - x0[0], slice.pts[n + q] - x0[1], slice.pts[2 * n + q] - x0[2]);
45 if (r <= rmin) continue;
46 const v = Math.hypot(u[2 * q], u[2 * q + 1]);
47 if (Number.isFinite(v)) vals.push(v);
48 }
49 }
50 if (!vals.length) return { umax: 1, ref: 1 };
51 vals.sort((x, y) => x - y);
52 return { umax: vals[Math.floor(0.99 * (vals.length - 1))], ref: vals[Math.floor(0.5 * (vals.length - 1))] };
53}
55/** Position in [0, 1] on the colormap of a complex value, or NaN. */
56export function colorParam(re: number, im: number, mode: Mode, c: number, s: number, umax: number, ref: number): number {
57 if (mode === 'pressure') return 0.5 + (0.5 * (re * c + im * s)) / umax;
58 const db = 20 * Math.log10(Math.hypot(re, im) / ref);
59 return (db - DB_MIN) / (DB_MAX - DB_MIN);
60}