Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Harmonize frequency range parameter #589

Merged
merged 10 commits into from Apr 26, 2024
Merged

Harmonize frequency range parameter #589

Show file tree
Hide file tree
Changes from 5 commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
9fe757f
Harmonize frequency range parameter
hoyer-a Apr 15, 2024
017b374
Fix deprecation warning
hoyer-a Apr 15, 2024
707d397
fix failed tests
hoyer-a Apr 15, 2024
7ad65e7
freq-parameter check in erb_frequencies
hoyer-a Apr 15, 2024
f44689a
Restructure deprecation warning and freq-range parameter check; add f…
hoyer-a Apr 16, 2024
4c2944e
New approach using decorator function
hoyer-a Apr 17, 2024
94c362f
restructure _rename_arg() for more general use
hoyer-a Apr 17, 2024
8600f30
changes according to f-brinkmann
hoyer-a Apr 18, 2024
73d54b2
remove decorator function from utils.py
hoyer-a Apr 18, 2024
7e61868
changes according to f-brinkmann
hoyer-a Apr 19, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
66 changes: 49 additions & 17 deletions pyfar/dsp/dsp.py
View file
Open in desktop
Expand Up @@ -718,9 +718,9 @@ def _time_window_symmetric_interval_four(interval, window):


def regularized_spectrum_inversion(
signal, freq_range,
signal, frequency_range=None,
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

IHm, its a bit annoying that we have to change the default behaviour to make the change. It would be better if we could keep

Suggested change
signal, frequency_range=None,
signal, frequency_range,

without the default None, but that would cause other problems. Any ideas anyone?

Copy link
Contributor Author

@hoyer-a hoyer-a Apr 17, 2024 •
edited

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What about using a decorator to intercept "freq_range" and rename the argument to "frequency_range" and throw a Deprecation warning? If that isn't too error-prone. And we wouldn't really have a "freq_range"-argument anymore.

In that case we wouldn't have to touch the functions at all (apart from changing all freq_range to frequency_range).

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think with fabians and my suggestions, it would already work as expected. No additional funtionallity and the old version would still work.

Copy link
Contributor Author

@hoyer-a hoyer-a Apr 17, 2024 •
edited

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Wouldn't that throw an error if "frequency_range" is a positional argument and only "freq_range" is passed?

e.g.

pf.dsp.regularized_spectrum_inversion(sweep, freq_range=(20, 20e3))

---------------------------------------------------------------------------
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
Cell In[6], [line 1](vscode-notebook-cell:?execution_count=6&line=1)
----> [1](vscode-notebook-cell:?execution_count=6&line=1) pf.dsp.regularized_spectrum_inversion(sweep, freq_range=(20, 20e3))

TypeError: regularized_spectrum_inversion() missing 1 required positional argument: 'frequency_range'

with

def regularized_spectrum_inversion(
        signal, frequency_range,
        regu_outside=1., regu_inside=10**(-200/20), regu_final=None,
        normalized=True, *, freq_range=None):

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, the TypeError is the problem I thought about. Should we maybe discuss this in the weekly. The decorator sounds interesting to me, but I'm not sure if Anne has a better soloution in mind.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

ahh i see, thank you. this might help: https://stackoverflow.com/questions/74564140/best-practice-to-rename-a-method-parameter-in-a-deployed-python-module

Copy link
Contributor Author

@hoyer-a hoyer-a Apr 17, 2024 •
edited

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yeah, that's where i got the idea with the decorator from. Should have mentioned it...
I already have a code snippet based on the stackoverflow question, that i want to try.

Maybe i can show you on friday.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In that case I think you could go straight ahead and implement it :)

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ok, perfect!

regu_outside=1., regu_inside=10**(-200/20), regu_final=None,
normalized=True):
normalized=True, *, freq_range=None):
r"""Invert the spectrum of a signal applying frequency dependent
regularization.

Expand All @@ -747,7 +747,7 @@ def regularized_spectrum_inversion(
----------
signal : Signal
The signals which spectra are to be inverted.
freq_range : tuple, array_like, double
frequency_range : tuple, array_like, double
The upper and lower frequency limits outside of which the
regularization factor is to be applied.
regu_outside : float, optional
Expand All @@ -763,6 +763,11 @@ def regularized_spectrum_inversion(
normalized : bool
Flag to indicate if the normalized spectrum (according to
`signal.fft_norm`) should be inverted. The default is ``True``.
freq_range: tuple, array_like, double
hoyer-a marked this conversation as resolved.
Show resolved Hide resolved
The upper and lower frequency limits outside of which the
regularization factor is to be applied.
``'freq_range'`` parameter will be deprecated in pyfar 0.8.0 in favor
of ``'frequency_range'``.


Returns
Expand All @@ -780,6 +785,17 @@ def regularized_spectrum_inversion(
numerical aspects of linear inversion. Philadelphia: SIAM, 1998.

"""
# Check frequency range parameter
if frequency_range is None and freq_range is None:
raise ValueError("Frequency range must be provided")
hoyer-a marked this conversation as resolved.
Show resolved Hide resolved
if freq_range is not None:
frequency_range = freq_range
# Deprecation warning for freq_range parameter
warnings.warn((
'freq_range parameter will be deprecated in pyfar 0.8.0 in favor'
' of frequency_range'),
PyfarDeprecationWarning)

if not isinstance(signal, pyfar.Signal):
raise ValueError("The input signal needs to be of type pyfar.Signal.")

Expand All @@ -791,9 +807,9 @@ def regularized_spectrum_inversion(
else:
data = signal.freq_raw

freq_range = np.asarray(freq_range)
frequency_range = np.asarray(frequency_range)

if freq_range.size < 2:
if frequency_range.size < 2:
raise ValueError(
"The frequency range needs to specify lower and upper limits.")

Expand All @@ -802,14 +818,15 @@ def regularized_spectrum_inversion(
regu_outside = np.ones(signal.n_bins, dtype=np.double) * regu_outside

idx_xfade_lower = signal.find_nearest_frequency(
[freq_range[0]/np.sqrt(2), freq_range[0]])
[frequency_range[0]/np.sqrt(2), frequency_range[0]])

regu_final = _cross_fade(regu_outside, regu_inside, idx_xfade_lower)

if freq_range[1] < signal.sampling_rate/2:
if frequency_range[1] < signal.sampling_rate/2:
idx_xfade_upper = signal.find_nearest_frequency([
freq_range[1],
np.min([freq_range[1]*np.sqrt(2), signal.sampling_rate/2])])
frequency_range[1],
np.min([frequency_range[1]*np.sqrt(2),
signal.sampling_rate/2])])

regu_final = _cross_fade(regu_final, regu_outside, idx_xfade_upper)

Expand Down Expand Up @@ -1418,8 +1435,8 @@ def find_impulse_response_start(
return np.squeeze(start_sample)


def deconvolve(system_output, system_input, fft_length=None, freq_range=None,
**kwargs):
def deconvolve(system_output, system_input, fft_length=None,
frequency_range=None, *, freq_range=None, **kwargs):
r"""Calculate transfer functions by spectral deconvolution of two signals.

The transfer function :math:`H(\omega)` is calculated by spectral
Expand Down Expand Up @@ -1453,7 +1470,7 @@ def deconvolve(system_output, system_input, fft_length=None, freq_range=None,
The system input signal (e.g., used to perform a measurement).
The system input signal is zero padded, if it is shorter than the
system output signal.
freq_range : tuple, array_like, double
frequency_range : tuple, array_like, double
The upper and lower frequency limits outside of which the
regularization factor is to be applied. The default ``None``
bypasses the regularization, which might cause numerical
Expand All @@ -1467,6 +1484,14 @@ def deconvolve(system_output, system_input, fft_length=None, freq_range=None,
kwargs : key value arguments
Key value arguments to control the inversion of :math:`H(\omega)` are
passed to to :py:func:`~pyfar.dsp.regularized_spectrum_inversion`.
freq_range: tuple, array_like, double
The upper and lower frequency limits outside of which the
regularization factor is to be applied. The default ``None``
bypasses the regularization, which might cause numerical
instabilities in case of band-limited `system_input`. Also see
:py:func:`~pyfar.dsp.regularized_spectrum_inversion`.
``'freq_range'`` parameter will be deprecated in pyfar 0.8.0 in favor
of ``'frequency_range'``.


Returns
Expand All @@ -1482,20 +1507,27 @@ def deconvolve(system_output, system_input, fft_length=None, freq_range=None,
sweeps. Directors cut." J. Audio Eng. Soc. 49(6):443-471,
(2001, June).
"""

# Check if system_output and system_input are both type Signal
if not isinstance(system_output, pyfar.Signal):
raise TypeError('system_output has to be of type pyfar.Signal')
if not isinstance(system_input, pyfar.Signal):
raise TypeError('system_input has to be of type pyfar.Signal')

# Check frequency range parameter
if frequency_range is None and freq_range is None:
frequency_range = (0, system_input.sampling_rate/2)
hoyer-a marked this conversation as resolved.
Show resolved Hide resolved
if freq_range is not None:
frequency_range = freq_range
# Deprecation warning for freq_range parameter
warnings.warn((
'freq_range parameter will be deprecated in pyfar 0.8.0 in favor'
' of frequency_range'),
PyfarDeprecationWarning)

# Check if both signals have the same sampling rate
if not system_output.sampling_rate == system_input.sampling_rate:
raise ValueError("The two signals have different sampling rates!")

if freq_range is None:
freq_range = (0, system_input.sampling_rate/2)

hoyer-a marked this conversation as resolved.
Show resolved Hide resolved
# Set fft_length to the max n_samples of both signals,
# if it is not explicitly set to a value
if fft_length is None:
Expand All @@ -1518,7 +1550,7 @@ def deconvolve(system_output, system_input, fft_length=None, freq_range=None,
# multiply system_output signal with regularized inversed system_input
# signal to get the system response
inverse_input = regularized_spectrum_inversion(
system_input, freq_range, **kwargs)
system_input, frequency_range, **kwargs)
system_response = system_output * inverse_input

# Check if the signals have any comments,
Expand Down
31 changes: 24 additions & 7 deletions pyfar/dsp/filter/fractional_octaves.py
View file
Open in desktop
Expand Up @@ -2,6 +2,7 @@
import numpy as np
import scipy.signal as spsignal
import pyfar as pf
from pyfar.classes.warnings import PyfarDeprecationWarning


def fractional_octave_frequencies(
Expand Down Expand Up @@ -157,8 +158,9 @@ def fractional_octave_bands(
signal,
num_fractions,
sampling_rate=None,
freq_range=(20.0, 20e3),
order=14):
frequency_range=(20.0, 20e3),
order=14,
*, freq_range=None):
"""Create and/or apply an energy preserving fractional octave filter bank.

The filters are designed using second order sections of Butterworth
Expand Down Expand Up @@ -189,11 +191,17 @@ def fractional_octave_bands(
sampling_rate : None, int
The sampling rate in Hz. Only required if signal is ``None``. The
default is ``None``.
freq_range : array, tuple, optional
frequency_range : array, tuple, optional
The lower and upper frequency limits. The default is
``frequency_range=(20, 20e3)``.
order : int, optional
Order of the Butterworth filter. The default is ``14``.
freq_range: array, tuple, optional
The lower and upper frequency limits. The default is
``frequency_range=(20, 20e3)``.
``'freq_range'`` parameter will be deprecated in pyfar 0.8.0 in favor
of ``'frequency_range'``.


Returns
-------
Expand All @@ -215,7 +223,7 @@ def fractional_octave_bands(
>>> # generate the data
>>> x = pf.signals.impulse(2**17)
>>> y = pf.dsp.filter.fractional_octave_bands(
... x, 1, freq_range=(20, 8e3))
... x, 1, frequency_range=(20, 8e3))
>>> # frequency domain plot
>>> y_sum = pf.FrequencyData(
... np.sum(np.abs(y.freq)**2, 0), y.frequencies)
Expand All @@ -225,6 +233,15 @@ def fractional_octave_bands(
... "Filter bands and the sum of their squared magnitudes")

"""
# Check frequency range parameter
if freq_range is not None:
frequency_range = freq_range
# Deprecation warning for freq_range parameter
warnings.warn((
'freq_range parameter will be deprecated in pyfar 0.8.0 in favor'
' of frequency_range'),
PyfarDeprecationWarning)

# check input
if (signal is None and sampling_rate is None) \
or (signal is not None and sampling_rate is not None):
Expand All @@ -234,7 +251,7 @@ def fractional_octave_bands(

sos = _coefficients_fractional_octave_bands(
sampling_rate=fs, num_fractions=num_fractions,
freq_range=freq_range, order=order)
frequency_range=frequency_range, order=order)

filt = pf.FilterSOS(sos, fs)
filt.comment = (
Expand All @@ -253,7 +270,7 @@ def fractional_octave_bands(

def _coefficients_fractional_octave_bands(
sampling_rate, num_fractions,
freq_range=(20.0, 20e3), order=14):
frequency_range=(20.0, 20e3), order=14):
"""Calculate the second order section filter coefficients of a fractional
octave band filter bank.

Expand Down Expand Up @@ -283,7 +300,7 @@ def _coefficients_fractional_octave_bands(
"""

f_crit = fractional_octave_frequencies(
num_fractions, freq_range, return_cutoff=True)[2]
num_fractions, frequency_range, return_cutoff=True)[2]

freqs_upper = f_crit[1]
freqs_lower = f_crit[0]
Expand Down