Repository: pyfar/pyfar
Pull request #734 (2024-12-14 12:24:30)

add air attenuation after iso
- [x] implement air attenuataion after ISO
- [x] separate vapour pressure calcualtion into private function (there is just a sligt difference between this and the magnus formula in #770) 
![grafik](https://github.com/user-attachments/assets/47e7ec54-5640-4aed-9d48-4628745301e3)


ahms5: In the standard, there is a table of results that I use for comparison during the tests. The Magnus formula in the ITA-Toolbox and AKtools was not well documented, as it referenced the manual of a humidity measuring device. Therefore, I searched for a more reliable reference.
ahms5: @CristobalAndrade pointed out that ISO Appendix B provides a defined method for calculating saturation vapor pressure. Therefore, I removed the old saturation_vapor_pressure function (still required for ideal gas in #770) and incorporated the calculation directly into the air attenuation function, following the standard's guidelines. 

Patch:
diff --git a/pyfar/classes/transmission_matrix.py b/pyfar/classes/transmission_matrix.py
index 7f086562..307718cc 100644
--- a/pyfar/classes/transmission_matrix.py
+++ b/pyfar/classes/transmission_matrix.py
@@ -32,6 +32,7 @@ from __future__ import annotations # required for Python <= 3.9
 import numpy as np
 import numpy.testing as npt
 from pyfar.classes.audio import FrequencyData
+from numbers import Number
 
 
 class TransmissionMatrix(FrequencyData):
@@ -739,6 +740,92 @@ class TransmissionMatrix(FrequencyData):
         frequencies = transducer_constant.frequencies
         return TransmissionMatrix.from_abcd(0, B, C, 0, frequencies)
 
+    @staticmethod
+    def create_transmission_line(
+            kl: FrequencyData,
+            characteristic_impedance: Number | FrequencyData,
+            ) -> TransmissionMatrix:
+        r"""Create a transmission matrix representing a transmission line.
+
+        The transmission matrix is calculated from the product :math:`kl` of
+        the wavenumber :math:`k` and length :math:`l`,
+        and characteristic impedance :math:`Z_0` following
+        Equation (5-10) of Reference [1]_:
+
+        .. math::
+            T = \begin{bmatrix}
+                \cos{kl} & j Z_0 \sin{kl} \\
+                j/Z_0 \sin{kl} & \cos{kl}
+                \end{bmatrix}
+
+        The inputs need to be broadcastable into one `cshape`.
+
+        Parameters
+        ----------
+        kl : FrequencyData
+            The product of wavenumber and length,
+            used as argument in the cosine and sine functions to calculate the
+            entries of the matrix.
+            If all imaginary parts of ``kl.freq`` are zero, the transmission
+            matrix represents a lossless transmission line.
+        characteristic_impedance : Number | FrequencyData
+            The characteristic impedance of the transmission line.
+
+        Returns
+        -------
+        tmat : TransmissionMatrix
+            A transmission matrix representing a transmission line.
+
+
+        Example
+        -------
+        Lossless acoustic duct with approximate ear canal dimensions.
+        Plot its input impedance with rigid termination.
+
+        .. plot::
+
+            >>> import pyfar as pf
+            >>> import numpy as np
+            >>> # Transmission line parameters
+            >>> frequencies = np.linspace(20, 20e3, 1000)
+            >>> omega = 2*np.pi*frequencies
+            >>> wavenumber = omega/pf.constants.reference_speed_of_sound
+            >>> length = 30e-3
+            >>> kl = pf.FrequencyData(wavenumber*length, frequencies)
+            >>> cross_section = np.pi*(3.5e-3)**2
+            >>> Z0 = pf.constants.reference_air_impedance/cross_section
+            >>> # Create transmission matrix
+            >>> T = pf.TransmissionMatrix.create_transmission_line(kl, Z0)
+            >>> # Plot input impedance with rigid termination
+            >>> ax = pf.plot.freq(T.input_impedance(np.inf))
+            >>> ax.set_ylim(100, 200)
+
+        """
+        if not isinstance(kl, FrequencyData):
+            raise ValueError("The input kl must be a FrequencyData object.")
+        if isinstance(characteristic_impedance, FrequencyData):
+            if not np.allclose(
+                    kl.frequencies,
+                    characteristic_impedance.frequencies, atol=1e-15):
+                raise ValueError(
+                    "The frequencies do not match.")
+            else:
+                characteristic_impedance = characteristic_impedance.freq
+        elif not isinstance(characteristic_impedance, Number):
+            raise ValueError(
+                "characteristic impedance must be a Number or " \
+                "FrequencyData object.")
+        # broadcastable shapes are checked by from_abcd
+
+        A = np.cos(kl.freq)
+        B = 1j*characteristic_impedance*np.sin(kl.freq)
+        C = 1j/characteristic_impedance*np.sin(kl.freq)
+        D = A.copy()
+
+        return TransmissionMatrix.from_abcd(
+            A, B, C, D, kl.frequencies)
+
+
     def __repr__(self):
         """String representation of TransmissionMatrix class."""
         repr_string = (
diff --git a/pyfar/constants/__init__.py b/pyfar/constants/__init__.py
index 3bdca282..cca81ad0 100644
--- a/pyfar/constants/__init__.py
+++ b/pyfar/constants/__init__.py
@@ -11,9 +11,14 @@ from .speed_of_sound import (
     speed_of_sound_ideal_gas,
 )
 
+from .constants import (
+    air_attenuation,
+)
+
 __all__ = [
     'saturation_vapor_pressure_magnus',
     'speed_of_sound_simple',
+    'air_attenuation',
     'speed_of_sound_ideal_gas',
 ]
 
diff --git a/pyfar/constants/constants.py b/pyfar/constants/constants.py
index 8db6fb3d..ff80ffc3 100644
--- a/pyfar/constants/constants.py
+++ b/pyfar/constants/constants.py
@@ -1,5 +1,315 @@
 """Constant calculation."""
 import numpy as np
+import pyfar as pf
+
+
+def air_attenuation(
+        temperature, frequencies, relative_humidity,
+        atmospheric_pressure=None):
+    r"""Calculate the pure tone attenuation of sound in air according to
+    ISO 9613-1.
+
+    Calculation is in accordance with ISO 9613-1 [#]_. The shape of the
+    outputs is broadcasted from the shapes of the ``temperature``,
+    ``relative_humidity``, and ``atmospheric_pressure``.
+
+    Parameters
+    ----------
+    temperature : float, array_like
+        Temperature in degree Celsius.
+        Must be in the range of -20°C to 50°C for accuracy of +/-10% or
+        must be greater than -70°C for accuracy of +/-50%.
+    frequencies : float, array_like
+        Frequency in Hz. Must be greater than 50 Hz.
+        Just one dimensional array is allowed.
+    relative_humidity : float, array_like
+        Relative humidity in the range from 0 to 1.
+    atmospheric_pressure : float, array_like, optional
+        Atmospheric pressure in Pascal, by default
+        :py:attr:`reference_atmospheric_pressure`.
+
+    Returns
+    -------
+    alpha : :py:class:`~pyfar.FrequencyData`
+        Pure tone air attenuation coefficient in decibels per meter for
+        atmospheric absorption.
+    m : :py:class:`~pyfar.FrequencyData`
+        Pure tone air attenuation coefficient per meter for
+        atmospheric absorption. The parameter ``m`` is calculated as
+        :math:`m = \alpha / (10 \cdot \log_{10}(e))`.
+    accuracy : :py:class:`~pyfar.FrequencyData`
+        accuracy of the results according to the standard:
+
+        ``10``, +/- 10% accuracy
+            - molar concentration of water vapour: 0.05% to 5 %.
+            - air temperature: 253.15 K to 323.15 (-20 °C to +50°C)
+            - atmospheric pressure: less than 200 000 Pa (2 atm)
+            - frequency-to-pressure ratio: 0.0004 Hz/Pa to 10 Hz/Pa.
+
+        ``20``, +/- 20% accuracy
+            - molar concentration of water vapour: 0.005 % to 0.05 %,
+              and greater than 5%
+            - air temperature: 253.15 K to 323.15 (-20 °C to +50°C)
+            - atmospheric pressure: less than 200 000 Pa (2 atm)
+            - frequency-to-pressure ratio: 0.0004 Hz/Pa to 10 Hz/Pa.
+
+        ``50``, +/- 50% accuracy
+            - molar concentration of water vapour: less than 0.005%
+            - air temperature: greater than 200 K (- 73 °C)
+            - atmospheric pressure: less than 200 000 Pa (2 atm)
+            - frequency-to-pressure ratio: 0.0004 Hz/Pa to 10 Hz/Pa.
+
+        ``-1``, no valid result
+            else.
+
+    References
+    ----------
+    .. [#] ISO 9613-1:1993, Acoustics -- Attenuation of sound during
+           propagation outdoors -- Part 1: Calculation of the absorption of
+           sound by the atmosphere.
+    """
+    if atmospheric_pressure is None:
+        atmospheric_pressure = pf.constants.reference_atmospheric_pressure
+    # check inputs
+    if not isinstance(temperature, (int, float, np.ndarray, list, tuple)):
+        raise TypeError(
+            'temperature must be a number or array of numbers')
+    if not isinstance(frequencies, (int, float, np.ndarray, list, tuple)):
+        raise TypeError(
+            'frequencies must be a number or array of numbers')
+    if not isinstance(
+            relative_humidity, (int, float, np.ndarray, list, tuple)):
+        raise TypeError(
+            'relative_humidity must be a number or array of numbers')
+    if np.array(frequencies).ndim > 1:
+        raise ValueError('frequencies must be one dimensional.')
+    if not isinstance(
+            atmospheric_pressure, (int, float, np.ndarray, list, tuple)):
+        raise TypeError(
+            'atmospheric_pressure must be a number or array of numbers')
+
+    # check if broadcastable
+    try:
+        _ = np.broadcast_shapes(
+            np.atleast_1d(temperature).shape,
+            np.atleast_1d(relative_humidity).shape,
+            np.atleast_1d(atmospheric_pressure).shape)
+    except ValueError as e:
+        raise ValueError(
+            'temperature, relative_humidity, and atmospheric_pressure must '
+            'have the same shape or be broadcastable.') from e
+
+    # check limits
+    if np.any(np.array(temperature) < -73):
+        raise ValueError("Temperature must be greater than -73°C.")
+    if np.any(np.array(frequencies) < 50):
+        raise ValueError("frequencies must be greater than 50 Hz.")
+    if np.any(np.array(relative_humidity) < 0) or np.any(
+            np.array(relative_humidity) > 1):
+        raise ValueError("Relative humidity must be between 0 and 1.")
+    if np.any(np.array(atmospheric_pressure) > 200000):
+        raise ValueError("Atmospheric pressure must be less than 200 kPa.")
+    # convert arrays
+    temperature = np.array(
+        temperature, dtype=float)[..., np.newaxis]
+    relative_humidity = np.array(
+        relative_humidity, dtype=float)[..., np.newaxis]
+    atmospheric_pressure = np.array(
+        atmospheric_pressure, dtype=float)[..., np.newaxis]
+    frequencies = np.array(frequencies, dtype=float)
+
+    # calculate air attenuation
+    p_atmospheric_ref = pf.constants.reference_atmospheric_pressure
+    t_degree_ref = pf.constants.reference_air_temperature_celsius
+
+    h_r = relative_humidity*100
+    p_a = atmospheric_pressure
+    p_r = p_atmospheric_ref
+    f = frequencies
+    T = temperature - pf.constants.absolute_zero_celsius
+    T_0 = t_degree_ref - pf.constants.absolute_zero_celsius
+
+    # saturation vapour pressure
+    p_sat = _saturation_vapour_pressure_iso(temperature)
+
+    # molar concentration of water vapor as a percentage (Equation B.1)
+    h = h_r * (p_sat / p_r) * (p_a / p_r)
+
+    # Oxygen relaxation frequency (Eq. 3)
+    f_rO = (p_a/p_r)*(24+4.04e4*h*(0.02+h)/(0.391+h))
+
+    # Nitrogen relaxation frequency (Eq. 4)
+    f_rN = (p_a/p_r)*(T/T_0)**(-1/2)*(9+280*h*np.exp(
+        -4.17*((T/T_0)**(-1/3)-1)))
+
+    # air attenuation (Eq. 5)
+    air_attenuation = 8.686*f**2*((1.84e-11*p_r/p_a*(T/T_0)**(1/2)) + \
+        (T/T_0)**(-5/2)*(0.01275*np.exp(-2239.1/T)*(f_rO + (f**2/f_rO))**(-1)
+        +0.1068*np.exp(-3352/T) * (f_rN + (f**2/f_rN))**(-1)))
+
+    alpha = pf.FrequencyData(
+        air_attenuation, frequencies=frequencies)
+
+    # Equation 3: ISO 17497-1:2004
+    m = alpha / (10*np.log10(np.exp(1)))
+
+    # calculate accuracy
+    accuracy = _air_attenuation_accuracy(
+        h, temperature, atmospheric_pressure, frequencies, alpha.freq.shape)
+
+    return alpha, m, accuracy
+
+
+def _air_attenuation_accuracy(
+        concentration_water_vapour, temperature, atmospheric_pressure,
+        frequencies, shape):
+    """Calculate the accuracy of the air attenuation calculation.
+
+    Calculation is in accordance with ISO 9613-1 [#]_. This method is used in
+    :py:func:`~pyfar.constants.air_attenuation` to calculate the accuracy of
+    the air attenuation calculation.
+
+    Parameters
+    ----------
+    concentration_water_vapour : float, array_like
+        Molar concentration of water vapor as a percentage.
+        Must be between 0% and 100%.
+    temperature : float, array_like
+        Temperature in degree Celsius.
+        Must be above -273.15°C.
+    atmospheric_pressure : float, array_like
+        Atmospheric pressure in Pascal.
+        Must be above 0Pa.
+    frequencies : float, array_like
+        Frequency in Hz.
+        Must be larger than 0 Hz.
+    shape : tuple
+        Shape of the output.
+
+    Returns
+    -------
+    accuracy : :py:class:`~pyfar.classes.FrequencyData`
+        accuracy of the results according to the standard:
+
+            ``10``, +/- 10% accuracy
+                - molar concentration of water vapour: 0.05% to 5 %.
+                - air temperature: 253.15 K to 323.15 (-20 °C to +50°C)
+                - atmospheric pressure: less than 200 000 Pa (2 atm)
+                - frequency-to-pressure ratio: 4 x 10-4 Hz/Pa to 10 Hz/Pa.
+
+            ``20``, +/- 20% accuracy
+                - molar concentration of water vapour: 0.005 % to 0.05 %,
+                  and greater than 5%
+                - air temperature: 253.15 K to 323.15 (-20 °C to +50°C)
+                - atmospheric pressure: less than 200 000 Pa (2 atm)
+                - frequency-to-pressure ratio: 4 x 10-4 Hz/Pa to 10 Hz/Pa.
+
+            ``50``, +/- 50% accuracy
+                - molar concentration of water vapour: less than 0.005%
+                - air temperature: greater than 200 K (- 73 °C)
+                - atmospheric pressure: less than 200 000 Pa (2 atm)
+                - frequency-to-pressure ratio: 4 x 10-4 Hz/Pa to 10 Hz/Pa.
+
+            ``-1``, no valid result
+                else.
+
+    References
+    ----------
+    .. [#] ISO 9613-1:1993, Acoustics -- Attenuation of sound during
+           propagation outdoors -- Part 1: Calculation of the absorption of
+           sound by the atmosphere.
+    """
+    if np.any(np.array(concentration_water_vapour) < 0) or np.any(
+            np.array(concentration_water_vapour) > 100):
+        raise ValueError(
+            r"Concentration of water vapour must be between 0% and 100%.")
+    if np.any(np.array(temperature) < -273.15):
+        raise ValueError(
+            "Temperature must be greater than -273.15°C.")
+    if np.any(np.array(atmospheric_pressure) < 0):
+        raise ValueError(
+            "Atmospheric pressure must be greater than 0 Pa.")
+    if np.any(np.array(frequencies) < 0):
+        raise ValueError(
+            "Frequencies must be positive.")
+
+    # broadcast inputs
+    atmospheric_pressure = np.broadcast_to(atmospheric_pressure, shape)
+    h_water_vapor = np.broadcast_to(concentration_water_vapour, shape)
+    frequency_pressure_ratio = frequencies/atmospheric_pressure
+    accuracy = np.zeros(shape) - 1
+
+    # atmospheric pressure < 200 kPa
+    atm_mask = atmospheric_pressure < 200000
+    # frequency-to-pressure ratio: 4 x 10-4 Hz/Pa to 10 Hz/Pa
+    frequency_pressure_ratio_mask = (4e-4 <= frequency_pressure_ratio) & (
+        frequency_pressure_ratio <= 10)
+    common_mask = atm_mask & frequency_pressure_ratio_mask
+
+    # molar concentration of water vapour: 0.05% to 5 %
+    vapor_10_mask = (0.05 <= h_water_vapor) & (h_water_vapor <= 5)
+    # molar concentration of water vapour: 0.005% to 0.05 % and greater than 5%
+    vapor_20_mask = (5 < h_water_vapor) | (
+        (0.005 <= h_water_vapor) & (h_water_vapor < 0.05))
+    # molar concentration of water vapour: less than 0.005%
+    vapor_50_mask = (0.005 > h_water_vapor)
+
+    # air temperature: 253,15 K to 323,15 (-20 °C to +50°C)
+    temp_20_mask = (-20 <= temperature) & (temperature <= 50)
+    # air temperature: greater than 200 K (- 73 °C)
+    temp_50_mask = (-73 <= temperature)
+
+    # apply masks
+    accuracy_50 = common_mask & temp_50_mask & (
+        vapor_10_mask | vapor_20_mask | vapor_50_mask)
+    accuracy[accuracy_50] = 50
+
+    accuracy_20 = common_mask & temp_20_mask & (
+        vapor_10_mask | vapor_20_mask)
+    accuracy[accuracy_20] = 20
+
+    accuracy_10 = vapor_10_mask & common_mask & temp_20_mask
+    accuracy[accuracy_10] = 10
+
+    # return FrequencyData object
+    return pf.FrequencyData(accuracy, frequencies=frequencies)
+
+
+def _saturation_vapour_pressure_iso(temperature):
+    """Calculates the saturation vapour pressure after ISO 9613-1:1993.
+
+    This method is used in the :py:func:`air_attenuation` function to calculate
+    the saturation vapour pressure of water vapour in air as a close
+    approximation to
+    those calculated by the World Meteorological Organization.
+    The method is described in Equation B.2 and B.3 in [#]_.
+
+    Parameters
+    ----------
+    temperature : float, array_like
+        Temperature in degree Celsius.
+        Must be in the range of -20°C to 50°C for accuracy of +/-10% or
+        must be greater than -70°C for accuracy of +/-50%.
+
+    Returns
+    -------
+    p_sat : float, array_like
+        Saturation vapour pressure in Pascal.
+
+    References
+    ----------
+    .. [#] ISO 9613-1:1993, Acoustics -- Attenuation of sound during
+           propagation outdoors -- Part 1: Calculation of the absorption of
+           sound by the atmosphere.
+    """
+    T = temperature - pf.constants.absolute_zero_celsius
+    p_r = pf.constants.reference_atmospheric_pressure
+    # triple-point isotherm temperature of 273.16 K
+    T_01 = 273.16
+
+    # saturation vapour pressure (Equation B.2 and B.3)
+    C = -6.8346*(T_01/T)**1.261+4.6151
+    return 10**C * p_r
 
 
 def saturation_vapor_pressure_magnus(temperature):


Test patch:
diff --git a/tests/test_constants_air_attenuation.py b/tests/test_constants_air_attenuation.py
new file mode 100644
index 00000000..835e153e
--- /dev/null
+++ b/tests/test_constants_air_attenuation.py
@@ -0,0 +1,247 @@
+import pytest
+import numpy as np
+import pyfar as pf
+import numpy.testing as npt
+import re
+
+
+@pytest.mark.parametrize((
+        "temperature", "frequency", "relative_humidity", "expected",
+        "expected_accuracy"), [
+    (10, 1000, .1, 2.16e1*1e-3, 10),
+    (10, 100, .1, 5.85e-1*1e-3, 10),
+    (-20, 800, .3, 4.92e-3, 20),
+    (35, 4000, .6, 2.58*10e-3, 10),
+    (20, 5000, .1, 1.33e2*1e-3, 10),
+    (20, 5000, .8, 3.06e1*1e-3, 10),
+    (50, 1000, .7, 8.03e-3, 20),
+])
+def test_air_attenuation_after_table_iso(
+        temperature, frequency, relative_humidity, expected,
+        expected_accuracy):
+    """
+    Check alpha after value table in ISO 9613-1. The accuracy is calculated
+    manually and not provided by the standard.
+    """
+    alpha, m, accuracy = pf.constants.air_attenuation(
+        temperature, frequency, relative_humidity)
+    npt.assert_allclose(alpha.freq, expected, atol=1e-3)
+    npt.assert_allclose(accuracy.freq, expected_accuracy)
+
+
+@pytest.mark.parametrize("temperature", [
+        np.array([[10, 10, 10, 10]]).T,
+        [10, 10],
+        10,
+])
+@pytest.mark.parametrize("frequency", [
+        np.array([1000, 20000, 40000]),
+        [1000, 20000, 40000],
+        1000,
+])
+@pytest.mark.parametrize("relative_humidity", [
+        np.array([[.1, .1, .1, .1]]).T,
+        [.1, .1],
+        .1,
+])
+@pytest.mark.parametrize("atmospheric_pressure", [
+        np.array([[101325, 101325, 101325, 101325]]).T,
+        [101325, 101325],
+        101325,
+])
+def test_air_attenuation_array(
+        temperature, frequency, relative_humidity, atmospheric_pressure):
+    alpha, m, accuracy = pf.constants.air_attenuation(
+        temperature, frequency, relative_humidity, atmospheric_pressure)
+    # test air attenuation alpha
+    expected = 2.16e1*1e-3 + np.zeros_like(alpha.freq[..., 0])
+    npt.assert_allclose(alpha.freq[..., 0], expected, atol=1e-3)
+    npt.assert_allclose(alpha.frequencies, frequency, atol=1e-3)
+    # test air attenuation factor m
+    expected = 2.16e1*1e-3/4.34 + np.zeros_like(alpha.freq[..., 0])
+    npt.assert_allclose(m.freq[..., 0], expected, atol=1e-3)
+    npt.assert_allclose(m.frequencies, frequency, atol=1e-3)
+    # test accuracy
+    expected_accuracy = 10 + np.zeros_like(accuracy.freq)
+    npt.assert_allclose(accuracy.freq, expected_accuracy, atol=1e-3)
+    npt.assert_allclose(accuracy.frequencies, frequency, atol=1e-3)
+
+
+def test_air_attenuation_inputs():
+    temperature = 10
+    frequency = 1000
+    relative_humidity = .1
+    with pytest.raises(TypeError, match='must be a number or'):
+        pf.constants.air_attenuation(
+            'test', frequency, relative_humidity)
+    with pytest.raises(TypeError, match='must be a number or'):
+        pf.constants.air_attenuation(
+            temperature, 'frequency', relative_humidity)
+    with pytest.raises(TypeError, match='must be a number or'):
+        pf.constants.air_attenuation(
+            temperature, frequency, 'relative_humidity')
+    with pytest.raises(
+            ValueError, match='frequencies must be one dimensional'):
+        pf.constants.air_attenuation(
+            temperature, [[1, 1]], relative_humidity)
+    with pytest.raises(
+            TypeError,
+            match='atmospheric_pressure must be a number or array of numbers'):
+        pf.constants.air_attenuation(
+            temperature, frequency, relative_humidity, 'test')
+
+
+def test_air_attenuation_broadcastable():
+    frequency = 1000
+    with pytest.raises(ValueError, match='same shape or be broadcastable'):
+        pf.constants.air_attenuation(
+            [10, 10], frequency, [.1, .1, .1])
+
+
+def test_air_attenuation_limits():
+    temperature = 10
+    frequency = 1000
+    relative_humidity = .1
+    atmospheric_pressure = 101325
+    match = 'Temperature must be greater than -73°C.'
+    with pytest.raises(ValueError,match=match):
+        pf.constants.air_attenuation(
+            -100, frequency, relative_humidity, atmospheric_pressure)
+    match = 'frequencies must be greater than 50 Hz'
+    with pytest.raises(ValueError,match=match):
+        pf.constants.air_attenuation(
+            temperature, 20, relative_humidity, atmospheric_pressure)
+    match = 'Relative humidity must be between 0 and 1.'
+    with pytest.raises(ValueError, match=match):
+        pf.constants.air_attenuation(
+            temperature, frequency, -.1, atmospheric_pressure)
+    with pytest.raises(ValueError, match=match):
+        pf.constants.air_attenuation(
+            temperature, frequency, 1.1, atmospheric_pressure)
+    match = 'Atmospheric pressure must be less than 200 kPa.'
+    with pytest.raises(ValueError, match=match):
+        pf.constants.air_attenuation(
+            temperature, frequency, relative_humidity, 200001)
+
+
+@pytest.mark.parametrize((
+        "concentration_water_vapour", "expected_accuracy"), [
+    (0.05, 10),
+    (5, 10),
+    (1, 10),
+    (10, 20),
+    (0.005, 20),
+    (0.006, 20),
+    (0.0006, 50),
+    (100, 20),
+])
+def test_accuracy_concentration_water_vapour(
+        concentration_water_vapour, expected_accuracy):
+    shape = (1,)
+    temperature = 10
+    atmospheric_pressure = 101325
+    frequencies = 1e3
+    accuracy = pf.constants.constants._air_attenuation_accuracy(
+        concentration_water_vapour, temperature, atmospheric_pressure,
+        frequencies, shape)
+    npt.assert_array_equal(accuracy.freq, expected_accuracy)
+
+
+@pytest.mark.parametrize((
+        "temperature", "expected_accuracy"), [
+    (-20, 10),
+    (50, 10),
+    (0, 10),
+    (-100, -1),
+    (-50, 50),
+    (-73, 50),
+    (100, 50),
+])
+def test_accuracy_temperature(
+        temperature, expected_accuracy):
+    shape = (1,)
+    concentration_water_vapour = 1
+    atmospheric_pressure = 101325
+    frequencies = 1e3
+    accuracy = pf.constants.constants._air_attenuation_accuracy(
+        concentration_water_vapour, temperature, atmospheric_pressure,
+        frequencies, shape)
+    npt.assert_array_equal(accuracy.freq, expected_accuracy)
+
+
+@pytest.mark.parametrize((
+        "atmospheric_pressure", "expected_accuracy"), [
+    (101325, 10),
+    (200000-1, 10),
+    (200000, -1),
+    (300000, -1),
+])
+def test_accuracy_atmospheric_pressure(
+        atmospheric_pressure, expected_accuracy):
+    shape = (1,)
+    temperature = 0
+    concentration_water_vapour = 1
+    frequencies = 1e3
+    accuracy = pf.constants.constants._air_attenuation_accuracy(
+        concentration_water_vapour, temperature, atmospheric_pressure,
+        frequencies, shape)
+    npt.assert_array_equal(accuracy.freq, expected_accuracy)
+
+
+
+@pytest.mark.parametrize((
+        "frequency_pressure_ratio", "expected_accuracy"), [
+    (4e-4, 10),
+    (10, 10),
+    (1, 10),
+    (3e-4, -1),
+    (11, -1),
+])
+def test_accuracy_frequency_pressure_ratio(
+        frequency_pressure_ratio, expected_accuracy):
+    shape = (1,)
+    temperature = 0
+    concentration_water_vapour = 1
+    atmospheric_pressure = 101325
+    frequencies = frequency_pressure_ratio * atmospheric_pressure
+    accuracy = pf.constants.constants._air_attenuation_accuracy(
+        concentration_water_vapour, temperature, atmospheric_pressure,
+        frequencies, shape)
+    npt.assert_array_equal(accuracy.freq, expected_accuracy)
+
+
+def test_accuracy_invalid_vapor():
+    match = r"Concentration of water vapour must be between 0% and 100%."
+    with pytest.raises(
+            ValueError, match=re.escape(match)):
+        pf.constants.constants._air_attenuation_accuracy(
+            101, 20, 101325, 1000, (1,))
+    with pytest.raises(
+            ValueError, match=re.escape(match)):
+        pf.constants.constants._air_attenuation_accuracy(
+            -1, 20, 101325, 1000, (1,))
+
+
+def test_accuracy_invalid_pressure():
+    match = "Atmospheric pressure must be greater than 0 Pa."
+    with pytest.raises(
+            ValueError, match=re.escape(match)):
+        pf.constants.constants._air_attenuation_accuracy(
+            .05, 20, -101325, 1000, (1,))
+
+
+def test_accuracy_invalid_temperature():
+    match = "Temperature must be greater than -273.15°C."
+    with pytest.raises(
+            ValueError, match=re.escape(match)):
+        pf.constants.constants._air_attenuation_accuracy(
+            .05, -280, 101325, 1000, (1,))
+
+
+def test_accuracy_invalid_frequency():
+    match = "Frequencies must be positive."
+    with pytest.raises(
+            ValueError, match=re.escape(match)):
+        pf.constants.constants._air_attenuation_accuracy(
+            .05, 20, 101325, -1, (1,))
+
diff --git a/tests/test_transmission_matrix.py b/tests/test_transmission_matrix.py
index 121c7d63..add12c72 100644
--- a/tests/test_transmission_matrix.py
+++ b/tests/test_transmission_matrix.py
@@ -269,6 +269,64 @@ def test_tmatrix_create_gyrator(transducer_constant, frequencies):
     Zin = tmat.input_impedance(Zl)
     npt.assert_allclose(Zin.freq, Zin_expected, atol = 1e-15)
 
+def test_create_transmission_line_number():
+    """Test `create_transmission_line` with impedance as number."""
+    kl = FrequencyData([1j, 2, 3j], [1, 2, 3])
+    Z = 4+0.1j
+    tmat = TransmissionMatrix.create_transmission_line(kl, Z)
+    assert isinstance(tmat, TransmissionMatrix)
+    _compare_tmat_vs_abcd(
+        tmat,
+        np.cos(kl.freq), 1j*Z*np.sin(kl.freq),
+        1j/Z*np.sin(kl.freq), np.cos(kl.freq))
+
+def test_create_transmission_line_frequency_data():
+    """Test `create_transmission_line` with impedance as FrequencyData."""
+    kl = FrequencyData([1j, 2, 3j], [1, 2, 3])
+    Z = FrequencyData([1+1j, 2+2j, 3+1j], [1, 2, 3])
+    tmat = TransmissionMatrix.create_transmission_line(kl, Z)
+    assert isinstance(tmat, TransmissionMatrix)
+    _compare_tmat_vs_abcd(
+        tmat,
+        np.cos(kl.freq), 1j*Z.freq*np.sin(kl.freq),
+        1j/Z.freq*np.sin(kl.freq), np.cos(kl.freq))
+
+def test_create_transmission_line_broadcasting():
+    """Test `create_transmission_line` broadcasting."""
+    kl = FrequencyData(np.array([[1j, 2, 3j],[4j, 5, 6j]]) , [1, 2, 3])
+    Z = FrequencyData([1+1j, 2+2j, 3+1j], [1, 2, 3])
+    tmat = TransmissionMatrix.create_transmission_line(kl, Z)
+    assert isinstance(tmat, TransmissionMatrix)
+    _compare_tmat_vs_abcd(
+        tmat,
+        np.cos(kl.freq), 1j*Z.freq*np.sin(kl.freq),
+        1j/Z.freq*np.sin(kl.freq), np.cos(kl.freq))
+
+@pytest.mark.parametrize(
+        "kl",
+        [2j, np.array([1, 2, 3]), "term"])
+def test_create_transmission_line_kl_error(kl):
+    """Test `create_transmission_line` error for kl input."""
+    Z = 4+0.1j
+    with pytest.raises(
+        ValueError, match="The input kl"):
+        TransmissionMatrix.create_transmission_line(kl, Z)
+
+@pytest.mark.parametrize("Z", [np.array([1, 2, 3]), "imp"])
+def test_create_transmission_line_characteristic_impedance_errors(Z):
+    """Test `create_transmission_line` errors for impedance parameter."""
+    kl = FrequencyData([1j, 2, 3j], [1, 2, 3])
+    with pytest.raises(
+        ValueError, match="characteristic impedance must be"):
+        TransmissionMatrix.create_transmission_line(kl, Z)
+
+def test_create_transmission_line_frequency_matching():
+    """Test `create_transmission_line` frequency matching."""
+    kl = FrequencyData([1j, 2, 3j], [1, 2, 3])
+    Z = FrequencyData([1+1j, 2+2j, 3+1j], [1, 2, 4])
+    with pytest.raises(
+        ValueError, match="The frequencies do not match"):
+        TransmissionMatrix.create_transmission_line(kl, Z)
 
 def test_tmatrix_slicing(frequencies):
     """Test whether slicing a T-Matrix object return T-Matrix or raises correct