From dda72739a9a444520cddbd6fb76b0b303c9aed2a Mon Sep 17 00:00:00 2001 From: MasterBabol Date: Sun, 30 Aug 2026 10:22:13 +0900 Subject: [PATCH] feat: add physical 3d receive pattern core --- packages/acoustics-core/beam-pattern.js | 407 ++++++++++++++++++++ tests/unit/beam-pattern.test.js | 479 ++++++++++++++++++++++++ 2 files changed, 886 insertions(+) diff --git a/packages/acoustics-core/beam-pattern.js b/packages/acoustics-core/beam-pattern.js index c573a50..abbba32 100644 --- a/packages/acoustics-core/beam-pattern.js +++ b/packages/acoustics-core/beam-pattern.js @@ -9,6 +9,11 @@ export const MAXIMUM_RECEIVE_BEAM_PATTERN_MICROPHONES = 1_024; export const MAXIMUM_RECEIVE_BEAM_PATTERN_FREQUENCY_HZ = 1_000_000; export const MAXIMUM_RECEIVE_BEAM_PATTERN_BASELINE_COMPARISONS = 600_000; export const RECEIVE_BEAM_PATTERN_METRIC_PHASE_STEP_LIMIT_RAD = Math.PI / 8; +export const DEFAULT_RECEIVE_BEAM_PATTERN_SPHERE_AZIMUTH_SAMPLES = 360; +export const DEFAULT_RECEIVE_BEAM_PATTERN_SPHERE_ELEVATION_SAMPLES = 181; +export const MAXIMUM_RECEIVE_BEAM_PATTERN_SPHERE_AZIMUTH_SAMPLES = 1_440; +export const MAXIMUM_RECEIVE_BEAM_PATTERN_SPHERE_ELEVATION_SAMPLES = 721; +export const MAXIMUM_RECEIVE_BEAM_PATTERN_SPHERE_SAMPLES = 262_144; const HALF_POWER = 0.5; const RADIANS_PER_DEGREE = Math.PI / 180; @@ -82,6 +87,29 @@ export function receiveBeamPatternDirection(cut, angleDeg) { return [0, sine, cosine]; } +/** + * Convert standard spherical azimuth/elevation to an array-to-source unit + * vector. Azimuth is measured in the XY plane from +X toward +Y. Elevation is + * measured from the XY plane toward +Z, with +90 degrees at the north pole. + */ +export function receiveBeamPatternSphereDirection(azimuthDeg, elevationDeg) { + const azimuthRadians = normalizeBeamPatternAngle(azimuthDeg) * + RADIANS_PER_DEGREE; + const elevation = requireFinite(elevationDeg, "elevationDeg"); + if (elevation < -90 || elevation > 90) { + throw new RangeError("elevationDeg must lie inside -90..90 degrees."); + } + if (elevation === 90) return [0, 0, 1]; + if (elevation === -90) return [0, 0, -1]; + const elevationRadians = elevation * RADIANS_PER_DEGREE; + const cosineElevation = Math.cos(elevationRadians); + return [ + cosineElevation * Math.cos(azimuthRadians), + cosineElevation * Math.sin(azimuthRadians), + Math.sin(elevationRadians), + ]; +} + function cutCoordinates(position, cut) { if (cut === "xy") return [position[0], position[1]]; if (cut === "xz") return [position[2], position[0]]; @@ -289,6 +317,32 @@ function maximumProjectedBaseline(coordinates) { return maximum; } +function maximumSpatialBaseline(positions) { + let maximum = 0; + for (let left = 0; left < positions.length; left += 1) { + for (let right = left + 1; right < positions.length; right += 1) { + maximum = Math.max( + maximum, + Math.hypot( + positions[left][0] - positions[right][0], + positions[left][1] - positions[right][1], + positions[left][2] - positions[right][2], + ), + ); + } + } + return maximum; +} + +function requireSphereSampleCount(value, label, minimum, maximum) { + if (!Number.isInteger(value) || value < minimum || value > maximum) { + throw new RangeError( + label + " must be an integer inside " + minimum + ".." + maximum + ".", + ); + } + return value; +} + function recommendedSampleCount(waveNumber, projectedApertureM) { if (projectedApertureM === 0) return 36; const ratio = RECEIVE_BEAM_PATTERN_METRIC_PHASE_STEP_LIMIT_RAD / @@ -623,3 +677,356 @@ export function computeReceiveBeamPattern(options) { }, }; } + +/** + * Sample the same ideal conventional receive array factor over a full sphere. + * + * The returned rectangular grid is equirectangular: azimuth samples span + * [0, 360) without a duplicate endpoint, while elevation samples include both + * poles. Pole rows therefore contain coincident directions for every azimuth. + * All response arrays are raw elevation-major samples; no interpolation, + * smoothing, frequency folding, or sampled-maximum normalization is applied. + */ +export function computeReceiveBeamPatternSphere(options) { + const array = options?.array; + const positions = requireMicrophones(array); + const frequencyHz = requirePositive(options?.frequencyHz, "frequencyHz"); + if (frequencyHz > MAXIMUM_RECEIVE_BEAM_PATTERN_FREQUENCY_HZ) { + throw new RangeError( + "frequencyHz must not exceed " + + MAXIMUM_RECEIVE_BEAM_PATTERN_FREQUENCY_HZ + " Hz.", + ); + } + const speedOfSound = requirePositive(options?.speedOfSound ?? 343, "speedOfSound"); + const steeringAzimuthDeg = normalizeBeamPatternAngle( + options?.steeringAzimuthDeg ?? 0, + ); + const steeringElevationDeg = requireFinite( + options?.steeringElevationDeg ?? 90, + "steeringElevationDeg", + ); + if (steeringElevationDeg < -90 || steeringElevationDeg > 90) { + throw new RangeError( + "steeringElevationDeg must lie inside -90..90 degrees.", + ); + } + const azimuthSampleCount = requireSphereSampleCount( + options?.azimuthSampleCount ?? + DEFAULT_RECEIVE_BEAM_PATTERN_SPHERE_AZIMUTH_SAMPLES, + "azimuthSampleCount", + 4, + MAXIMUM_RECEIVE_BEAM_PATTERN_SPHERE_AZIMUTH_SAMPLES, + ); + const elevationSampleCount = requireSphereSampleCount( + options?.elevationSampleCount ?? + DEFAULT_RECEIVE_BEAM_PATTERN_SPHERE_ELEVATION_SAMPLES, + "elevationSampleCount", + 3, + MAXIMUM_RECEIVE_BEAM_PATTERN_SPHERE_ELEVATION_SAMPLES, + ); + const sampleCount = azimuthSampleCount * elevationSampleCount; + if ( + !Number.isSafeInteger(sampleCount) || + sampleCount > MAXIMUM_RECEIVE_BEAM_PATTERN_SPHERE_SAMPLES + ) { + throw new RangeError( + "sphere sample product exceeds the receive beam-pattern limit of " + + MAXIMUM_RECEIVE_BEAM_PATTERN_SPHERE_SAMPLES + ".", + ); + } + const dbFloor = requireFinite( + options?.dbFloor ?? DEFAULT_RECEIVE_BEAM_PATTERN_DB_FLOOR, + "dbFloor", + ); + if (dbFloor < -300 || dbFloor > -1) { + throw new RangeError("dbFloor must lie inside -300..-1 dB."); + } + const maximumEvaluations = options?.maximumEvaluations ?? + MAXIMUM_RECEIVE_BEAM_PATTERN_EVALUATIONS; + if (!Number.isSafeInteger(maximumEvaluations) || maximumEvaluations < 1) { + throw new RangeError("maximumEvaluations must be a positive safe integer."); + } + const microphoneCount = positions.length; + const sampleEvaluationCount = microphoneCount * sampleCount; + if ( + !Number.isSafeInteger(sampleEvaluationCount) || + sampleEvaluationCount > maximumEvaluations + ) { + throw new RangeError( + "receive beam-pattern sphere requires " + sampleEvaluationCount + + " microphone-direction evaluations; limit is " + + maximumEvaluations + ".", + ); + } + const maximumBaselineComparisons = options?.maximumBaselineComparisons ?? + MAXIMUM_RECEIVE_BEAM_PATTERN_BASELINE_COMPARISONS; + if ( + !Number.isSafeInteger(maximumBaselineComparisons) || + maximumBaselineComparisons < 0 + ) { + throw new RangeError( + "maximumBaselineComparisons must be a non-negative safe integer.", + ); + } + const spatialBaselineComparisons = + microphoneCount * (microphoneCount - 1) / 2; + if (spatialBaselineComparisons > maximumBaselineComparisons) { + throw new RangeError( + "receive beam-pattern sphere requires " + spatialBaselineComparisons + + " spatial-baseline comparisons; limit is " + + maximumBaselineComparisons + ".", + ); + } + + const centroid = [0, 0, 0]; + for (const position of positions) { + centroid[0] += position[0] / microphoneCount; + centroid[1] += position[1] / microphoneCount; + centroid[2] += position[2] / microphoneCount; + } + const centeredPositions = positions.map((position) => [ + position[0] - centroid[0], + position[1] - centroid[1], + position[2] - centroid[2], + ]); + const spatialApertureM = maximumSpatialBaseline(centeredPositions); + const waveNumber = TAU * frequencyHz / speedOfSound; + const maximumPhaseMagnitudeRadBound = 2 * waveNumber * spatialApertureM; + if ( + !Number.isFinite(waveNumber) || + !Number.isFinite(maximumPhaseMagnitudeRadBound) + ) { + throw new RangeError( + "frequency, speed of sound, and array coordinates overflow phase.", + ); + } + + const azimuthStepDeg = 360 / azimuthSampleCount; + const elevationStepDeg = 180 / (elevationSampleCount - 1); + const azimuthStepRad = TAU / azimuthSampleCount; + const elevationStepRad = Math.PI / (elevationSampleCount - 1); + const azimuthAdjacentDirectionChordBound = + 2 * Math.sin(azimuthStepRad / 2); + const elevationAdjacentDirectionChordBound = + 2 * Math.sin(elevationStepRad / 2); + const maximumCellDiagonalDirectionChordBound = 2 * Math.hypot( + Math.sin(elevationStepRad / 2), + Math.cos(elevationStepRad / 2) * Math.sin(azimuthStepRad / 2), + ); + const maximumAzimuthAdjacentPhaseStepRadBound = + waveNumber * spatialApertureM * azimuthAdjacentDirectionChordBound; + const maximumElevationAdjacentPhaseStepRadBound = + waveNumber * spatialApertureM * elevationAdjacentDirectionChordBound; + const maximumCellDiagonalPhaseStepRadBound = + waveNumber * spatialApertureM * + maximumCellDiagonalDirectionChordBound; + const resolutionCriterionSatisfied = + maximumCellDiagonalPhaseStepRadBound <= + RECEIVE_BEAM_PATTERN_METRIC_PHASE_STEP_LIMIT_RAD * + (1 + 64 * Number.EPSILON); + + const azimuthAxisDeg = new Float64Array(azimuthSampleCount); + const elevationAxisDeg = new Float64Array(elevationSampleCount); + const azimuthCosines = new Float64Array(azimuthSampleCount); + const azimuthSines = new Float64Array(azimuthSampleCount); + const elevationCosines = new Float64Array(elevationSampleCount); + const elevationSines = new Float64Array(elevationSampleCount); + for (let azimuthIndex = 0; azimuthIndex < azimuthSampleCount; azimuthIndex += 1) { + const azimuthDeg = azimuthIndex * azimuthStepDeg; + const azimuthRadians = azimuthIndex * azimuthStepRad; + azimuthAxisDeg[azimuthIndex] = azimuthDeg; + azimuthCosines[azimuthIndex] = Math.cos(azimuthRadians); + azimuthSines[azimuthIndex] = Math.sin(azimuthRadians); + } + for ( + let elevationIndex = 0; + elevationIndex < elevationSampleCount; + elevationIndex += 1 + ) { + const elevationDeg = elevationIndex === elevationSampleCount - 1 + ? 90 + : -90 + elevationIndex * elevationStepDeg; + const elevationRadians = elevationDeg * RADIANS_PER_DEGREE; + elevationAxisDeg[elevationIndex] = elevationDeg; + if (elevationIndex === 0 || elevationIndex === elevationSampleCount - 1) { + elevationCosines[elevationIndex] = 0; + elevationSines[elevationIndex] = elevationIndex === 0 ? -1 : 1; + } else { + elevationCosines[elevationIndex] = Math.cos(elevationRadians); + elevationSines[elevationIndex] = Math.sin(elevationRadians); + } + } + + const steeringDirection = receiveBeamPatternSphereDirection( + steeringAzimuthDeg, + steeringElevationDeg, + ); + const steeringProjections = new Float64Array(microphoneCount); + for (let microphoneIndex = 0; microphoneIndex < microphoneCount; microphoneIndex += 1) { + const position = centeredPositions[microphoneIndex]; + steeringProjections[microphoneIndex] = + position[0] * steeringDirection[0] + + position[1] * steeringDirection[1] + + position[2] * steeringDirection[2]; + } + + const responseReal = new Float64Array(sampleCount); + const responseImaginary = new Float64Array(sampleCount); + const responsePower = new Float64Array(sampleCount); + const responseLinear = new Float64Array(sampleCount); + const responseDb = new Float64Array(sampleCount); + const displayDb = new Float64Array(sampleCount); + let sampledMaximumLinear = 0; + let sampledMaximumPower = 0; + for ( + let elevationIndex = 0; + elevationIndex < elevationSampleCount; + elevationIndex += 1 + ) { + const cosineElevation = elevationCosines[elevationIndex]; + const directionZ = elevationSines[elevationIndex]; + for ( + let azimuthIndex = 0; + azimuthIndex < azimuthSampleCount; + azimuthIndex += 1 + ) { + const directionX = cosineElevation * azimuthCosines[azimuthIndex]; + const directionY = cosineElevation * azimuthSines[azimuthIndex]; + let summedReal = 0; + let summedImaginary = 0; + for ( + let microphoneIndex = 0; + microphoneIndex < microphoneCount; + microphoneIndex += 1 + ) { + const position = centeredPositions[microphoneIndex]; + const arrivalProjection = + position[0] * directionX + + position[1] * directionY + + position[2] * directionZ; + const phase = waveNumber * + (arrivalProjection - steeringProjections[microphoneIndex]); + summedReal += Math.cos(phase); + summedImaginary += Math.sin(phase); + } + const flatIndex = + elevationIndex * azimuthSampleCount + azimuthIndex; + const real = summedReal / microphoneCount; + const imaginary = summedImaginary / microphoneCount; + const power = real * real + imaginary * imaginary; + if ( + !Number.isFinite(real) || + !Number.isFinite(imaginary) || + !Number.isFinite(power) + ) { + throw new RangeError( + "receive beam-pattern sphere phasor accumulation became non-finite.", + ); + } + const magnitude = Math.sqrt(power); + const decibels = power > 0 ? 10 * Math.log10(power) : -Infinity; + responseReal[flatIndex] = real; + responseImaginary[flatIndex] = imaginary; + responsePower[flatIndex] = power; + responseLinear[flatIndex] = magnitude; + responseDb[flatIndex] = decibels; + displayDb[flatIndex] = Math.max(dbFloor, decibels); + sampledMaximumLinear = Math.max(sampledMaximumLinear, magnitude); + sampledMaximumPower = Math.max(sampledMaximumPower, power); + } + } + + const returnedFloat64Count = + sampleCount * 6 + azimuthSampleCount + elevationSampleCount; + return { + model: "far-field-narrowband-conventional-delay-and-sum-receive-array-factor", + surface: "full-sphere-equirectangular-samples", + assumptions: Object.freeze({ + wavefront: "far-field-plane-wave", + bandwidth: "narrowband", + elements: "ideal-calibrated-omnidirectional-point-microphones", + environment: "free-field", + coordinates: "meters", + }), + coordinateSystem: Object.freeze({ + handedness: "right-handed", + azimuth: "XY-plane-from-+X-toward-+Y", + elevation: "from-XY-plane-toward-+Z", + directionEquation: + "[cos(elevation)*cos(azimuth), cos(elevation)*sin(azimuth), sin(elevation)]", + }), + frequencyHz, + speedOfSound, + wavelengthM: speedOfSound / frequencyHz, + waveNumberRadPerM: waveNumber, + steeringAzimuthDeg, + steeringElevationDeg, + steeringDirection, + microphoneCount, + azimuthSampleCount, + elevationSampleCount, + sampleCount, + gridShape: [elevationSampleCount, azimuthSampleCount], + flatIndexing: "elevation-major-row-major", + evaluationCount: sampleEvaluationCount, + dbFloor, + grid: { + topology: "equirectangular-azimuth-elevation", + shape: [elevationSampleCount, azimuthSampleCount], + flatIndex: "elevationIndex * azimuthSampleCount + azimuthIndex", + azimuthEndpointDuplicated: false, + elevationEndpointsIncluded: true, + poleAzimuthSamplesCoincident: true, + uniqueDirectionCount: + azimuthSampleCount * (elevationSampleCount - 2) + 2, + interpolation: "none", + smoothing: "none", + sampledMaximumNormalization: false, + }, + workload: { + spatialBaselineComparisons, + maximumBaselineComparisons, + sampleEvaluationCount, + maximumEvaluationCount: sampleEvaluationCount, + maximumEvaluations, + returnedFloat64Count, + returnedFloat64Bytes: + returnedFloat64Count * Float64Array.BYTES_PER_ELEMENT, + maximumSurfaceSamples: MAXIMUM_RECEIVE_BEAM_PATTERN_SPHERE_SAMPLES, + }, + angularSampling: { + mode: "fixed-equirectangular-grid", + criterion: "maximum-cell-diagonal-pair-phase-step", + azimuthStepDeg, + elevationStepDeg, + spatialApertureM, + maximumPhaseMagnitudeRadBound, + azimuthAdjacentDirectionChordBound, + elevationAdjacentDirectionChordBound, + maximumCellDiagonalDirectionChordBound, + maximumAzimuthAdjacentPhaseStepRadBound, + maximumElevationAdjacentPhaseStepRadBound, + maximumCellDiagonalPhaseStepRadBound, + phaseStepLimitRad: + RECEIVE_BEAM_PATTERN_METRIC_PHASE_STEP_LIMIT_RAD, + resolutionCriterionSatisfied, + interpretation: + "sampling-bound-only-not-proof-of-complete-lobe-detection", + }, + azimuthAxisDeg, + elevationAxisDeg, + responseReal, + responseImaginary, + responsePower, + responseLinear, + responseDb, + displayDb, + sampledSummary: { + maximumLinearOnGrid: sampledMaximumLinear, + maximumPowerOnGrid: sampledMaximumPower, + maximumDbOnGrid: sampledMaximumPower > 0 + ? 10 * Math.log10(sampledMaximumPower) + : -Infinity, + }, + }; +} diff --git a/tests/unit/beam-pattern.test.js b/tests/unit/beam-pattern.test.js index 1ac43fd..3d71dfd 100644 --- a/tests/unit/beam-pattern.test.js +++ b/tests/unit/beam-pattern.test.js @@ -1,6 +1,8 @@ import { computeReceiveBeamPattern, + computeReceiveBeamPatternSphere, receiveBeamPatternDirection, + receiveBeamPatternSphereDirection, } from "/packages/acoustics-core/beam-pattern.js"; function assertNear(actual, expected, tolerance, label) { @@ -413,3 +415,480 @@ test("co-located pattern accepts an exact sample-only evaluation cap", () => { } assert(error instanceof RangeError, "strict cap rejects one fewer evaluation"); }); + +function sphereFlatIndex(pattern, azimuthDeg, elevationDeg) { + let azimuthIndex = -1; + let elevationIndex = -1; + for (let index = 0; index < pattern.azimuthAxisDeg.length; index += 1) { + if (Math.abs(pattern.azimuthAxisDeg[index] - azimuthDeg) < 1e-10) { + azimuthIndex = index; + break; + } + } + for (let index = 0; index < pattern.elevationAxisDeg.length; index += 1) { + if (Math.abs(pattern.elevationAxisDeg[index] - elevationDeg) < 1e-10) { + elevationIndex = index; + break; + } + } + assert(azimuthIndex >= 0, "requested azimuth is present on sphere grid"); + assert(elevationIndex >= 0, "requested elevation is present on sphere grid"); + return elevationIndex * pattern.azimuthSampleCount + azimuthIndex; +} + +test("sphere directions use standard azimuth and elevation coordinates", () => { + const cases = [ + [0, 0, [1, 0, 0]], + [90, 0, [0, 1, 0]], + [180, 0, [-1, 0, 0]], + [270, 0, [0, -1, 0]], + [37, 90, [0, 0, 1]], + [241, -90, [0, 0, -1]], + ]; + for (const [azimuthDeg, elevationDeg, expected] of cases) { + const direction = receiveBeamPatternSphereDirection( + azimuthDeg, + elevationDeg, + ); + direction.forEach((value, axis) => + assertNear( + value, + expected[axis], + 1e-12, + "sphere direction axis " + axis, + )); + } + const diagonal = receiveBeamPatternSphereDirection(45, 30); + assertNear(diagonal[0], Math.sqrt(6) / 4, 1e-12, "diagonal x"); + assertNear(diagonal[1], Math.sqrt(6) / 4, 1e-12, "diagonal y"); + assertNear(diagonal[2], 0.5, 1e-12, "diagonal z"); +}); + +test("one microphone returns an unshaped raw full-sphere grid", () => { + const pattern = computeReceiveBeamPatternSphere({ + array: arrayFromPositions([[12.5, -7, 3]]), + frequencyHz: 49_125, + speedOfSound: 343, + azimuthSampleCount: 12, + elevationSampleCount: 7, + maximumEvaluations: 84, + }); + assertEqual(pattern.steeringAzimuthDeg, 0, "default sphere steering azimuth"); + assertEqual(pattern.steeringElevationDeg, 90, "default sphere steering is +Z"); + assertEqual(pattern.steeringDirection.join(","), "0,0,1", "default direction"); + assertEqual(pattern.gridShape.join(","), "7,12", "top-level grid shape"); + assertEqual( + pattern.flatIndexing, + "elevation-major-row-major", + "top-level flat indexing", + ); + assertEqual(pattern.grid.shape.join(","), "7,12", "elevation-major shape"); + assertEqual(pattern.microphoneCount, 1, "sphere microphone provenance"); + assertEqual(pattern.sampleCount, 84, "sphere sample product"); + assertEqual(pattern.evaluationCount, 84, "exact sphere evaluation count"); + assertEqual(pattern.responsePower.length, 84, "raw output sample count"); + assertEqual(pattern.azimuthAxisDeg.at(-1), 330, "360 endpoint is not duplicated"); + assertEqual(pattern.elevationAxisDeg[0], -90, "south pole is included"); + assertEqual(pattern.elevationAxisDeg.at(-1), 90, "north pole is included"); + assertEqual(pattern.grid.interpolation, "none", "no interpolation"); + assertEqual(pattern.grid.smoothing, "none", "no smoothing"); + assert( + !pattern.grid.sampledMaximumNormalization, + "sampled maximum normalization is disabled", + ); + assert(pattern.angularSampling.resolutionCriterionSatisfied, "zero aperture resolves"); + for (let index = 0; index < pattern.sampleCount; index += 1) { + assertNear(pattern.responseReal[index], 1, 0, "single-mic raw real"); + assertNear(pattern.responseImaginary[index], 0, 0, "single-mic raw imaginary"); + assertNear(pattern.responsePower[index], 1, 0, "single-mic raw power"); + assertNear(pattern.responseLinear[index], 1, 0, "single-mic raw magnitude"); + assertNear(pattern.responseDb[index], 0, 0, "single-mic raw dB"); + assertNear(pattern.displayDb[index], 0, 0, "single-mic display dB"); + } +}); + +test("half-wavelength X baseline has the expected three-dimensional null surface", () => { + const speedOfSound = 343; + const frequencyHz = 1_000; + const spacing = speedOfSound / frequencyHz / 2; + const pattern = computeReceiveBeamPatternSphere({ + array: arrayFromPositions([[-spacing / 2, 0, 0], [spacing / 2, 0, 0]]), + frequencyHz, + speedOfSound, + azimuthSampleCount: 72, + elevationSampleCount: 37, + }); + const positiveX = sphereFlatIndex(pattern, 0, 0); + const positiveY = sphereFlatIndex(pattern, 90, 0); + const negativeX = sphereFlatIndex(pattern, 180, 0); + const northPole = sphereFlatIndex(pattern, 0, 90); + const southPole = sphereFlatIndex(pattern, 0, -90); + assertNear(pattern.responseLinear[positiveX], 0, 1e-12, "+X endfire null"); + assertNear(pattern.responseLinear[negativeX], 0, 1e-12, "-X endfire null"); + assertNear(pattern.responseLinear[positiveY], 1, 1e-12, "+Y broadside"); + assertNear(pattern.responseLinear[northPole], 1, 1e-12, "+Z steering response"); + assertNear(pattern.responseLinear[southPole], 1, 1e-12, "-Z front-back tie"); + assert( + pattern.angularSampling.resolutionCriterionSatisfied, + "five-degree grid satisfies the disclosed cell-diagonal phase criterion", + ); + const expectedDiagonalChord = 2 * Math.hypot( + Math.sin(2.5 * Math.PI / 180), + Math.cos(2.5 * Math.PI / 180) * Math.sin(2.5 * Math.PI / 180), + ); + assertNear( + pattern.angularSampling.maximumCellDiagonalPhaseStepRadBound, + Math.PI * expectedDiagonalChord, + 1e-12, + "full-aperture cell-diagonal phase bound", + ); + for (let azimuthIndex = 1; azimuthIndex < pattern.azimuthSampleCount; azimuthIndex += 1) { + assertNear( + pattern.responsePower[northPole + azimuthIndex], + pattern.responsePower[northPole], + 0, + "north-pole azimuth samples are physically coincident", + ); + assertNear( + pattern.responsePower[southPole + azimuthIndex], + pattern.responsePower[southPole], + 0, + "south-pole azimuth samples are physically coincident", + ); + } +}); + +test("sphere raw complex response uses all three coordinates and the receive sign", () => { + const speedOfSound = 343; + const frequencyHz = 1_000; + const wavelengthM = speedOfSound / frequencyHz; + const positions = [ + [0, 0, 0], + [wavelengthM / 4, 0, wavelengthM / 16], + [0, wavelengthM / 8, 3 * wavelengthM / 16], + ]; + const pattern = computeReceiveBeamPatternSphere({ + array: arrayFromPositions(positions), + frequencyHz, + speedOfSound, + steeringAzimuthDeg: 0, + steeringElevationDeg: 90, + azimuthSampleCount: 8, + elevationSampleCount: 7, + }); + const targetDirection = receiveBeamPatternSphereDirection(45, 30); + const steeringDirection = [0, 0, 1]; + const centroid = [0, 0, 0]; + for (const position of positions) { + for (let axis = 0; axis < 3; axis += 1) { + centroid[axis] += position[axis] / positions.length; + } + } + let expectedReal = 0; + let expectedImaginary = 0; + for (const position of positions) { + const centered = position.map((coordinate, axis) => coordinate - centroid[axis]); + const phase = 2 * Math.PI / wavelengthM * ( + centered[0] * (targetDirection[0] - steeringDirection[0]) + + centered[1] * (targetDirection[1] - steeringDirection[1]) + + centered[2] * (targetDirection[2] - steeringDirection[2]) + ); + expectedReal += Math.cos(phase) / positions.length; + expectedImaginary += Math.sin(phase) / positions.length; + } + const targetIndex = sphereFlatIndex(pattern, 45, 30); + assertNear(pattern.responseReal[targetIndex], expectedReal, 1e-12, "3D raw real"); + assertNear( + pattern.responseImaginary[targetIndex], + expectedImaginary, + 1e-12, + "3D raw imaginary sign", + ); + assertNear( + pattern.responsePower[targetIndex], + expectedReal * expectedReal + expectedImaginary * expectedImaginary, + 1e-12, + "3D raw power", + ); +}); + +test("sphere samples agree with the existing XZ great-circle contract", () => { + const array = arrayFromPositions([ + [-0.12, 0.03, 0.04], + [0.04, -0.08, 0.02], + [0.15, 0.07, -0.11], + ]); + const shared = { + array, + frequencyHz: 3_200, + speedOfSound: 343, + }; + const cut = computeReceiveBeamPattern({ + ...shared, + cut: "xz", + steeringAngleDeg: 0, + sampleCount: 360, + }); + const sphere = computeReceiveBeamPatternSphere({ + ...shared, + steeringAzimuthDeg: 0, + steeringElevationDeg: 90, + azimuthSampleCount: 4, + elevationSampleCount: 3, + }); + const mappings = [ + [0, 0, 90], + [90, 0, 0], + [180, 0, -90], + [270, 180, 0], + ]; + for (const [cutAngleDeg, azimuthDeg, elevationDeg] of mappings) { + const cutIndex = [...cut.anglesDeg].findIndex( + (angleDeg) => angleDeg === cutAngleDeg, + ); + const sphereIndex = sphereFlatIndex(sphere, azimuthDeg, elevationDeg); + assert(cutIndex >= 0, "mapped XZ angle is sampled"); + assertNear( + sphere.responseReal[sphereIndex], + cut.responseReal[cutIndex], + 1e-12, + "sphere/XZ raw real at " + cutAngleDeg, + ); + assertNear( + sphere.responseImaginary[sphereIndex], + cut.responseImaginary[cutIndex], + 1e-12, + "sphere/XZ raw imaginary at " + cutAngleDeg, + ); + assertNear( + sphere.responsePower[sphereIndex], + cut.responsePower[cutIndex], + 1e-12, + "sphere/XZ raw power at " + cutAngleDeg, + ); + } +}); + +test("sphere array factor remains invariant to translating a three-dimensional array", () => { + const basePositions = [ + [-0.12, 0.03, 0.04], + [0.04, -0.08, 0.02], + [0.15, 0.07, -0.11], + ]; + const translation = [10, -20, 5]; + const translatedPositions = basePositions.map((position) => + position.map((coordinate, axis) => coordinate + translation[axis])); + const options = { + frequencyHz: 7_500, + speedOfSound: 343, + steeringAzimuthDeg: 45, + steeringElevationDeg: 30, + azimuthSampleCount: 8, + elevationSampleCount: 7, + }; + const base = computeReceiveBeamPatternSphere({ + ...options, + array: arrayFromPositions(basePositions), + }); + const translated = computeReceiveBeamPatternSphere({ + ...options, + array: arrayFromPositions(translatedPositions), + }); + for (let index = 0; index < base.sampleCount; index += 1) { + assertNear( + translated.responseReal[index], + base.responseReal[index], + 1e-11, + "translation-invariant sphere real " + index, + ); + assertNear( + translated.responseImaginary[index], + base.responseImaginary[index], + 1e-11, + "translation-invariant sphere imaginary " + index, + ); + } +}); + +test("sphere response preserves analog frequency and display floor only affects display dB", () => { + const highFrequencyHz = 49_125; + const speedOfSound = 343; + const halfHighWavelength = speedOfSound / (2 * highFrequencyHz); + const array = arrayFromPositions([ + [-halfHighWavelength / 2, 0, 0], + [halfHighWavelength / 2, 0, 0], + ]); + const shared = { + array, + speedOfSound, + steeringElevationDeg: 90, + azimuthSampleCount: 4, + elevationSampleCount: 3, + }; + const highAnalog = computeReceiveBeamPatternSphere({ + ...shared, + frequencyHz: highFrequencyHz, + dbFloor: -30, + }); + const highAnalogDeepFloor = computeReceiveBeamPatternSphere({ + ...shared, + frequencyHz: highFrequencyHz, + dbFloor: -100, + }); + const foldedTemporalAlias = computeReceiveBeamPatternSphere({ + ...shared, + frequencyHz: 1_125, + }); + const positiveX = sphereFlatIndex(highAnalog, 0, 0); + assertEqual(highAnalog.frequencyHz, highFrequencyHz, "analog frequency is preserved"); + assertNear( + highAnalog.wavelengthM, + speedOfSound / highFrequencyHz, + 1e-15, + "sphere uses analog wavelength", + ); + assert( + Math.abs( + highAnalog.responseLinear[positiveX] - + foldedTemporalAlias.responseLinear[positiveX] + ) > 0.1, + "sphere response is not computed from a folded ADC alias", + ); + for (let index = 0; index < highAnalog.sampleCount; index += 1) { + assertNear( + highAnalog.responsePower[index], + highAnalogDeepFloor.responsePower[index], + 0, + "sphere raw power is display-floor independent", + ); + } + assertEqual(highAnalog.displayDb[positiveX], -30, "sphere shallow display floor"); + assertEqual( + highAnalogDeepFloor.displayDb[positiveX], + -100, + "sphere deep display floor", + ); + assert( + highAnalog.responseDb[positiveX] < -300, + "sphere raw null is retained below display floors", + ); +}); + +test("coarse sphere grids expose failed resolution provenance without shaping raw samples", () => { + const pattern = computeReceiveBeamPatternSphere({ + array: arrayFromPositions([[-0.5, 0, 0], [0.5, 0, 0]]), + frequencyHz: 100_000, + speedOfSound: 343, + azimuthSampleCount: 12, + elevationSampleCount: 7, + }); + assert( + !pattern.angularSampling.resolutionCriterionSatisfied, + "coarse full-sphere phase-step criterion fails", + ); + assert( + pattern.angularSampling.maximumCellDiagonalPhaseStepRadBound > + pattern.angularSampling.phaseStepLimitRad, + "failed cell-diagonal phase bound is exposed", + ); + assertEqual(pattern.responsePower.length, 84, "unresolved raw sphere is retained"); + assertEqual( + pattern.angularSampling.interpretation, + "sampling-bound-only-not-proof-of-complete-lobe-detection", + "resolution claim remains limited", + ); +}); + +test("sphere preflight checks axis, product, baseline, and evaluation caps before output allocation", () => { + const single = arrayFromPositions([[0, 0, 0]]); + const exact = computeReceiveBeamPatternSphere({ + array: single, + frequencyHz: 1_000, + azimuthSampleCount: 4, + elevationSampleCount: 3, + maximumEvaluations: 12, + }); + assertEqual(exact.evaluationCount, 12, "exact sphere evaluation cap is accepted"); + assertEqual( + exact.workload.returnedFloat64Bytes, + (12 * 6 + 4 + 3) * Float64Array.BYTES_PER_ELEMENT, + "returned allocation provenance", + ); + + let evaluationError = null; + try { + computeReceiveBeamPatternSphere({ + array: single, + frequencyHz: 1_000, + azimuthSampleCount: 4, + elevationSampleCount: 3, + maximumEvaluations: 11, + }); + } catch (error) { + evaluationError = error; + } + assert(evaluationError instanceof RangeError, "sphere evaluation cap rejects work"); + assert(evaluationError.message.includes("12"), "evaluation requirement is explicit"); + + let productError = null; + try { + computeReceiveBeamPatternSphere({ + array: single, + frequencyHz: 1_000, + azimuthSampleCount: 720, + elevationSampleCount: 721, + }); + } catch (error) { + productError = error; + } + assert(productError instanceof RangeError, "sphere sample product is bounded"); + assert( + productError.message.includes("262144"), + "sphere sample-product limit is explicit", + ); + + let axisError = null; + try { + computeReceiveBeamPatternSphere({ + array: single, + frequencyHz: 1_000, + azimuthSampleCount: 1_441, + elevationSampleCount: 3, + }); + } catch (error) { + axisError = error; + } + assert(axisError instanceof RangeError, "sphere azimuth axis is bounded"); + + let baselineError = null; + try { + computeReceiveBeamPatternSphere({ + array: arrayFromPositions([[0, 0, 0], [0.1, 0, 0]]), + frequencyHz: 1_000, + azimuthSampleCount: 4, + elevationSampleCount: 3, + maximumBaselineComparisons: 0, + }); + } catch (error) { + baselineError = error; + } + assert(baselineError instanceof RangeError, "sphere baseline work is bounded"); + assert( + baselineError.message.includes("1 spatial-baseline"), + "sphere baseline requirement is explicit", + ); + + let steeringError = null; + try { + computeReceiveBeamPatternSphere({ + array: single, + frequencyHz: 1_000, + steeringElevationDeg: 91, + azimuthSampleCount: 4, + elevationSampleCount: 3, + }); + } catch (error) { + steeringError = error; + } + assert(steeringError instanceof RangeError, "invalid sphere elevation is rejected"); +});