Comment rechercher des sinus de fréquences différentes dans une table de conversion de taille fixe?

Je suis en train d’échantillonner une onde sinusoïdale à 48 kHz. La gamme de fréquences de mon onde sinusoïdale peut varier de 0 à 20000 Hz avec un pas d’environ 100 Hz. J’utilise une approche de table de consultation. Je génère donc 4096 échantillons pour une onde sinusoïdale pour 4096 phases différentes. Je pense que l’idée générale derrière cela est d’augmenter la taille des pas et d’utiliser différentes tailles de pas pour différentes fréquences. Donc, je fais ce qui suit (pseudo code). Mais je ne suis pas sûr de savoir comment la taille de pas va être liée à la fréquence que je veux générer les échantillons de l’onde sinusoïdale de? Par exemple, si ma fréquence est de 15 000 Hz, quelle est la taille de pas que je dois traverser? La taille de mon échantillon (4096) est-elle trop basse pour cela?

// Pseudocode uint16_t audio_sample[4096] = {...}; NSTEP = freq; //???How is the step size going to be related to the freq here for(int i = 0; i < 4096; i = i+NSTEP) { sine_f(i) = audio_sample[i]; } 

Merci d’avance.

Vous êtes sur la bonne voie – nous devons d’abord générer une LUT à onde sinusoïdale:

 const int Fs = 48000; // sample rate (Hz) const int LUT_SIZE = 1024; // lookup table size int16_t LUT[LUT_SIZE]; // our sine wave LUT for (int i = 0; i < LUT_SIZE; ++i) { LUT[i] = (int16_t)roundf(SHRT_MAX * sinf(2.0f * M_PI * (float)i / LUT_SIZE)); } // fill LUT with 16 bit sine wave sample values 

Notez que nous n'avons besoin de générer cette LUT qu'une seule fois, par exemple lors de l'initialisation.

Maintenant que nous avons une table sinusoïdale, nous pouvons l'utiliser pour générer toute fréquence désirée avec un accumulateur de phase :

 const int BUFF_SIZE = 4096; // size of output buffer (samples) int16_t buff[BUFF_SIZE]; // output buffer const int f = 1000; // frequency we want to generate (Hz) const float delta_phi = (float) f / Fs * LUT_SIZE; // phase increment float phase = 0.0f; // phase accumulator // generate buffer of output for (int i = 0; i < BUFF_SIZE; ++i) { int phase_i = (int)phase; // get integer part of our phase buff[i] = LUT[phase_i]; // get sample value from LUT phase += delta_phi; // increment phase if (phase >= (float)LUT_SIZE) // handle wraparound phase -= (float)LUT_SIZE; } 

Remarque: pour une sortie de meilleure qualité, vous pouvez utiliser une interpolation linéaire entre les valeurs de LUT à phase_i et phase_i + 1 , mais l'approche ci-dessus est suffisante pour la plupart des applications audio.

Avec l’approche de la table de consultation, on peut souhaiter utiliser efficacement la mémoire et ne stocker que le premier quadrant de l’onde sinusoïdale.

 Then second quadrant = sin[180 - angle] ; // for angles 90..180 degrees third quadrant = -sin[angle-180]; // for 180..270 fourth quadrant = -sin[-angle+360] // for 270..360 

Je recommanderais une interpolation linéaire, mais il existe également une approche basée sur la rotation vectorielle pour la génération de sinus (produisant simultanément le sinus et le cosinus)

 x_0 = 1, y_0 = 0; x_n+1 = x_n * cos(alpha) - y_n * sin(alpha) y_n+1 = x_n * sin(alpha) + y_n * cos(alpha), 

où alpha = différence de phase de la fréquence == 2pi * fHz / Ts, avec fHz la fréquence à produire (en hertz) et Ts le temps d’échantillonnage (ou 1 / Ts = fréquence d’échantillonnage, par exemple 44100 Hz).

ce qui conduit à une approche de filtre de retour résonant , dont la fonction de transfert f (z) a deux pôles conjugués au cercle unitaire (z = e ^ jomegaT).

 y[n] = y[n-1] - 2*cos(alpha) * y[n-2], with y[0] = 0, y[1] = sin(alpha) 

La partie amusante est que l’on peut modifier l’alpha (cos (alpha)) à la volée. L’inconvénient de cette approche de filtre IIR est qu’il est instable par définition. Les imprécisions de virgule flottante et surtout de virgule fixe s’accumulent et mènent soit à une décroissance exponentielle, soit à une croissance exponentielle de la magnitude. Cela peut toutefois être guéri en permettant une légère distorsion de phase.

Au lieu de comme dans CORDIC, rotation ayant un facteur d’amplification connu par itération:

 K = sqrt(1+sin(alpha)^2); x' = x - y*sin(alpha); y' = y + x*sin(alpha); one can do x' = x - y * sin(alpha); y'' = y + x' * sin(alpha); 

qui ne produit pas de cercles parfaits pour (x ‘, y’ ‘) mais des ellipses stables même avec une arithmétique en virgule fixe. (Notez que cela suppose des valeurs alpha relativement faibles, ce qui signifie également des fréquences relativement basses.)

Si vous voulez une haute précision, vous pouvez utiliser les identités sortinggonomésortingques pour avoir à la fois une petite LUT et des ondes sinusoïdales nettes.

 static float gTrigSinHi[256], gTrigSinLo[256]; static float gTrigCosHi[256], gTrigCosLo[256]; //////////////////////////////////////// // Sets up the fast sortingg functions void FastTrigInit() { unsigned i; for(i = 0; i < 256; ++i) { gTrigSinHi[i] = sin(2.0 * M_PI / 0xFFFF * (i << 8)); gTrigSinLo[i] = sin(2.0 * M_PI / 0xFFFF * (i << 0)); gTrigCosHi[i] = cos(2.0 * M_PI / 0xFFFF * (i << 8)); gTrigCosLo[i] = cos(2.0 * M_PI / 0xFFFF * (i << 0)); } } //////////////////////////////////////// // Implements sin as // sin(u+v) = sin(u)*cos(v) + cos(u)*sin(v) float FastSin(unsigned short val) { unsigned char hi = (val >> 8); unsigned char lo = (val & 0xFF); return gTrigSinHi[hi] * gTrigCosLo[lo] + gTrigCosHi[hi] * gTrigSinLo[lo]; } //////////////////////////////////////// // Implements cos as // cos(u+v) = cos(u)*cos(v) - sin(u)*sin(v) float FastCos(unsigned short val) { unsigned char hi = (val >> 8); unsigned char lo = (val & 0xFF); return gTrigCosHi[hi] * gTrigCosLo[lo] - gTrigSinHi[hi] * gTrigSinLo[lo]; } 

Très bonne réponse, il s’agit d’un logiciel classique DDS. Face au même problème ces jours-ci. Il n’y a pas besoin d’utiliser des flotteurs

  UInt16 f = 400; // Hz, up to 16384 :) UInt16 delta = (UInt16)(((UInt32)f * LUT_SIZE ) / fmt_SamplesPerSec); UInt16 phase = 0; for (int i = 0; i < BUFF_SIZE; ++i) { buff[i] = LUT[phase]; phase += delta; phase &= LUT_SIZE-1; } 

Laisser la phase enveloppante LUT taille comme masque. Et ne vous souciez pas de l’utilisation des quadrants car j’ai déjà un énorme MIPS pour répondre à ces besoins.

Selon “http://en.wikipedia.org/wiki/X86_instruction_listings” si vous avez x80387 ou une version plus récente, il existe une instruction sine, appelez-la directement. Vous devez juste comprendre comment append un langage d’assemblage en ligne à votre programme. De cette façon, vous n’avez pas à vous inquiéter si votre valeur d’entrée ne correspond pas exactement à ce qui se trouve dans votre table.