# How to create a filter that is a single period of a time-warped raised cosine?

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

Here is Figure 9:

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:

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

• 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 Commented Jun 26, 2023 at 13:15
• Here is a link to the original paper. @RodrigodeAzevedo.
– Bryn
Commented Jun 26, 2023 at 15:33
• @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.
– Bryn
Commented Jun 26, 2023 at 15:36
• Please provide a full reference. In the question itself, not in the comment section Commented Jun 26, 2023 at 15:43
• Maybe not related but there's an erratum published. Commented Jun 26, 2023 at 20:59

According to the erratum

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

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:

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.

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:

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()