2313
|
1 /*
|
|
2 * PCM time-domain equalizer
|
|
3 *
|
|
4 * Copyright (C) 2002-2005 Felipe Rivera <liebremx at users sourceforge net>
|
|
5 *
|
|
6 * This program is free software; you can redistribute it and/or modify
|
|
7 * it under the terms of the GNU General Public License as published by
|
|
8 * the Free Software Foundation; either version 2 of the License, or
|
|
9 * (at your option) any later version.
|
|
10 *
|
|
11 * This program is distributed in the hope that it will be useful,
|
|
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
14 * GNU General Public License for more details.
|
|
15 *
|
|
16 * You should have received a copy of the GNU General Public License
|
|
17 * along with this program; if not, write to the Free Software
|
|
18 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
|
|
19 *
|
|
20 * $Id: iir_fpu.c,v 1.3 2005/11/13 20:02:58 lisanet Exp $
|
|
21 */
|
|
22
|
|
23 #include <strings.h>
|
|
24 #include <stdlib.h>
|
|
25 #include <glib.h>
|
|
26 #include "iir.h"
|
|
27 #include "iir_fpu.h"
|
|
28
|
|
29 static sXYData data_history[EQ_MAX_BANDS][EQ_CHANNELS] __attribute__((aligned));
|
|
30 static sXYData data_history2[EQ_MAX_BANDS][EQ_CHANNELS] __attribute__((aligned));
|
|
31 float gain[EQ_MAX_BANDS][EQ_CHANNELS] __attribute__((aligned));
|
|
32 /* random noise */
|
|
33 sample_t dither[256];
|
|
34 gint di;
|
|
35
|
|
36 void set_gain(gint index, gint chn, float val)
|
|
37 {
|
|
38 gain[index][chn] = val;
|
|
39 }
|
|
40
|
|
41 void clean_history()
|
|
42 {
|
|
43 gint n;
|
|
44 /* Zero the history arrays */
|
|
45 bzero(data_history, sizeof(sXYData) * EQ_MAX_BANDS * EQ_CHANNELS);
|
|
46 bzero(data_history2, sizeof(sXYData) * EQ_MAX_BANDS * EQ_CHANNELS);
|
|
47 /* this is only needed if we use fpu code and there's no other place for
|
|
48 the moment to init the dither array*/
|
|
49 for (n = 0; n < 256; n++) {
|
|
50 dither[n] = (rand() % 4) - 2;
|
|
51 }
|
|
52 di = 0;
|
|
53 }
|
|
54
|
|
55 __inline__ int iir(gpointer * d, gint length, gint nch)
|
|
56 {
|
|
57 // FTZ_ON;
|
|
58 gint16 *data = (gint16 *) * d;
|
|
59 /* Indexes for the history arrays
|
|
60 * These have to be kept between calls to this function
|
|
61 * hence they are static */
|
|
62 static gint i = 2, j = 1, k = 0;
|
|
63
|
|
64 gint index, band, channel;
|
|
65 gint tempgint, halflength;
|
|
66 sample_t out[EQ_CHANNELS], pcm[EQ_CHANNELS];
|
|
67
|
|
68 #if 0
|
|
69 // Load the correct filter table according to the sampling rate if needed
|
|
70 if (srate != rate)
|
|
71 {
|
|
72 band_count = eqcfg.band_num;
|
|
73 rate = srate;
|
|
74 iir_cf = get_coeffs(&band_count, rate, eqcfg.use_xmms_original_freqs);
|
|
75 clean_history();
|
|
76 }
|
|
77 #endif
|
|
78
|
|
79 #ifdef BENCHMARK
|
|
80 start_counter();
|
|
81 #endif //BENCHMARK
|
|
82
|
|
83 /**
|
|
84 * IIR filter equation is
|
|
85 * y[n] = 2 * (alpha*(x[n]-x[n-2]) + gamma*y[n-1] - beta*y[n-2])
|
|
86 *
|
|
87 * NOTE: The 2 factor was introduced in the coefficients to save
|
|
88 * a multiplication
|
|
89 *
|
|
90 * This algorithm cascades two filters to get nice filtering
|
|
91 * at the expense of extra CPU cycles
|
|
92 */
|
|
93 /* 16bit, 2 bytes per sample, so divide by two the length of
|
|
94 * the buffer (length is in bytes)
|
|
95 */
|
|
96 halflength = (length >> 1);
|
|
97 for (index = 0; index < halflength; index+=nch)
|
|
98 {
|
|
99 /* For each channel */
|
|
100 for (channel = 0; channel < nch; channel++)
|
|
101 {
|
|
102 pcm[channel] = data[index+channel] * 4;
|
|
103 /* Preamp gain */
|
|
104 pcm[channel] *= (preamp[channel] / 2);
|
|
105
|
|
106 /* add random noise */
|
|
107 pcm[channel] += dither[di];
|
|
108
|
|
109 out[channel] = 0.0;
|
|
110 /* For each band */
|
|
111 for (band = 0; band < band_count; band++)
|
|
112 {
|
|
113 /* Store Xi(n) */
|
|
114 data_history[band][channel].x[i] = pcm[channel];
|
|
115 /* Calculate and store Yi(n) */
|
|
116 data_history[band][channel].y[i] =
|
|
117 (
|
|
118 /* = alpha * [x(n)-x(n-2)] */
|
|
119 iir_cf[band].alpha * ( data_history[band][channel].x[i]
|
|
120 - data_history[band][channel].x[k])
|
|
121 /* + gamma * y(n-1) */
|
|
122 + iir_cf[band].gamma * data_history[band][channel].y[j]
|
|
123 /* - beta * y(n-2) */
|
|
124 - iir_cf[band].beta * data_history[band][channel].y[k]
|
|
125 );
|
|
126 /*
|
|
127 * The multiplication by 2.0 was 'moved' into the coefficients to save
|
|
128 * CPU cycles here */
|
|
129 /* Apply the gain */
|
|
130 out[channel] += data_history[band][channel].y[i]*gain[band][channel]; // * 2.0;
|
|
131 } /* For each band */
|
|
132
|
|
133 if (cfg.eq_extra_filtering)
|
|
134 {
|
|
135 /* Filter the sample again */
|
|
136 for (band = 0; band < band_count; band++)
|
|
137 {
|
|
138 /* Store Xi(n) */
|
|
139 data_history2[band][channel].x[i] = out[channel];
|
|
140 /* Calculate and store Yi(n) */
|
|
141 data_history2[band][channel].y[i] =
|
|
142 (
|
|
143 /* y(n) = alpha * [x(n)-x(n-2)] */
|
|
144 iir_cf[band].alpha * (data_history2[band][channel].x[i]
|
|
145 - data_history2[band][channel].x[k])
|
|
146 /* + gamma * y(n-1) */
|
|
147 + iir_cf[band].gamma * data_history2[band][channel].y[j]
|
|
148 /* - beta * y(n-2) */
|
|
149 - iir_cf[band].beta * data_history2[band][channel].y[k]
|
|
150 );
|
|
151 /* Apply the gain */
|
|
152 out[channel] += data_history2[band][channel].y[i]*gain[band][channel];
|
|
153 } /* For each band */
|
|
154 }
|
|
155
|
|
156 /* Volume stuff
|
|
157 Scale down original PCM sample and add it to the filters
|
|
158 output. This substitutes the multiplication by 0.25
|
|
159 Go back to use the floating point multiplication before the
|
|
160 conversion to give more dynamic range
|
|
161 */
|
|
162 out[channel] += pcm[channel]*0.25;
|
|
163
|
|
164 /* remove random noise */
|
|
165 out[channel] -= dither[di]*0.25;
|
|
166
|
|
167 /* Round and convert to integer */
|
|
168 #ifdef ARCH_PPC
|
|
169 tempgint = round_ppc(out[channel]);
|
|
170 #else
|
|
171 #ifdef ARCH_X86
|
|
172 tempgint = round_trick(out[channel]);
|
|
173 #else
|
|
174 tempgint = (int)out[channel];
|
|
175 #endif
|
|
176 #endif
|
|
177
|
|
178 /* Limit the output */
|
|
179 if (tempgint < -32768)
|
|
180 data[index+channel] = -32768;
|
|
181 else if (tempgint > 32767)
|
|
182 data[index+channel] = 32767;
|
|
183 else
|
|
184 data[index+channel] = tempgint;
|
|
185 } /* For each channel */
|
|
186
|
|
187 /* Wrap around the indexes */
|
|
188 i = (i+1)%3;
|
|
189 j = (j+1)%3;
|
|
190 k = (k+1)%3;
|
|
191 /* random noise index */
|
|
192 di = (di + 1) % 256;
|
|
193
|
|
194 }/* For each pair of samples */
|
|
195
|
|
196 #ifdef BENCHMARK
|
|
197 timex += get_counter();
|
|
198 blength += length;
|
|
199 if (count++ == 1024)
|
|
200 {
|
|
201 printf("FLOATING POINT: %f %d\n",timex/1024.0, blength/1024);
|
|
202 blength = 0;
|
|
203 timex = 0.;
|
|
204 count = 0;
|
|
205 }
|
|
206 #endif // BENCHMARK
|
|
207
|
|
208 // FTZ_OFF;
|
|
209 return length;
|
|
210 }
|