ObieWebApp

1.1 What is a FRF?

We want to express the violin mathematically as a transfer function (also called the frequency response function) — it converts some input (in this case a force from the hammer) to some output (in this case sound or a pressure fluctuation):

Taking the Fourier Transform, we can find the transfer function for the violin as:

where \(P(f)\) and \(F(f)\) are the Fourier transforms of pressure and force and have real and imaginary parts. The Fourier Transform is defined as:

We can simplify — get the denominator to be real — by multiplying numerator and denominator by \(F^*\):

1.2 Least Squares Fit — H1 and H2

In reality, we take 5–10 hits per position so that we can get a "best fit" reading. We know that:

So the error of the best fit for measured values \(F_i\) and \(P_i\) is:

Realizing that \(H(f)\) is a complex number \((a+ib)\) and that \(F, P\) are complex as well:

Summing over all readings and finding the minimum by differentiating with respect to \(a\) and \(b\) and setting equal to zero:

Multiplying the second equation by \(i\) and then adding:

Adding them together:

This gives the H1 estimator (minimises output noise):

Alternatively, we can minimise the effect of noise on the input by using as the error:

Expanding the squared error:

Differentiating and setting to zero:

Multiplying the second by \(i\) and adding:

Adding and rearranging:

This gives the H2 estimator (minimises input noise):

If the data is clean and the noise is truly random, H1 and H2 should be very close and their ratio — the coherence — should be 1. George used the geometric mean as the working estimate:

1.3 Coherence

The coherence \(C\) measures how consistently the system behaves across repeated measurements:

The numerator is the cross-spectrum squared — how much \(F\) and \(P\) move together. The denominator is the product of the individual power spectra — how much each signal varies on its own. It is the frequency-domain analogue of \(R^2\) from linear regression, and answers: "how much of the output can be explained by the input?"

• \(C \approx 1\): output is well explained by input — clean measurement, linear system.
• \(C \approx 0\): output has little to do with input at that frequency — noise dominated, or the system simply doesn't respond there.
• \(0 < C < 1\): partial explanation — could be noise, nonlinearity, or multiple inputs contributing to the output.

We can look at the average distribution of energy by breaking the FRF into frequency bands and averaging the magnitudes over each band:

We also calculate the spectral centroid of the band — the first moment — which tells you where the energy is concentrated within the band:

In Python (Python/processing/bands.py):

f_lo  = float(band['start'])
f_hi  = float(band['end'])
mask     = (freq >= f_lo) & (freq <= f_hi)
H_band   = H_mag[mask]          # linear magnitude, NOT in dB
f_band   = freq[mask]
avg_db   = 20.0 * math.log10(max(float(np.mean(H_band)), 1e-12))
centroid = float(np.sum(f_band * H_band) / np.sum(H_band))

We can estimate what an instrument will sound like by convolving its FRF with a recording — ideally captured with a force sensor at the bridge rather than a microphone. The convolution of two functions \(f(t)\) and \(g(t)\) is:

What it does

A violin has resonances — certain frequencies ring longer and louder than others. The FRF is a fingerprint of exactly how the violin shapes sound. Convolution applies that fingerprint to any audio signal. Feed in a Tchaikovsky recording; get out what it would sound like coloured by that violin's resonance structure.

FFT acceleration

Computing the integral directly is slow. The convolution theorem says we can instead multiply the two signals in the frequency domain:

In Python (Python/processing/convolution.py):

y = np.fft.irfft(
    np.fft.rfft(wav, n_fft) * np.fft.rfft(ir.astype(np.float64), n_fft),
    n_fft,
)[:N].real

The output is the input signal with the violin's character stamped onto it — its resonant peaks amplified, its dead spots suppressed, its natural decay times imposed on every note.

IR length

The impulse response (IR) length determines frequency resolution. ObieApp uses an adaptive length: the next power of 2 above \(\text{sample\_rate} / \Delta f\), capped at 65536:

def _adaptive_ir_length(freqs, sample_rate):
    df = freqs[1] - freqs[0]
    natural = int(sample_rate / df)
    next_pow2 = 1 << (natural - 1).bit_length()
    return min(next_pow2, 65536)

Minimum-phase reconstruction

For real-only FRF files (AvR format, no phase information), the phase can be estimated by treating the system as minimum-phase — for such systems the log-magnitude and phase are a Hilbert transform pair, so knowing one determines the other. In practice ObieApp stores all FRFs as complex (AvC or TRF), so minimum-phase reconstruction is not required.

All multi-byte fields are little-endian. Fields prefixed f in the MtxVec TSingleArray descriptor are 32-bit single-precision floats; all others are as specified. Sources: george.pas (Dew Research / DSP Master) and the LabVIEW writer code.

4.1 TRF / TRV — Transfer Function Format

Binary format used by Rational Acoustics / MtxVec. The header is 110 bytes, followed immediately by the data payload.

TFRFHeader (110 bytes, packed — no padding)

4.2 AvR / AvC — Averaged Spectrum Format

AvR files store a real-valued averaged spectrum; AvC files store a complex-valued averaged spectrum. The file layout is identical except for the data section, which is \(N \times 8\) bytes for AvR and \(N \times 16\) bytes for AvC. \(N\) is given by fLength in the TSingleArray descriptor.

This section derives the linearised acoustic wave equation from first principles and applies it to the mechanics of violin plates and the sound hole.

5.1 Wave Equation Derivation

Assume that the pressure is \(p\), the speed of a fluid particle is \(u\) with a mean of zero, and the density everywhere is \(\rho_0 + \rho\). The 1D continuity equation gives:

Expanding and dropping the nonlinear term \(\partial(\rho u)/\partial x\):

Differentiating again with respect to \(t\):

The momentum equation (1D, no gravity or friction) linearised for small fluctuations:

Differentiating with respect to \(x\) and dropping the middle (nonlinear) term:

Combining the continuity and momentum equations:

Adiabatic closure

For an adiabatic process \(p = K\rho^\gamma\), so:

Linearising the second time derivative:

Substituting into the combined equation gives the acoustic wave equation:

where the speed of sound \(c = \sqrt{\gamma RT}\).

5.2 Acoustic Impedance

Specific acoustic impedance \(z\) is the ratio of pressure to acoustic particle velocity: \(p = zu\), analogous to voltage driving a current.

For a plate with mass \(m\), stiffness \(k\), and damping \(\beta\), area \(A\), the equation of motion is:

Solving for impedance:

Assuming a harmonic solution \(u = \dot{x} = Ce^{i\omega t}\):

where the lumped parameters are:

This is the direct mechanical analogue of an electrical RLC circuit: resistance \(r\), inductance (mass) \(M\), and compliance (capacitance) \(C\).

For the sound hole, modelled with frictional resistance \(R_h\), and the top plate with restoring force:

5.3 Sound Intensity

Sound intensity is defined as the sound power per unit surface area. For a point source of power \(P\):

Sound Power Level and Sound Intensity Level are referenced to \(P_0 = 10^{-12}\,\mathrm{W}\) and \(I_0 = 10^{-12}\,\mathrm{W/m^2}\) respectively:

Since \(I = pu = p^2/z\), the Sound Pressure Level is:

referenced to \(p_0 = 20\,\mu\mathrm{Pa}\).

5.4 Quality Factor

The quality factor \(Q\) characterises the sharpness of a resonance peak. For the damped oscillator \(m\ddot{x} + R\dot{x} + kx = 0\):

where \(\Delta f\) is the half-power bandwidth (the −3 dB width of the peak), \(r\) is the exponential decay rate, and \(R\) is the damping coefficient. Higher \(Q\) means a sharper, more persistent resonance.