home *** CD-ROM | disk | FTP | other *** search
/ Magazyn Amiga 14 / MA_Cover_14.iso / source / c / pegase_src / fastfftc.cpp < prev    next >
Encoding:
C/C++ Source or Header  |  1999-06-14  |  21.3 KB  |  793 lines
  1. /*
  2. **
  3. ** FastFFT.cpp
  4. **
  5. ** $VER: FastFFT.cpp 1.0 (11.07.98)
  6. **
  7. ** Calcul de la transformée de Fourier rapide
  8. **
  9. ** $Revision: 1.4 $
  10. ** $State: Exp $
  11. ** $Date: 1998/08/16 19:03:44 $
  12. **
  13. ** $Log: FastFFTC.cpp $
  14. ** Revision 1.4  1998/08/16 19:03:44  kakace
  15. ** Version Beta3+
  16. **
  17. ** Revision 1.3  1998/08/09 19:21:36  kakace
  18. ** Correction d'un bug dans la routine FFT (avec tableà
  19. **
  20. ** Revision 1.2  1998/08/02 15:17:10  kakace
  21. ** Fonctionnement OK (Layer 1/2 + Stereo/JStereo)
  22. **
  23. */
  24.  
  25.  
  26. #ifndef  _INCLUDE_MATH_H
  27. #include <math.h>
  28. #endif
  29.  
  30. #ifndef  _MEMORY_HPP
  31. #include "Memory.hpp"
  32. #endif
  33.  
  34. #ifndef  _FASTFFT_CLASS_HPP
  35. #include "FastFFTC.hpp"
  36. #endif
  37.  
  38. #ifndef  _PEGASECOND_HPP
  39. #include "PegaseCond.hpp"
  40. #endif
  41.  
  42.  
  43. //----------------------------------------------------------------------------------------------------
  44.  
  45. // Fonctions "inline"
  46.  
  47. /// FastFFTC::GetShuffleIndex_?()
  48. /****** Class FastFFTC/GetShuffleIndex_10 ***********************************
  49. *
  50. *   NAME
  51. *   FastFFTC::GetShuffleIndex_10 -- Obtenir un index sur 10 bits.
  52. *   FastFFTC::GetShuffleIndex_8  -- Obtenir un index sur 8 bits.
  53. *
  54. *   SYNOPSIS
  55. *   Index = FastFFTC::GetShuffleIndex_10(input)
  56. *   Index = FastFFTC::GetShuffleIndex_8(input)
  57. *
  58. *   inline ULONG FastFFTC::GetShuffleIndex_10(ULONG);
  59. *   inline ULONG FastFFTC::GetShuffleIndex_8(ULONG);
  60. *
  61. *   FUNCTION
  62. *   Retourne la valeur transmise en argument après avoir inversé le champs
  63. *   binaire (bitreverse). La longueur du champs de bits doit correspondre à
  64. *   longueur indiquée dans le nom de la fonction.
  65. *
  66. *   INPUTS
  67. *   input - Valeur à traiter.
  68. *
  69. *   RESULT
  70. *   Index - Index obtenu après inversion du masque binaire.
  71. *
  72. *****************************************************************************
  73. *
  74. */
  75.  
  76.  
  77. static UBYTE ShuffleArray[] = {0, 8, 4, 12, 2, 10, 6, 14, 1, 9, 5, 13, 3, 11, 7, 15};
  78.  
  79. inline ULONG FastFFTC::GetShuffleIndex_10(ULONG input)
  80. {
  81.     ULONG k;
  82.     ULONG bitrev = 0;
  83.  
  84.     for (k = 0; k < 3; k++)
  85.     {
  86.         bitrev = (bitrev * 16) | (::ShuffleArray[input & 15]);
  87.         input >>= 4;
  88.     }
  89.  
  90.     return (bitrev >> 2);
  91. }
  92.  
  93.  
  94. inline ULONG FastFFTC::GetShuffleIndex_8(ULONG input)
  95. {
  96.     return ( (::ShuffleArray[input & 15] * 16) | (::ShuffleArray[input >> 4]) );
  97. }
  98. ///
  99.  
  100. // Constructeur et destructeur.
  101.  
  102. /// FastFFTC::FastFFTC
  103. FastFFTC::FastFFTC()
  104. {
  105.     int i;
  106.  
  107. #if _USE_SINCOS_TABLE_ == 1
  108.  
  109.     // Initialiser la table SinCos (68020/68030)
  110.  
  111.     REAL cos_alpha;
  112.     REAL f = - (PI / 512.0);
  113.  
  114.     for (i = 0; i < SINCOSSIZE / 3; i++)
  115.     {
  116.         cos_alpha = cos(f * i);
  117.         FFT_dSinCos[i]                      =   cos_alpha;
  118.         FFT_dSinCos[SINCOSSIZE * 2 / 3 - i] = - cos_alpha;
  119.         FFT_dSinCos[SINCOSSIZE * 2 / 3 + i] = - cos_alpha;
  120.     }
  121.  
  122.     FFT_dSinCos[SINCOSSIZE / 3] = 0;
  123.  
  124. #else
  125.  
  126.     // Initialiser la table SinCos réduite (68040/68060)
  127.  
  128.     REAL f = - (PI / 512.0);
  129.  
  130.     for (i = 9; i >= 0; i--)
  131.     {
  132.         FFT_dSinCos[0][i] = cos(f);
  133.         FFT_dSinCos[1][i] = sin(f);
  134.         f += f;
  135.     }
  136.  
  137. #endif
  138.  
  139.     FFT_MemoryError = FALSE;
  140.  
  141.     FFT_MemoryError |= !(FFT_pdWinSmplShuffled = (T_FFTBUFF *) AllocAlignedBuffer(sizeof(T_FFTBUFF)));
  142.     FFT_MemoryError |= !(FFT_pdHannWindowShuffled = (T_HANN *) AllocAlignedBuffer(sizeof(T_HANN)));
  143.  
  144.     // Initialiser la table HannWindow.
  145.     // Cette table est réordonnée pour faciliter l'initialisation des tampons d'entrée
  146.     // des FFTs.
  147.  
  148.     {
  149.         REAL cos_alpha, sin_alpha;
  150.         REAL cos_delta, sin_delta;
  151.         REAL f;
  152.  
  153.         f = - (PI / 1024.0);
  154.         cos_alpha = cos(f);
  155.         sin_alpha = sin(f);
  156.  
  157.         f = cos_alpha + cos_alpha;                  // = 2.cos(a)
  158.         cos_delta = f * cos_alpha - 1.0;            // cos(2a) = 2.cos²(a) - 1
  159.         sin_delta = - f * sin_alpha;                // sin(2a) = 2.sin(a).cos(a)
  160.  
  161.         for (i = 0; i < BLKSIZE; i++)
  162.         {
  163.             f = 0.5 * (1.0 - cos_alpha);
  164.             (*FFT_pdHannWindowShuffled)[GetShuffleIndex_10(i)] = f;
  165.  
  166.             f = cos_alpha;
  167.             cos_alpha = f * cos_delta - sin_alpha * sin_delta;
  168.             sin_alpha = f * sin_delta + sin_alpha * cos_delta;
  169.         }
  170.     }
  171. }
  172. ///
  173. /// FastFFTC::~FastFFTC()
  174. FastFFTC::~FastFFTC()
  175. {
  176.     if (FFT_pdWinSmplShuffled)  FreeAlignedBuffer(FFT_pdWinSmplShuffled);
  177.     if (FFT_pdHannWindowShuffled) FreeAlignedBuffer(FFT_pdHannWindowShuffled);
  178. }
  179. ///
  180.  
  181. // Fonctions privées.
  182.  
  183. /// FastFFTC::FastFFT()
  184. /****** Class FastFFTC/FastFFT **********************************************
  185. *
  186. *   NAME
  187. *   FastFFTC::FastFFT -- Calculer la transformée de Fourier.
  188. *
  189. *   SYNOPSIS
  190. *   FastFFTC::FastFFT(pdFFTIn, iFFTSize)
  191. *
  192. *   void FastFFTC::FastFFT(REAL *, UWORD);
  193. *
  194. *   FUNCTION
  195. *   Calcule la transformée de Fourier sur un ensemble d'échantillons.
  196. *
  197. *   INPUTS
  198. *   pdFFTIn  - Tableau contenant les échantillons (série de nombres complexes
  199. *       dans l'ordre : [partie réelle][partie imaginaire].
  200. *   iFFTSize - Nombre d'échantillons à traiter.
  201. *
  202. *   NOTES
  203. *   Cette fonction a été fortement optimisée en tirant partie de quelques
  204. *   angles particuliers. Elle est capable de calculer une FFT sur un ensemble
  205. *   de 1024 échantillons au maximum.
  206. *
  207. *   BUGS
  208. *   Fonction contrôlée conforme le 20.06.98
  209. *
  210. *   SEE ALSO
  211. *   FastFFTC::ComputeFastFFT_Stereo()
  212. *
  213. *****************************************************************************
  214. *
  215. *   Benchmark (3148 cycles par seconde (global)) :
  216. *       21.06.98 : 1'45.53
  217. *
  218. *   06.07.98 :
  219. *       - Utilisation d'une table SinCos[] en lieu et place des calculs
  220. *         trigonométriques (si _USE_SINCOS_TABLE_ == 1).
  221. */
  222.  
  223. #if _USE_ASM_FCT_ == 0
  224.  
  225. void FastFFTC::FastFFT(T_FFTDATAS *pdFFTIn, UWORD iFFTSize)
  226. {
  227. #define MAX_FFT_BLOCK 256
  228.  
  229.     ULONG i, j, delta;
  230.     LONG  k;
  231.     ULONG BlockSize, BlockEnd, n;
  232.     ULONG iBandSize = (iFFTSize < MAX_FFT_BLOCK) ? iFFTSize : MAX_FFT_BLOCK;
  233.  
  234.     double tr, ti;
  235.     T_FFTDATAS *pd_J, *pd_K, *pd_J2, *pd_K2;
  236.  
  237. #if _USE_SINCOS_TABLE_ == 1
  238.     ULONG step;
  239.     float *Cosin = &FFT_dSinCos[0];
  240.     float *Sinus = &FFT_dSinCos[SINCOSSIZE / 3];
  241. #endif
  242.  
  243.     // Appliquer la FFT sur le tampon d'entrée.
  244.     // Les passes 0 à 7 sont calculées ici, par colonnes de "iBandSize" échantillons
  245.     // (pour tirer un meilleur parti des caches).
  246.     // Les passes 8 et 9 sont quant à elles calculées à part.
  247.  
  248.     // On commence par la fin pour profiter du fait que le tampon d'entrée vient d'être
  249.     // initialisé, et donc que les derniers échantillons sont toujours dans le cache.
  250.  
  251.     for (k = iFFTSize - iBandSize; k >= 0; k -= iBandSize)
  252.     {
  253.         BlockEnd = 1;
  254.  
  255. #if _USE_SINCOS_TABLE_ == 0
  256.         delta = 0;
  257. #else
  258.         step = 512;
  259. #endif
  260.  
  261.         for (BlockSize = 2; BlockSize <= iBandSize; BlockSize += BlockSize)
  262.         {
  263. #if _USE_SINCOS_TABLE_ == 1
  264.             delta = 0;
  265. #endif
  266.  
  267.             n = BlockSize + BlockSize;
  268.  
  269.             pd_J = pd_J2 = &pdFFTIn[k + k];
  270.             pd_K = pd_K2 = &pd_J2[BlockSize];
  271.  
  272.             // Première partie : ar = 1, ai = 0 (commune à toutes les passes)
  273.  
  274.             for (i = 0; i < iBandSize; i += BlockSize)      // Boucle exécutée 1020 fois
  275.             {
  276.                 tr = pd_K[0]; ti = pd_J[0];
  277.                 pd_K[0] = ti - tr;
  278.                 pd_J[0] = ti + tr;
  279.  
  280.                 tr = pd_J[1]; ti = pd_K[1];
  281.                 pd_K[1] = tr - ti;
  282.                 pd_J[1] = tr + ti;
  283.  
  284.                 pd_K = &pd_K[n];                            // <==> pd_J += n
  285.                 pd_J = &pd_J[n];                            // <==> pd_J += n
  286.             }
  287.             pd_J2 += 2;
  288.             pd_K2 += 2;
  289.  
  290.             // Deuxième partie : ar et ai quelconques.
  291.  
  292.             if (BlockEnd > 4)
  293.             {
  294.                 double ar, ai;
  295.  
  296. #if _USE_SINCOS_TABLE_ == 0
  297.                 double cos_alpha, sin_alpha;
  298.  
  299.                 ar = cos_alpha = FFT_dSinCos[0][delta];     // = cos(2.PI/BlockSize)
  300.                 ai = sin_alpha = FFT_dSinCos[1][delta];     // = sin(2.PI/BlockSize)
  301. #endif
  302.  
  303.                 for (j = 1; j < BlockEnd; j++)
  304.                 {
  305.                     pd_J = pd_J2;
  306.                     pd_K = pd_K2;
  307.  
  308. #if _USE_SINCOS_TABLE_ == 1
  309.                     delta += step;
  310.                     ar = Cosin[delta];
  311.                     ai = Sinus[delta];
  312. #endif
  313.                     for (i = 0; i < iBandSize; i += BlockSize)
  314.                     {
  315.                         ti = pd_K[0];
  316.                         tr = ar * ti - ai * pd_K[1];
  317.                         ti = ar * pd_K[1] + ai * ti;
  318.  
  319.                         pd_K[0] = pd_J[0] - tr;
  320.                         pd_J[0] = pd_J[0] + tr;
  321.  
  322.                         tr = pd_J[1];
  323.                         pd_K[1] = tr - ti;
  324.                         pd_J[1] = tr + ti;
  325.  
  326.                         pd_K = &pd_K[n];            // <==> pd_J += n
  327.                         pd_J = &pd_J[n];            // <==> pd_J += n
  328.                     }
  329.                     pd_J2 += 2;
  330.                     pd_K2 += 2;
  331.  
  332. #if _USE_SINCOS_TABLE_ == 0
  333.                     tr = ar;
  334.                     ar = ar * cos_alpha - ai * sin_alpha;
  335.                     ai = tr * sin_alpha + ai * cos_alpha;
  336. #endif
  337.                 }
  338.             }
  339.             else
  340.             {
  341.                 // Troisième partie : ar = 0, ai = -1.
  342.  
  343.                 pd_J = pd_J2;
  344.                 pd_K = pd_K2;
  345.  
  346.                 if (BlockEnd == 2)
  347.                 {
  348.                     for (i = 0; i < iBandSize; i += 4)
  349.                     {
  350.                         double tz;
  351.  
  352.                         ti = pd_K[1];                   // ai = -1
  353.                         tr = pd_K[0];                   // ar = 0
  354.  
  355.                         tz = pd_J[0];
  356.                         pd_K[0] = tz - ti;
  357.                         pd_J[0] = tz + ti;
  358.  
  359.                         tz = pd_J[1];
  360.                         pd_K[1] = tz + tr;
  361.                         pd_J[1] = tz - tr;
  362.  
  363.                         pd_K    += 8;
  364.                         pd_J    += 8;
  365.                     }
  366.                 }
  367.                 // Quatrième partie : ar = {1; 0.707; 0; -0.707}, ai = {0; -0.707; -1; -0.707}
  368.  
  369.                 if (BlockEnd == 4)
  370.                 {
  371.                     double ar, tz;
  372.  
  373. #if _USE_SINCOS_TABLE_ == 0
  374.                     ar = FFT_dSinCos[0][delta];         // = cos(2.PI/8) = sqrt(2) / 2
  375. #else
  376.                     ar = Cosin[SINCOSSIZE / 6];
  377. #endif
  378.  
  379.                     for (i = 0; i < iBandSize; i += 8)
  380.                     {
  381.                         ti = pd_K[0];
  382.                         tz = pd_K[1];
  383.                         tr = ar * (ti + tz);            // ar = -ai = sqrt(2) / 2
  384.                         ti = ar * (tz - ti);
  385.  
  386.                         tz = pd_J[0];
  387.                         pd_K[0] = tz - tr;
  388.                         pd_J[0] = tz + tr;
  389.  
  390.                         tz = pd_J[1];
  391.                         pd_K[1] = tz - ti;
  392.                         pd_J[1] = tz + ti;
  393.  
  394.                         tr = pd_K[2];                   // ar = 0
  395.                         ti = pd_K[3];                   // ai = -1
  396.  
  397.                         tz = pd_J[2];
  398.                         pd_K[2] = tz - ti;
  399.                         pd_J[2] = tz + ti;
  400.  
  401.                         tz = pd_J[3];
  402.                         pd_K[3] = tz + tr;
  403.                         pd_J[3] = tz - tr;
  404.  
  405.                         ti = pd_K[4];
  406.                         tz = pd_K[5];
  407.                         tr = ar * (ti - tz);            // ar = ai = - sqrt(2) / 2
  408.                         ti = ar * (tz + ti);
  409.  
  410.                         tz = pd_J[4];
  411.                         pd_K[4] = tz + tr;
  412.                         pd_J[4] = tz - tr;
  413.  
  414.                         tz = pd_J[5];
  415.                         pd_K[5] = tz + ti;
  416.                         pd_J[5] = tz - ti;
  417.  
  418.                         pd_K    += 16;
  419.                         pd_J    += 16;
  420.                     }
  421.                 }
  422.             }
  423.  
  424.             BlockEnd = BlockSize;
  425.  
  426. #if _USE_SINCOS_TABLE_ == 0
  427.             ++delta;
  428. #else
  429.             step >>= 1;
  430. #endif
  431.  
  432.         }
  433.     }
  434.  
  435.     // ================================================================
  436.     // Dernière passe : BlockSize = 512, BlockEnd = 256 (passes 8 et 9)
  437.     // ================================================================
  438.  
  439.     for ( ; BlockSize <= iFFTSize; BlockSize += BlockSize)
  440.     {
  441. #if _USE_SINCOS_TABLE_ == 0
  442.         double cos_alpha, sin_alpha, ar, ai;
  443.  
  444.         cos_alpha = FFT_dSinCos[0][delta];
  445.         sin_alpha = FFT_dSinCos[1][delta];
  446. #else
  447.         double ar, ai;
  448.         delta = 0;
  449. #endif
  450.  
  451.         pd_J = &pdFFTIn[0];
  452.         pd_K = &pdFFTIn[BlockSize];
  453.  
  454.         ar = 1.0;
  455.         ai = 0.0;
  456.  
  457.         for (j = 0; j < BlockEnd; j++)
  458.         {
  459.             if (BlockEnd == 256 && iFFTSize == 1024)
  460.             {
  461.                 ti = pd_K[1024];
  462.                 tr = ar * ti - ai * pd_K[1025];
  463.                 ti = ar * pd_K[1025] + ai * ti;
  464.  
  465.                 pd_K[1025] = pd_J[1025] - ti;
  466.                 pd_J[1025] = pd_J[1025] + ti;
  467.  
  468.                 ti = pd_J[1024];
  469.                 pd_K[1024] = ti - tr;
  470.                 pd_J[1024] = ti + tr;
  471.             }
  472.  
  473.             ti = pd_K[0];
  474.             tr = ar * ti - ai * pd_K[1];
  475.             ti = ar * pd_K[1] + ai * ti;
  476.  
  477.             pd_K[0] = pd_J[0] - tr;
  478.             pd_J[0] = pd_J[0] + tr;
  479.  
  480.             tr = pd_J[1];
  481.             pd_K[1] = tr - ti;
  482.             pd_J[1] = tr + ti;
  483.  
  484.             pd_J += 2; pd_K += 2;
  485.  
  486. #if _USE_SINCOS_TABLE_ == 0
  487.             tr = ar;
  488.             ar = ar * cos_alpha - ai * sin_alpha;
  489.             ai = tr * sin_alpha + ai * cos_alpha;
  490. #else
  491.             delta += step;
  492.             ar = Cosin[delta];
  493.             ai = Sinus[delta];
  494. #endif
  495.         }
  496.  
  497.         BlockEnd = BlockSize;
  498.  
  499. #if _USE_SINCOS_TABLE_ == 0
  500.         ++delta;
  501. #else
  502.         step >>= 1;
  503. #endif
  504.     }
  505. }
  506.  
  507. #endif  // _USE_ASM_FCT_
  508. ///
  509. /// FastFFTC::FixFFT_Stereo()
  510. /****** Class FastFFTC/FixFFT_Stereo ****************************************
  511. *
  512. *   NAME
  513. *   FastFFTC::FixFFT_Stereo -- Calcule une FFT.
  514. *
  515. *   SYNOPSIS
  516. *   FastFFTC::FixFFT_Stereo(pd_FFTIn, iFFTSize)
  517. *
  518. *   void FastFFTC::FixFFT_Stereo(REAL *, UWORD);
  519. *
  520. *   FUNCTION
  521. *   Calcule la FFT d'un ensemble d'échantillons stéréo.
  522. *
  523. *   INPUTS
  524. *   pd_FFTIn - Tableau de sortie de la FFT.
  525. *   iFFTSize - Nombre d'échantillons à traiter.
  526. *
  527. *   NOTES
  528. *   La routine FFT effectue son calcul sur les deux canaux simultanément (les
  529. *   échantillons ayant été placés dans l'ordre correct auparavant).
  530. *   Dès lors, il suffit de réordonner les valeurs fournies par la FFT en
  531. *   profitant des symétries pour obtenir la véritable FFT de chacun des deux
  532. *   canaux.
  533. *
  534. *   BUGS
  535. *   Fonction contrôlée conforme le 20.06.98
  536. *
  537. *   SEE ALSO
  538. *   FastFFTC::FastFFT()
  539. *
  540. *****************************************************************************
  541. *
  542. *   Benchmark (3 cycles par seconde) :
  543. *       21.06.98 : 0'19.61
  544. */
  545.  
  546. #if _USE_ASM_FCT_ == 0
  547.  
  548. void FastFFTC::FixFFT_Stereo(T_FFTDATAS *pd_FFTIn, UWORD iFFTSize)
  549. {
  550.     T_FFTDATAS *pd_FFTOut = &pd_FFTIn[BLKSIZE * 2];       // &(*FFT_pdWinSmplShuffled)[1][0][0];
  551.     T_FFTDATAS *pd_FFTIn2 = &pd_FFTIn[iFFTSize * 2];      // &(*FFT_pdWinSmplShuffled)[0][iFFTSize][0];
  552.     double dOneHalf = 0.5;
  553.  
  554.     ULONG j;
  555.     ULONG iFFTLen = iFFTSize >> 1;
  556.  
  557.     // Décoder les données en sortie de FFT pour obtenir la table finale.
  558.  
  559.     pd_FFTOut[0] = pd_FFTIn[1];
  560.     pd_FFTIn[1]  = pd_FFTOut[1] = 0.0;
  561.  
  562.     pd_FFTIn  += 2;
  563.     pd_FFTOut += 2;
  564.  
  565.     for (j = iFFTLen; j > 0; j--)
  566.     {
  567.         pd_FFTIn2 += (-2);
  568.  
  569.         *pd_FFTOut++ = dOneHalf * (pd_FFTIn[1]  + pd_FFTIn2[1]);
  570.         *pd_FFTOut++ = dOneHalf * (pd_FFTIn2[0] - pd_FFTIn[0]);
  571.  
  572.         *pd_FFTIn++ = dOneHalf * (*pd_FFTIn + pd_FFTIn2[0]);
  573.         *pd_FFTIn++ = dOneHalf * (*pd_FFTIn - pd_FFTIn2[1]);
  574.     }
  575. }
  576.  
  577. #endif  // _USE_ASM_FCT_
  578. ///
  579. /// FastFFTC::FixFFT_Mono()
  580. /****** Class FastFFTC/FixFFT_Mono ******************************************
  581. *
  582. *   NAME
  583. *   FastFFTC::FixFFT_Mono -- Calcule une FFT
  584. *
  585. *   SYNOPSIS
  586. *   FastFFTC::FixFFT_Mono(pd_FFTIn, iFFTSize)
  587. *
  588. *   void FastFFTC::FixFFT_Mono(REAL *, UWORD);
  589. *
  590. *   FUNCTION
  591. *   Calcule une FFT sur un groupe d'échantillons (fichier mono).
  592. *
  593. *   INPUTS
  594. *   pd_FFTIn - Tableau de sortie de la FFT.
  595. *   iFFTSize - Taille de la fenêtre FFT (est égale à la moitié du nombre d'
  596. *       échantillons à analyser).
  597. *
  598. *   NOTES
  599. *   Les échantillons sont placés dans un ordre particulier afin de limiter
  600. *   le nombre de calculs (échantillons d'index pair comme partie réelle,
  601. *   échantillons d'index impair comme partie imaginaire).
  602. *   Cette solution permet de limiter la taille de la fenêtre FFT à la moitié
  603. *   de sa taille normale.
  604. *   Il suffit de réordonner les nombres complexes en sortie de FFT pour
  605. *   obtenir le résultat recherché.
  606. *
  607. *   BUGS
  608. *   Fonction contrôlée conforme le 05/07/98
  609. *
  610. *   SEE ALSO
  611. *   FastFFTC::FastFFT()
  612. *
  613. *****************************************************************************
  614. *
  615. */
  616.  
  617. #if _USE_ASM_FCT_ == 0
  618.  
  619. void FastFFTC::FixFFT_Mono(T_FFTDATAS *pd_FFTIn, UWORD iFFTSize)
  620. {
  621.     T_FFTDATAS *pd_FFTIn2 = &pd_FFTIn[iFFTSize * 2];     // &(*FFT_pdWinSmplShuffled)[0][iFFTSize][0];
  622.     T_FFTDATAS *pd_FFTOrg = pd_FFTIn;
  623.  
  624.     ULONG i;
  625.     ULONG iFFTLen = iFFTSize >> 1;
  626.  
  627.     double wr, wi;
  628.  
  629. #if _USE_SINCOS_TABLE_ == 0
  630.     double  wpr, wpi;
  631. #else
  632.     ULONG delta = 0;
  633.     float *Cosin = &FFT_dSinCos[0];
  634.     float *Sinus = &FFT_dSinCos[SINCOSSIZE / 3];
  635. #endif
  636.  
  637.     double h1r, h1i, h2r, h2i;
  638.     double dOneHalf = 0.5;
  639.  
  640.     // Décoder les données en sortie de FFT pour obtenir la table finale.
  641.  
  642. #if _USE_SINCOS_TABLE_ == 0
  643.     wpr   = FFT_dSinCos[0][9];                  // cos(PI/iFFTSize);
  644.     wpi   = FFT_dSinCos[1][9];                  // sin(PI/iFFTSize);
  645.  
  646.     wr    = wpr;                                // = cos(theta)
  647.     wi    = wpi;                                // = sin(theta)
  648. #endif
  649.  
  650.     pd_FFTIn += 2;
  651.  
  652.     for (i = 1; i <= iFFTLen; i++)
  653.     {
  654. #if _USE_SINCOS_TABLE_ == 1
  655.         ++delta;
  656.         wr = Cosin[delta];
  657.         wi = Sinus[delta];
  658. #endif
  659.  
  660.         pd_FFTIn2 += (-2);
  661.  
  662.         h2r = dOneHalf * (pd_FFTIn[1]  + pd_FFTIn2[1]);
  663.         h2i = dOneHalf * (pd_FFTIn2[0] - pd_FFTIn[0]);
  664.  
  665.         h1r = h2r;
  666.         h2r = wr * h1r - wi * h2i;
  667.         h2i = wr * h2i + wi * h1r;
  668.  
  669.         h1r = dOneHalf * (pd_FFTIn[0]  + pd_FFTIn2[0]);
  670.         h1i = dOneHalf * (pd_FFTIn[1]  - pd_FFTIn2[1]);
  671.  
  672.         *pd_FFTIn++  = h2r + h1r;               // + wr * h2r - wi * h2i;
  673.         *pd_FFTIn++  = h2i + h1i;               // + wr * h2i + wi * h2r;
  674.         pd_FFTIn2[0] = h1r - h2r;               // - wr * h2r + wi * h2i;
  675.         pd_FFTIn2[1] = h2i - h1i;               // + wr * h2i + wi * h2r;
  676.  
  677. #if _USE_SINCOS_TABLE_ == 0
  678.         h1r = wr;
  679.         wr = h1r * wpr - wi * wpi;              // = cos(x).cos(theta) - sin(x).sin(theta) = cos(x+theta)
  680.         wi = wi * wpr + h1r * wpi;              // = sin(x).cos(theta) + cos(x).sin(theta) = sin(x+theta)
  681. #endif
  682.     }
  683.  
  684.     iFFTSize *= 2;
  685.  
  686.     pd_FFTOrg[iFFTSize]  = pd_FFTOrg[0] - pd_FFTOrg[1];
  687.     pd_FFTOrg[0]        += pd_FFTOrg[1];
  688.  
  689.     ++pd_FFTOrg;
  690.     pd_FFTOrg[0] = pd_FFTOrg[iFFTSize] = 0.0;
  691. }
  692.  
  693. #endif  // _USE_ASM_FCT_
  694. ///
  695. /// FastFFTC::SetupWinSamples_1024()
  696. /****** Class FastFFTC/SetupWinSamples_1024 *********************************
  697. *
  698. *   NAME
  699. *   FastFFTC::SetupWinSamples_1024 -- Préparer les données pour la FFT.
  700. *
  701. *   SYNOPSIS
  702. *   FastFFTC::SetupWinSamples_1024(pSamples, wNumChannels)
  703. *
  704. *   void FastFFTC::SetupWinSamples_1024(WORD *, UWORD);
  705. *
  706. *   FUNCTION
  707. *   Initialise le tampon d'entrée de la FFT à partir d'un bloc contenant 1024
  708. *   échantillons.
  709. *   Lorsque le fichier audio traité est stéréo, les deux canaux sont préparés
  710. *   simultanément.
  711. *
  712. *   INPUTS
  713. *   pSamples     - Pointeur vers les échantillons.
  714. *   wNumChannels - Nombre de canaux.
  715. *
  716. *   NOTES
  717. *   Afin de tirer partie des symétries de la FFT, les deux canaux sont
  718. *   préparés de sorte que l'un corresponde à la partie réelle, et l'autre
  719. *   à la partie imaginaire. De fait, la FFT effectue son calcul sur les deux
  720. *   canaux simultanément.
  721. *
  722. *   Pour un fichier mono, les échantillons sont alternativement placés dans
  723. *   la partie réelle (échantillons d'indices pairs) ou la partie imaginaire
  724. *   (échantillons d'indices impairs) du tampon d'entrée. De cette façon, la
  725. *   fonction FFT peut être fortement optimisée.
  726. *
  727. *   BUGS
  728. *
  729. *   SEE ALSO
  730. *   FastFFTC::ComputeFastFFT_Stereo()
  731. *
  732. *****************************************************************************
  733. *
  734. *   Benchmarks (12 cycles par seconde) :
  735. *       21.06.98 : 0'25.42
  736. *
  737. */
  738.  
  739. #if _USE_ASM_FCT_ == 0
  740.  
  741. void FastFFTC::SetupWinSamples_1024(WORD *pSamples, UWORD wNumChannels)
  742. {
  743.     const T_FFTDATAS *pdHANNWindow = &(*FFT_pdHannWindowShuffled)[0];
  744.           T_FFTDATAS *pdWSmpls     = &(*FFT_pdWinSmplShuffled)[0][0][0];
  745.           double f;
  746.  
  747.     WORD  *pwSmpl;
  748.     ULONG i, k;
  749.  
  750.     if (wNumChannels == 2)
  751.     {
  752.         for (i = 0; i < 1024 / 4; i++)
  753.         {
  754.             k = GetShuffleIndex_8(i) * 2;
  755.             pwSmpl = &pSamples[k];
  756.  
  757.             f = *pdHANNWindow++;
  758.             *pdWSmpls++ = pwSmpl[0] * f;              // Canal gauche.
  759.             *pdWSmpls++ = pwSmpl[1] * f;              // Canal droit.
  760.  
  761.             f = *pdHANNWindow++;
  762.             *pdWSmpls++ = pwSmpl[512 * 2]     * f;    // Canal gauche.
  763.             *pdWSmpls++ = pwSmpl[512 * 2 + 1] * f;    // Canal droit.
  764.  
  765.             f = *pdHANNWindow++;
  766.             *pdWSmpls++ = pwSmpl[256 * 2]     * f;    // Canal gauche.
  767.             *pdWSmpls++ = pwSmpl[256 * 2 + 1] * f;    // Canal droit.
  768.  
  769.             f = *pdHANNWindow++;
  770.             *pdWSmpls++ = pwSmpl[768 * 2]     * f;    // Canal gauche.
  771.             *pdWSmpls++ = pwSmpl[768 * 2 + 1] * f;    // Canal droit.
  772.         }
  773.     }
  774.     else
  775.     {
  776.         for (i = 0; i < 1024 / 4; i++)
  777.         {
  778.             k = GetShuffleIndex_8(i) * 2;
  779.             pwSmpl = &pSamples[k];
  780.  
  781.             *pdWSmpls++ = pwSmpl[0]   * *pdHANNWindow++;
  782.             *pdWSmpls++ = pwSmpl[1]   *  pdHANNWindow[511];
  783.             *pdWSmpls++ = pwSmpl[512] * *pdHANNWindow++;
  784.             *pdWSmpls++ = pwSmpl[513] *  pdHANNWindow[511];
  785.         }
  786.     }
  787. }
  788.  
  789. #endif  // _USE_ASM_FCT_
  790. ///
  791.  
  792.  
  793.