.. only:: html
.. note::
:class: sphx-glr-download-link-note
Click :ref:`here ` to download the full example code
.. rst-class:: sphx-glr-example-title
.. _sphx_glr_auto_examples_plot_higherOrderModes.py:
The empty circular duct with flow and higher-order modes
========================================================
In this example we extract higher-order modes from measurement data in a circular duct with moderate mean flow.
The data is part of `this study `_\,
which is referred here for further details.
.. image:: ../../image/wave2.JPG
:width: 800
1. Initialization
-----------------
First, we import the packages needed for this example.
.. code-block:: default
import matplotlib.pyplot as plt
import numpy
import acdecom
The duct has a circular cross section and is filled with air for which we use the standard properties.
The duct radius is 0.075 m. The distance from the microphones to the reference cross section is 0.55 m.
The highest frequency of interest is 3200 Hz.
.. code-block:: default
section = "circular"
radius = 0.075 # m
distance = 0.55 # [m]
f_max = 3200 # Hz
We create objects for the upstream and the downstream section of the duct.
.. code-block:: default
td_us = acdecom.WaveGuide(dimensions=(radius,), cross_section=section, f_max=f_max, damping="dokumaci",
distance=distance, flip_flow=True)
td_ds = acdecom.WaveGuide(dimensions=(radius,), cross_section=section, f_max=f_max, damping="dokumaci",
distance=distance)
.. note::
The standard flow direction is in :math:`P_+` direction. On the upstream side, the Mach-number therefore must be
either set negative or the argument *flip_flow* must be set to *True*.
2. Sensor Positions
-------------------
The microphone coordinates are saved in the
`emptyUS.mic `_ and
`emptyDS.mic `_ file.
.. code-block:: default
td_us.read_microphonefile("data/emptyUS.mic", cylindrical_coordinates=True)
td_ds.read_microphonefile("data/emptyDS.mic", cylindrical_coordinates=True)
.. note::
In this case, the microphone coordinates are defined in a cylindrical coordinate system with the circumferential
position in deg. Therefore, we set the argument *cylindrical_coordinates* to *True*. This will transform the
circumferential position from deg. to radians.
3. Decomposition
----------------
The measurement must be pre-processed in a format that is understood by the *WaveGuide* object. Generally,
this must be a numpy.ndarray, wherein the columns contain the measurement data, such as the measured frequency, the
pressure values for that frequency, the bulk Mach-number, and the temperature. The rows can be different frequencies
or different sound excitations (cases). In this example, the
measurement was post-processed into the
`higherOrderModes.txt `_ file and
can be loaded with the `numpy.loadtxt `_
function.
.. note::
The pressure used for the decomposition must be pre-processed, for example to account for microphone.
.. code-block:: default
pressure = numpy.loadtxt("data/higherOrderModes.txt",dtype=complex, delimiter=",", skiprows=1)
We examine the file header to understand how the data is stored in our input file.
.. code-block:: default
with open("data/higherOrderModes.txt") as pressurefile:
print(pressurefile.readline().split(","))
.. rst-class:: sphx-glr-script-out
Out:
.. code-block:: none
['Mach-Number', 'temperature', 'f', 'Mic1', 'Mic2', 'Mic3', 'Mic4', 'Mic5', 'Mic6', 'Mic7', 'Mic8', 'Mic9', 'Mic10', 'Mic11', 'Mic12', 'Mic13', 'Mic14', 'Mic15', 'Mic16', 'Mic17', 'Mic18', 'Mic19', 'Mic20', 'Mic21', 'Mic22', 'Mic23', 'Mic24', 'case\n']
Mach-number, temperature, and frequency are stored in columns 0, 1, and 2. The upstream microphones 1-12 are in
columns 3 - 14, the downstream microphones 13-24 are in columns 15-26, and the case number is in the last column.
.. code-block:: default
Mach_number = 0
temperature = 1
f = 2
mics_us = range(3, 15)
Mics_ds = range(15, 27)
case = -1
Now, we can decompose the sound-fields into the propagating modes. We decompose the sound-fields on the upstream
and downstream side of the duct, using the two *WaveGuide* objects defined earlier.
.. code-block:: default
decomp_us, headers_us = td_us.decompose(pressure, f, mics_us, temperature_col=temperature, case_col=case,
Mach_col=Mach_number)
decomp_ds, headers_ds = td_ds.decompose(pressure, f, Mics_ds, temperature_col=temperature, case_col=case,
Mach_col=Mach_number)
.. rst-class:: sphx-glr-script-out
Out:
.. code-block:: none
/home/docs/checkouts/readthedocs.org/user_builds/acdecom/checkouts/latest/src/acdecom.py:1322: UserWarning: The Modal analysis is ill-conditioned for some of the frequencies.
warnings.warn("The Modal analysis is ill-conditioned for some of the frequencies.")
.. note ::
The decomposition may show warnings for ill-conditioned modal matrices. This typically happens for frequencies close
to the cut-on of a mode. However, it can also indicate that the microphone array is insufficient to separate the
modes. The condition number of the wave decomposition is stored in the data returned by *decompose* and
should be checked in case a warning is triggered.
4. Further Post-processing
--------------------------
We can print the *headers_ds* to see the names of the columns of the arrays that store the decomposed sound fields.
.. code-block:: default
print(headers_ds)
.. rst-class:: sphx-glr-script-out
Out:
.. code-block:: none
['(0,0) plus Direction', '(-1,0) plus Direction', '(1,0) plus Direction', '(-2,0) plus Direction', '(2,0) plus Direction', '(0,1) plus Direction', '(0,0) minus Direction', '(-1,0) minus Direction', '(1,0) minus Direction', '(-2,0) minus Direction', '(2,0) minus Direction', '(0,1) minus Direction', 'f', 'Mach_number', 'temperature', 'Ps', 'condition number', 'case']
We use that information to extract the modal data for the different frequencies and cases.
.. code-block:: default
plusmodes = [0,1,2,3,4,5]
minusmodes = [6,7,8,9,10,11]
Furthermore, we can get the unique decomposed frequency points.
.. code-block:: default
frequs = numpy.abs(numpy.unique(decomp_us[:,headers_us.index("f")]))
nof = frequs.shape[0]
For each of the frequencies we can compute the scattering matrix by solving a linear system of equations
:math:`S = p_+ p_-^{-1}`\, where :math:`S` is the scattering matrix and :math:`p_{\pm}` are matrices containing the
acoustic modes palded in rows and the different test cases placed in columns.
.. note::
Details for the computation of the Scattering Matrix and the procedure to measure the different test-cases can be
found in `this study `_\.
.. code-block:: default
S = numpy.zeros((12,12,nof), dtype=complex)
for fIndx, f in enumerate(frequs):
frequ_rows = numpy.where(decomp_us[:,headers_us.index("f")] == f)
ppm_us = decomp_us[frequ_rows]
ppm_ds = decomp_ds[frequ_rows]
pp = numpy.concatenate((ppm_us[:,plusmodes].T, ppm_ds[:,plusmodes].T))
pm = numpy.concatenate((ppm_us[:,minusmodes].T, ppm_ds[:,minusmodes].T))
S[:,:,fIndx] = numpy.dot(pp,numpy.linalg.pinv(pm))
5. Plot
-------
Finally, we can plot the transmission and reflection coefficients of the 6 propagating modes.
.. code-block:: default
mode_names = td_us.mode_vector
fig,axs=plt.subplots(6,1,figsize=(10,10))
axs[0].set_title("Empty Circular Duct with Higher Order Modes")
for mode in range (6):
axs[mode].plot(frequs, numpy.abs(S[mode,mode+6,:]),
color="#67A3C1", label = str(mode_names[mode]) + "-Mode Transmission")
axs[mode].plot(frequs, numpy.abs(S[mode,mode,:]), ls="--",
color="#D38D7B", label = str(mode_names[mode]) + "-Mode Reflection")
axs[mode].set_xlim([0, 3200])
axs[mode].set_ylim([0, 1.1])
axs[mode].set_xticks([])
axs[mode].legend()
axs[2].set_ylabel("Scattering Magnitude")
axs[5].set_xticks(range(0,3200,250))
plt.xlabel("Frequency [Hz]")
plt.show()
.. image:: /auto_examples/images/sphx_glr_plot_higherOrderModes_001.png
:alt: Empty Circular Duct with Higher Order Modes
:class: sphx-glr-single-img
.. rst-class:: sphx-glr-timing
**Total running time of the script:** ( 1 minutes 3.698 seconds)
.. _sphx_glr_download_auto_examples_plot_higherOrderModes.py:
.. only :: html
.. container:: sphx-glr-footer
:class: sphx-glr-footer-example
.. container:: sphx-glr-download sphx-glr-download-python
:download:`Download Python source code: plot_higherOrderModes.py `
.. container:: sphx-glr-download sphx-glr-download-jupyter
:download:`Download Jupyter notebook: plot_higherOrderModes.ipynb `
.. only:: html
.. rst-class:: sphx-glr-signature
`Gallery generated by Sphinx-Gallery `_