Расчет коэффициента корреляции FFT

Question

Я хотел бы рассчитать коэффициент корреляции для 2 образцов звука, используя AForge 2.2.5.

Я читал от here формулу расчета кросс-корреляции.
И here Я читал о формуле для вычисления коэффициента корреляции.

В настоящее время у меня есть:
До вызова CrossCorrelation() выполнено FFT.

static Complex[] CrossCorrelation(Complex[] ffta, Complex[] fftb) 
{ 
    var conj = ffta.Select(i => new Complex(i.Re, -i.Im)).ToArray(); 

    for (int a = 0; a < conj.Length; a++) 
     conj[a] = Complex.Multiply(conj[a], fftb[a]); 

    FourierTransform.FFT(conj, FourierTransform.Direction.Backward); 

    return conj; 
} 

static double CorrelationCoefficient(Complex[] ffta, Complex[] fftb) 
{ 
    var correlation = CrossCorrelation(ffta, fftb); 
    var a = CrossCorrelation(ffta, ffta); 
    var b = CrossCorrelation(fftb, fftb); 

    // Not sure if this part is correct.. 
    var numerator = correlation.Select(i => i.SquaredMagnitude).Max(); 
    var denominatora = a.Select(i => i.Magnitude).Max(); 
    var denominatorb = b.Select(i => i.Magnitude).Max(); 

    return numerator/(denominatora * denominatorb); 
} 

Я не уверен, что это правильный путь для реализации функции (или обработки данных), как я новичок в ЦОС. Будет очень признателен, если кто-то сможет указать мне в правильном направлении.

Answer 1

Чтобы сделать кросс-корреляцию с FFTs и Afrog:

• Padding сигналов с нулями: Согласно AForge документации FFT: Метод принимает массив данных 2n размера только, где п может изменяться в [1, 14].

Так что вам нужно будет убедиться, что размер входного правильно проложенный до длины, которая является степенью 2, и в указанном диапазоне: принимая во внимание, по крайней мере половина волны «пустой» (нули)

исх:

https://dsp.stackexchange.com/questions/741/why-should-i-zero-pad-a-signal-before-taking-the-fourier-transform

https://dsp.stackexchange.com/questions/1919/efficiently-calculating-autocorrelation-using-ffts

• принять БПФ обоих сигналов
• многократно один с конъюгатом другого (поэлементного умножения)
• сделать обратное FFT
• получить максимальное значение в качестве коэффициента corellation и его индекса в качестве задержки (лаг сигнала)

Почему максимальное значение результирующего IFFT:

от wikipedia, кросс-корреляции являются полезными для определения временной задержки между двумя сигналами, например, для определения временных задержек для распространения акустических сигналов через микрофонную решетку. 2 [3] [разъяснение необходимо] После вычисления взаимной корреляции между двумя сигналами максимальный (или минимальный, если сигналы отрицательно коррелированы) функции взаимной корреляции указывает момент времени, когда сигналы наилучшим образом выровнены, т.е.задержка времени между двумя сигналами определяется аргументом максимума, или агд макс взаимной корреляции, как и в

исх: https://math.stackexchange.com/questions/1080709/why-is-the-maximum-value-of-cross-correlation-achieved-at-similar-section

На основании приведенных выше пунктов, Крест Расчет рассчитывается используя следующий код:

//same as OP 
    public Complex[] CrossCorrelation(Complex[] ffta, Complex[] fftb) 
    { 
    var conj = ffta.Select(i => new Complex(i.Re, -i.Im)).ToArray(); 

    conj = conj.Zip(fftb, (v1, v2) => Complex.Multiply(v1, v2)).ToArray(); 
    FourierTransform.FFT(conj, FourierTransform.Direction.Backward); 

    //get that data and plot in Excel, to show the max peak 
    Console.WriteLine("---------rr[i]---------"); 
    double[] rr = new double[conj.Length]; 
    rr = conj.Select(i => i.Magnitude).ToArray(); 

    for (int i = 0; i < conj.Length; i++) 
     Console.WriteLine(rr[i]); 

    Console.WriteLine("----end-----"); 
    return conj; 
    } 

//tuble.Item1: Cor. Coefficient 
//tuble.Item2: Delay of signal (Lag) 
public Tuple<double, int> CorrelationCoefficient(Complex[] ffta, Complex[] fftb) 
{ 
    Tuple<double, int> tuble; 
    var correlation = CrossCorrelation(ffta, fftb); 
    var seq = correlation.Select(i => i.Magnitude); 
    var maxCoeff = seq.Max(); 

    int maxIndex = seq.ToList().IndexOf(maxCoeff); 
    Console.WriteLine("max: {0}", maxIndex); 
    tuble = new Tuple<double, int>(maxCoeff, maxIndex); 
    return tuble; 
} 
    // Pad signal with zeros up to 2^n and convert to complex 
    public List<Complex> ToComplexWithPadding(List<double> sample, int padding = 1) 
    { 
     //As per AForge documentation: 
     // The method accepts data array of 2n size only, where n may vary in the [1, 14] range 
     //So you would need to make sure the input size is correctly padded to a length that is a power of 2, and in the specified range: 

     double logLength = Math.Ceiling(Math.Log(sample.Count * padding, 2.0)); 
     int paddedLength = (int)Math.Pow(2.0, Math.Min(Math.Max(1.0, logLength), 14.0)); 
     Complex[] complex = new Complex[paddedLength]; 
     var samples = sample.ToArray(); 
     // copy all input samples 
     int i = 0; 
     for (; i < sample.Count; i++) 
     { 
      complex[i] = new Complex(samples[i], 0); 
      Console.WriteLine(complex[i].Re); 

     } 
     // pad with zeros 
     for (; i < paddedLength; i++) 
     { 
      complex[i] = new Complex(0, 0); 
      Console.WriteLine(complex[i].Re); 
     } 
     return complex.ToList(); 

    } 

    // extra for signal generation for testing. You can find in the link of the life demo. 

Вы можете запустить жизнь демо с образцом два сигнала, генерируемого с задержкой 11 и результат соответствует фактической задержки сигнала

Life demo with two signals generated

Выходной результат:

correlation Coef: 0.33867796353274 | Delay: 11| actual delay: 11