Forums

Remez Exchange Algorithm

Started by Unknown June 22, 2006
Hi,

I want to optimize a 32 band EQ plugin. Currently it uses inverse FFT 
for the creation of the FIR filter kernel. But this requires rather long 
FIR kernels to process low frequencies properly. While computing power 
is nowadays not an issue the latency is.

I heard about the "Remez Exchange Algorithmn". Unfortunatly I did not 
find a decription how to implement this in the web. Only the keyword is 
everywhere.

Does anybody know where to get information about this algorithmn?


Marcel
Marcel M=FCller wrote:
> Hi, > > I want to optimize a 32 band EQ plugin. Currently it uses inverse FFT > for the creation of the FIR filter kernel. But this requires rather long > FIR kernels to process low frequencies properly. While computing power > is nowadays not an issue the latency is. > > I heard about the "Remez Exchange Algorithmn". Unfortunatly I did not > find a decription how to implement this in the web. Only the keyword is > everywhere. > > Does anybody know where to get information about this algorithmn?
The Remez algorithm is a somewhat "hazy" tool in DSP, hovering out there somewhere in the mist, halfway between sheer wizardry and black craft. It works, and is a classical tool for what it does, but comprehensable, detailed information is hard to find. Comprehensable information is easy to find, it's available in almost every text on DSP. But the texts I know of are not detailed enough to be used as recipes for implementation. On the other hand, the descriptions that are detailed enough to be used as program recipes are all but incomprehensable. Having said that, I am not sure the Remez is the correct solution for your problem. The Remez is a tool to *design* a filter, i.e. find the coefficients that ultimately constitute the filter.
>From your question, it seems to me as if you need an efficient way
of *implementing* a long filter, get the run-time down to a minimum. These are two different questions that will have different answers. Could you elaborate on exactly what you try to do, and what the problem is? Rune
Rune Allnor schrieb:
> The Remez algorithm is a somewhat "hazy" tool in DSP, hovering > out there somewhere in the mist, halfway between sheer wizardry > and black craft. It works, and is a classical tool for what it does, > but > comprehensable, detailed information is hard to find. Comprehensable > information is easy to find, it's available in almost every text on > DSP. > But the texts I know of are not detailed enough to be used as recipes > for implementation. > > On the other hand, the descriptions that are detailed enough to be > used as program recipes are all but incomprehensable. > > Having said that, I am not sure the Remez is the correct solution > for your problem. The Remez is a tool to *design* a filter, i.e. > find the coefficients that ultimately constitute the filter. > >>From your question, it seems to me as if you need an efficient way > of *implementing* a long filter, get the run-time down to a minimum.
The filter implementation is already done and the run-time is acceptable too (FFT convolution). But the FIR order to achieve an accuracy in the order of +-1dB is about 32768. This introduces a latency in the order of 1s. That is the reason why I look for a shorter kernel. So the Idea is to design a shorter FIR kernel that matches a given function with some tolerance. In case of an EQ this is only a Frequency response. The filter is linear phase. The only idea I have so far is, that knowing the frequency response of the window funtion (currently a slightly modified Hamming) it should be possible to tune the coefficients before the IFFT to compensate for that. But this is a deconvolution task and I have no idea of doing this efficiently. Marcel
Marcel M=FCller wrote:
> Rune Allnor schrieb:
.=2E.
> > From your question, it seems to me as if you need an efficient way > > of *implementing* a long filter, get the run-time down to a minimum. > > The filter implementation is already done and the run-time is acceptable > too (FFT convolution). But the FIR order to achieve an accuracy in the > order of +-1dB is about 32768. This introduces a latency in the order of > 1s. That is the reason why I look for a shorter kernel. > > So the Idea is to design a shorter FIR kernel that matches a given > function with some tolerance. In case of an EQ this is only a Frequency > response. The filter is linear phase.
You could use a minimum-phase kernel and an inhomogeneous frame division of the impulse response to simulatenously reduce the filter and processing latency. Check out http://www.music.miami.edu/programs/mue/Research/jvandekieft/jvchapter2.htm on how to do this. However, for your application (EQ), I would suggest a series of IIR biquads instead of a FIR filter. Even at 32 bands, this is more efficient and far simpler to implement than (low-latency) frequency domain filtering. Another upside is that the order of the filter does not depend on the range of EQ frequencies. Regards, Andor
Marcel M�ller wrote:
> Hi, > > I want to optimize a 32 band EQ plugin. Currently it uses inverse FFT > for the creation of the FIR filter kernel. But this requires rather long > FIR kernels to process low frequencies properly. While computing power > is nowadays not an issue the latency is. > > I heard about the "Remez Exchange Algorithmn". Unfortunatly I did not > find a decription how to implement this in the web. Only the keyword is > everywhere. > > Does anybody know where to get information about this algorithmn?
The Remez exchange algorithm is a procedure used to approximate arbitrary rational functions. It is used in DSP to primarily as the core of a program that optimizes the coefficients of FIR filters. Its use for that purpose was introduced by Parks and McClellan. I doubt that the algorithm itself (as opposed to programs that use it) would be of much use to you, but if you want to know more than you're likely to learn via Google, it is described briefly in "An Introduction to the Approximation of Functions? by Theodore Rivlin, available in inexpensive reprint from Dover Publications (ISBN 0-486-64069-8). The original paper is Remez, E. Ya. "General Computational Methods of Chebychev Approximation. The Problems with Linear Real Parameters", AEC-tr-4491 (US Atomic Energy Commission, 1962) Translated from Russian in two volumes. You don't want to read it! "Ralston, A. "Rational Chebychev Approximation by Remes' [sic] Algorithm" Numer. Math. 7 (1965) pp. 322-330 is a somewhat more accessible paper. Jerry -- Engineering is the art of making what you want from things you can get. �����������������������������������������������������������������������
To rephrase what you are doing:

You are designing a system where the user can define a filter by moving
32 knobs, with each knob representing a frequency area.
Now, you want to take the information from these knobs, and compute the
coefficients for an FIR filter.  Then you want to use these
coefficients to filter the actual incoming data.

You would like to try the Remez algorithm instead of doing an inverse
FFT to get the coefficients, hoping to reduce the number of
coefficients.

I don't think you will get the efficiency you are looking for with the
Remez algorithm.  I would build several prototype models, and use
Matlab to compute the coefficients with different algorithms.
While I haven't tried it, I think the Remez algorithm will reduce the
number of coefficients you need by 10% or so, while you would probably
like to reduce the number of coefficients by 90%

"Marcel M&#2013266172;ller" <news.5.maazl@spamgourmet.com> wrote in message 
news:449a595e$0$11079$9b4e6d93@newsread4.arcor-online.net...
> Hi, > > I want to optimize a 32 band EQ plugin. Currently it uses inverse FFT for > the creation of the FIR filter kernel. But this requires rather long FIR > kernels to process low frequencies properly. While computing power is > nowadays not an issue the latency is. > > I heard about the "Remez Exchange Algorithmn". Unfortunatly I did not find > a decription how to implement this in the web. Only the keyword is > everywhere. > > Does anybody know where to get information about this algorithmn?
Jerry has given you good text references. Here is a quick description of the algorithm which can be applied to FIR filters among other functional approximation applications: [This algorithm is particularly suited to minimax error criterion] 1) Define a desired function D(w) (as a "perfect" filter would be a desired function) 2) Choose a set of basis functions that will be used to add up, suitably weighted by individual coefficients, to an approximation of the desired function. Harmonically related cosines are the typical basis functions for symmetrical/linear phase FIR filters. Call any sum the approximant A(w) 3) Define the error E(w) as D(w) - A(w). 4) Define an error value of +/-|e| 4) Select an arbitrary set of N+1 points that spans the frequency range of interest. 5) Solve a linear system of N+1 equations that equalizes the error |e| with alternating sign of the error (as a function of frequency): E(k) = D(k) - A(k) = e*1^-k Now the error is guaranteed to have alternating signs. 5) Find a new set of N+1 points that have the largest *peaks* of E(w) AND have alternating signs. 6) Solve a new linear system of N+1 equations as above. Continue 5 and 6 until the peaks and the equalized error values are equal. The equalized error value will increase at each step and the absolute maximum peak will reduce at each step until they are equal. That said, I doubt it will yield a shorter filter. Fred
"Fred Marshall" <fmarshallx@remove_the_x.acm.org> wrote in message 
news:nMadnfLZSMz-4gbZnZ2dnUVZ_tSdnZ2d@centurytel.net...
> > "Marcel M&#2013266172;ller" <news.5.maazl@spamgourmet.com> wrote in message > news:449a595e$0$11079$9b4e6d93@newsread4.arcor-online.net... >> Hi, >> >> I want to optimize a 32 band EQ plugin. Currently it uses inverse FFT for >> the creation of the FIR filter kernel. But this requires rather long FIR >> kernels to process low frequencies properly. While computing power is >> nowadays not an issue the latency is. >> >> I heard about the "Remez Exchange Algorithmn". Unfortunatly I did not >> find a decription how to implement this in the web. Only the keyword is >> everywhere. >> >> Does anybody know where to get information about this algorithmn? > > Jerry has given you good text references. > > Here is a quick description of the algorithm which can be applied to FIR > filters among other functional approximation applications: > > [This algorithm is particularly suited to minimax error criterion] > > > 1) Define a desired function D(w) (as a "perfect" filter would be a > desired function) > > 2) Choose a set of basis functions that will be used to add up, suitably > weighted by individual coefficients, to an approximation of the desired > function. Harmonically related cosines are the typical basis functions > for symmetrical/linear phase FIR filters. Call any sum the approximant > A(w) > > 3) Define the error E(w) as D(w) - A(w). > > 4) Define an error value of +/-|e| > > 4) Select an arbitrary set of N+1 points that spans the frequency range of > interest. > > 5) Solve a linear system of N+1 equations that equalizes the error |e| > with alternating sign of the error (as a function of frequency): > > E(k) = D(k) - A(k) = e*(-1)^-k > > Now the error is guaranteed to have alternating signs. > > 5) Find a new set of N+1 points that have the largest *peaks* of E(w) AND > have alternating signs. > > 6) Solve a new linear system of N+1 equations as above. > > Continue 5 and 6 until the peaks and the equalized error values are equal. > The equalized error value will increase at each step and the absolute > maximum peak will reduce at each step until they are equal. > > That said, I doubt it will yield a shorter filter. > > Fred
It was late and I was in a hurry. It should have said, says now, e*(-1)^k Here's a bit more information: The N+1 equations to be solved have N basis function coefficients and "e" as the N+1 variables. So, the equations look more like this: A(k) + e*(-1)^k = D(k) so: a11*r1(1') + a12*r2(1') + .... + a1N*rN(1') + e*(-1)^1 = D(1') a21*r1(2') + a22*r2(2') + .... + a2N*rN(2') + e*(-1)^2 = D(2') .... aN1*r1(N+1') + aN2*r2(N+1') + .... + aNN*rN(N+1') + e*(-1)^(N+1) = D(N+1') where 1', 2', 3' ...., N+1' refer to the ordered values of frequency at which the error will be equalized at each step and the r's are the basis functions. One can also continuously weight the approximation error using D(w) - A(w) = W(w)*E(w) So, an otherwise equiripple approximant can have growing ripple at the edges of passbands, smaller stopband gain in regions of stopbands, etc. Large W yields smaller E. Fred
I have written an adaptable filter using a modification of the standard 
Parks-McClellan code available on the web.  Someone else wrote a front-end in 
X-Windows to adjust the sliders.  It isn't hard to do this though you won't 
use more than a hundred taps of so for the equalizer.  You can also use the 
Inverse DFT for this calculation.  I did both the the IFFT was faster and more 
stable.  32K taps is very noisy computationally and you wouldn't to use that 
many for any realtime processing unless you were (very, very) desperate.  

What you may want to do is some analysis on the minimum frequency and 
binwidth to determine the number of taps and their location.  Don't forget 
that if you use multiple filters with different delay characteristics, you 
may introduce differential delay distortion that may cause you some problems..

In article <449a595e$0$11079$9b4e6d93@newsread4.arcor-online.net>, 
=?ISO-8859-15?Q?Marcel_M=FCller?= <news.5.maazl@spamgourmet.com> wrote:
>Hi, > >I want to optimize a 32 band EQ plugin. Currently it uses inverse FFT >for the creation of the FIR filter kernel. But this requires rather long >FIR kernels to process low frequencies properly. While computing power >is nowadays not an issue the latency is. > >I heard about the "Remez Exchange Algorithmn". Unfortunatly I did not >find a decription how to implement this in the web. Only the keyword is >everywhere. > >Does anybody know where to get information about this algorithmn? > > >Marcel
Many thanks for all the answers.
At least they brought me much closer to a solution.

Fred Marshall schrieb:
>>Here is a quick description of the algorithm which can be applied to FIR >>filters among other functional approximation applications:
[...] The algorithm is really funny and I did not really understand it. Especially the trick with the alternation errors. However, I did not implement it so far either because maybe you are right and it is far too complicated for the yield.
>>That said, I doubt it will yield a shorter filter.
I think there is potential anyway. I did use a much easier approach now which gives me reasonable results without as much effort. Maybe it is suitable. I used a very aggessive window function (.75+.25*cos...) and compensated for their frequency respone by deconvoluting the desired filter respose by the frequency response of the window function before the IDCT. Because the natural frequency scale in audio applications is logarithmic, I have to apply this correction only to the first 100 points (or so). The roll-off is significantly improoved. Now I get +-3dB with the FIR length of 12288. Before the optimization it was about +-8dB. (Hopefully the simulation was right.) Marcel