Éviter l'arithmétique en virgule flottante

Question

J'ai écrit un petit logiciel de synthèse pour iPhone.
Pour améliorer les performances, j'ai mesuré mon application avec Shark et constaté que je perdais beaucoup de temps dans les conversions float / SInt16.
J'ai donc réécrit certaines parties afin de contourner les conversions en pré-calculant les tables de consultation qui renvoient "Prêt à l'emploi". SInt16 échantillons. Cela fonctionne très bien jusqu'à présent.
Actuellement, j'essaie de réécrire certains filtres et l'implémentation de mon enveloppe ADSR pour utiliser uniquement l'arithmétique entière, mais je pourrais utiliser quelques conseils pour effectuer des multiplications / divisions sans flottants.
Je cible le format canonique de l'iPhone :

• LPCM
• échantillons de nombre entier 16 bits

Quelles sont les bonnes approches pour appliquer une amplitude à mon échantillon final sans utiliser de float?

Modifier:
La seule chose que j’ai compris jusqu’à présent, c’est que je peux diviser par deux par deux en décalant vers la droite mon échantillon actuel.

inBuffer[frame] = wavetable[i % cycleLengthInSamples] >> 4;

Mais je ne vois pas de moyen élégant de créer une enveloppe ADSR lisse avec cela.

Modifier2: Merci pour toutes vos bonnes réponses!
Mon approche actuelle:

• apporter toutes mes valeurs d'enveloppe ADSR dans la plage positive SInt16
• multiplier par la valeur actuelle de la table d'ondes (stocker les intermédiaires en tant que SInt32)
• décale le résultat de 16 vers la droite

cela semble fonctionner:)

Answer 1

Le point fixe est bon, car dans ce cas, vous utilisez 16 bits. Le moyen le plus simple est de multiplier par une puissance de 10 en fonction de la précision dont vous avez besoin. Si vous pouvez utiliser les ints 32 bits comme intermédiaire, vous devriez pouvoir obtenir une précision décente. À la fin, vous pouvez reconvertir en 16 bits, arrondir ou tronquer à votre guise.

Modifier: Vous voulez décaler à gauche, pour rendre les valeurs plus grandes. Stockez le résultat du décalage dans un type avec plus de précision (32 ou 64 bits, selon vos besoins). Le décalage simple ne fonctionnera pas si vous utilisez des types signés

Faites attention si vous multipliez ou divisez deux nombres à virgule fixe. Multiplier est (a * n) * (b n) et vous obtiendrez un b n ^ 2 au lieu de b n. La division est (a n) / (b n) qui est (a / b) au lieu de ((a n) / b). C'est pourquoi j'ai suggéré d'utiliser des puissances de 10, il est facile de trouver vos erreurs si vous n'êtes pas familier avec le point fixe.

Lorsque vous avez terminé vos calculs, vous revenez à droite pour revenir à un int de 16 bits. Si vous voulez avoir l’imagination, vous pouvez également arrondir avant de changer de poste.

Je vous suggère de lire un peu si vous êtes vraiment intéressé par la mise en œuvre de points fixes efficaces. http://www.digitalsignallabs.com/fp.pdf

Answer 2

Les réponses à cette question SO sont assez complet en termes de mise en œuvre. Voici un peu plus d'explication que ce que j'ai vu là-bas:

Une approche consiste à forcer tous vos nombres dans une plage, disons [-1.0,1.0). Ensuite, vous mappez ces nombres dans la plage [-2 ^ 15, (2 ^ 15) -1]. Par exemple,

Half = round(0.5*32768); //16384
Third = round((1.0/3.0)*32768); //10923

Lorsque vous multipliez ces deux nombres, vous obtenez

Temp = Half*Third; //178962432
Result = Temp/32768; //5461 = round(1.0/6.0)*32768

La dernière ligne est divisée par 32768 comme patros à propos de multiplie nécessitant une étape supplémentaire de mise à l'échelle. Cela a plus de sens si vous écrivez explicitement la mise à l'échelle 2 ^ N:

x1 = x1Float*(2^15);
x2 = x2Float*(2^15);
Temp = x1Float*x2Float*(2^15)*(2^15);
Result = Temp/(2^15); //get back to 2^N scaling

Voilà donc l'arithmétique. Pour l'implémentation, notez que la multiplication de deux entiers 16 bits nécessite un résultat de 32 bits, donc Temp doit être de 32 bits. De plus, 32768 n'est pas représentable dans une variable de 16 bits, alors sachez que le compilateur créera une version immédiate de 32 bits. Et comme vous l'avez déjà noté, vous pouvez passer à multiplier / diviser par 2 pour pouvoir écrire

N = 15;
SInt16 x1 = round(x1Float * (1 << N));
SInt16 x2 = round(x2Float * (1 << N));
SInt32 Temp = x1*x2;
Result = (SInt16)(Temp >> N);
FloatResult = ((double)Result)/(1 << N);

Mais supposons que [-1,1) ne soit pas la bonne plage? Si vous préférez limiter vos nombres à, par exemple, [-4.0,4.0), vous pouvez utiliser N = 13. Ensuite, vous avez 1 bit de signe, deux bits avant le point binaire et 13 après. Ceux-ci sont appelés respectivement 1.15 et 3.13 types fractionnaires à point fixe. Vous échangez la précision dans la fraction de la marge disponible.

L'ajout et la soustraction de types fractionnaires fonctionnent correctement tant que vous recherchez la saturation. Pour diviser, comme dit Patros, la mise à l'échelle s'annule. Donc, vous devez faire

Quotient = (x1/x2) << N;

ou, pour préserver la précision

Quotient = (SInt16)(((SInt32)x1 << N)/x2); //x1 << N needs wide storage

Multiplier et diviser par des nombres entiers fonctionne normalement. Par exemple, pour diviser par 6, vous pouvez simplement écrire

Quotient = x1/6; //equivalent to x1Float*(2^15)/6, stays scaled

Et en cas de division par une puissance de 2,

Quotient = x1 >> 3; //divides by 8, can't do x1 << -3 as Patros pointed out

L'ajout et la soustraction de nombres entiers ne fonctionnent toutefois pas naïvement. Vous devez d’abord voir si le nombre entier s’intègre dans votre type x.y, créer le type fractionnaire équivalent et continuer.

J'espère que cela vous aidera avec l'idée. Regardez le code de cette autre question pour des implémentations propres.

Answer 3

Consultez cette page qui décrit les algorithmes de multiplication rapide.

http://www.newton.dep.anl.gov /askasci/math99/math99199.htm

Answer 4

En général, supposons que vous utilisiez une représentation 16,16 signée en points fixes. Ainsi, un entier de 32 bits aura une partie entière signée de 16 bits et une partie fractionnaire de 16 bits. Ensuite, je ne sais pas quelle langue est utilisée dans le développement d'un iPhone (Objective-C peut-être?), Mais cet exemple est en C:

#include <stdint.h>

typedef fixed16q16_t int32_t ;
#define FIXED16Q16_SCALE 1 << 16 ;

fixed16q16_t mult16q16( fixed16q16_t a, fixed16q16_t b )
{
    return (a * b) / FIXED16Q16_SCALE ;
}

fixed16q16_t div16q16( fixed16q16_t a, fixed16q16_t b )
{
    return (a * FIXED16Q16_SCALE) / b ;
}

Notez que ce qui précède est une implémentation simpliste et n’offre aucune protection contre les débordements arithmétiques. Par exemple, dans div16q16 (), je multiple avant la division pour conserver la précision, mais l'opération peut déborder en fonction des opérandes. Vous pouvez utiliser un intermédiaire 64 bits pour résoudre ce problème. De plus, la division est toujours arrondie car elle utilise une division entière. Cela donne les meilleures performances, mais peut affecter la précision des calculs itératifs. Les correctifs sont simples mais s’ajoutent à la surcharge.

Notez que lors de la multiplication ou de la division par une puissance de deux constante, la plupart des compilateurs détectent l'optimisation triviale et utilisent un décalage. Cependant, C ne définit pas le comportement pour un décalage à droite d'un entier signé négatif. Je laisse donc au compilateur le soin de le résoudre pour des raisons de sécurité et de portabilité. YMV sur la langue que vous utilisez.

Dans un langage OO, fixed16q16_t serait naturellement un candidat pour une classe avec surcharge d’opérateur afin que vous puissiez l’utiliser comme un type arithmétique normal.

Vous pouvez trouver utile de convertir entre les types:

double fixed16q16_to_double( fixed16q16_t fix )
{
    return (double)fix / FIXED16Q16_SCALE ;
}

int fixed16q16_to_int( fixed16q16_t fix )
{
    // Note this rounds to nearest rather than truncates
    return ((fix + FIXED16Q16_SCALE/2)) / FIXED16Q16_SCALE ;
}

fixed16q16_t int_to_fixed16q16( int i )
{
    return i * FIXED16Q16_SCALE ;
}

fixed16q16_t double_to_fixed16q16( double d )
{
    return (int)(d * FIXED16Q16_SCALE) ;
}

C’est la base, il est possible de devenir plus sophistiqué et d’ajouter trig et d’autres fonctions mathématiques.

L'addition et la soustraction fixes fonctionnent avec les opérateurs + et - intégrés et leurs variantes.