[FEATURE] add pytta.iso3382_3 with the parameters of open-plan offices

docs/pyttadoc/index.rst

@@ -64,4 +64,5 @@ Modules:
    functions
    generate
    rooms
+   iso3382_3
    properties

docs/pyttadoc/iso3382_3.rst

@@ -0,0 +1,7 @@
+.. currentmodule:: pytta
+
+ISO 3382-3
+==========
+
+.. automodule:: pytta.iso3382_3
+   :members:

pytta/__init__.py

@@ -72,6 +72,9 @@
     * pytta.iso3741:
         calculations according to the standard;
 
+    * pytta.iso3382_3:
+        open-plan offices (ISO 3382-3);
+
 
 You can find usage examples in the toolbox's examples folder.
 
@@ -108,6 +111,7 @@
 from . import utils
 from . import rooms
 from . import iso3741
+from . import iso3382_3
 
 from .functions import read_wav, write_wav, merge, split, \
     list_devices, print_devices, get_device_from_user,\
@@ -130,6 +134,7 @@
            'generate',
            'utils',
            'iso3741',
+           'iso3382_3',
            'rooms',
 
            # Functions

pytta/iso3382_3.py

@@ -0,0 +1,71 @@
+# -*- coding: utf-8 -*-
+"""
+Room acoustic parameters of open-plan offices, ISO 3382-3:2012 (read in its
+Brazilian adoption, ABNT NBR ISO 3382-3): the A-weighted speech level at
+each position from measurements with a calibrated omnidirectional source
+(Eq. 1 to 4), the spatial decay rate of speech D2,S and the level at 4 m
+Lp,A,S,4m (Eq. 5), and the distraction and privacy distances rD and rP from
+the speech transmission index along the measurement line.
+
+The STI at each position is computed separately, by the full method of IEC
+60268-16, e.g. with pytta.iec60268.sti.
+
+>>> LpAS = pytta.iso3382_3.speech_levels(Lw_source, Lp_source)
+>>> D2S, Lp4m = pytta.iso3382_3.spatial_decay(r, LpAS)
+>>> rD = pytta.iso3382_3.distance_at_sti(r, sti, 0.5)
+
+"""
+
+import numpy as np
+
+__all__ = ['OCTAVES', 'SPEECH_LW', 'SPEECH_LP_1M', 'A_WEIGHTING', 'source_level_1m',
+           'speech_levels', 'spatial_decay', 'distance_at_sti']
+
+# Table 1: normal unisex speech in the octave bands 125 Hz to 8 kHz.
+OCTAVES = np.array([125, 250, 500, 1000, 2000, 4000, 8000])
+SPEECH_LW = np.array([60.9, 65.3, 69.0, 63.0, 55.8, 49.8, 44.5])        # dB re 1 pW
+SPEECH_LP_1M = np.array([49.9, 54.3, 58.0, 52.0, 44.8, 38.8, 33.5])     # omnidirectional, dB
+A_WEIGHTING = np.array([-16.1, -8.6, -3.2, 0.0, 1.2, 1.0, -1.1])
+
+
+def source_level_1m(Lw):
+    """Free field level at 1 m of a source of octave band power levels Lw,
+    Lw + 10 lg(1 / (4 pi 1^2)), about Lw - 11 dB (Eq. 1)."""
+    return np.asarray(Lw, dtype=float) + 10 * np.log10(1 / (4 * np.pi))
+
+
+def speech_levels(Lw_source, Lp_source, speech=SPEECH_LP_1M):
+    """
+    A-weighted level of normal speech at each measurement position,
+    Lp,A,S,n (Eq. 2 to 4). The attenuation of each octave band from 1 m in
+    free field to position n, Dn,i = Lp,Ls,1m,i - Lp,Ls,n,i, measured with
+    the source of power levels Lw_source (7 octave bands), is applied to the
+    speech levels at 1 m of Table 1, and the bands are A-weighted and summed.
+    Lp_source holds the levels of the source at the positions, shape (N, 7).
+    """
+    D = source_level_1m(Lw_source) - np.atleast_2d(np.asarray(Lp_source, dtype=float))
+    LpS = np.asarray(speech, dtype=float) - D
+    return 10 * np.log10(np.sum(10 ** ((LpS + A_WEIGHTING) / 10), axis=1))
+
+
+def spatial_decay(r, LpAS, rmin=2., rmax=16.):
+    """
+    Spatial decay rate of speech per distance doubling D2,S and the speech
+    level at 4 m Lp,A,S,4m from the regression line of Lp,A,S,n on lg(r / 1 m)
+    over the positions between rmin and rmax (Eq. 5, 3.3 and 6.2).
+    """
+    r, L = np.asarray(r, dtype=float), np.asarray(LpAS, dtype=float)
+    keep = (r >= rmin) & (r <= rmax)
+    slope, intercept = np.polyfit(np.log10(r[keep]), L[keep], 1)
+    return -np.log10(2) * slope, intercept + slope * np.log10(4.)
+
+
+def distance_at_sti(r, sti, limit=0.5):
+    """
+    Distance at which the regression line of the STI on the distance (linear
+    axis) falls to limit: the distraction distance rD for 0.50 and the
+    privacy distance rP for 0.20 (3.6, 3.7 and 6.3). nan when the line does
+    not fall.
+    """
+    slope, intercept = np.polyfit(np.asarray(r, dtype=float), np.asarray(sti, dtype=float), 1)
+    return (limit - intercept) / slope if slope < -1e-12 else np.nan

tests/test_iso3382_3.py

@@ -0,0 +1,52 @@
+"""
+Tests for pytta.iso3382_3 against ISO 3382-3:2012: the values of Table 1 and
+analytic sound fields along a measurement line.
+"""
+import unittest
+
+import numpy as np
+
+from pytta import iso3382_3 as iso
+
+R = np.array([2., 3., 4., 6., 8., 11., 16.])
+
+
+class TestISO338203(unittest.TestCase):
+
+    def test_table_1_and_equation_1(self):
+        # The omnidirectional speech of Table 1 sums to 57.4 dB(A), and the
+        # level at 1 m is about Lw - 11 dB.
+        self.assertAlmostEqual(10 * np.log10(np.sum(10 ** ((iso.SPEECH_LP_1M + iso.A_WEIGHTING) / 10))),
+                               57.4, places=1)
+        self.assertAlmostEqual(float(iso.source_level_1m(100.)), 100 - 10.99, places=2)
+
+    def test_free_field_decays_6_db_per_doubling(self):
+        # A source of 90 dB power in every band in free field: the speech
+        # level falls by 20 lg 2 = 6.02 dB per doubling of the distance and
+        # is 57.4 - 20 lg 4 = 45.36 dB(A) at 4 m.
+        Lp = iso.source_level_1m(np.full(7, 90.)) - 20 * np.log10(R)[:, None]
+        LpAS = iso.speech_levels(np.full(7, 90.), Lp)
+        D2S, L4 = iso.spatial_decay(R, LpAS)
+        L1 = 10 * np.log10(np.sum(10 ** ((iso.SPEECH_LP_1M + iso.A_WEIGHTING) / 10)))
+        print(f'\n  free field: D2,S {D2S:.3f} dB, Lp,A,S,4m {L4:.3f} dB', end='')
+        self.assertAlmostEqual(D2S, 20 * np.log10(2), places=10)
+        self.assertAlmostEqual(L4, L1 - 20 * np.log10(4), places=10)
+
+    def test_positions_outside_2_to_16_m_are_left_out(self):
+        # A wall reflection that raises the level at 20 m does not change
+        # D2,S, which uses 2 m to 16 m only.
+        r = np.r_[R, 20.]
+        L = 60 - 5 * np.log2(r)
+        L[-1] += 10
+        self.assertAlmostEqual(iso.spatial_decay(r, L)[0], 5., places=10)
+
+    def test_distraction_and_privacy_distances(self):
+        # STI = 0.8 - 0.05 r: 0.50 at 6 m and 0.20 at 12 m.
+        sti = 0.8 - 0.05 * R
+        self.assertAlmostEqual(iso.distance_at_sti(R, sti, 0.5), 6., places=10)
+        self.assertAlmostEqual(iso.distance_at_sti(R, sti, 0.2), 12., places=10)
+        self.assertTrue(np.isnan(iso.distance_at_sti(R, np.full(R.size, 0.6), 0.5)))
+
+
+if __name__ == '__main__':
+    unittest.main()