~ubuntu-branches/ubuntu/lucid/mpg123/lucid

« back to all changes in this revision

Viewing changes to src/dct64.c

Tags: upstream-0.60
ImportĀ upstreamĀ versionĀ 0.60

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/*
 
2
        dct64.c: DCT64, the plain C version
 
3
 
 
4
        copyright ?-2006 by the mpg123 project - free software under the terms of the LGPL 2.1
 
5
        see COPYING and AUTHORS files in distribution or http://mpg123.de
 
6
        initially written by Michael Hipp
 
7
*/
 
8
 
 
9
/*
 
10
 * Discrete Cosine Tansform (DCT) for subband synthesis
 
11
 *
 
12
 * -funroll-loops (for gcc) will remove the loops for better performance
 
13
 * using loops in the source-code enhances readabillity
 
14
 *
 
15
 *
 
16
 * TODO: write an optimized version for the down-sampling modes
 
17
 *       (in these modes the bands 16-31 (2:1) or 8-31 (4:1) are zero 
 
18
 */
 
19
 
 
20
#include "config.h"
 
21
#include "mpg123.h"
 
22
 
 
23
void dct64(real *out0,real *out1,real *samples)
 
24
{
 
25
  real bufs[64];
 
26
 
 
27
 {
 
28
  register int i,j;
 
29
  register real *b1,*b2,*bs,*costab;
 
30
 
 
31
  b1 = samples;
 
32
  bs = bufs;
 
33
  costab = pnts[0]+16;
 
34
  b2 = b1 + 32;
 
35
 
 
36
  for(i=15;i>=0;i--)
 
37
    *bs++ = (*b1++ + *--b2); 
 
38
  for(i=15;i>=0;i--)
 
39
    *bs++ = REAL_MUL((*--b2 - *b1++), *--costab);
 
40
 
 
41
  b1 = bufs;
 
42
  costab = pnts[1]+8;
 
43
  b2 = b1 + 16;
 
44
 
 
45
  {
 
46
    for(i=7;i>=0;i--)
 
47
      *bs++ = (*b1++ + *--b2); 
 
48
    for(i=7;i>=0;i--)
 
49
      *bs++ = REAL_MUL((*--b2 - *b1++), *--costab);
 
50
    b2 += 32;
 
51
    costab += 8;
 
52
    for(i=7;i>=0;i--)
 
53
      *bs++ = (*b1++ + *--b2); 
 
54
    for(i=7;i>=0;i--)
 
55
      *bs++ = REAL_MUL((*b1++ - *--b2), *--costab);
 
56
    b2 += 32;
 
57
  }
 
58
 
 
59
  bs = bufs;
 
60
  costab = pnts[2];
 
61
  b2 = b1 + 8;
 
62
 
 
63
  for(j=2;j;j--)
 
64
  {
 
65
    for(i=3;i>=0;i--)
 
66
      *bs++ = (*b1++ + *--b2); 
 
67
    for(i=3;i>=0;i--)
 
68
      *bs++ = REAL_MUL((*--b2 - *b1++), costab[i]);
 
69
    b2 += 16;
 
70
    for(i=3;i>=0;i--)
 
71
      *bs++ = (*b1++ + *--b2); 
 
72
    for(i=3;i>=0;i--)
 
73
      *bs++ = REAL_MUL((*b1++ - *--b2), costab[i]);
 
74
    b2 += 16;
 
75
  }
 
76
 
 
77
  b1 = bufs;
 
78
  costab = pnts[3];
 
79
  b2 = b1 + 4;
 
80
 
 
81
  for(j=4;j;j--)
 
82
  {
 
83
    *bs++ = (*b1++ + *--b2); 
 
84
    *bs++ = (*b1++ + *--b2);
 
85
    *bs++ = REAL_MUL((*--b2 - *b1++), costab[1]);
 
86
    *bs++ = REAL_MUL((*--b2 - *b1++), costab[0]);
 
87
    b2 += 8;
 
88
    *bs++ = (*b1++ + *--b2); 
 
89
    *bs++ = (*b1++ + *--b2);
 
90
    *bs++ = REAL_MUL((*b1++ - *--b2), costab[1]);
 
91
    *bs++ = REAL_MUL((*b1++ - *--b2), costab[0]);
 
92
    b2 += 8;
 
93
  }
 
94
  bs = bufs;
 
95
  costab = pnts[4];
 
96
 
 
97
  for(j=8;j;j--)
 
98
  {
 
99
    real v0,v1;
 
100
    v0=*b1++; v1 = *b1++;
 
101
    *bs++ = (v0 + v1);
 
102
    *bs++ = REAL_MUL((v0 - v1), (*costab));
 
103
    v0=*b1++; v1 = *b1++;
 
104
    *bs++ = (v0 + v1);
 
105
    *bs++ = REAL_MUL((v1 - v0), (*costab));
 
106
  }
 
107
 
 
108
 }
 
109
 
 
110
 
 
111
 {
 
112
  register real *b1;
 
113
  register int i;
 
114
 
 
115
  for(b1=bufs,i=8;i;i--,b1+=4)
 
116
    b1[2] += b1[3];
 
117
 
 
118
  for(b1=bufs,i=4;i;i--,b1+=8)
 
119
  {
 
120
    b1[4] += b1[6];
 
121
    b1[6] += b1[5];
 
122
    b1[5] += b1[7];
 
123
  }
 
124
 
 
125
  for(b1=bufs,i=2;i;i--,b1+=16)
 
126
  {
 
127
    b1[8]  += b1[12];
 
128
    b1[12] += b1[10];
 
129
    b1[10] += b1[14];
 
130
    b1[14] += b1[9];
 
131
    b1[9]  += b1[13];
 
132
    b1[13] += b1[11];
 
133
    b1[11] += b1[15];
 
134
  }
 
135
 }
 
136
 
 
137
 
 
138
  out0[0x10*16] = bufs[0];
 
139
  out0[0x10*15] = bufs[16+0]  + bufs[16+8];
 
140
  out0[0x10*14] = bufs[8];
 
141
  out0[0x10*13] = bufs[16+8]  + bufs[16+4];
 
142
  out0[0x10*12] = bufs[4];
 
143
  out0[0x10*11] = bufs[16+4]  + bufs[16+12];
 
144
  out0[0x10*10] = bufs[12];
 
145
  out0[0x10* 9] = bufs[16+12] + bufs[16+2];
 
146
  out0[0x10* 8] = bufs[2];
 
147
  out0[0x10* 7] = bufs[16+2]  + bufs[16+10];
 
148
  out0[0x10* 6] = bufs[10];
 
149
  out0[0x10* 5] = bufs[16+10] + bufs[16+6];
 
150
  out0[0x10* 4] = bufs[6];
 
151
  out0[0x10* 3] = bufs[16+6]  + bufs[16+14];
 
152
  out0[0x10* 2] = bufs[14];
 
153
  out0[0x10* 1] = bufs[16+14] + bufs[16+1];
 
154
  out0[0x10* 0] = bufs[1];
 
155
 
 
156
  out1[0x10* 0] = bufs[1];
 
157
  out1[0x10* 1] = bufs[16+1]  + bufs[16+9];
 
158
  out1[0x10* 2] = bufs[9];
 
159
  out1[0x10* 3] = bufs[16+9]  + bufs[16+5];
 
160
  out1[0x10* 4] = bufs[5];
 
161
  out1[0x10* 5] = bufs[16+5]  + bufs[16+13];
 
162
  out1[0x10* 6] = bufs[13];
 
163
  out1[0x10* 7] = bufs[16+13] + bufs[16+3];
 
164
  out1[0x10* 8] = bufs[3];
 
165
  out1[0x10* 9] = bufs[16+3]  + bufs[16+11];
 
166
  out1[0x10*10] = bufs[11];
 
167
  out1[0x10*11] = bufs[16+11] + bufs[16+7];
 
168
  out1[0x10*12] = bufs[7];
 
169
  out1[0x10*13] = bufs[16+7]  + bufs[16+15];
 
170
  out1[0x10*14] = bufs[15];
 
171
  out1[0x10*15] = bufs[16+15];
 
172
 
 
173
}
 
174
 
 
175