Comment puis-je vectoriser les calculs pour un filtre récursif de premier ordre?

9

J'ai un simple filtre passe-bas unipolaire (pour le lissage des paramètres) qui peut être expliqué par la formule suivante:

y[n]=(1a)y[n1]+ax[n]

L'architecture que j'utilise a accès à des instructions à instruction unique et à données multiples (SIMD) qui peuvent effectuer plusieurs calculs vectorisés en parallèle. Je voudrais profiter de cette capacité, mais je ne sais pas comment faire cela pour un filtre récursif comme celui-ci. Le problème est que chaque calcul a besoin d'un résultat précédent.

user1132968
la source
Quelqu'un pourrait-il préciser pourquoi cela a été fermé comme "hors sujet"?
Paul R
La question chevauche ici et Stack Overflow. La question d'origine demandait spécifiquement comment l'implémenter à l'aide des extensions ARM NEON. La question vivra sur les deux sites; il a été édité ici pour en faire davantage une discussion théorique sur la structuration du problème pour tirer parti du parallélisme.
Jason R
@PaulR J'ai demandé à ce qu'il soit migré ici hier, mais phonon était fermement convaincu qu'il appartenait à SO, et non ici. J'avoue, je ne connais pas ARM NEON en particulier, et il a probablement raison et je respecte son jugement. Voir cette méta-question . Ma réponse supprimée disait essentiellement que je suis d'accord avec Jason R et que les utilisateurs devraient signaler ces questions pour la migration
Lorem Ipsum
1
Merci à tous les deux pour la clarification - j'ai fait de mon mieux sur SO pour faire migrer ici les questions liées au DSP afin de rehausser notre profil, donc je me demandais si c'était la bonne chose à faire lorsque cette question a été fermée - heureux de le voir à nouveau ouvert sous une forme plus générale.
Paul R

Réponses:

7

En supposant que vous effectuez des opérations vectorielles M éléments à la fois, vous pouvez dérouler l'équation de la différence par un facteur de Massez facilement pour le simple filtre unipolaire. Supposons que vous avez déjà calculé toutes les sorties jusqu'ày[n]. Ensuite, vous pouvez calculer les suivants comme suit:

y[n+1]=(1a)y[n]+ax[n+1]y[n+2]=(1a)y[n+1]+ax[n+2]=(1a)((1a)y[n]+ax[n+1])+ax[n+2]=(1a)2y[n]+a(1a)x[n+1]+ax[n+2]

En général, vous pouvez écrire y[n+k] comme:

y[n+k]=(1a)ky[n]+i=1ka(1a)kix[n+i]

Pour chaque indice d'échantillon n+k, cela ressemble à un filtre FIR avec k+1 taps: un tap multiplie la dernière sortie du filtre y[n], et l'autre k les robinets multiplient les entrées du filtre x[n+1],,x[n+k]. La bonne chose est que les coefficients utilisés pour tous ces taps peuvent être précalculés, vous permettant de dérouler le filtre récursif enM M+1-appuyez sur des filtres non récursifs parallèles (ces filtres calculent les échantillons de sortie y[n+1],,y[n+M]), en mettant à jour le terme de rétroaction tous les Méchantillons de sortie. Donc, étant donné une condition initialey[n] (qui est supposé être la dernière sortie calculée sur la précédente itération vectorisée), vous pouvez calculer la prochaine M sorties en parallèle.

Il y a quelques mises en garde à cette approche:

  • Si Mdevient grand, vous finissez par multiplier un ensemble de nombres afin d'obtenir les coefficients FIR effectifs pour les filtres déroulés. Selon le format de votre numéro et la valeur dea, cela pourrait avoir des implications de précision numérique.

  • De plus, vous ne recevez pas de Mmultipliez par deux cette approche: vous finissez par calculer y[n+k] avec ce qui équivaut à un k-filtre le filtre FIR. Bien que vous calculiezM sorties en parallèle, le fait que vous devez faire kles opérations de multiplication-accumulation (MAC) au lieu de la simple mise en œuvre récursive de premier ordre diminuent certains des avantages de la vectorisation. L'approche non vectorisée utilise 2 MAC par sortie, vous avez donc besoin2M opérations pour calculer Mles sorties. Le schéma vectorisé calcule laM sorties à la fois, nécessitant M+1MAC en cours. Ainsi, la réduction des opérations peut être exprimée en fonction deM comme:

    R=M+12M=12(1+1M)

    Donc, même avec Mtrès grand, vous pouvez obtenir une réduction maximale de 50% du nombre de calculs en utilisant cette méthode. Pour les valeurs communes deM=4 et M=8, la réduction est respectivement de 37,5% et 43,75%. Cependant, en plus de la réduction pure du nombre d'opérations, vous pouvez bénéficier d'avantages de performances supplémentaires en raison du modèle d'accès à la mémoire différent utilisé par l'approche non déroulée; pour une implémentation plus simple, vous pouvez rencontrer des retards en raison de la dépendance de chaque échantillon de sortiey[n] sur l'échantillon précédent y[n1]. Cet effet est évidemment très dépendant de la plateforme.

Jason R
la source
3

En général, vous ne pouvez vectoriser que des ensembles de calculs complètement indépendants. Mais dans votre passe-bas IIR, chaque sortie dépend d'une autre (sauf la 1ère), donc la vectorisation n'est pas possible.

Si votre variable "a" est suffisamment grande pour que (1-a) ^ n décroisse rapidement en dessous du plancher de bruit souhaité ou de l'erreur autorisée, vous pouvez remplacer une courte approximation de filtre FIR pour votre IIR et vectoriser cette convolution à la place. Mais ce n'est pas susceptible d'être plus rapide.

hotpaw2
la source
2

Vous ne pouvez vraiment vraiment le vectoriser efficacement que si vous avez plus d'un signal auquel vous souhaitez appliquer le même filtre, par exemple s'il s'agit d'un signal audio stéréo, vous pouvez traiter les canaux gauche et droit en parallèle. Quatre ou huit canaux en parallèle seraient évidemment encore meilleurs.

Paul R
la source
0

Que diriez-vous d'étendre les équations à 4 étapes et d'utiliser la multiplication matricielle? a est constant de sorte qu'une matrice peut être précalculée


la source
0

Il est plus efficace de vectoriser en parallèle le filtrage de plusieurs flux indépendants.

Si vous n'avez qu'un seul long flux, vous pouvez également diviser le flux en plusieurs sections qui se chevauchent et les filtrer comme s'il s'agissait de flux indépendants.

Vous voulez un chevauchement car les premiers échantillons de chaque flux seront incorrects. La quantité d'échantillons que vous devez jeter dépendra de la valeur de a et de la précision souhaitée. Vous voudrez éliminer environ n échantillons où (1 / a) (1-a) ** n <eps afin que le résultat soit précis à eps * max (| x [i] |).

Par exemple, si a = 0,1, alors dans 128 échantillons, l'erreur provoquée par (1 / a) (1-a) ^ n serait inférieure au LSB dans un entier de 16 bits, tandis que pour a = 0,9, vous n'auriez qu'à supprimer environ 6 échantillons.

Peter de Rivaz
la source