Evitar aritmética de coma flotante

Question

Escribí un pequeño sintetizador de software para el iPhone.
Para ajustar aún más el rendimiento, midí mi aplicación con Shark y descubrí que estoy perdiendo mucho tiempo en las conversiones float / SInt16.
Así que reescribí algunas partes para evitar las conversiones calculando previamente las tablas de búsqueda que devuelven "listo para usar". SInt16 muestras. Esto funciona bien hasta ahora.
Actualmente estoy tratando de reescribir algunos filtros y mi implementación de envolvente ADSR para usar solo aritmética de enteros, pero podría usar algunos consejos para realizar multiplicaciones / divisiones sin flotantes.
Estoy apuntando a formato canónico para iPhone :

• LPCM
• Muestras enteras de 16 bits

¿Cuáles son los buenos enfoques para aplicar una amplitud a mi muestra final sin usar un flotador?

Editar:
Lo único que he descubierto hasta ahora es que puedo dividir por potencias de 2 desplazando a la derecha mi muestra actual.

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

Pero no se me ocurre ninguna forma elegante de crear un sobre ADSR sin problemas con eso.

Editar2: ¡Gracias por todas sus excelentes respuestas!
Mi enfoque actual:

• traer todos mis valores de sobre ADSR en el rango positivo SInt16
• multiplicar con el valor actual de la tabla de ondas (almacenar intermedios como SInt32)
• cambia el resultado en 16 a la derecha

esto parece funcionar :)

Answer 1

El punto fijo es bueno, ya que en este caso estás usando 16 bits. La forma más simple es multiplicar por una potencia de 10 dependiendo de la precisión que necesite. Si puede usar entradas de 32 bits como intermediario, debería poder obtener una precisión decente. Al final, puede volver a convertir a 16 bits int, redondeando o truncando como prefiera.

Editar: Desea desplazarse hacia la izquierda para aumentar los valores. Almacene el resultado del cambio en un tipo con más precisión (32 o 64 bits dependiendo de lo que necesite). el cambio simple no funcionará si está utilizando tipos con signo

Tenga cuidado si está multiplicando o dividiendo dos números de punto fijo. Multiplicar termina siendo (a * n) * (b n) y terminarás con un b n ^ 2 en lugar de a b n. La división es (a n) / (b n) que es (a / b) en lugar de ((a n) / b). Es por eso que sugerí usar potencias de 10, hace que sea fácil encontrar sus errores si no está familiarizado con el punto fijo.

Cuando haya terminado sus cálculos, cambia de nuevo a la derecha para volver a un int de 16 bits. Si quieres ponerte elegante, también puedes redondear antes de cambiar.

Le sugiero que lea un poco si realmente está interesado en implementar un punto fijo eficiente. http://www.digitalsignallabs.com/fp.pdf

Answer 2

Las respuestas a esta pregunta SO son bastante completo en términos de implementación. Aquí hay una explicación un poco más de lo que vi allí:

Un enfoque es forzar todos sus números en un rango, digamos [-1.0,1.0). Luego mapeas esos números en el rango [-2 ^ 15, (2 ^ 15) -1]. Por ejemplo,

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

Cuando multiplicas estos dos números obtienes

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

Dividiendo por 32768 en la última línea está el punto Patros sobre multiplica y necesita un paso de escala adicional. Esto tiene más sentido si escribe la escala 2 ^ N explícitamente:

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

Entonces esa es la aritmética. Para la implementación, tenga en cuenta que la multiplicación de dos enteros de 16 bits necesita un resultado de 32 bits, por lo que Temp debe ser de 32 bits. Además, 32768 no es representable en una variable de 16 bits, así que tenga en cuenta que el compilador hará inmediatos de 32 bits. Y como ya ha notado, puede cambiar a multiplicar / dividir por potencias de 2 para poder escribir

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);

Pero supongamos que [-1,1) no es el rango correcto? Si prefiere limitar sus números a, digamos, [-4.0,4.0), puede usar N = 13. Luego tiene 1 bit de signo, dos bits antes del punto binario y 13 después. Estos se denominan 1.15 y 3.13 tipos fraccionales de punto fijo respectivamente. Cambia la precisión en la fracción por altura libre.

Sumar y restar tipos fraccionarios funciona bien siempre que tenga en cuenta la saturación. Para dividir, como dijo Patros, la escala en realidad se cancela. Entonces tienes que hacer

Quotient = (x1/x2) << N;

o, para preservar la precisión

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

Multiplicar y dividir entre números enteros funciona normalmente. Por ejemplo, para dividir entre 6 simplemente puede escribir

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

Y en el caso de dividir por una potencia de 2,

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

Sin embargo, sumar y restar números enteros no funciona ingenuamente. Primero debe ver si el número entero encaja en su tipo x.y, hacer el tipo fraccional equivalente y continuar.

Espero que esto ayude con la idea, mire el código en esa otra pregunta para implementaciones limpias.

Answer 3

Eche un vistazo a esta página que describe algoritmos de multiplicación rápida.

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

Answer 4

En general, digamos que usará una representación de punto fijo 16.16 con signo. Para que un entero de 32 bits tenga una parte entera de 16 bits con signo y una parte fraccional de 16 bits. Entonces no sé qué lenguaje se usa en el desarrollo de iPhone (¿quizás Objective-C?), Pero este ejemplo 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 ;
}

Tenga en cuenta que lo anterior es una implementación simplista y no proporciona protección contra el desbordamiento aritmético. Por ejemplo, en div16q16 (), multiplico antes de la división para mantener la precisión, pero dependiendo de los operandos, la operación puede desbordarse. Puede usar un intermedio de 64 bits para superar esto. Además, la división siempre se redondea hacia abajo porque usa la división de enteros. Esto proporciona el mejor rendimiento, pero puede afectar la precisión de los cálculos iterativos. Las soluciones son simples pero se agregan a los gastos generales.

Tenga en cuenta que al multiplicar o dividir por un poder constante de dos, la mayoría de los compiladores detectarán la optimización trivial y usarán un cambio. Sin embargo, C no define el comportamiento para un desplazamiento a la derecha de un entero con signo negativo, por lo que se lo he dejado al compilador para que lo resuelva por seguridad y portabilidad. YMV en cualquier idioma que esté utilizando.

En un lenguaje OO, fixed16q16_t sería, naturalmente, un candidato para una clase con sobrecarga de operadores, por lo que puede usarlo como un tipo aritmético normal.

Puede resultarle útil convertir entre tipos:

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) ;
}

Eso es lo básico, es posible ser más sofisticado y agregar trigonometría y otras funciones matemáticas.

La suma y resta fija funciona con los operadores + y - incorporados y sus variantes.