transfo.c
上传用户:sun1608
上传日期:2007-02-02
资源大小:6116k
文件大小:14k
源码类别:

流媒体/Mpeg4/MP4

开发平台:

Visual C++

  1. /************************* MPEG-2 NBC Audio Decoder **************************
  2.  *                                                                           *
  3. "This software module was originally developed by
  4. AT&T, Dolby Laboratories, Fraunhofer Gesellschaft IIS in the course of
  5. development of the MPEG-2 NBC/MPEG-4 Audio standard ISO/IEC 13818-7,
  6. 14496-1,2 and 3. This software module is an implementation of a part of one or more
  7. MPEG-2 NBC/MPEG-4 Audio tools as specified by the MPEG-2 NBC/MPEG-4
  8. Audio standard. ISO/IEC  gives users of the MPEG-2 NBC/MPEG-4 Audio
  9. standards free license to this software module or modifications thereof for use in
  10. hardware or software products claiming conformance to the MPEG-2 NBC/MPEG-4
  11. Audio  standards. Those intending to use this software module in hardware or
  12. software products are advised that this use may infringe existing patents.
  13. The original developer of this software module and his/her company, the subsequent
  14. editors and their companies, and ISO/IEC have no liability for use of this software
  15. module or modifications thereof in an implementation. Copyright is not released for
  16. non MPEG-2 NBC/MPEG-4 Audio conforming products.The original developer
  17. retains full right to use the code for his/her  own purpose, assign or donate the
  18. code to a third party and to inhibit third party from using the code for non
  19. MPEG-2 NBC/MPEG-4 Audio conforming products. This copyright notice must
  20. be included in all copies or derivative works."
  21. Copyright(c)1996.
  22.  *                                                                           *
  23.  ****************************************************************************/
  24. /*
  25.  * $Id: transfo.c,v 1.5 2002/01/09 22:25:41 wmay Exp $
  26.  */
  27. #include <math.h>
  28. #include "all.h"
  29. #include "transfo.h"
  30. #include "fastfft.h"
  31. /*
  32. #define INTEL_SPL
  33. */
  34. #ifdef INTEL_SPL
  35. #define nsp_UsesAll
  36. #include <nsp.h>
  37. #include <nspwin.h>
  38. /* choose the library that best fits your processor */
  39. #pragma comment (lib, "nsppxl.lib")
  40. #endif
  41. #ifndef M_PI
  42. #define M_PI        3.14159265358979323846
  43. #endif
  44. #ifndef M_PI_2
  45. #define M_PI_2      1.57079632679489661923
  46. #endif
  47. void MDCT_Long(faacDecHandle hDecoder, fftw_real *data)
  48. {
  49. #ifdef INTEL_SPL
  50.     SCplx FFT_data[512];
  51. #else
  52.     fftw_complex FFTarray[512];    /* the array for in-place FFT */
  53. #endif
  54.     fftw_real tempr, tempi, c, s, cold, cfreq, sfreq; /* temps for pre and post twiddle */
  55. /*  fftw_real freq = 0.0030679616611450911f; */
  56.     fftw_real fac,cosfreq8,sinfreq8;
  57.     int i, n;
  58.     int b = 2048 >> 1;
  59.     int N4 = 2048 >> 2;
  60.     int N2 = 2048 >> 1;
  61.     int a = 2048 - b;
  62.     int a2 = a >> 1;
  63.     int a4 = a >> 2;
  64.     int b4 = b >> 2;
  65.     int unscambled;
  66.     /* Choosing to allocate 2/N factor to Inverse Xform! */
  67.     fac = 2.; /* 2 from MDCT inverse  to forward */
  68.     /* prepare for recurrence relation in pre-twiddle */
  69.     cfreq = 0.99999529123306274f;
  70.     sfreq = 0.0030679567717015743f;
  71.     cosfreq8 = 0.99999994039535522f;
  72.     sinfreq8 = 0.00038349519824312089f;
  73.     c = cosfreq8;
  74.     s = sinfreq8;
  75.     for (i = 0; i < N4; i++)
  76.     {
  77.         /* calculate real and imaginary parts of g(n) or G(p) */
  78.     n = 2048 / 2 - 1 - 2 * i;
  79.         if (i < b4) {
  80.       tempr = data [a2 + n] + data [2048 + a2 - 1 - n]; /* use second form of e(n) for n = N / 2 - 1 - 2i */
  81.         } else {
  82.             tempr = data [a2 + n] - data [a2 - 1 - n]; /* use first form of e(n) for n = N / 2 - 1 - 2i */
  83.         }
  84.         n = 2 * i;
  85.         if (i < a4) {
  86.             tempi = data [a2 + n] - data [a2 - 1 - n]; /* use first form of e(n) for n=2i */
  87.         } else {
  88.       tempi = data [a2 + n] + data [2048 + a2 - 1 - n]; /* use second form of e(n) for n=2i*/
  89.         }
  90.         /* calculate pre-twiddled FFT input */
  91. #ifdef INTEL_SPL
  92.     FFT_data[i].re = tempr * c + tempi * s;
  93.     FFT_data[i].im = tempi * c - tempr * s;
  94. #else
  95.         FFTarray [i].re = tempr * c + tempi * s;
  96.         FFTarray [i].im = tempi * c - tempr * s;
  97. #endif
  98.         /* use recurrence to prepare cosine and sine for next value of i */
  99.         cold = c;
  100.         c = c * cfreq - s * sfreq;
  101.         s = s * cfreq + cold * sfreq;
  102.     }
  103.     /* Perform in-place complex FFT of length N/4 */
  104. #ifdef INTEL_SPL
  105.     nspcFft(FFT_data, 9, NSP_Forw);
  106. #else
  107.     pfftw_512(FFTarray);
  108. #endif
  109.     /* prepare for recurrence relations in post-twiddle */
  110.     c = cosfreq8;
  111.     s = sinfreq8;
  112.     /* post-twiddle FFT output and then get output data */
  113.     for (i = 0; i < N4; i++)
  114.     {
  115.       /* get post-twiddled FFT output  */
  116.       /* Note: fac allocates 4/N factor from IFFT to forward and inverse */
  117. #ifdef INTEL_SPL
  118.       tempr = fac * (FFT_data[i].re * c + FFT_data[i].im * s);
  119.       tempi = fac * (FFT_data[i].im * c - FFT_data[i].re * s);
  120. #else
  121.       unscambled = hDecoder->unscambled512[i];
  122.       tempr = fac * (FFTarray [unscambled].re * c + FFTarray [unscambled].im * s);
  123.       tempi = fac * (FFTarray [unscambled].im * c - FFTarray [unscambled].re * s);
  124. #endif
  125.       /* fill in output values */
  126.       data [2 * i] = -tempr;   /* first half even */
  127.       data [N2 - 1 - 2 * i] = tempi;  /* first half odd */
  128.       data [N2 + 2 * i] = -tempi;  /* second half even */
  129.       data [2048 - 1 - 2 * i] = tempr;  /* second half odd */
  130.       /* use recurrence to prepare cosine and sine for next value of i */
  131.       cold = c;
  132.       c = c * cfreq - s * sfreq;
  133.       s = s * cfreq + cold * sfreq;
  134.     }
  135. }
  136. void MDCT_Short(faacDecHandle hDecoder, fftw_real *data)
  137. {
  138. #ifdef INTEL_SPL
  139.     SCplx FFT_data[64];
  140. #else
  141.     fftw_complex FFTarray[64];    /* the array for in-place FFT */
  142. #endif
  143.     fftw_real tempr, tempi, c, s, cold, cfreq, sfreq; /* temps for pre and post twiddle */
  144. /*  fftw_real freq = 0.024543693289160728f; */
  145.     fftw_real fac,cosfreq8,sinfreq8;
  146.     int i, n;
  147.     int b = 256 >> 1;
  148.     int N4 = 256 >> 2;
  149.     int N2 = 256 >> 1;
  150.     int a = 256 - b;
  151.     int a2 = a >> 1;
  152.     int a4 = a >> 2;
  153.     int b4 = b >> 2;
  154.     int unscambled;
  155.     /* Choosing to allocate 2/N factor to Inverse Xform! */
  156.     fac = 2.; /* 2 from MDCT inverse  to forward */
  157.     /* prepare for recurrence relation in pre-twiddle */
  158.     cfreq = 0.99969881772994995f;
  159.     sfreq = 0.024541229009628296f;
  160.     cosfreq8 = 0.99999529123306274f;
  161.     sinfreq8 = 0.0030679568483393833f;
  162.     c = cosfreq8;
  163.     s = sinfreq8;
  164.     for (i = 0; i < N4; i++)
  165.     {
  166.         /* calculate real and imaginary parts of g(n) or G(p) */
  167.     n = 256 / 2 - 1 - 2 * i;
  168.         if (i < b4) {
  169.       tempr = data [a2 + n] + data [256 + a2 - 1 - n]; /* use second form of e(n) for n = N / 2 - 1 - 2i */
  170.         } else {
  171.             tempr = data [a2 + n] - data [a2 - 1 - n]; /* use first form of e(n) for n = N / 2 - 1 - 2i */
  172.         }
  173.         n = 2 * i;
  174.         if (i < a4) {
  175.             tempi = data [a2 + n] - data [a2 - 1 - n]; /* use first form of e(n) for n=2i */
  176.         } else {
  177.       tempi = data [a2 + n] + data [256 + a2 - 1 - n]; /* use second form of e(n) for n=2i*/
  178.         }
  179.         /* calculate pre-twiddled FFT input */
  180. #ifdef INTEL_SPL
  181.     FFT_data[i].re = tempr * c + tempi * s;
  182.     FFT_data[i].im = tempi * c - tempr * s;
  183. #else
  184.         FFTarray [i].re = tempr * c + tempi * s;
  185.         FFTarray [i].im = tempi * c - tempr * s;
  186. #endif
  187.         /* use recurrence to prepare cosine and sine for next value of i */
  188.         cold = c;
  189.         c = c * cfreq - s * sfreq;
  190.         s = s * cfreq + cold * sfreq;
  191.     }
  192.     /* Perform in-place complex FFT of length N/4 */
  193. #ifdef INTEL_SPL
  194.     nspcFft(FFT_data, 6, NSP_Forw);
  195. #else
  196.     pfftw_64(FFTarray);
  197. #endif
  198.     /* prepare for recurrence relations in post-twiddle */
  199.     c = cosfreq8;
  200.     s = sinfreq8;
  201.     /* post-twiddle FFT output and then get output data */
  202.     for (i = 0; i < N4; i++) {
  203.         /* get post-twiddled FFT output  */
  204.         /* Note: fac allocates 4/N factor from IFFT to forward and inverse */
  205. #ifdef INTEL_SPL
  206.     tempr = fac * (FFT_data[i].re * c + FFT_data[i].im * s);
  207.     tempi = fac * (FFT_data[i].im * c - FFT_data[i].re * s);
  208. #else
  209.         unscambled = hDecoder->unscambled64[i];
  210.     tempr = fac * (FFTarray [unscambled].re * c + FFTarray [unscambled].im * s);
  211.     tempi = fac * (FFTarray [unscambled].im * c - FFTarray [unscambled].re * s);
  212. #endif
  213.         /* fill in output values */
  214.         data [2 * i] = -tempr;   /* first half even */
  215.         data [N2 - 1 - 2 * i] = tempi;  /* first half odd */
  216.         data [N2 + 2 * i] = -tempi;  /* second half even */
  217.     data [256 - 1 - 2 * i] = tempr;  /* second half odd */
  218.         /* use recurrence to prepare cosine and sine for next value of i */
  219.         cold = c;
  220.         c = c * cfreq - s * sfreq;
  221.         s = s * cfreq + cold * sfreq;
  222.     }
  223. }
  224. void IMDCT_Long(faacDecHandle hDecoder, fftw_real *data)
  225. {
  226. #ifdef INTEL_SPL
  227.   SCplx FFT_data[512];    /* the array for in-place FFT */
  228. #else
  229.     fftw_complex FFTarray[512];    /* the array for in-place FFT */
  230. #endif
  231.     fftw_real tempr, tempi, c, s, cold, cfreq, sfreq; /* temps for pre and post twiddle */
  232. /*  fftw_real freq = 0.0030679616611450911f; */
  233.     fftw_real fac, cosfreq8, sinfreq8;
  234.     int i;
  235.   int Nd2 = 2048 >> 1;
  236.   int Nd4 = 2048 >> 2;
  237.   int Nd8 = 2048 >> 3;
  238.     int unscambled;
  239.     /* Choosing to allocate 2/N factor to Inverse Xform! */
  240.   fac = 0.0009765625f;
  241.     /* prepare for recurrence relation in pre-twiddle */
  242.   cfreq = 0.99999529123306274f;
  243.   sfreq = 0.0030679567717015743f;
  244.   cosfreq8 = 0.99999994039535522f;
  245.   sinfreq8 = 0.00038349519824312089f;
  246.     c = cosfreq8;
  247.     s = sinfreq8;
  248.     for (i = 0; i < Nd4; i++) {
  249.         /* calculate real and imaginary parts of g(n) or G(p) */
  250.         tempr = -data [2 * i];
  251.         tempi = data [Nd2 - 1 - 2 * i];
  252.         /* calculate pre-twiddled FFT input */
  253. #ifdef INTEL_SPL
  254.     FFT_data[i].re = tempr * c - tempi * s;
  255.     FFT_data[i].im = tempi * c + tempr * s;
  256. #else
  257.         unscambled = hDecoder->unscambled512[i];
  258.         FFTarray [unscambled].re = tempr * c - tempi * s;
  259.         FFTarray [unscambled].im = tempi * c + tempr * s;
  260. #endif
  261.         /* use recurrence to prepare cosine and sine for next value of i */
  262.         cold = c;
  263.         c = c * cfreq - s * sfreq;
  264.         s = s * cfreq + cold * sfreq;
  265.     }
  266.     /* Perform in-place complex IFFT of length N/4 */
  267. #ifdef INTEL_SPL
  268.   nspcFft(FFT_data, 9, NSP_Inv | NSP_NoScale);
  269. #else
  270.   pfftwi_512(FFTarray);
  271. #endif
  272.     /* prepare for recurrence relations in post-twiddle */
  273.     c = cosfreq8;
  274.     s = sinfreq8;
  275.     /* post-twiddle FFT output and then get output data */
  276.     for (i = 0; i < Nd4; i++) {
  277.         /* get post-twiddled FFT output  */
  278. #ifdef INTEL_SPL
  279.   tempr = fac * (FFT_data[i].re * c - FFT_data[i].im * s);
  280.       tempi = fac * (FFT_data[i].im * c + FFT_data[i].re * s);
  281. #else
  282.     tempr = fac * (FFTarray[i].re * c - FFTarray[i].im * s);
  283.         tempi = fac * (FFTarray[i].im * c + FFTarray[i].re * s);
  284. #endif
  285.         /* fill in output values */
  286.         data [Nd2 + Nd4 - 1 - 2 * i] = tempr;
  287.         if (i < Nd8) {
  288.             data [Nd2 + Nd4 + 2 * i] = tempr;
  289.         } else {
  290.             data [2 * i - Nd4] = -tempr;
  291.         }
  292.         data [Nd4 + 2 * i] = tempi;
  293.         if (i < Nd8) {
  294.             data [Nd4 - 1 - 2 * i] = -tempi;
  295.         } else {
  296.       data [Nd4 + 2048 - 1 - 2*i] = tempi;
  297.         }
  298.         /* use recurrence to prepare cosine and sine for next value of i */
  299.         cold = c;
  300.         c = c * cfreq - s * sfreq;
  301.         s = s * cfreq + cold * sfreq;
  302.   }
  303. }
  304. void IMDCT_Short(faacDecHandle hDecoder, fftw_real *data)
  305. {
  306. #ifdef INTEL_SPL
  307.   SCplx FFT_data[64];    /* the array for in-place FFT */
  308. #else
  309.     fftw_complex FFTarray[64];    /* the array for in-place FFT */
  310. #endif
  311.     fftw_real tempr, tempi, c, s, cold, cfreq, sfreq; /* temps for pre and post twiddle */
  312. /*  fftw_real freq = 0.024543693289160728f; */
  313.     fftw_real fac, cosfreq8, sinfreq8;
  314.     int i;
  315.     int Nd2 = 256 >> 1;
  316.     int Nd4 = 256 >> 2;
  317.     int Nd8 = 256 >> 3;
  318.     int unscambled;
  319.     /* Choosing to allocate 2/N factor to Inverse Xform! */
  320.   fac = 0.0078125f; /* remaining 2/N from 4/N IFFT factor */
  321.     /* prepare for recurrence relation in pre-twiddle */
  322.   cfreq = 0.99969881772994995f;
  323.   sfreq = 0.024541229009628296f;
  324.   cosfreq8 = 0.99999529123306274f;
  325.   sinfreq8 = 0.0030679568483393833f;
  326.     c = cosfreq8;
  327.     s = sinfreq8;
  328.     for (i = 0; i < Nd4; i++)
  329.   {
  330.         /* calculate real and imaginary parts of g(n) or G(p) */
  331.         tempr = -data [2 * i];
  332.         tempi = data [Nd2 - 1 - 2 * i];
  333.         /* calculate pre-twiddled FFT input */
  334. #ifdef INTEL_SPL
  335.     FFT_data[i].re = tempr * c - tempi * s;
  336.     FFT_data[i].im = tempi * c + tempr * s;
  337. #else
  338.         unscambled = hDecoder->unscambled64[i];
  339.         FFTarray [unscambled].re = tempr * c - tempi * s;
  340.         FFTarray [unscambled].im = tempi * c + tempr * s;
  341. #endif
  342.         /* use recurrence to prepare cosine and sine for next value of i */
  343.         cold = c;
  344.         c = c * cfreq - s * sfreq;
  345.         s = s * cfreq + cold * sfreq;
  346.     }
  347.     /* Perform in-place complex IFFT of length N/4 */
  348. #ifdef INTEL_SPL
  349.   nspcFft(FFT_data, 6, NSP_Inv | NSP_NoScale);
  350. #else
  351.     pfftwi_64(FFTarray);
  352. #endif
  353.     /* prepare for recurrence relations in post-twiddle */
  354.     c = cosfreq8;
  355.     s = sinfreq8;
  356.     /* post-twiddle FFT output and then get output data */
  357.     for (i = 0; i < Nd4; i++) {
  358.         /* get post-twiddled FFT output  */
  359. #ifdef INTEL_SPL
  360.       tempr = fac * (FFT_data[i].re * c - FFT_data[i].im * s);
  361.       tempi = fac * (FFT_data[i].im * c + FFT_data[i].re * s);
  362. #else
  363.     tempr = fac * (FFTarray[i].re * c - FFTarray[i].im * s);
  364.         tempi = fac * (FFTarray[i].im * c + FFTarray[i].re * s);
  365. #endif
  366.         /* fill in output values */
  367.         data [Nd2 + Nd4 - 1 - 2 * i] = tempr;
  368.         if (i < Nd8) {
  369.             data [Nd2 + Nd4 + 2 * i] = tempr;
  370.         } else {
  371.             data [2 * i - Nd4] = -tempr;
  372.         }
  373.         data [Nd4 + 2 * i] = tempi;
  374.         if (i < Nd8) {
  375.             data [Nd4 - 1 - 2 * i] = -tempi;
  376.         } else {
  377.       data [Nd4 + 256 - 1 - 2*i] = tempi;
  378.         }
  379.         /* use recurrence to prepare cosine and sine for next value of i */
  380.         cold = c;
  381.         c = c * cfreq - s * sfreq;
  382.         s = s * cfreq + cold * sfreq;
  383.     }
  384. }