MUSIC Spectrum — Control Systems/Spectral Measurements

MUSIC Spectrum

Control Systems / Spectral Measurements

The MUSIC pseudospectrum of a signal frame – MATLAB's pmusic. The frame's correlation matrix R of order m is decomposed as R = V·diag(λ)·V′, eigenvalues in descending order; the p largest span the signal subspace and the other m−p eigenvectors v_k the noise subspace orthogonal to it. A sinusoid in the frame lies in the signal subspace, so its frequency is where the noise subspace has almost nothing:

S(f) = 1 ÷ Σ_k>p |V_k(f)|²,   V_k(f) = Σ_n v_k[n]·e^−j2πfn,   f = k/nfft cycles per sample.

R is built by the covariance method: with L = N−m+1, row i of the data matrix X holds samples i+m−1 down to i, and R = X′X/L, exactly corrmtx(x, m−1, 'cov') multiplied by its own transpose. Peaks mark the sinusoids' frequencies; their heights are not powers.

Ports

• x – the signal frame, an [N,1] column with the oldest sample first. N is read off the wire and must be at least 2m−1 – as many rows in the data matrix as it has columns – and at most 128. A Tapped Delay with Delay Order Oldest turns a stream into such a frame.
• S – the pseudospectrum, [floor(nfft/2)+1, 1] one-sided (bins 0 to nfft/2, DC first) or [nfft, 1] two-sided; bin k is at frequency k/nfft of the sampling rate.

Parameters

• Signal Subspace Dimension – p, the number of complex exponentials in the frame: two per real sinusoid. A whole number from 1 to m−1. Defaults to 4, two real sinusoids.
• Correlation Matrix Order – m, the size of the m×m correlation matrix and the length of each eigenvector, a whole number from 2 to 16. pmusic uses 2p, and so does the default, 8.
• Dimension Rule – how the signal subspace dimension is decided:
  • Fixed – it is p.
  • Threshold – pmusic's two-element form [p thr]: it is the number of eigenvalues greater than Threshold times the smallest one, but never more than p, and p when none is.
• Threshold – thr, the multiple of the smallest eigenvalue the Threshold rule compares against; a number of 0 or more. Defaults to 10. Ignored by the Fixed rule.
• FFT Length – nfft, the number of frequency points around the whole circle, a whole number from 2 to 1024. Defaults to 256, pmusic's. An nfft below m samples each eigenvector's transform at nfft points all the same, which is what pmusic's own transform does.
• Spectrum Range – One-sided (the default, pmusic's for a real signal): bins 0 to floor(nfft/2), with nothing doubled, because a pseudospectrum is not a density. Two-sided: all nfft bins.
• Sampling Time (s) – zero or less inherits the solver's rate; a positive value runs the block at that period.

Code export

All ten targets: Python, MATLAB, Java, Rust, C, C++, VHDL, Verilog, SystemVerilog and PLC Structured Text, all from one program, so every target does the same arithmetic in the same order as the simulation: the correlation matrix, a cyclic Jacobi eigenvalue solve, the sort, and one transform per eigenvector. All six settings and the frame length are structural; re-export after changing any of them.

The three HDL targets run that program in real arithmetic and quantize only at the ports: simulation-only, not offered as synthesizable. An eigenvalue solve takes square roots and divides by the data, which a Q16.16 datapath does not carry.

Simulink bridge

None (Support::None). pmusic is a Signal Processing Toolbox function and that toolbox ships no Simulink library; a search of every DSP System Toolbox and Simulink library on R2026a found no MUSIC, subspace or pseudospectrum block, so there is no path a diagram could name. The bridge reports this block rather than dropping it silently, and it has no parity testbench; code export verification covers all ten languages.

Notes

• Algebraic: the output depends on this step's frame alone.
• Verified against R2026a: the same program agrees with pmusic to 2.2e−11 relative on every bin over 1400 frames of tones in noise, the threshold form and both ranges included. With m = 2p the block is pmusic(x, p, nfft); for any other m it is pmusic(corrmtx(x, m−1, 'cov'), p, nfft). Add 'twosided' for the two-sided range and write p as [p thr] for the threshold rule.
• When the signal subspace is the smaller one (p ≤ m−p) the block does what pmusic does: it sums the p signal eigenvectors and answers 1 ÷ |m − Σ|, the same number by orthonormality.
• ⚠ There is no rootmusic mode, and the reason is measured. rootmusic roots the noise-subspace polynomial and reports the p frequencies whose roots sit closest to the unit circle, with their powers. Over 500 perturbations of a 24-sample frame the FREQUENCIES were steady (1.5e−14 under a 1e−15 perturbation, 1e−6 rad under one Q16.16 quantum) but the POWERS were not: a 1e−15 perturbation moved them by 1.4e−10 relative and one Q16.16 quantum moved them by 9.8 %, ten times the band the three HDL targets are held to. A mode whose answer no rig can certify is not shipped; the peaks of this pseudospectrum locate the same frequencies.
• An eigenvector's sign is arbitrary, and this block's solver and MATLAB's SVD choose differently; the pseudospectrum reads each one only through |V_k(f)|², which no sign changes.
• A frame whose signal and noise eigenvalues are nearly equal – no sinusoid, or more sinusoids than p/2 – has no well-defined split and the answer then depends on rounding. Drive it with tones.
• The three HDL targets' output port is Q16.16, whose largest value is 32767; a frame with almost no noise gives a peak far beyond it, which those three targets cannot carry while the other seven do.