Find the equivalent of this python remez specs in C++ remez or Matlab firpm

I want to start with that I am quite new in this community, so if this question does not belong here, please point me to the right place. Thank you!

Background

I am learning and writing a half-band filter in C++ using the famous Remez exchange algorithm. And I am implementing the filter following A “TRICK” for the Design of FIR Half-Band Filters. I am also using this blog as a reference. This is the source remez.cpp.

In the blog, section Half-Band Filter Design with NumPy, filter coefficients are generated using the python signal package from scipy module.

g = signal.remez(N//2+1, [0., 2*Fpb/Fs, .5, .5], [1, 0], [1, 1]) # The Trick


with numbers plugged in:

g = signal.remez(6, [0., 0.4, .5, .5], [1, 0], [1, 1]) # The Trick


And it outputs

g = [ 0.10753661 -0.18306401  0.62633938  0.62633938 -0.18306401  0.10753661]


Problem

Calling remez() in remez.cpp using the same specs, however, does not return the same coefficients.

/********************
* remez
*=======
* Calculates the optimal (in the Chebyshev/minimax sense)
* FIR filter impulse response given a set of band edges,
* the desired response on those bands, and the weight given to
* the error in those bands.
*
* INPUT:
* ------
* int     numtaps     - Number of filter coefficients
* int     numband     - Number of bands in filter specification
* double  bands[]     - User-specified band edges [2 * numband]
* double  des[]       - User-specified band responses [numband]
* double  weight[]    - User-specified error weights [numband]
* int     type        - Type of filter
*
* OUTPUT:
* -------
* double h[]      - Impulse response of final filter [numtaps]
* returns         - true on success, false on failure to converge
********************/

int numOrder   = 5;
std::vector<double> h(numOrder + 1);
std::vector<double> bandsEdges  = {0, 0.4, 0.5, 0.5};
std::vector<double> desiredAmps = {1, 0};
std::vector<double> weights     = {1, 1};
remez(h.data(),6,2,bandsEdges.data(),desiredAmps.data(),weights.data(),1,16)

// Pseudo code
// h = [0.0154523, 0.115447, 0.343533, 0.343533, 0.115447, 0.0154523];


I've also tried in Matlab using the following command, and I also couldn't get the same coefficients as Python ones.

h = firpm(5,[0 0.8 1 1],[1,1,0,0],[1,1]);
h = [0.089673, -0.176451, 0.624020, 0.624020, -0.176451, 0.089673];


Question

I am really baffled by the three completely different sets of results. Am I using the wrong specs in either Matlab or C++ remez?

After some trial and error, I figured out the correct specs for C++ remez and Matlab.

C++

Treat it as a one band response. Thus numBand should be set to 1.

int numOrder   = 5;
std::vector<double> h(numOrder + 1);
std::vector<double> bandsEdges  = {0, 0.4, 0.5, 0.5};
std::vector<double> desiredAmps = {1, 1, 0, 0};
std::vector<double> weights     = {1, 1};
remez(h.data(),6,1,bandsEdges.data(),desiredAmps.data(),weights.data(),1,16);


Matlab

Similarly, since we only care about the passband response, we can treat the rest of bandwidth as don't cares.

h = firpm(5,[0 0.8],[1 1]);

• Yes, that looks correct!
– Peter K.
Oct 21 '21 at 19:37

To start with: the Python and Matlab versions are not that different.

The plot below shows all three of them. The red and green ones are the Python and Matlab versions, respectively. These overlap each other, with just a small error. The black one is the C++ result. It is significantly different from the others.

The thing to know about the Remez algorithm is that it is iterative and needs to converge. The difference between the three implementations is probably the choice of how it's determined if the result has converged.

It looks like what Yihan says is correct.

If I use

h_matlab_2 = firpm(5,[0 0.8],[1 1]);


instead, then the Python and Matlab versions are much closer.

• Thanks for the info and for trying it out. I think the way each implementation interprets arguments is different from each others. Also, it seems that the C++ implementation has a mistake in the description of des argument. It should be the same length as bands, just like how firpm in Matlab is. Fixing that gives me the correst result. Oct 21 '21 at 20:50
• @YihanHu Good to hear! You figured it out yourself. :-) Oh well, it was a fun thing to try.
– Peter K.
Oct 21 '21 at 21:46