Minimum phase transformation

[Aleksey Vaneev]

Why minimum phase FIR?

Because it offers a different sound without pre-ringing, and FIR because you can define filters of any complexity with it. It can be hard to get the same thing with a serie of biquads.

[takoha]

It's not really possible. But you can get around that by defining presets with generic filenames like 'impulse1.wav', 'impulse2.wav', etc and then replacing these WAV files. Not the best way, of course.

Oki, that is a way.
storing predefined config files and then restore them (as data file) and generate appropriate impulse file.

Another way maybe, generate config file of the vst host (like Console, cause it memorise all vst plugin's parametre) but still not the best way

Regards
Avi

[takoha]

we did it using the hilbert transformate.

Could you show the algo for that?

maybe see my first post for a link on the hilbert transform, but this one didn't work fine in all the case for me (certainly i forgot something in the algo)
(it s the same that i found in wikipedia, which gives me the same result)

Avi

[mystran]

Why minimum phase FIR?

Because it offers a different sound without pre-ringing, and FIR because you can define filters of any complexity with it. It can be hard to get the same thing with a serie of biquads.

Hehe, actually only after having asked did I realize one might need a really complex filter and prefer implementing it with FFT convolution or something.

Couldn't think of a domain where I'd prefer minimum phase yet would be happy with a static filter design, but I guess your pre-ringing point is a valid one. Now that I think of that though, I'm really started to wonder how to balance the ringing so that it optimally falls into the pre- and post-masking regions (that would be "almost minimum phase I guess"). Guess it's another night in bed with mathematics.

Anyway, I said I'm gonna ask a stupid question 'cos I had the feeling it must have been a stupid one, and I now realize I actually can think of more and more domains where that would be useful. Just shows I've spent too much time tinkering with various IIR structures (saturation behaviors, transient responses, all that crap) and probably should be spending more time with FIR tricks in the future.

Well.. such is life. Thanks anyway.

[Zaphod (giancarlo)]

we did it using the hilbert transformate.

Could you show the algo for that?

... do you mean for the hilbert transformate? we are using a commercial library coming from intel for that (ipps), so it is pretty fast and we kept the implementation really simple. I don't remember it very well, but you should just perform the hilber transformate of the log of the module, and take the imaginary part of it. It is the new phase. But my memory is weak, and my code is not simple to understand even for me ( ), it's better if you google it.

[Aleksey Vaneev]

we did it using the hilbert transformate.

Could you show the algo for that?

... do you mean for the hilbert transformate? we are using a commercial library coming from intel for that (ipps), so it is pretty fast and we kept the implementation really simple. I don't remember it very well, but you should just perform the hilber transformate of the log of the module, and take the imaginary part of it. It is the new phase. But my memory is weak, and my code is not simple to understand even for me ( ), it's better if you google it.

OK, understood. I dislike Hilbert transform (HT) myself for such tasks because, well, you can't have an ideal HT implemented in finite math because response of HT is infinite, and performance of such operation will be limited by window's length (the most low and most high frequencies will be missing from such transform).

[Henke]

I've successfully used real cepstrum for this task.

[kasaudio]

As I found Aleksey Vaneev's solution to work perfect, I enclose source code that I made according to his algorithm.

The code is quite self explanatory, it is in c# and is only for real 1D transformation. For FFT I use Intel IPP's routines around which I made a .net wrapper.

Code: Select all

        /// <summary>
        /// Transforms input filter coefficients to minimum phase coefficients.
        /// </summary>
        /// <param name="input">Input filter coefficients.</param>
        /// <returns>Minimum phase transformed filter coefficients.</returns>
        public static double[] Transform(double[] input)
        {
            int n = input.Length;
            int power = (int)Math.Ceiling(Math.Log(n, 2d)) + 3; // Larger fft size instead double min required to prevent aliasing...
            int fftSize = (int)Math.Pow(2, power);
            int halfFttSize = fftSize / 2;
            int halfFftSizeP1 = halfFttSize + 1;
            double[] kernel = new double[fftSize];
            double[] kernelCss = new double[fftSize + 2];
            double[] magnitudes = new double[halfFftSizeP1];
            Array.Copy(input, kernel, n);
            // Double precision FFT transformer.
            FFTDouble fft = new FFTDouble(fftSize);

            // Forward transformation. Result is halfFttSize+1 real-imaginary pairs.
            fft.FftForward(kernel, 0, kernelCss, 0);

            int offset=0, offsetP1;
            double real, imaginary, magnitude;

            for (int i = 0; i < halfFftSizeP1; i++)
            {
                offsetP1 = offset + 1;
                real = kernelCss[offset];
                imaginary = kernelCss[offsetP1];
                magnitude = Math.Sqrt(real * real + imaginary * imaginary);
                magnitudes[i] = magnitude;
                kernelCss[offset] = Math.Log(magnitude);
                kernelCss[offsetP1] = 0d;
                offset += 2;
            }

            // Backward transformation.
            fft.FftBackward(kernelCss, 0, kernel, 0);

            for (int i = halfFttSize; i < fftSize; i++)
                kernel[i] *= -1d;

            fft.FftForward(kernel, 0, kernelCss, 0);

            offset = 0;

            for (int i = 0; i < halfFftSizeP1; i++)
            {
                offsetP1 = offset + 1;
                imaginary = kernelCss[offsetP1];
                real = Math.Cos(imaginary) * magnitudes[i];
                imaginary = Math.Sin(imaginary) * magnitudes[i];
                kernelCss[offset] = real;
                kernelCss[offsetP1] = imaginary;
                offset += 2;
            }

            fft.FftBackward(kernelCss, 0, kernel, 0);
            fft.Free();
            fft = null;
            double[] result = new double[n];
            Array.Copy(kernel, result, n);
            return result;
        }

Best regards, kasaudio