Merge pull request #71 from wantysal/sii

Speech intelligibility index (ANSI S3.5)

Author: Martin Glesser

CITATION.cff

@@ -4,9 +4,9 @@ authors:
   - family-names: Green Forge Coop
     given-names:
 title: MOSQITO
-version: 1.0.8
+version: 1.1
 doi: 10.5281/zenodo.6675733
-date-released: 2022-06-21
+date-released: 2023-12-21
 keywords:
   - audio
   - python

mosqito/__init__.py

@@ -1,4 +1,6 @@
 from mosqito.sound_level_meter.noct_spectrum.noct_spectrum import noct_spectrum
+from mosqito.sound_level_meter.noct_spectrum.noct_synthesis import noct_synthesis
+from mosqito.sound_level_meter.spectrum import spectrum
 
 from mosqito.sq_metrics.loudness.loudness_ecma.loudness_ecma import loudness_ecma
 from mosqito.sq_metrics.loudness.loudness_zwst.loudness_zwst import loudness_zwst
@@ -27,6 +29,10 @@
 from mosqito.sq_metrics.sharpness.sharpness_din.sharpness_din_perseg import sharpness_din_perseg
 from mosqito.sq_metrics.sharpness.sharpness_din.sharpness_din_freq import sharpness_din_freq
 
+from mosqito.sq_metrics.speech_intelligibility.sii_ansi.sii_ansi import sii_ansi
+from mosqito.sq_metrics.speech_intelligibility.sii_ansi.sii_ansi_freq import sii_ansi_freq
+from mosqito.sq_metrics.speech_intelligibility.sii_ansi.sii_ansi_level import sii_ansi_level
+
 from mosqito.sq_metrics.loudness.utils.sone_to_phon import sone_to_phon
 from mosqito.utils.isoclose import isoclose
 from mosqito.utils.load import load
@@ -36,8 +42,9 @@
 # Colors and linestyles
 COLORS = [
     "#69c3c5",
-    "#9969c4",
-    "#c46b69",
-    "#95c469",
-    "#2a6c6e",
+    "#ffd788",
+    "#ff8b88",
+    "#7894cf",
+    "#228080",
+    "#a8e2e2"
 ]

mosqito/sound_level_meter/__init__.py

@@ -1,3 +1,4 @@
 from mosqito.sound_level_meter.noct_spectrum.noct_spectrum import noct_spectrum
 from mosqito.sound_level_meter.noct_spectrum.noct_synthesis import noct_synthesis
-
+from mosqito.sound_level_meter.spectrum import spectrum
+from mosqito.sound_level_meter.freq_band_synthesis import freq_band_synthesis

mosqito/sound_level_meter/freq_band_synthesis.py

@@ -0,0 +1,74 @@
+# -*- coding: utf-8 -*-
+
+# Standard library import
+from numpy import array, concatenate, zeros, log10, power, argmin, split, arange, interp, iscomplex
+
+def freq_band_synthesis(spectrum, freqs, fmin, fmax):
+    """Adapt input spectrum to frequency band levels
+    
+    Convert the input spectrum to frequency band spectrum
+    between "fmin" and "fmax".
+    
+    Parameters
+    ----------
+    spectrum : numpy.ndarray
+        One-sided spectrum of the signal in [dB], size (nperseg, nseg).
+    freqs : list
+        List of input frequency , size (nperseg) or (nperseg, nseg).
+    fmin : float
+        Min frequency band [Hz].
+    fmax : float
+        Max frequency band [Hz].
+    n : int
+        Number of bands pr octave.
+    G : int
+        System for specifying the exact geometric mean frequencies.
+        Can be base 2 or base 10.
+    fr : int
+        Reference frequency. Shall be set to 1 kHz for audible frequency
+        range, to 1 Hz for infrasonic range (f < 20 Hz) and to 1 MHz for
+        ultrasonic range (f > 31.5 kHz).
+    Outputs
+    -------
+    spec : numpy.ndarray
+        Third octave band spectrum of signal sig [dB re.2e-5 Pa], size (nbands, nseg).
+    fpref : numpy.ndarray
+        Corresponding preferred third octave band center frequencies, size (nbands).
+    """
+    if iscomplex(spectrum).any():
+        raise ValueError('Input spectrum must be in dB, no complex value allowed.')
+    
+    if (fmin.min() < freqs.min()):
+        print("[WARNING] Input spectrum minimum frequency if higher than fmin. Empty values will be filled with 0.")
+        df = freqs[1] - freqs[0]
+        spectrum = interp(arange(fmin.min(),fmax.max()+df, df), freqs, spectrum)
+        freqs = arange(fmin.min(),fmax.max()+df, df)
+
+    if (fmax.max() > freqs.max()):
+        print("[WARNING] Input spectrum maximum frequency if lower than fmax. Empty values will be filled with 0.")
+        df = freqs[1] - freqs[0]
+        spectrum = interp(arange(fmin.min(),fmax.max()+df, df), freqs, spectrum)
+        freqs = arange(fmin.min(),fmax.max()+df, df)
+    
+    # Find the lower and upper index of each band
+    idx_low = argmin(abs(freqs[:,None] - fmin), axis=0)
+    idx_up = argmin(abs(freqs[:,None] - fmax), axis=0)
+    idx = concatenate((idx_low, [idx_up[-1]]))
+    
+    # Split the given spectrum in several bands
+    bands = array(split(spectrum, idx), dtype=object)[1:-1]
+    
+    # Compute the bands level
+    band_spectrum = zeros((len(bands)))
+    i = 0
+    for s in bands:
+        band_spectrum[i] = 10*log10(sum(power(10,s/10)))      
+        i += 1
+
+    return band_spectrum, (fmin+fmax)/2
+
+
+        
+
+
+

mosqito/sound_level_meter/noct_spectrum/noct_synthesis.py

@@ -11,6 +11,7 @@
 
 def noct_synthesis(spectrum, freqs, fmin, fmax, n=3, G=10, fr=1000):
     """Adapt input spectrum to nth-octave band spectrum
+    
     Convert the input spectrum to third-octave band spectrum
     between "fc_min" and "fc_max".
     Parameters

mosqito/sq_metrics/__init__.py

@@ -32,3 +32,7 @@
 from mosqito.sq_metrics.sharpness.sharpness_din.sharpness_din_freq import sharpness_din_freq
 
 from mosqito.sq_metrics.loudness.utils.sone_to_phon import sone_to_phon
+
+from mosqito.sq_metrics.speech_intelligibility.sii_ansi.sii_ansi import sii_ansi
+from mosqito.sq_metrics.speech_intelligibility.sii_ansi.sii_ansi_freq import sii_ansi_freq
+from mosqito.sq_metrics.speech_intelligibility.sii_ansi.sii_ansi_level import sii_ansi_level

mosqito/sq_metrics/loudness/loudness_zwst/_main_loudness.py

@@ -226,7 +226,7 @@ def _main_loudness(spec_third, field_type):
     # taking into account the dependance of absolute threshold
     # within this critical band
     korry = 0.4 + 0.32 * nm[0] ** 0.2
-    nm[0, korry <= 1] *= korry
+    nm[0, korry <= 1] *= korry[korry <= 1]
     # nm[:, -1] = 0
     # if spec_third.ndim == 1 or spec_third.shape[1] == 1:
     #     # This is only for test_loudness_zwicker_3oct because only one array of one col is given and this routine needs 2 or more