From 7e560af4a11e6c9c53f71d28e746808826853e67 Mon Sep 17 00:00:00 2001 From: Sergio Aguirre Date: Sat, 31 Oct 2026 21:36:40 +0100 Subject: [PATCH] [FEATURE] add pytta.iso10534 for impedance tube measurements Two-microphone transfer function method of ISO 10534-2: transfer function from the two microphone positions, reflection factor (Equation 17, with the calibration factor of 7.5.2), absorption coefficient, specific impedance ratio, speed of sound from temperature and the working frequency range of a circular tube. Tests synthesize H12 with Equation (D.7) for a rigid, an anechoic and a complex termination and recover r exactly. Co-Authored-By: Claude Fable 5.1 --- pytta/__init__.py | 2 ++ pytta/iso10534.py | 68 ++++++++++++++++++++++++++++++++++++++++++ tests/test_iso10534.py | 67 +++++++++++++++++++++++++++++++++++++++++ 3 files changed, 137 insertions(+) create mode 100644 pytta/iso10534.py create mode 100644 tests/test_iso10534.py diff --git a/pytta/__init__.py b/pytta/__init__.py index d8740aa..31433fb 100644 --- a/pytta/__init__.py +++ b/pytta/__init__.py @@ -108,6 +108,7 @@ from . import utils from . import rooms from . import iso3741 +from . import iso10534 from .functions import read_wav, write_wav, merge, split, \ list_devices, print_devices, get_device_from_user,\ @@ -130,6 +131,7 @@ 'generate', 'utils', 'iso3741', + 'iso10534', 'rooms', # Functions diff --git a/pytta/iso10534.py b/pytta/iso10534.py new file mode 100644 index 0000000..aaa7fbe --- /dev/null +++ b/pytta/iso10534.py @@ -0,0 +1,68 @@ +""" +PyTTa ISO 10534-2 impedance tube: +-------------------------------- +Normal-incidence reflection factor, absorption coefficient and impedance +of a sample in an impedance tube, by the two-microphone transfer function +method of ISO 10534-2:1998. + + >>> f, H12 = pytta.iso10534.transfer_function(mic1, mic2) + >>> r = pytta.iso10534.reflection_factor(H12, f, x1, s, c0) + >>> alpha = pytta.iso10534.absorption(r) + +Microphone 1 is the one further from the sample, at x1 from it; microphone +2 is s closer. The time convention is that of numpy.fft, exp(+j w t), the +same as the standard's. +""" +import numpy as np + + +def speed_of_sound(temperature): + """Speed of sound in m/s at `temperature` in degrees Celsius, + Equation (5): c0 = 343.2 sqrt(T / 293), T in kelvin.""" + return 343.2 * np.sqrt((temperature + 273.15) / 293) + + +def frequency_range(d, s, c0): + """ + Working frequency range (f_low, f_up) in Hz of a circular tube of inner + diameter d with microphones s apart, both in m: f_up from Equations (2) + and (4) (f_up d < 0.58 c0, f_up s < 0.45 c0), and f_low where s is 5 % + of the wavelength, the general guide of 4.2. + """ + return 0.05 * c0 / s, min(0.58 * c0 / d, 0.45 * c0 / s) + + +def transfer_function(mic1, mic2): + """Frequencies and transfer function H12 = P2 / P1 from the recordings + or impulse responses of the two microphone positions (first channel of + each SignalObj).""" + p1 = np.asarray(mic1.timeSignal[:, 0], dtype=float) + p2 = np.asarray(mic2.timeSignal[:, 0], dtype=float) + return np.fft.rfftfreq(p1.size, 1 / mic1.samplingRate)[1:], \ + (np.fft.rfft(p2) / np.fft.rfft(p1))[1:] + + +def reflection_factor(H12, f, x1, s, c0, calibration=None): + """ + Normal-incidence reflection factor, Equation (17): + r = (H12 - HI) / (HR - H12) exp(2 j k x1), with HI = exp(-j k s) and + HR = exp(j k s) (D.5, D.6). `calibration` is the factor Hc of 7.5.2, + which divides H12, if the microphones were not interchanged. + """ + H12 = np.asarray(H12) if calibration is None \ + else np.asarray(H12) / calibration + k = 2 * np.pi * np.asarray(f) / c0 + return (H12 - np.exp(-1j * k * s)) / (np.exp(1j * k * s) - H12) \ + * np.exp(2j * k * x1) + + +def absorption(r): + """Normal-incidence absorption coefficient, Equation (18): + 1 - |r|^2.""" + return 1 - np.abs(r) ** 2 + + +def impedance_ratio(r): + """Specific acoustic impedance ratio Z / (rho c0), Equation (19): + (1 + r) / (1 - r).""" + return (1 + np.asarray(r)) / (1 - np.asarray(r)) diff --git a/tests/test_iso10534.py b/tests/test_iso10534.py new file mode 100644 index 0000000..d4f32ab --- /dev/null +++ b/tests/test_iso10534.py @@ -0,0 +1,67 @@ +""" +Tests for pytta.iso10534, the two-microphone transfer function method of +ISO 10534-2 for impedance tubes. The transfer functions are synthesized +with Equation (D.7) for known reflection factors, so the expected values +are exact. Each test pins one rule. +""" +import os +import unittest + +os.environ.setdefault('MPLBACKEND', 'Agg') + +import numpy as np + +import pytta +from pytta import iso10534 + +F = np.linspace(100, 1800, 200) +C = iso10534.speed_of_sound(20.0) +X1, S = 0.10, 0.05 # m: sample to far microphone, spacing + + +def h12(r): + """Equation (D.7): H12 = p2 / p1 for a reflection factor r at x = 0.""" + k = 2 * np.pi * F / C + x2 = X1 - S + return (np.exp(1j * k * x2) + r * np.exp(-1j * k * x2)) \ + / (np.exp(1j * k * X1) + r * np.exp(-1j * k * X1)) + + +class TestISO10534(unittest.TestCase): + + def test_speed_of_sound(self): + # Equation (5): c0 = 343.2 sqrt(T / 293), T in kelvin. + self.assertAlmostEqual(iso10534.speed_of_sound(19.85), 343.2, 6) + + def test_known_reflection_factors_come_back(self): + # Rigid termination (r = 1, alpha = 0), anechoic (r = 0, alpha = 1, + # Z = rho c) and a complex r with frequency-dependent phase. + cases = (np.ones(F.size), np.zeros(F.size), + 0.5 * np.exp(1j * (np.pi / 4 + F / 1000))) + for r in cases: + got = iso10534.reflection_factor(h12(r), F, X1, S, C) + self.assertTrue(np.allclose(got, r, atol=1e-9)) + alpha = iso10534.absorption(got) + self.assertTrue(np.allclose(alpha, 1 - np.abs(r) ** 2, atol=1e-9)) + z = iso10534.impedance_ratio(np.zeros(F.size)) + self.assertTrue(np.allclose(z, 1)) + + def test_transfer_function_from_two_responses(self): + # p2 is p1 delayed by d samples and halved: H12 = 0.5 exp(-j w d). + fs, d = 48000, 7 + p1 = pytta.SignalObj(np.eye(1, 4096, 100)[0], 'time', fs) + p2 = pytta.SignalObj(0.5 * np.eye(1, 4096, 100 + d)[0], 'time', fs) + f, H = iso10534.transfer_function(p1, p2) + self.assertTrue(np.allclose(H, 0.5 * np.exp(-2j * np.pi * f * d / fs))) + + def test_working_frequency_range(self): + # Equations (2) and (4) and the 5 % rule: a 100 mm circular tube + # with microphones 50 mm apart. + c = 343.2 + low, high = iso10534.frequency_range(0.1, 0.05, c) + self.assertAlmostEqual(high, min(0.58 * c / 0.1, 0.45 * c / 0.05)) + self.assertAlmostEqual(low, 0.05 * c / 0.05) + + +if __name__ == '__main__': + unittest.main()