FM synthesis using a phase battery

Question

I am trying to implement a FM synthesis operator with feedback using a phase accumulator in C. In Tomisawa, an original patent , the phase accumulator getting into the adder for both negative and positive indices, from -2 ^ (n-1} in the phase of the sine wave from -pi to 2 ^ (n-1) in the pi phase. For simplicity, for example, I would use a phase accumulator that takes into account only positive values, using the upper bytes of an unused 32-bit integer as an index in the search for the sine table.

I experimented with this, and unfortunately, I cannot get the algorithm to generate the expected results when using feedback. Adding the sine wave output to the phase battery should result in a sawtooth waveform, but I cannot figure out how to properly add the sine wave output (which is a 16-bit int signature) for an unsigned phase battery to create it. Any suggestions would be appreciated.

Edit:

Some refinements, probably in order. Here are some diagrams from the original Tomisawa:

When both a phase battery and a sinusoidal signal are signed, the algorithm is simple enough to implement. The phase battery starts at -1 and runs up to 1, and the sine wave output is also between -1 and 1. In Python, the algorithm looks something like this to generate 1000 samples:

table = [] feedback = 0.25 accumulator = -1 for i in xrange(1000): output = math.sin(math.pi*(accumulator + feedback*output) table[i] = output accumulator += 0.005 if accumulator > 1: accumulator = -1 

Which creates an output that looks like this:

I am trying to adapt this algorithm to C. In C for computational efficiency, I would like the phase battery to be a 32-bit unsigned integer rather than a signed integer. Thus, I can use the first two bits of the high byte of the battery as a quadrant index, and the second high byte as an index in an array of 256 16-bit sine values ​​for a sine table of size 1024. How:

 XXXXXXQQ.IIIIIIII.XXXXXXXX.XXXXXXXX ^^ ^^^^^^^^ quadrant index 

My problem is that it is difficult for me to adapt the FM algorithm as indicated for an unsigned phase battery. If the phase storage device is an unsigned 32-bit int, and the output of the sinusoidal table is a 16-bit integer (signed or unsigned), then how can I adapt the algorithm, as shown in the patent, and the above Python code to work with this format, and make the same exit?

Answer 1

First of all, we can try to write you pironic code in C

 #include <stdio.h> #include <math.h> void main() { double table[1000]; double feedback = 0.25; double accumulator = -1; int i; for (i=0;i<1000;i++) { double output = sin(M_PI*(accumulator + feedback*output)); table[i]=output; accumulator += 0.005; if (accumulator > 1) accumulator = -1; printf("%f\n",output); } } 

The next step is to use the calculated values ​​for sin

 #include <stdio.h> #include <math.h> void main() { double table[1000]; double feedback = 0.25; double accumulator = 1; int i; double sinvalue[1024]; for (i=0;i<1024;i++) { sinvalue[i]=sin(M_PI*i/512); } for (i=0;i<1000;i++) { double output = sinvalue[(int)(512*(accumulator + feedback*output))%1024]; printf("%0.6f %0.6f %0.6f\t",accumulator,feedback,output); table[i]=output; accumulator += 0.005; if (accumulator > 2) accumulator = 0; printf("%f\n",output); } } 

The next step is to use the 16-bit values ​​of sin and output. In this version, the value in "output", for example XXXXXXQQ.IIIIIIII.XXXXXXXX.XXXXXXXX In addition, we have lost some accuracy.

 #include <stdio.h> #include <math.h> #define ONE ((int)(2*256*256*256/M_PI)) void main() { double table[1000]; double feedback = 0.25; double accumulator = 1; double accumulatorDelta = 0.005; unsigned int feedback_i = ONE*feedback/32768; unsigned int accumulator_i = ONE*accumulator; unsigned int accumulatorDelta_i = ONE*accumulatorDelta; int i; double sinvalue[1025]; short int sinvalue_i[1025]; for (i=0;i<1025;i++) { sinvalue[i]=sin(M_PI*i/512); sinvalue_i[i]=32786*sinvalue[i]; if (sinvalue[i]*32768>32768) sinvalue_i[i]=32768; if (sinvalue[i]*32768<-32767) sinvalue_i[i]=-32767; } for (i=0;i<1000;i++) { double output = sin(M_PI*(accumulator + feedback*output)); short int output_i = sinvalue_i[ ((unsigned int) ((accumulator_i + feedback_i*output_i)*M_PI)>>16)%1024 ]; table[i]=output; accumulator += 0.005; if (accumulator > 2) accumulator = 0; accumulator_i += accumulatorDelta_i; if (accumulator_i > 2*ONE) accumulator_i = 0; printf("%f %f %04X\n",output,(float)output_i/32768,(unsigned short int)output_i); } } 

But we lost some time for the conversion int-> double-> int If we change the ONE constant, we will lose the ability to quickly get the quadrant, but get rid of the conversion

 #include <stdio.h> #include <math.h> #define ONE ((int)(2*256*256*256)) void main() { short int table[1000]; unsigned int feedback_i = ONE*0.25/32768; unsigned int accumulator_i = ONE*1; unsigned int accumulatorDelta_i = ONE*0.005; int i; short int sinvalue_i[1025]; for (i=0;i<1025;i++) { double sinvalue=sin(M_PI*i/512); sinvalue_i[i]=32786*sinvalue; if (sinvalue*32768>32768) sinvalue_i[i]=32768; if (sinvalue*32768<-32767) sinvalue_i[i]=-32767; } for (i=0;i<1000;i++) { short int output_i = sinvalue_i[ ( (accumulator_i + feedback_i*output_i)>>16)%1024 ]; table[i]=output_i; accumulator_i += accumulatorDelta_i; if (accumulator_i > 2*ONE) accumulator_i = 0; printf("%f %04X\n",(float)output_i/32768,(unsigned short int)output_i); } }