diff --git a/pytta/__init__.py b/pytta/__init__.py index d8740aa..0d81a8a 100644 --- a/pytta/__init__.py +++ b/pytta/__init__.py @@ -108,6 +108,7 @@ from . import utils from . import rooms from . import iso3741 +from . import iso717 from .functions import read_wav, write_wav, merge, split, \ list_devices, print_devices, get_device_from_user,\ @@ -130,6 +131,7 @@ 'generate', 'utils', 'iso3741', + 'iso717', 'rooms', # Functions diff --git a/pytta/iso717.py b/pytta/iso717.py new file mode 100644 index 0000000..ba5d1b1 --- /dev/null +++ b/pytta/iso717.py @@ -0,0 +1,137 @@ +""" +PyTTa ISO 717 rating: +-------------------- +Single-number quantities of airborne sound insulation, ISO 717-1:2013, and +of impact sound insulation, ISO 717-2. + + >>> Rw, C, Ctr = pytta.iso717.airborne(R) + >>> Lnw, CI = pytta.iso717.impact(Ln) + +""" +import numpy as np + +# ISO 717-1:2013, Table 3 (reference values) and Table 4 (spectra No. 1 and +# No. 2 for C and Ctr), from 100 Hz to 3 150 Hz in thirds and from 125 Hz +# to 2 000 Hz in octaves. +_TABLES = { + 16: dict(reference=[33, 36, 39, 42, 45, 48, 51, 52, 53, 54, 55, 56, 56, + 56, 56, 56], + spectra=([-29, -26, -23, -21, -19, -17, -15, -13, -12, -11, -10, + -9, -9, -9, -9, -9], + [-20, -20, -18, -16, -15, -14, -13, -12, -11, -9, -8, + -9, -10, -11, -13, -15]), + limit=320, at500=7), + 5: dict(reference=[36, 45, 52, 55, 56], + spectra=([-21, -14, -8, -5, -4], [-14, -10, -7, -4, -6]), + limit=100, at500=2), +} + + +def _half_up(x, step): + """Round to a multiple of step, halves up (ISO 80000-1), clearing the + binary error of the input first.""" + return np.floor(np.round(np.asarray(x, dtype=float) / step, 6) + 0.5) + + +def _shift(tenths, reference, limit): + """Largest shift s, in dB, of the reference curve (in tenths of a dB) + such that the deviations of reference + s above tenths sum to no more + than limit.""" + shift = int(np.floor((tenths - reference).min() / 10)) # no deviation + while np.maximum(reference + (shift + 1) * 10 - tenths, 0).sum() \ + <= limit: + shift += 1 + return shift + + +def airborne(R): + """ + Weighted sound reduction index and spectrum adaptation terms of + ISO 717-1:2013, from R, R', Dn or DnT by band. + + Parameters: + ----------- + * R (array-like): + 16 one-third-octave bands from 100 Hz to 3 150 Hz or 5 octave + bands from 125 Hz to 2 000 Hz, in dB. + + Return: + ------- + * (Xw, C, Ctr), integers in dB: Rw (or R'w, Dn,w, DnT,w) and the + adaptation terms for pink noise and for urban traffic noise. + + R is rounded to 0.1 dB, halves up (4.4). The reference curve is shifted + in steps of 1 dB until the sum of unfavourable deviations is as large as + possible but not more than 32.0 dB in thirds or 10.0 dB in octaves, and + Xw is its value at 500 Hz. C and Ctr follow Equations (1) and (2). + Values are kept in tenths of a decibel as integers, so the limit is + compared exactly. + """ + tables = _TABLES.get(len(R)) + if tables is None: + raise ValueError('R must have 16 one-third-octave or 5 octave bands.') + tenths = _half_up(R, 0.1) + reference = np.array(tables['reference']) * 10 + shift = _shift(tenths, reference, tables['limit']) + Xw = tables['reference'][tables['at500']] + shift + C = [int(_half_up(-10 * np.log10(np.sum(10 ** ((np.array(L) - tenths / 10) + / 10))), 1)) - Xw + for L in tables['spectra']] + return Xw, C[0], C[1] + + +# ISO 717-2 reference values for impact sound, from 100 Hz to 3 150 Hz in +# thirds and from 125 Hz to 2 000 Hz in octaves, as restated by Crocker +# (ed.), Handbook of Noise and Vibration Control (2007), chapter 106, Table 13 +# and Figure 4. The octave rating is 5 dB below the curve at 500 Hz, and CI +# sums the bands up to 2 500 Hz (the first 15 thirds, all 5 octaves). +_IMPACT = { + 16: dict(reference=[62, 62, 62, 62, 62, 62, 61, 60, 59, 58, 57, 54, 51, + 48, 45, 42], + limit=320, at500=7, offset=0, upto2500=15), + 5: dict(reference=[67, 67, 65, 62, 49], + limit=100, at500=2, offset=-5, upto2500=5), +} + + +def impact(Ln): + """ + Weighted normalized impact sound pressure level and spectrum adaptation + term of ISO 717-2, from Ln, L'n or L'nT by band. + + Parameters: + ----------- + * Ln (array-like): + 16 one-third-octave bands from 100 Hz to 3 150 Hz, 5 octave + bands from 125 Hz to 2 000 Hz, or 19 one-third-octave bands from + 50 Hz to 3 150 Hz, in dB. + + Return: + ------- + * (Lnw, CI), integers in dB: Ln,w (or L'n,w, L'nT,w) and the + adaptation term for walking noise; (Lnw, CI, CI50) with the + 19 bands, CI50 being CI,50-2500. + + Ln is rounded to 0.1 dB, halves up. An unfavourable deviation is Ln above + the reference curve. The curve is shifted in steps of 1 dB until the sum + of unfavourable deviations is as large as possible but not more than + 32.0 dB in thirds or 10.0 dB in octaves; Ln,w is its value at 500 Hz, + less 5 dB in octaves. CI = Ln,sum - 15 - Ln,w, with Ln,sum the energetic + sum of the bands from 100 Hz (125 Hz in octaves) to 2 500 Hz; + CI,50-2500 sums from 50 Hz. Ln,w is rated on 100 Hz to 3 150 Hz only. + """ + tenths = _half_up(Ln, 0.1) + extended = len(tenths) == 19 # 50, 63 and 80 Hz come first + low, tenths = (tenths[:3], tenths[3:]) if extended else ([], tenths) + tables = _IMPACT.get(len(tenths)) + if tables is None: + raise ValueError('Ln must have 16 or 19 one-third-octave or 5 octave ' + 'bands.') + reference = np.array(tables['reference']) * 10 + # an impact deviation is the airborne one with both curves negated + shift = -_shift(-tenths, -reference, tables['limit']) + Lnw = tables['reference'][tables['at500']] + shift + tables['offset'] + bands = tenths[:tables['upto2500']] + ranges = (bands, np.r_[low, bands]) if extended else (bands,) + return (Lnw, *[int(_half_up(10 * np.log10(np.sum(10 ** (L / 100))), 1)) + - 15 - Lnw for L in ranges]) diff --git a/tests/test_iso717.py b/tests/test_iso717.py new file mode 100644 index 0000000..0b9ce4c --- /dev/null +++ b/tests/test_iso717.py @@ -0,0 +1,50 @@ +""" +Tests for pytta.iso717.airborne, the single-number rating of airborne +sound insulation of ISO 717-1:2013. Each test pins one rule. +""" +import os +import unittest + +os.environ.setdefault('MPLBACKEND', 'Agg') + +import pytta + +# ISO 717-1:2013, Annex C, Table C.1: sound reduction index of a building +# element, one-third-octave bands from 100 Hz to 3 150 Hz. +R_ANNEX_C = [20.4, 16.3, 17.7, 22.6, 22.4, 22.7, 24.8, 26.6, 28.0, 30.5, 31.8, + 32.5, 33.4, 33.0, 31.0, 25.5] + + +class TestISO717(unittest.TestCase): + + def test_example_of_annex_c(self): + # Rw(C;Ctr) = 30 (-2; -3) dB, with 31.8 dB of unfavourable deviation + self.assertEqual(pytta.iso717.airborne(R_ANNEX_C), (30, -2, -3)) + + def test_deviation_limit_is_inclusive(self): + # Octave bands: R 5 dB above the reference curve. Shifting it up by + # 7 dB gives 2 dB below R in each of the 5 bands, exactly 10.0 dB, + # which is allowed; 8 dB would give 15 dB. + reference = [36, 45, 52, 55, 56] + Rw, C, Ctr = pytta.iso717.airborne([r + 5 for r in reference]) + self.assertEqual(Rw, 52 + 7) + + def test_input_is_rounded_half_up_to_one_decimal(self): + # The standard rounds R to 0.1 dB before rating, halves up. These + # five bands end in 5, and round() of Python, which rounds the + # binary value, takes each of them 0.1 dB down: the deviations would + # sum 32.3 dB and Rw would drop to 29. Rounded half up, the curve is + # that of Annex C. + R = list(R_ANNEX_C) + for band, value in ((5, 22.65), (8, 27.95), (9, 30.45), (14, 30.95), + (15, 25.45)): + R[band] = value + self.assertEqual(pytta.iso717.airborne(R), (30, -2, -3)) + + def test_rejects_other_band_counts(self): + with self.assertRaises(ValueError): + pytta.iso717.airborne(R_ANNEX_C[:-1]) + + +if __name__ == '__main__': + unittest.main() diff --git a/tests/test_iso717_2.py b/tests/test_iso717_2.py new file mode 100644 index 0000000..7fa24a0 --- /dev/null +++ b/tests/test_iso717_2.py @@ -0,0 +1,100 @@ +""" +Tests for pytta.iso717.impact, the single-number rating of impact sound +insulation of ISO 717-2. Each test pins one rule. + +The text of ISO 717-2 is not at hand; the worked example is the one of +Crocker (ed.), Handbook of Noise and Vibration Control (2007), chapter 106, +Table 13, which follows ISO 717-2. +""" +import os +import unittest + +import numpy as np + +os.environ.setdefault('MPLBACKEND', 'Agg') + +import pytta + +REFERENCE_THIRDS = [62, 62, 62, 62, 62, 62, 61, 60, 59, 58, 57, 54, 51, 48, + 45, 42] +REFERENCE_OCTAVES = [67, 67, 65, 62, 49] + + +def unfavourable(Ln, reference, Lnw, at500): + """Sum of the deviations of Ln above the reference curve shifted to Lnw + at 500 Hz.""" + shifted = np.array(reference) - reference[at500] + Lnw + return np.maximum(np.array(Ln) - shifted, 0).sum() + + +class TestISO717Impact(unittest.TestCase): + + def test_example_of_crocker_table_13(self): + # Bare floor: Ln,w(CI) = 77 (-11) dB, 29.4 dB of unfavourable + # deviation. Floor with covering: Ln,w = 40 dB with 29.2 dB; the + # table prints CI = -3 dB, but its own Ln,sum = 57.9 dB gives + # 57.9 - 15 - 40 = +2.9 dB, and the floor covering term of the same + # table (-11 - Cr = -14, Cr = +3) confirms the sign is a misprint. + bare = [65.5, 66.5, 67, 67.4, 67.8, 68.2, 68.6, 69, 69.3, 69.6, 70, + 70.3, 70.5, 70.7, 71, 71.2] + covered = [56.2, 50, 47.4, 43.6, 39.1, 36.4, 33.2, 31.7, 31, 29.6, + 24.7, 23.5, 22.6, 24.3, 19.1, 14.3] + for Ln, expected, deviation in ((bare, (77, -11), 29.4), + (covered, (40, 3), 29.2)): + result = pytta.iso717.impact(Ln) + sum_ = unfavourable(Ln, REFERENCE_THIRDS, result[0], 7) + print(f'\n{result} vs {expected}; deviations {sum_:.1f} dB vs ' + f'{deviation} dB') + self.assertEqual(result, expected) + self.assertAlmostEqual(sum_, deviation, places=6) + + def test_octave_bands_rate_5_dB_below_the_curve_at_500_Hz(self): + # Ln 5 dB above the octave reference curve. Shifting it up by 3 dB + # leaves 2 dB in each of the 5 bands, exactly 10.0 dB, which is + # allowed; 2 dB would give 15 dB. The curve then reads 65 + 3 dB at + # 500 Hz and Ln,w is 5 dB less. + Ln = [r + 5 for r in REFERENCE_OCTAVES] + Lnw, CI = pytta.iso717.impact(Ln) + print(f'\nLn,w = {Lnw} dB, deviations ' + f'{unfavourable(Ln, REFERENCE_OCTAVES, Lnw + 5, 2):.1f} dB') + self.assertEqual(Lnw, 63) + + def test_octave_rating_within_1_dB_of_the_third_octave_rating(self): + # Crocker, chapter 106: the octave reference curve keeps the ratings + # from thirds and from octaves within 1 dB. The floors of Table 13, + # summed to octaves, rate 77 and 41 dB against 77 and 40 dB in + # thirds. This pins the octave curve, 49 dB at 2 kHz included: with + # 54 dB there the bare floor would rate 73 dB. + floors = ([65.5, 66.5, 67, 67.4, 67.8, 68.2, 68.6, 69, 69.3, 69.6, 70, + 70.3, 70.5, 70.7, 71, 71.2], + [56.2, 50, 47.4, 43.6, 39.1, 36.4, 33.2, 31.7, 31, 29.6, + 24.7, 23.5, 22.6, 24.3, 19.1, 14.3]) + for Ln in floors: + octaves = [round(10 * np.log10(np.sum(10**(np.array(Ln[i:i + 3]) + / 10))), 1) + for i in range(0, 15, 3)] + thirds, octave = (pytta.iso717.impact(Ln)[0], + pytta.iso717.impact(octaves)[0]) + print(f'\nLn,w {thirds} dB in thirds, {octave} dB in octaves') + self.assertLessEqual(abs(octave - thirds), 1) + + def test_adaptation_terms_sum_up_to_2500_Hz(self): + # Flat 59.6 dB from 50 Hz to 3 150 Hz (19 thirds). Shifting the curve + # to 66 dB leaves 28.4 dB of deviation (65 dB would give 33.0 dB). + # CI sums the 15 bands from 100 Hz to 2 500 Hz: + # 59.6 + 10 lg 15 - 15 - 66 = -9.64, so -10; with the 3 150 Hz band + # it would be -9.36, so -9. CI,50-2500 adds 50, 63 and 80 Hz: + # 59.6 + 10 lg 18 - 15 - 66 = -8.85, so -9. + result = pytta.iso717.impact([59.6] * 19) + CI = 59.6 + 10 * np.log10(15) - 15 - 66 + CI50 = 59.6 + 10 * np.log10(18) - 15 - 66 + print(f'\n{result}; CI {CI:.2f}, CI,50-2500 {CI50:.2f}') + self.assertEqual(result, (66, -10, -9)) + + def test_rejects_other_band_counts(self): + with self.assertRaises(ValueError): + pytta.iso717.impact([60.0] * 15) + + +if __name__ == '__main__': + unittest.main()