Comment fonctionne un filtre passe-bas par programmation?

9

Je travaille sur un simple filtre passe-bas pour une mesure <100 Hz dans mon application. Mais jusqu'à présent, je lutte avec la théorie derrière tout cela. C'est cool que je l'ai fait fonctionner, mais je l'apprécierais vraiment si je savais comment / pourquoi cela fonctionne.

J'ai trouvé le code suivant:

void getLPCoefficientsButterworth2Pole(const int samplerate, const double cutoff, double* const ax, double* const by)
{
    double PI = M_PI;
    double sqrt2 = sqrt(2);

    double QcRaw  = (2 * PI * cutoff) / samplerate; // Find cutoff frequency in [0..PI]
    double QcWarp = tan(QcRaw); // Warp cutoff frequency

    double gain = 1 / ( 1 + sqrt2 / QcWarp + 2 / ( QcWarp * QcWarp ) );

    by[2] = ( 1 - sqrt2 / QcWarp + 2 / ( QcWarp * QcWarp ) ) * gain;
    by[1] = ( 2 - 2 * 2 / ( QcWarp * QcWarp ) ) * gain;
    by[0] = 1;

    ax[0] = 1 * gain;
    ax[1] = 2 * gain;
    ax[2] = 1 * gain;
}

Pour calculer les coefficients. Ensuite, dans les échantillons audio, je les «passe-bas» de cette façon:

        xv[2] = xv[1];
        xv[1] = xv[0];

        xv[0] = pData[j];
        yv[2] = yv[1];
        yv[1] = yv[0];

        yv[0] = (ax[0] * xv[0] + ax[1] * xv[1] + ax[2] * xv[2]
                   - by[1] * yv[0]
                   - by[2] * yv[1]);

        pData[j] = yv[0];

Pour obtenir une conception passe-bas.

Je me pose quelques questions:

  1. Je reçois les échantillons audio dans un simple tableau flottant *. Quel est ce nombre flottant? La seule chose que je vois est un nombre, comment est-ce un morceau de son?
  2. Le code utilise des calculs antérieurs (trois d'entre eux) dans le nouveau calcul par échantillon. Est-ce à dire que les 2 premiers échantillons de données ne sont pas filtrés correctement? (ce n'est pas important car ce ne sont que 2 échantillons, mais je me demande simplement)
  3. En essayant de tout apprendre, j'ai trouvé quelques formules pour le filtre Butterworth (2e pôle). Comment ces formules sont-elles reflétées dans ce code? Aucune des formules que j'ai trouvées n'a ces calculs que vous pouvez voir dans la fonction 'getLPCoefficientsButterworth2Pole ()'.
Niek van der Steen
la source
1
Je n'essaie pas d'être irrespectueux ici, mais votre déclaration "comment est-ce un morceau de son?" semble indiquer que vous ne comprenez pas les principes de base du traitement discret dans le temps. Il sera très difficile de résoudre tout problème DSP sans comprendre les bases telles que le théorème d'échantillonnage, la quantification, le système LTI, etc. Je recommanderais un peu de temps avec un bon manuel. Celui-ci est gratuit dspguide.com/pdfbook.htm
Hilmar

Réponses:

6
  1. Le nombre float * tableau est un pointeur sur le tableau. Il s'agit d'un nombre unique qui contient l'adresse du premier élément du tableau de valeurs flottantes.

  2. Habituellement, la condition initiale (c'est-à-dire les éléments «passés» initiaux de x et y) est 0, mais si leurs valeurs ne sont pas égales à 0 ce n'est pas un gros problème non plus, car après un certain temps, les conditions initiales n'ont aucun effet sur la sortie signal pour tout filtre stable. Et votre filtre est évidemment stable.

  3. Une fonction de transfert passe-bas du deuxième ordre avec caractéristique de Butterworth et fréquence de coupure ωune=2πFune (avec Fune en Hertz) est donné par

(1)H(s)=ωune2s2+2ωunes+ωune2
Il s'agit d'un résultat standard que vous pouvez facilement trouver sur le Web. Afin d'obtenir un filtre à temps discret, vous pouvez appliquer une soi-disant transformation bilinéaire:

(2)s=2Fsz-1z+1

Fsest la fréquence d'échantillonnage. Cela est nécessaire car l'axe de fréquence des signaux analogiques (0F) doit être mappé sur la plage de fréquences autorisée pour les signaux à temps discret (0FFs/2). Comme cette transformation déforme les fréquences, nous devons calculer la fréquence de coupure analogique souhaitée à partir de la fréquence de coupure donnéeF dans le domaine du temps discret:

ωune=2Fsbronzer(ω2) avec ω=2πF/Fs

Si vous insérez (2) dans (1) vous obtenez

(3)H(z)=kz2+b1z+b2z2+une1z+une2

avec

k=α21+2α+α2 (Gain)
une1=2(α2-1)1+2α+α2,une2=1-2α+α21+2α+α2,b1=2,b2=1

où j'ai utilisé α=bronzer(ω2). La fonction de transfert (3) correspond à l'équation du domaine temporel

y(n)=kX(n)+kb1X(n-1)+kb2X(n-2)-une1y(n-1)-une2y(n-2)

Comme vous pouvez le voir, il y a une certaine ressemblance avec votre code. Cependant, il existe également des différences et vous devez vérifier la réponse en fréquence de votre filtre. Je pense que les équations ci-dessus sont correctes, mais c'est à vous de les vérifier et de décider quelle version est la bonne.

Matt L.
la source
Merci. Le premier, je comprends. Mais c'est toujours une série de flotteurs. Je veux dire: que représente une seule valeur flottante?
Niek van der Steen
2
Une seule valeur flottante est l'un de vos échantillons d'entrée.
Matt L.
3
En bref: le son est une variation de la pression atmosphérique, convertie par un transducteur (microphone) en un signal électrique dont les variations suivent le même schéma que la pression atmosphérique. L'amplitude du signal électrique est mesurée N fois (= échantillonnage) chaque seconde pour produire la séquence de nombres que vous voyez dans votre tableau flottant *. Il existe des dizaines de présentations en ligne sur le traitement (audio) numérique expliquant ce processus de manière plus détaillée.
pichenettes
J'ai modifié ma réponse pour répondre à la dernière partie de votre question.
Matt L.
Vous utilisez x (n - quelque chose), qu'est-ce que x? Je suppose que «n» est l'échantillon d'entrée? Bonne réponse!
Niek van der Steen
11

Vous avez demandé comment fonctionne un filtre passe-bas et avez mentionné que le filtre utilise les valeurs passées de vos données. Il s'agit d'une discussion non technique de ce qui se passe dans un filtre passe-bas.

Le filtre passe-bas prend des vues différentes (décalées dans le temps) de votre signal, les met à l'échelle et les ajoute ensemble. Vous pouvez imaginer dessiner votre signal 3 fois, l'un étant courant, le second étant décalé d'un temps d'échantillonnage, le troisième étant décalé de 2 temps d'échantillonnage.

Aux basses fréquences, toutes les vues sont très similaires (le décalage d'un seul échantillon change à peine où vous êtes sur le signal à tout instant). Dans ce cas, les trois versions s'ajouteront de manière constructive (ou du moins non destructive), de sorte que le signal passe à travers le filtre.

Passant maintenant à des fréquences plus élevées, chaque version décalée du signal devient plus distincte à un instant donné (point d'échantillonnage) et peut même en fait s'inverser en signe. À ces fréquences plus élevées, les trois versions de votre signal ont tendance à s'annuler (ajoutées de manière destructive) afin que le signal soit atténué.

Différents types de filtres organisent ce type d'interférence constructive / destructrice sur les bandes de fréquences appropriées pour créer des filtres passe-bas, passe-bande ou passe-haut.

user2718
la source
1

Cela dépend de la façon dont vous le souhaitez. Pour moi, j'ai mis en œuvre en C. Donc, ce que j'ai fait, c'est que j'ai généré des coefficients de filtre en matlabutilisant firls. Ensuite, j'ai stocké ces coefficients de filtre dans un tableau, puis je passe ces coefficients de filtre à la fonction de filtre dans C. Il matlabest assez facile de générer les coefficients de filtre, puis d'utiliser l'opération de convolution ou d'utiliser un filtre spécifique. Mais Cvous devez également implémenter le filtre. Pour ma part, je dois le faire pour les DSP, donc implémenté dans Cet obtenu les résultats en matlabutilisant la fonction mex. Pour calculer les coefficients, utilisez cette commande dans matlab:

n=17                  %filter order
f=[0 0.167 0.333 1]   %Frequency band edges
a=[1 1 0 0]           %Desired amplitudes
fir= firls (n, f,a )
DX
la source