2
0

dct32.c 10 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280
  1. /* ***** BEGIN LICENSE BLOCK *****
  2. * Version: RCSL 1.0/RPSL 1.0
  3. *
  4. * Portions Copyright (c) 1995-2002 RealNetworks, Inc. All Rights Reserved.
  5. *
  6. * The contents of this file, and the files included with this file, are
  7. * subject to the current version of the RealNetworks Public Source License
  8. * Version 1.0 (the "RPSL") available at
  9. * http://www.helixcommunity.org/content/rpsl unless you have licensed
  10. * the file under the RealNetworks Community Source License Version 1.0
  11. * (the "RCSL") available at http://www.helixcommunity.org/content/rcsl,
  12. * in which case the RCSL will apply. You may also obtain the license terms
  13. * directly from RealNetworks. You may not use this file except in
  14. * compliance with the RPSL or, if you have a valid RCSL with RealNetworks
  15. * applicable to this file, the RCSL. Please see the applicable RPSL or
  16. * RCSL for the rights, obligations and limitations governing use of the
  17. * contents of the file.
  18. *
  19. * This file is part of the Helix DNA Technology. RealNetworks is the
  20. * developer of the Original Code and owns the copyrights in the portions
  21. * it created.
  22. *
  23. * This file, and the files included with this file, is distributed and made
  24. * available on an 'AS IS' basis, WITHOUT WARRANTY OF ANY KIND, EITHER
  25. * EXPRESS OR IMPLIED, AND REALNETWORKS HEREBY DISCLAIMS ALL SUCH WARRANTIES,
  26. * INCLUDING WITHOUT LIMITATION, ANY WARRANTIES OF MERCHANTABILITY, FITNESS
  27. * FOR A PARTICULAR PURPOSE, QUIET ENJOYMENT OR NON-INFRINGEMENT.
  28. *
  29. * Technology Compatibility Kit Test Suite(s) Location:
  30. * http://www.helixcommunity.org/content/tck
  31. *
  32. * Contributor(s):
  33. *
  34. * ***** END LICENSE BLOCK ***** */
  35. /**************************************************************************************
  36. * Fixed-point MP3 decoder
  37. * Jon Recker (jrecker@real.com), Ken Cooke (kenc@real.com)
  38. * June 2003
  39. *
  40. * dct32.c - optimized implementations of 32-point DCT for matrixing stage of
  41. * polyphase filter
  42. **************************************************************************************/
  43. #include "coder.h"
  44. #include "assembly.h"
  45. #define COS0_0 0x4013c251 /* Q31 */
  46. #define COS0_1 0x40b345bd /* Q31 */
  47. #define COS0_2 0x41fa2d6d /* Q31 */
  48. #define COS0_3 0x43f93421 /* Q31 */
  49. #define COS0_4 0x46cc1bc4 /* Q31 */
  50. #define COS0_5 0x4a9d9cf0 /* Q31 */
  51. #define COS0_6 0x4fae3711 /* Q31 */
  52. #define COS0_7 0x56601ea7 /* Q31 */
  53. #define COS0_8 0x5f4cf6eb /* Q31 */
  54. #define COS0_9 0x6b6fcf26 /* Q31 */
  55. #define COS0_10 0x7c7d1db3 /* Q31 */
  56. #define COS0_11 0x4ad81a97 /* Q30 */
  57. #define COS0_12 0x5efc8d96 /* Q30 */
  58. #define COS0_13 0x41d95790 /* Q29 */
  59. #define COS0_14 0x6d0b20cf /* Q29 */
  60. #define COS0_15 0x518522fb /* Q27 */
  61. #define COS1_0 0x404f4672 /* Q31 */
  62. #define COS1_1 0x42e13c10 /* Q31 */
  63. #define COS1_2 0x48919f44 /* Q31 */
  64. #define COS1_3 0x52cb0e63 /* Q31 */
  65. #define COS1_4 0x64e2402e /* Q31 */
  66. #define COS1_5 0x43e224a9 /* Q30 */
  67. #define COS1_6 0x6e3c92c1 /* Q30 */
  68. #define COS1_7 0x519e4e04 /* Q28 */
  69. #define COS2_0 0x4140fb46 /* Q31 */
  70. #define COS2_1 0x4cf8de88 /* Q31 */
  71. #define COS2_2 0x73326bbf /* Q31 */
  72. #define COS2_3 0x52036742 /* Q29 */
  73. #define COS3_0 0x4545e9ef /* Q31 */
  74. #define COS3_1 0x539eba45 /* Q30 */
  75. #define COS4_0 0x5a82799a /* Q31 */
  76. // faster in ROM
  77. static const int dcttab[48] PROGMEM = {
  78. /* first pass */
  79. COS0_0, COS0_15, COS1_0, /* 31, 27, 31 */
  80. COS0_1, COS0_14, COS1_1, /* 31, 29, 31 */
  81. COS0_2, COS0_13, COS1_2, /* 31, 29, 31 */
  82. COS0_3, COS0_12, COS1_3, /* 31, 30, 31 */
  83. COS0_4, COS0_11, COS1_4, /* 31, 30, 31 */
  84. COS0_5, COS0_10, COS1_5, /* 31, 31, 30 */
  85. COS0_6, COS0_9, COS1_6, /* 31, 31, 30 */
  86. COS0_7, COS0_8, COS1_7, /* 31, 31, 28 */
  87. /* second pass */
  88. COS2_0, COS2_3, COS3_0, /* 31, 29, 31 */
  89. COS2_1, COS2_2, COS3_1, /* 31, 31, 30 */
  90. -COS2_0, -COS2_3, COS3_0, /* 31, 29, 31 */
  91. -COS2_1, -COS2_2, COS3_1, /* 31, 31, 30 */
  92. COS2_0, COS2_3, COS3_0, /* 31, 29, 31 */
  93. COS2_1, COS2_2, COS3_1, /* 31, 31, 30 */
  94. -COS2_0, -COS2_3, COS3_0, /* 31, 29, 31 */
  95. -COS2_1, -COS2_2, COS3_1, /* 31, 31, 30 */
  96. };
  97. #define D32FP(i, s0, s1, s2) { \
  98. a0 = buf[i]; a3 = buf[31-i]; \
  99. a1 = buf[15-i]; a2 = buf[16+i]; \
  100. b0 = a0 + a3; b3 = MULSHIFT32(*cptr++, a0 - a3) << (s0); \
  101. b1 = a1 + a2; b2 = MULSHIFT32(*cptr++, a1 - a2) << (s1); \
  102. buf[i] = b0 + b1; buf[15-i] = MULSHIFT32(*cptr, b0 - b1) << (s2); \
  103. buf[16+i] = b2 + b3; buf[31-i] = MULSHIFT32(*cptr++, b3 - b2) << (s2); \
  104. }
  105. /**************************************************************************************
  106. * Function: FDCT32
  107. *
  108. * Description: Ken's highly-optimized 32-point DCT (radix-4 + radix-8)
  109. *
  110. * Inputs: input buffer, length = 32 samples
  111. * require at least 6 guard bits in input vector x to avoid possibility
  112. * of overflow in internal calculations (see bbtest_imdct test app)
  113. * buffer offset and oddblock flag for polyphase filter input buffer
  114. * number of guard bits in input
  115. *
  116. * Outputs: output buffer, data copied and interleaved for polyphase filter
  117. * no guarantees about number of guard bits in output
  118. *
  119. * Return: none
  120. *
  121. * Notes: number of muls = 4*8 + 12*4 = 80
  122. * final stage of DCT is hardcoded to shuffle data into the proper order
  123. * for the polyphase filterbank
  124. * fully unrolled stage 1, for max precision (scale the 1/cos() factors
  125. * differently, depending on magnitude)
  126. * guard bit analysis verified by exhaustive testing of all 2^32
  127. * combinations of max pos/max neg values in x[]
  128. *
  129. * TODO: code organization and optimization for ARM
  130. * possibly interleave stereo (cut # of coef loads in half - may not have
  131. * enough registers)
  132. **************************************************************************************/
  133. // about 1ms faster in RAM
  134. /* attribute__ ((section (".data"))) */ void FDCT32(int *buf, int *dest, int offset, int oddBlock, int gb)
  135. {
  136. int i, s, tmp, es;
  137. const int *cptr = dcttab;
  138. int a0, a1, a2, a3, a4, a5, a6, a7;
  139. int b0, b1, b2, b3, b4, b5, b6, b7;
  140. int *d;
  141. /* scaling - ensure at least 6 guard bits for DCT
  142. * (in practice this is already true 99% of time, so this code is
  143. * almost never triggered)
  144. */
  145. es = 0;
  146. if (gb < 6) {
  147. es = 6 - gb;
  148. for (i = 0; i < 32; i++)
  149. buf[i] >>= es;
  150. }
  151. /* first pass */
  152. D32FP(0, 1, 5, 1);
  153. D32FP(1, 1, 3, 1);
  154. D32FP(2, 1, 3, 1);
  155. D32FP(3, 1, 2, 1);
  156. D32FP(4, 1, 2, 1);
  157. D32FP(5, 1, 1, 2);
  158. D32FP(6, 1, 1, 2);
  159. D32FP(7, 1, 1, 4);
  160. /* second pass */
  161. for (i = 4; i > 0; i--) {
  162. a0 = buf[0]; a7 = buf[7]; a3 = buf[3]; a4 = buf[4];
  163. b0 = a0 + a7; b7 = MULSHIFT32(*cptr++, a0 - a7) << 1;
  164. b3 = a3 + a4; b4 = MULSHIFT32(*cptr++, a3 - a4) << 3;
  165. a0 = b0 + b3; a3 = MULSHIFT32(*cptr, b0 - b3) << 1;
  166. a4 = b4 + b7; a7 = MULSHIFT32(*cptr++, b7 - b4) << 1;
  167. a1 = buf[1]; a6 = buf[6]; a2 = buf[2]; a5 = buf[5];
  168. b1 = a1 + a6; b6 = MULSHIFT32(*cptr++, a1 - a6) << 1;
  169. b2 = a2 + a5; b5 = MULSHIFT32(*cptr++, a2 - a5) << 1;
  170. a1 = b1 + b2; a2 = MULSHIFT32(*cptr, b1 - b2) << 2;
  171. a5 = b5 + b6; a6 = MULSHIFT32(*cptr++, b6 - b5) << 2;
  172. b0 = a0 + a1; b1 = MULSHIFT32(COS4_0, a0 - a1) << 1;
  173. b2 = a2 + a3; b3 = MULSHIFT32(COS4_0, a3 - a2) << 1;
  174. buf[0] = b0; buf[1] = b1;
  175. buf[2] = b2 + b3; buf[3] = b3;
  176. b4 = a4 + a5; b5 = MULSHIFT32(COS4_0, a4 - a5) << 1;
  177. b6 = a6 + a7; b7 = MULSHIFT32(COS4_0, a7 - a6) << 1;
  178. b6 += b7;
  179. buf[4] = b4 + b6; buf[5] = b5 + b7;
  180. buf[6] = b5 + b6; buf[7] = b7;
  181. buf += 8;
  182. }
  183. buf -= 32; /* reset */
  184. /* sample 0 - always delayed one block */
  185. d = dest + 64*16 + ((offset - oddBlock) & 7) + (oddBlock ? 0 : VBUF_LENGTH);
  186. s = buf[ 0]; d[0] = d[8] = s;
  187. /* samples 16 to 31 */
  188. d = dest + offset + (oddBlock ? VBUF_LENGTH : 0);
  189. s = buf[ 1]; d[0] = d[8] = s; d += 64;
  190. tmp = buf[25] + buf[29];
  191. s = buf[17] + tmp; d[0] = d[8] = s; d += 64;
  192. s = buf[ 9] + buf[13]; d[0] = d[8] = s; d += 64;
  193. s = buf[21] + tmp; d[0] = d[8] = s; d += 64;
  194. tmp = buf[29] + buf[27];
  195. s = buf[ 5]; d[0] = d[8] = s; d += 64;
  196. s = buf[21] + tmp; d[0] = d[8] = s; d += 64;
  197. s = buf[13] + buf[11]; d[0] = d[8] = s; d += 64;
  198. s = buf[19] + tmp; d[0] = d[8] = s; d += 64;
  199. tmp = buf[27] + buf[31];
  200. s = buf[ 3]; d[0] = d[8] = s; d += 64;
  201. s = buf[19] + tmp; d[0] = d[8] = s; d += 64;
  202. s = buf[11] + buf[15]; d[0] = d[8] = s; d += 64;
  203. s = buf[23] + tmp; d[0] = d[8] = s; d += 64;
  204. tmp = buf[31];
  205. s = buf[ 7]; d[0] = d[8] = s; d += 64;
  206. s = buf[23] + tmp; d[0] = d[8] = s; d += 64;
  207. s = buf[15]; d[0] = d[8] = s; d += 64;
  208. s = tmp; d[0] = d[8] = s;
  209. /* samples 16 to 1 (sample 16 used again) */
  210. d = dest + 16 + ((offset - oddBlock) & 7) + (oddBlock ? 0 : VBUF_LENGTH);
  211. s = buf[ 1]; d[0] = d[8] = s; d += 64;
  212. tmp = buf[30] + buf[25];
  213. s = buf[17] + tmp; d[0] = d[8] = s; d += 64;
  214. s = buf[14] + buf[ 9]; d[0] = d[8] = s; d += 64;
  215. s = buf[22] + tmp; d[0] = d[8] = s; d += 64;
  216. s = buf[ 6]; d[0] = d[8] = s; d += 64;
  217. tmp = buf[26] + buf[30];
  218. s = buf[22] + tmp; d[0] = d[8] = s; d += 64;
  219. s = buf[10] + buf[14]; d[0] = d[8] = s; d += 64;
  220. s = buf[18] + tmp; d[0] = d[8] = s; d += 64;
  221. s = buf[ 2]; d[0] = d[8] = s; d += 64;
  222. tmp = buf[28] + buf[26];
  223. s = buf[18] + tmp; d[0] = d[8] = s; d += 64;
  224. s = buf[12] + buf[10]; d[0] = d[8] = s; d += 64;
  225. s = buf[20] + tmp; d[0] = d[8] = s; d += 64;
  226. s = buf[ 4]; d[0] = d[8] = s; d += 64;
  227. tmp = buf[24] + buf[28];
  228. s = buf[20] + tmp; d[0] = d[8] = s; d += 64;
  229. s = buf[ 8] + buf[12]; d[0] = d[8] = s; d += 64;
  230. s = buf[16] + tmp; d[0] = d[8] = s;
  231. /* this is so rarely invoked that it's not worth making two versions of the output
  232. * shuffle code (one for no shift, one for clip + variable shift) like in IMDCT
  233. * here we just load, clip, shift, and store on the rare instances that es != 0
  234. */
  235. if (es) {
  236. d = dest + 64*16 + ((offset - oddBlock) & 7) + (oddBlock ? 0 : VBUF_LENGTH);
  237. s = d[0]; CLIP_2N(s, 31 - es); d[0] = d[8] = (s << es);
  238. d = dest + offset + (oddBlock ? VBUF_LENGTH : 0);
  239. for (i = 16; i <= 31; i++) {
  240. s = d[0]; CLIP_2N(s, 31 - es); d[0] = d[8] = (s << es); d += 64;
  241. }
  242. d = dest + 16 + ((offset - oddBlock) & 7) + (oddBlock ? 0 : VBUF_LENGTH);
  243. for (i = 15; i >= 0; i--) {
  244. s = d[0]; CLIP_2N(s, 31 - es); d[0] = d[8] = (s << es); d += 64;
  245. }
  246. }
  247. }