1d linear convolution in ANSI-C-code?

Anstatt das Rad neu zu erfinden, ich Frage mich, wenn jemand könnte sich mich auf eine 1D linear convolution code-snippet in ANSI C? Ich habe eine Suche auf google und in-stack-overflow, aber konnte nicht finden, dass irgendetwas in C, die ich verwenden könnte.

Beispielsweise für Arrays A, B und C, alle mit doppelter Genauigkeit, wobei A und B sind die Eingänge und C ausgegeben wird, mit Längen len_A, len_B, und len_C = len_A + len_B - 1 bzw.

Mein array so klein sind und so jeder Erhöhung der Geschwindigkeit bei der Umsetzung der schnellen Faltung durch FFT ist nicht notwendig. Suchen für einfache Berechnungen.

Auf welcher Plattform sind Sie targeting? Es ist durchaus möglich, dass eine solche Funktion existiert bereits und Sie können es benutzen.
Ich bin mit gcc4.4.4 und centos 5.7 auf linux 64 bit-server.
Nach viel suchen, fand ich den folgenden code für LinearConvolution(), die in meiner Implementierung funktioniert es schnell und produziert die gleichen Ergebnisse wie Matlab, obwohl der code nicht ganz so einfach zu Lesen und zu verstehen, wie Alex ' s unten. Nicht sicher, ob es irgendwelche Unterschiede gibt ansonsten. Anfangs dachte ich, das war C++ - code, aber es scheint zu laufen gut, wenn kompiliert mit -ansi-Schalter in meinem C-Programm. Ich habe verbunden es hier im Fall ist es nützlich, um andere. dsprelated.com/showmessage/71405/1.php

InformationsquelleAutor ggkmath | 2011-12-07

  1. 20

    Hier ist, wie:

    #include <stddef.h>
#include <stdio.h>

void convolve(const double Signal[/* SignalLen */], size_t SignalLen,
              const double Kernel[/* KernelLen */], size_t KernelLen,
              double Result[/* SignalLen + KernelLen - 1 */])
{
  size_t n;

  for (n = 0; n < SignalLen + KernelLen - 1; n++)
  {
    size_t kmin, kmax, k;

    Result[n] = 0;

    kmin = (n >= KernelLen - 1) ? n - (KernelLen - 1) : 0;
    kmax = (n < SignalLen - 1) ? n : SignalLen - 1;

    for (k = kmin; k <= kmax; k++)
    {
      Result[n] += Signal[k] * Kernel[n - k];
    }
  }
}

void printSignal(const char* Name,
                 double Signal[/* SignalLen */], size_t SignalLen)
{
  size_t i;

  for (i = 0; i < SignalLen; i++)
  {
    printf("%s[%zu] = %f\n", Name, i, Signal[i]);
  }
  printf("\n");
}

#define ELEMENT_COUNT(X) (sizeof(X) / sizeof((X)[0]))

int main(void)
{
  double signal[] = { 1, 1, 1, 1, 1 };
  double kernel[] = { 1, 1, 1, 1, 1 };
  double result[ELEMENT_COUNT(signal) + ELEMENT_COUNT(kernel) - 1];

  convolve(signal, ELEMENT_COUNT(signal),
           kernel, ELEMENT_COUNT(kernel),
           result);

  printSignal("signal", signal, ELEMENT_COUNT(signal));
  printSignal("kernel", kernel, ELEMENT_COUNT(kernel));
  printSignal("result", result, ELEMENT_COUNT(result));

  return 0;
}

    Ausgabe:

    signal[0] = 1.000000
signal[1] = 1.000000
signal[2] = 1.000000
signal[3] = 1.000000
signal[4] = 1.000000

kernel[0] = 1.000000
kernel[1] = 1.000000
kernel[2] = 1.000000
kernel[3] = 1.000000
kernel[4] = 1.000000

result[0] = 1.000000
result[1] = 2.000000
result[2] = 3.000000
result[3] = 4.000000
result[4] = 5.000000
result[5] = 4.000000
result[6] = 3.000000
result[7] = 2.000000
result[8] = 1.000000
    Vielen Dank Alex, ich habe gerade kopiert/umgesetzt, die convolve () - Funktion oben in meinem Programm und es scheint zu schnell zu laufen und produzieren im wesentlichen die genau gleichen Ergebnisse als Matlab - Funktion conv, die tun sollten, die gleiche Sache (es ist z.B. nicht FFT-basiert). Der code ist auch Recht einfach zu Lesen, vielen Dank. Gibt es irgendwelche Annahmen oder anderen Dinge zu beachten (limits, etc.) hier?
    Längen und Hinweise scheinbar muss nicht 0 sein. Auch SignalLen + KernelLen - 1 <= (size_t)-1. Floating-point-Werte gelten als gut. Nichts ungewöhnliches.

    InformationsquelleAutor Alexey Frunze

  2. 3

    Nicht getestet, aber es scheint, wie es funktionieren würde...

    void conv(const double v1[], size_t n1, const double v2[], size_t n2, double r[])
{
    for (size_t n = 0; n < n1 + n2 - 1; n++)
        for (size_t k = 0; k < max(n1, n2); k++)
            r[n] += (k < n1 ? v1[k] : 0) * (n - k < n2 ? v2[n - k] : 0);
}

    Tipp: Wenn es dauert weniger Zeit, um ein Rad neu erfinden, als eine zu finden, halten die ehemalige.

    Danke Mehrdad, ich verstehe den beiden Schleifen, aber die bedingten Anweisungen in der r[n] Ausdruck-sind Sie einfach verhindert Multiplikation array-Werte mit negativen Indizes? Oder können Sie erklären Sie diese kurz? Außerdem sollte r[n] initialisiert zu 0 nach dem ersten loop und vor dem zweiten loop, oder sind Sie unter der Annahme r ist bereits mit null initialisiert, für alle Elemente?
    Ich gehe davon aus, die Ergebnisse sind zunächst null. Und ja, die Bedingungen sind nur zu verhindern, dass out-of-bound-Indizes.
    Ich glaube nicht, dass der max () - Funktion von ansi-C-kompatibel, so dass ich gerade eine Funktion zu tun, statt. Mit einem kleinen test-array funktionierte es gut. Wenn ich meinen real-arrays (Größe etwa 20K Elemente), es irgendwie enthalten ein ghost-image von einem der Eingabe-arrays in der Ausgabe. Entweder habe ich einen Fehler gemacht bei der Umsetzung der obigen Funktion, oder die Datensätze zwischen diesem test und den ersten test unterscheidet sich von dem Algorithmus s point of view (nicht sicher).
    Ich kann wirklich nichts sagen, ohne zu sehen, Ihr Beispiel.
    Sie nicht null r[n] bevor die innere Schleife wird zu erhöhen, IOW, conv() erwartet r[] initialisiert alle 0.0's. Das war vielleicht das problem.

    InformationsquelleAutor Mehrdad

  3. 0

    Da, wir sind unter Faltung von 2 finite-Länge-Sequenzen, damit der gewünschte Frequenzgang erzielt wird, wenn die kreisförmige Faltung durchgeführt wird anstatt einer linearen Faltung. Eine sehr einfache Umsetzung der zirkulären Faltung wird das gleiche Ergebnis zu erzielen als der Algorithmus, gegeben durch Alex.

    #define MOD(n, N) ((n<0)? N+n : n)
......
......

for(n=0; n < signal_Length + Kernel_Length - 1; n++)
{
    out[n] = 0;
    for(m=0; m < Kernel_Length; m++)
    {
        out[n] = h[m] * x[MOD(n-m, N)];
    }
}

    InformationsquelleAutor Siddhant Raman

  4. 0

    Ich verwendet @Mehrdad Ansatz, und erstellt die folgende Antwort:

    void conv(const double v1[], size_t n1, const double v2[], size_t n2, double r[])
{
    for (size_t n = 0; n < n1 + n2 - 1; n++)
        for (size_t k = 0; k < max(n1, n2) && n >= k; k++)
            r[n] += (k < n1 ? v1[k] : 0) * (n - k < n2 ? v2[n - k] : 0);
}

    Gibt es problem mit index-überschreitung der unteren Grenze, wenn in der zweiten loops k größer wird als n, also denke mal es sollten zusätzliche Bedingung zu verhindern.

    InformationsquelleAutor Imro Stocky

Schreibe einen Kommentar