3
$\begingroup$

I am attempting to implement a filter described in the paper Functional Count-Comparison Model for Binaural Decoding by Pulkki et al 2009:

enter image description here

Here is Figure 9:

enter image description here

However, I am quite certain that the equation for n is incorrect since the cosine is squared and there is a fourth root of n/fs. This equation does not produce a single period for every omega value (in my case between 100 Hz and 12 kHz).

I am also unsure why fs is in the equation for f at all, since it should be handled by the length of n. Am I correct on either or both of these points?

Another issue I am having is that the BF of the highest frequency band is 12 kHz, while fs = 20 kHz, which I suspect will cause some aliasing. Should I remove all frequency bands above 10 kHz?

The final issue I have is that given the equation, there does not seem to be a high enough fs to represent the filter for the highest frequencies as the length of n quickly goes to 0 well before 10 kHz.

Here is my attempt at implementing it in python:

import numpy as np
import matplotlib.pyplot as plt


# Define the function
def raised_cosine_filter(omega, fs=20000):
    # Calculate the period based on omega
    period = (2 * np.pi) / omega

    # Calculate the number of samples
    num_samples = int(np.ceil(period * fs))

    # Generate x values for one period
    n = np.linspace(0, period, num_samples)

    return ((np.cos(omega * (n / fs) ** 0.25 - np.pi) + 1) ** 2)


# Example usage:
f = 100  # Frequency in Hz
omega = 2 * np.pi * f  # Angular frequency in rad/s
y = raised_cosine_filter(omega)

# Plot the raised cosine pulse
plt.plot(y)
plt.xlabel('Sample number')
plt.ylabel('y')
plt.title(f'One Period of the Raised Cosine Filter (f = {f} Hz)')
plt.grid()
plt.show()

And here is the resulting plot, which looks incorrect as you can see:

enter image description here

Any help would be appreciated in fixing either my implementation or clarifying whether the paper has a mistake. Thank you.

$\endgroup$
7
  • $\begingroup$ What confuses me, is that the x-Axis are different in the both figures. It is not really clear to me what exactly is meant by amplitude in fig.9 $\endgroup$ Commented Jun 26, 2023 at 13:15
  • $\begingroup$ Here is a link to the original paper. @RodrigodeAzevedo. $\endgroup$
    – Bryn
    Commented Jun 26, 2023 at 15:33
  • $\begingroup$ @Irreducible. I was also confused by the axes of figure 9 from the paper. They plotted the phase response, but I'm not sure why. Based on the equations, the function should look similar with samples on the x-axis as well. $\endgroup$
    – Bryn
    Commented Jun 26, 2023 at 15:36
  • $\begingroup$ Please provide a full reference. In the question itself, not in the comment section $\endgroup$ Commented Jun 26, 2023 at 15:43
  • $\begingroup$ Maybe not related but there's an erratum published. $\endgroup$
    – Juha P
    Commented Jun 26, 2023 at 20:59

2 Answers 2

1
$\begingroup$

According to the erratum

Equation 5 is an error. Here is the actual equation:

enter image description here

And the code after the correction:

import numpy as np
import matplotlib.pyplot as plt


# Define the function
def raised_cosine_filter(fbf, fs=20000):

    # Calculate the number of samples
    num_samples = int(np.ceil(fs / fbf))

    # Generate x values for one period
    n = np.arange(num_samples)

    return 0.25 * (np.cos(2 * np.pi * (n * fbf / fs) ** 0.25 - np.pi) + 1) ** 2


# Example usage:
fbf = 750  # Frequency in Hz
y = raised_cosine_filter(fbf)

# Plot the raised cosine pulse
plt.plot(y)
plt.xlabel('Sample number')
plt.ylabel('y')
plt.title(f'One Period of the Raised Cosine Filter (f = {fbf} Hz)')
plt.grid()
plt.show()

Here is what the plot looks like:

enter image description here

I'm still concerned about the sampling rate, but not sure how to address that without making the filter inconsistent with the data I'm using.

$\endgroup$
0
$\begingroup$

Not a definite answer, but a tentative hint that's too long to put in a comment


I could be completely off mark but I believe there's either a typo in the paper or no typo but rather something off with your implementation. Regardless, if you take $f_s$ out of the equation for $f$ like you hinted at, and use *normalized* angular frequency, i.e. $\omega = 2\pi f/f_s$, you get something that looks very close (up to normalization) to what you're looking for.

I used $750 \,\texttt{Hz}$ to match what the authors do in the paper linked:

enter image description here

import numpy as np
import matplotlib.pyplot as plt


# Define the function
def raised_cosine_filter(omega, fs):
    # Calculate the period based on omega
    period = (2 * np.pi * fs) / omega

    n = np.linspace(0, period, 360)

    return (np.cos(omega * n ** 0.25 - np.pi) + 1) ** 2


# Example usage:
f = 750  # Frequency in Hz
fs = 20000 # sampling frequency
omega = 2 * np.pi * f / fs  # normalized angular frequency 
y = raised_cosine_filter(omega, fs)

# Plot the raised cosine pulse
plt.plot(y)
plt.xlabel('Sample number')
plt.ylabel('y')
plt.title(f'One Period of the Raised Cosine Filter (f = {f} Hz)')
plt.grid()
plt.show()
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.