I am trying to find a Python replacement for the Matlab function cfirpm (complex FIR Parks-McClellan). I use cfirpm in the design of a complex FIR magnitude and group delay compensation filter for a crystal IF filter. Since the frequency response of the crystal IF filter is not symmetric, the compensating filter isn't either. Therefore, I cannot simply use a lowpass filter design and multiply it by e^(jwft). I have searched for hints on how to do this on this site and many others. I have looked at 23147/design-filter-with-arbitrary-magnitude-and-phase-response-in-matlab where the answer was to find the inverse FFT of the two-sided complex frequency response, etc. but what I have from a frequency sweep using an RF signal generator and processing the I and Q samples is a vector representing the desired magnitude at evenly spaced frequencies and I'm not sure how to represent that as a two-sided complex frequency response. Any ideas on how this type of filter can be designed in Python would be greatly appreciated.
Here is an example I already had in Python demonstrating the Wiener-Hopf equations detailed here Compensating Loudspeaker frequency response in an audio signal to solve for the channel equalizer. This application by the OP is ideal for this given the static channel condition as long as deep spectral nulls don't exist, over a iterative LMS or RLS equalizer solution such as those detailed here (which are great for dynamic channels that don't exhibit deep frequency nulls in band): Recursive Least Square Adaptive Linear Equalizer For channels with deep spectral nulls in band (frequency selective fading), a decision-feedback equalizer would be the better approach due to degradations from noise enhancement with the linear equalizer approaches.
Eye diagram of undistorted waveform with no channel distortion, as it would be received after the 2nd RRC filter in the receiver:
Tx Spectrum of undistorted waveform prior to channel distortion:
Example complex channel distortion that we will equalize:
cir = [.019*np.exp(1j*.2), .021*np.exp(-1j*.24), .08*np.exp(1j*2.4), .032*np.exp(-1j*.14), -.032*np.exp(-1j*.14)]
This is the same eye diagram and Rx spectrum above after passing through our channel (yuck, what a mess! How do we possibly demodulate this????). Notice the assymetrical passband, a tell-tale sign that the compensating filter MUST be complex.
Compute Equalizer Solution Using Tx and Rx Waveforms
(Note this is Python 3.8, but nothing jumps out at me that wouldn't be compatible with earlier versions)
import scipy.linalg as linalg import numpy as np def convmtx(h,n): # creates the convolution (Toeplitz) matrix, which transforms convolution to matrix multiplication # h is input array, # n is length of array to convolve return linalg.toeplitz(np.hstack([h, np.zeros(n-1)]), np.hstack([h, np.zeros(n-1)])) # Compute equalizer given known tx and rx waveforms omit = 100 #initial samples to exclude shift = 0 # number of samples to shift dominant equalizer tap to the right ntaps = 20 depth = ntaps * 30 delay = ntaps//2 A = convmtx(rx[omit+shift: omit+shift+depth], ntaps) R = np.dot(np.conj(A).T, A) X = np.concatenate([np.zeros(delay), tx[omit: omit+depth], np.zeros(np.int(np.ceil(ntaps/2)-1))]) ro = np.dot(np.conj(A).T, X) equal_coeff = np.dot(linalg.inv(R),ro)
Below is the resulting frequency response for the equalizer and the equalized waveform. The details for the parameters used such as equalizer number of taps and depth are further detailed in the earlier links given. Note how the equalizer can only determine the compensation where signal energy exists (otherwise it has no way to measure the channel), therefore it is important to use a test (sounding) waveform that fully covers the bandwidth of the channel of interest. Pseudo-random sources or any white noise source is ideal given the broad spectral coverage.