MCUdude / SigmaDSP

A versatile Arduino library for interfacing with the ADAU1401, ADAU1701 and ADAU1702 audio DSPs
GNU Lesser General Public License v3.0
171 stars 36 forks source link

Convert some of the functions from floating point to integer math #4

Closed MCUdude closed 4 years ago

MCUdude commented 5 years ago

Since most low-end microcontrollers don't have a floating point unit built-in, working with floats is both slow and space hungry. With some clever programming, many of the functions provided in this library should be able to be converted into integer based functions.

@Jackjan4 you did a wonderful job with the millis() function in MCUdude_corefiles. Maybe you could help me out with a few examples? 🙂

Let's take a very simple function, a sine wave generator. At the moment it looks like this. Don't worry about the 0xff and 1.0 parameters. These are just needed parameters.

void SigmaDSP::sineSource(uint16_t startMemoryAddress, float frequency)
{
  float value = (1.00/24000.00)*frequency;

  safeload_write(startMemoryAddress, 0xff, value, 1.0);
}

The thing is that this DSP doesn't really need floats at all! before everything is sent to the DSP the data (value in this case) is multiplied by a bitshift, converted to an int32 and chopped up in 4 bytes.

The function above can be rewritten into this:

void SigmaDSP::sineSource(uint16_t startMemoryAddress, float frequency)
{
  float value = (1.00/24000.00)*frequency;
  int32_t converted_value = floatToInt(value);
  safeload_write(startMemoryAddress, 0xff, converted_value, 1.0);
}

int32_t SigmaDSP::floatToInt(float value)
{
  // Convert float 5.23 to int 28.0
  return (value * ((int32_t)1 << 23));
}

And, the 1.0 parameter could actually be replaced by the number 0x800000. In case you're wondering what kind of function safeload_write is, it is simply a template based function that calls one of these functions based on the parameters:

    void safeload_writeRegister(uint16_t memoryAddress, int32_t data, bool finished);
    void safeload_writeRegister(uint16_t memoryAddress, int16_t data, bool finished);
    void safeload_writeRegister(uint16_t memoryAddress, uint8_t data, bool finished);
    void safeload_writeRegister(uint16_t memoryAddress, float data, bool finished);

So the big question is: how can the function above be replaced with integer based math? Can the same bitshift magic be used here as well @JackJan4?

EDIT: If you don't know what a Sigma DSP is it is an audio DSP family from Analog Devices. The one this library is targeted at, the ADAU1701, is small (48-pin), cheap, versatile and easy to program through Analog Devices' SigmaStudio software. However they never really provided a good way of controlling this DSP using an external microcontroller. No API, only some exported arrays in a c file. IMO this is a terrible solution, so this is why this library exists.

Jackjan4 commented 5 years ago

Hi,

I looked through your case here and if I understand it right we need to address two multiplications, right?

This one return (value * ((int32_t)1 << 23));

and this one (1.00/24000.00)*frequency;

right?

EDIT: Oh, and I have another question. What type of numbers will "frequency" be? Also just integers?

MCUdude commented 5 years ago

Yes, you're correct. Or, if we do some simplifications:

1 << 23 = 8388608 (or 0x800000 in hex)

int32_val = frequency * (8388608/24000)

int32_val = frequency * 349.523333...

frequency will be an integer number, between 1 and 20000

Jackjan4 commented 5 years ago

Ah okay, I already got something here for you: Since you can only use bitshifting on integers you need to rewrite "frequency" to an integer datatype first.

After that you can use this for a multiplication. int32_val = (freq << 9) - (freq << 7) - (freq << 5) - (freq << 1) - freq + (freq >> 2) + (freq >> 5) This line will multiply by 349.5315 so it is only an approximation (since bitshifting can only approximate). If you need less accuracy you can cut the end off. If you need more, the bitshifting will be less efficient but theoretically it can of course be done.

MCUdude commented 5 years ago

Wow, amazing!

The next one is very similar. instead of 1/24000, it's 0.5/24000 (or 1/48000). Would you like to provide an integer variant of this one as well?

void SigmaDSP::sawtoothSource(uint16_t startMemoryAddress, float frequency)
{
  float value = (0.50/24000.00)*frequency;
  int32_t converted_value = floatToInt(value);
  safeload_write(startMemoryAddress, converted_value, 0x800000);
}

If you'll want a real puzzle you could have a look at the volume_slew function. In the SigmaStudio program, this function is represented as a simple volume slider. Don't worry about the slew parameter. The dB value ranges from -80 to +4, where -80 is the lowest possible volume, and +4 is the highest. In this case, an 8-bit signed integer value is a valid alternative for a float.

void SigmaDSP::volume_slew(uint16_t startMemoryAddress, float dB, uint8_t slew)
{
  float volume = pow(10, dB / 20); // 10^(dB / 20)
  int32_t converted_volume = floatToInt(volume);
  int32_t slewrate = 0x400000 / (1 << (slew - 1)); // 0x400000/2^(slew - 1))

  safeload_write(startMemoryAddress, converted_volume, slewrate);
}
Jackjan4 commented 5 years ago

Okay the second one will be the following: int32_val = (freq << 8) - (freq << 6) - (freq << 4) - freq - (freq>>1) + (freq >> 2) This one will multiply by 174.75 (the exact value should be 174.762666...), is it accurate enough for your purpose?

For the third one I will have to look into it tomorrow :)

MCUdude commented 5 years ago

thank you! It is definitely accurate enough!

MCUdude commented 5 years ago

Have you had the time to look at the third formula? And would you like to provide a simple explanation on how you do calculations based on bit shifts? Thanks!

MCUdude commented 5 years ago

I don't have any actual hardware to test the new multiplication on, so I wrote a simple C program to run on my computer:

The 174.75 multiplication is basically dead on, but the 349.5315 is a tiny bit off (about 0.25):

#include <stdio.h>
#include <stdint.h>

// Different "frequencies" to test
uint32_t mult_array[] = {5, 10, 20, 30, 50, 100, 200, 234, 500, 514, 1000, 1250, 2000, 5000, 10000, 12000, 15000, 18000, 20000};

int32_t multiply_349(int32_t frequency)
{
  // Multiply frequency by 349.5315 (ideally 349.525333333)
  int32_t value = (frequency << 9) - (frequency << 7) - (frequency << 5) - (frequency << 1) - frequency + (frequency >> 2) + (frequency >> 5);

  return value;
}

int32_t multiply_174(int32_t frequency)
{
  // Multiply frequency by 174.75 (ideally 174.762666667)
  int32_t value = (frequency << 8) - (frequency << 6) - (frequency << 4) - frequency - (frequency >> 1) + (frequency >> 2);

  return value;
}

int main(void)    
{

  printf("349.5315 multiplier\nFreq: \tOutput: \tFactor:\n");
  for(uint8_t i = 0; i < (sizeof(mult_array)/sizeof(uint32_t)); i++)
    printf("%d \t%d \t%.6f\n", mult_array[i], multiply_349(mult_array[i]), ((float)multiply_349(mult_array[i])/(float)mult_array[i]));

  printf("\n\n174.75 multiplier\nFreq: \tOutput: \tFactor:\n");  
  for(uint8_t i = 0; i < (sizeof(mult_array)/sizeof(uint32_t)); i++)
    printf("%d \t%d \t%.6f\n", mult_array[i], multiply_174(mult_array[i]), ((float)multiply_174(mult_array[i])/(float)mult_array[i]));  

  return 0;
}

Output:

349.5315 multiplier
Freq:   Output:     Factor:
5   1746    349.200012
10  3492    349.200012
20  6985    349.250000
30  10477   349.233337
50  17463   349.260010
100     34928   349.279999
200     69856   349.279999
234     81731   349.277771
500     174640  349.279999
514     179530  349.280151
1000    349281  349.281006
1250    436601  349.280792
2000    698562  349.281006
5000    1746406     349.281189
10000   3492812     349.281189
12000   4191375     349.281250
15000   5239218     349.281189
18000   6287062     349.281219
20000   6985625     349.281250

174.75 multiplier
Freq:   Output:     Factor:
5   874     174.800003
10  1747    174.699997
20  3495    174.750000
30  5242    174.733337
50  8737    174.740005
100     17475   174.750000
200     34950   174.750000
234     40891   174.747864
500     87375   174.750000
514     89821   174.749023
1000    174750  174.750000
1250    218437  174.749603
2000    349500  174.750000
5000    873750  174.750000
10000   1747500     174.750000
12000   2097000     174.750000
15000   2621250     174.750000
18000   3145500     174.750000
20000   3495000     174.750000
Jackjan4 commented 5 years ago

Yes you are right. I miscalculated my 349.5315 multiplication. It should be the following: int32_val = (freq << 9) - (freq << 7) - (freq << 5) - (freq << 1) - (freq >> 1) + (frequency >> 5)

Here is a little explanation on how this "magic" actually works:

Since we are in binary, bitshifting only means dividing or multiplying by 2: So x << 1 will multiply x with 2 and x >> 2 will divide x by 2.

x << 2 will multiply x with 4 x << 3 will multply x with 8 ...and so on, just the powers of two.

So with bitshifting we can precisly multiply and divide with the powers of 2 but nearly all other numbers we have to approximate. So how do we do this?

Lets take a simple example first - we want to multiply x with 7.5: So what is the nearest precise multiplication we can do with power of two? Of course we can multiply x with 8, so we do: x << 3 // Multiply x by 8

But now we multiplied too much. We need to get rid of 0.5x. But actually thats quite easy. We just subtract x >> 1 which is 0.5x or x/2. So in the end we get: x = x << 3 - x >> 1 // Multiply x by 7.5

This time our real world example - we want to multiply x by 349.525333... . Since this is an irregular number we can of course only approximate it. But lets look how precise we can go:

First, whats the best power of two that comes near to our number? Well, its 512, or in different words (x << 9). Okay, but now we are too large, so we subtract again as good as we can. We can subtract 128 which is (x << 7). So now we got (x << 9) - (x << 7) or in different words we got 384. We can subtract by 32 (which is x << 5) and we get 352. Further we subtract 2 and finally we got the following: (x << 9) - (x << 7) - (x << 5) - (x << 1) which results in 350 which is quite good at the moment. To come even closer we now use rightshifts which help us to get fractions. If we subtract by x >> 1 now we get 349.5 because we subtracted 0.5. (like in the example with 7.5)

And the last number thats left to get 349.5315 I will leave to you. ;)

So you see, bitshifting is playing with math and to explore how you can approximate a number by multiplication/division of 2. I hope this little explanation helps to understand the technique behind it.

MCUdude commented 5 years ago

Yes you are right. I miscalculated my 349.5315 multiplication. It should be the following: int32_val = (freq << 9) - (freq << 7) - (freq << 5) - (freq << 1) - (freq >> 1) + (freq >> 5)

Thanks! It's pretty much dead on now 🙂 But we're still theoretically multiplying by 349.5315, right? I'm adding some comments to the main code you see.

Freq:   Output:     Factor:
5   1748    349.600006
10  3495    349.500000
20  6990    349.500000
30  10485   349.500000
50  17476   349.519989
100     34953   349.529999
200     69906   349.529999
234     81790   349.529907
500     174765  349.529999
514     179659  349.531128
1000    349531  349.531006
1250    436914  349.531189
2000    699062  349.531006
5000    1747656     349.531189
10000   3495312     349.531189
12000   4194375     349.531250
15000   5242968     349.531189
18000   6291562     349.531219
20000   6990625     349.531250

So you see, bitshifting is playing with math and to explore how you can approximate a number by multiplication/division of 2. I hope this little explanation helps to understand the technique behind it.

Wow, this is pure gold! I've been googling around to see if I can find a simple guide for this. This is just clear and simple; exactly what I've been looking for. However, I need to sit down with my trusty old Texas calculator and do some number crunching to get the proper feel for it.

But this got me thinking, this is "just" a simple multiplication. How is the volume (dB) calculating formula handled? This involves a power of 10 and a division of the exponent. Is it even possible?

Jackjan4 commented 5 years ago

But we're still theoretically multiplying by 349.5315, right?

Yes, do you need it more precise, or is this okay?

But this got me thinking, this is "just" a simple multiplication. How is the volume (dB) calculating formula handled? This involves a power of 10 and a division of the exponent. Is it even possible?

It is possible, but this involves multiplying in a for-loop. I don't know if this approach results in less code and size in the end, but you could optimize some stuff, for example "slew / 20".

MCUdude commented 5 years ago

Yes, do you need it more precise, or is this okay?

It's OK. I just asked because I wanted to get the comment right. I'll push the new code as soon as I've tested on real hardware.

It is possible, but this involves multiplying in a for-loop. I don't know if this approach results in less code and size in the end, but you could optimize some stuff, for example, "slew / 20".

If it involves loops then we're just solving the equation in a numerical way. I don't think it's worth it really. How would you improve the slewrate formula? get rid of the division?

Jackjan4 commented 5 years ago

You could replace float volume = pow(10, dB / 20); with float volume = pow(10, dB >> 4 - dB >> 7 - dB >> 8);

1/20 equals 0.05 and this bitshift equals multiplying with 0.05078... This should maybe be a bit quicker and maybe also save a couple bytes. However this would be the only thing that I would change.

MaxPayne86 commented 5 years ago

Hello :)

in original AidaDSP.cpp I've decided to keep the user level simple and make use of floating point. 1) The code using Aida DSP apis is not supposed to run at audio sampling rate, but rather at much slower speed (~60Hz?) so converting float math to integer math is a waste of time imho 2) The slow part is the cast between float to integer, not the bit-shift 3) When I developed the original Aida DSP library (5 years ago) the overall trend was to switch from 8-bit microcontrollers to 32-bit, which notably have much higher cpu clock speeds, and some of which can do a little audio processing on their own (Axoloti?) so for sure can handle calculation done in low speed to manage a user hmi (i.e. track down a pot value)

Just a few thoughts from me

MCUdude commented 4 years ago

I've decided I'll put this dead for now. Most users probably are fine with fairly slow calculations. If not, why no just use a 32-bit microcontroller instead?