Mercurial > audlegacy

comparison src/libSAD/dither.c @ 4232:704607c1f858

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
1st attempt to integrate dithering and RG engine
author Eugene Zagidullin <e.asphyx@gmail.com>
date Wed, 30 Jan 2008 01:22:37 +0300
parents
children 74c6f3d3cf1d
comparison
equal deleted inserted replaced
4225:f6625617716b 4232:704607c1f858
1 /* Scale & Dither library (libSAD)
2 * High-precision bit depth converter with ReplayGain support
3 *
4 * (c)2007 by Eugene Zagidullin (e.asphyx@gmail.com)
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 * You should have received a copy of the GNU General Public License
16 * along with this program; if not, write to the Free Software
17 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18 */
19
20
21 #include "common.h"
22 #include "buffer.h"
23 #include "noisegen.h"
24
25 #include <assert.h>
26 #include <math.h>
27
28 /*
29 * Supported conversions:
30 *
31 * O U T P U T
32 * ,------------------+-----------------------------------------------.
33 * | |S8 U8 S16 U16 S24 U24 S32 U32 FLOAT FIXED-POINT|
34 * +------------------+-----------------------------------------------+
35 * | S8 |X X X X X X X X - - |
36 * | U8 |X X X X X X X X - - |
37 * I | S16 |X X X X X X X X - - |
38 * N | U16 |X X X X X X X X - - |
39 * P | S24 |X X X X X X X X - - |
40 * U | U24 |X X X X X X X X - - |
41 * T | S32 |X X X X X X X X - - |
42 * | U32 |X X X X X X X X - - |
43 * | FLOAT |X X X X X X X X X - |
44 * | FIXED-POINT |X X X X X X X X X - |
45 * `------------------+-----------------------------------------------'
46 */
47
48 #define SCALE(x,s) (s != 1.0 ? x * s : x)
49 #define MAXINT(a) (1L << ((a)-1))
50 #define CLIP(x,m) (x > m-1 ? m-1 : (x < -m ? -m : x))
51
52 /* private object */
53 typedef struct {
54 SAD_sample_format input_sample_format;
55 SAD_sample_format output_sample_format;
56 int input_bits;
57 int input_fracbits;
58 int output_bits;
59 int output_fracbits;
60 int channels;
61 SAD_channels_order input_chorder;
62 SAD_channels_order output_chorder;
63 SAD_get_sample_proc get_sample;
64 SAD_put_sample_proc put_sample;
65 int dither;
66 int hardlimit;
67 float scale;
68 float rg_scale;
69 } SAD_state_priv;
70
71 /* error code */
72
73 //static SAD_error SAD_last_error = SAD_ERROR_OK;
74
75 static inline double compute_hardlimit (double sample, float scale) {
76 sample *= scale;
77 const double k = 0.5; /* -6dBFS */
78 if (sample > k) {
79 return tanh((sample - k) / (1 - k)) * (1 - k) + k;
80 }
81 else if (sample < -k) {
82 return tanh((sample + k) / (1 - k)) * (1 - k) - k;
83 }
84 return sample;
85 }
86
87 /*
88 * Dither fixed-point normalized or integer sample to n-bits integer
89 * samples < -1 and > 1 will be clipped
90 */
91
92 static inline sad_sint32 __dither_sample_fixed_to_int (sad_sint32 sample, int inbits, int fracbits, int outbits, float scale, int dither,
93 int hardlimit)
94 {
95 int n_bits_to_loose, bitwidth, precision_loss;
96 sad_sint32 maxint = MAXINT(outbits);
97
98 n_bits_to_loose = 0;
99 bitwidth = inbits;
100 precision_loss = FALSE;
101
102 /*#ifdef DEEP_DEBUG
103 printf("f: __dither_sample_fixed_to_int\n");
104 #endif*/
105
106 if (fracbits == 0) {
107 if (inbits<29) {
108 /* convert to 4.28 fixed-point */
109 n_bits_to_loose = 29 - inbits;
110 sample <<= n_bits_to_loose;
111 bitwidth += n_bits_to_loose;
112 }
113
114 n_bits_to_loose += inbits - outbits;
115
116 if (inbits > outbits) {
117 precision_loss = TRUE;
118 #ifdef PRECISION_DEBUG
119 printf("Precision loss, reason: bitwidth loss %d --> %d\n", inbits, outbits);
120 #endif
121 }
122 } else {
123 n_bits_to_loose = fracbits + 1 - outbits;
124 bitwidth = fracbits;
125 precision_loss = TRUE;
126 #ifdef PRECISION_DEBUG
127 printf("Precision loss, reason: fixed-point input\n", inbits, outbits);
128 #endif
129 }
130
131 assert(n_bits_to_loose >=0 );
132
133 if (hardlimit) {
134 sample = (sad_sint32)(compute_hardlimit((double)sample/(double)MAXINT(bitwidth), scale) * (double)MAXINT(bitwidth));
135 #ifdef PRECISION_DEBUG
136 printf("Precision loss, reason: hard limiter\n", inbits, outbits);
137 #endif
138 precision_loss = TRUE;
139 } else {
140 sample = SCALE(sample, scale);
141 }
142
143 if (scale != 1.0){
144 precision_loss = TRUE;
145 #ifdef PRECISION_DEBUG
146 printf("Precision loss, reason: scale\n", inbits, outbits);
147 #endif
148 }
149
150 if (precision_loss && (n_bits_to_loose >= 1)) sample += (1L << (n_bits_to_loose - 1));
151
152 #ifdef DITHER_DEBUG
153 sad_sint32 val_wo_dither = sample, >> n_bits_to_loose;
154 val_wo_dither = CLIP(val_wo_dither, maxint);
155 #endif
156 if (dither && precision_loss && (n_bits_to_loose >= 1)) {
157 sad_sint32 dither_num = triangular_dither_noise(n_bits_to_loose + 1);
158 sample += dither_num;
159 }
160
161 sample >>= n_bits_to_loose;
162
163 /* Clipping */
164 #ifdef CLIPPING_DEBUG
165 sad_sint32 val_wo_clip = sample;
166 #endif
167 sample = CLIP(sample, maxint);
168 #ifdef CLIPPING_DEBUG
169 if (val_wo_clip != sample) {
170 printf("Clipping: %d --> %d\n", val_wo_clip, sample);
171 }
172 #endif
173 #ifdef DITHER_DEBUG
174 if (dither && precision_loss && (n_bits_to_loose >= 1)) printf("%d --> %d, noise: %d\n", val_wo_dither, sample, sample - val_wo_dither);
175 #endif
176 return sample;
177 }
178
179 /*
180 * Dither floating-point normalized sample to n-bits integer
181 * samples < -1 and > 1 will be clipped
182 */
183 static inline sad_sint32 __dither_sample_float_to_int (float sample, int nbits, float scale, int dither, int hardlimit) {
184
185 #ifdef DEEP_DEBUG
186 printf("f: __dither_sample_float_to_int\n");
187 #endif
188
189 sad_sint32 maxint = MAXINT(nbits);
190
191 if (hardlimit) {
192 sample = compute_hardlimit((double)sample, scale);
193 } else {
194 sample = SCALE(sample, scale);
195 }
196
197 sample *= maxint;
198 /* we want to round precisely */
199 sample = (sample < 0 ? sample - 0.5 : sample + 0.5);
200
201 #ifdef DITHER_DEBUG
202 sad_sint32 val_wo_dither = (sad_sint32) sample;
203 val_wo_dither = CLIP(val_wo_dither, maxint);
204 #endif
205 if (dither) {
206 float dither_num = triangular_dither_noise_f();
207 sample += dither_num;
208 }
209
210 /* Round and clipping */
211 sad_sint32 value = (sad_sint32) sample;
212 #ifdef CLIPPING_DEBUG
213 sad_sint32 val_wo_clip = value;
214 #endif
215 value = CLIP(value, maxint);
216 #ifdef CLIPPING_DEBUG
217 if (val_wo_clip != value) {
218 printf("Clipping: %d --> %d\n", val_wo_clip, value);
219 }
220 #endif
221
222 #ifdef DITHER_DEBUG
223 printf("%d --> %d, noise: %d\n", val_wo_dither, value, value - val_wo_dither);
224 #endif
225 return value;
226 }
227
228 static inline float __dither_sample_float_to_float (float sample, float scale, int hardlimit) {
229 #ifdef DEEP_DEBUG
230 printf("f: __dither_sample_float_to_float\n");
231 #endif
232 if (hardlimit) {
233 sample = compute_hardlimit((double)sample, scale);
234 } else {
235 sample = SCALE(sample, scale);
236 }
237 return sample;
238 }
239
240 static inline float __dither_sample_fixed_to_float (sad_sint32 sample, int inbits, int fracbits, float scale, int hardlimit) {
241 float fsample;
242
243 #ifdef DEEP_DEBUG
244 printf("f: __dither_sample_fixed_to_float\n");
245 #endif
246 if (fracbits == 0) {
247 fsample = (float)sample / (float)MAXINT(inbits);
248 } else {
249 fsample = (float)sample / (float)MAXINT(fracbits+1);
250 }
251 return __dither_sample_float_to_float (fsample, scale, hardlimit);
252 }
253
254
255
256
257
258 SAD_dither_t* SAD_dither_init(SAD_buffer_format *inbuf_format, SAD_buffer_format *outbuf_format, int *error) {
259 SAD_state_priv *priv;
260
261 DEBUG_MSG("f: SAD_dither_init\n",0);
262
263 priv = calloc(sizeof(SAD_state_priv), 1);
264
265 /* Check buffer formats and assign buffer ops */
266 SAD_buffer_ops* inops = SAD_assign_buf_ops(inbuf_format);
267
268 if (inbuf_format->sample_format != SAD_SAMPLE_FLOAT) {
269 if (inops != NULL) {
270 priv->get_sample = inops->get_sample;
271 } else {
272 free(priv);
273 *error = SAD_ERROR_INCORRECT_INPUT_SAMPLEFORMAT;
274 return NULL;
275 }
276 }
277
278 SAD_buffer_ops* outops = SAD_assign_buf_ops(outbuf_format);
279
280 if (outbuf_format->sample_format != SAD_SAMPLE_FLOAT) {
281 if (outops != NULL) {
282 priv->put_sample = outops->put_sample;
283 } else {
284 free(priv);
285 *error = SAD_ERROR_INCORRECT_OUTPUT_SAMPLEFORMAT;
286 return NULL;
287 }
288 }
289
290 priv->input_fracbits = 0;
291 priv->output_fracbits = 0;
292 priv->input_sample_format = inbuf_format->sample_format;
293 priv->output_sample_format = outbuf_format->sample_format;
294 priv->input_chorder = inbuf_format->channels_order;
295 priv->output_chorder = outbuf_format->channels_order;
296 priv->channels = inbuf_format->channels;
297 priv->scale = 1.0;
298 priv->rg_scale = 1.0;
299 priv->dither = TRUE;
300 priv->hardlimit = FALSE;
301
302 switch(outbuf_format->sample_format){
303 case SAD_SAMPLE_S8:
304 case SAD_SAMPLE_U8: priv->output_bits = 8; break;
305 case SAD_SAMPLE_S16:
306 case SAD_SAMPLE_U16: priv->output_bits = 16; break;
307 case SAD_SAMPLE_S24:
308 case SAD_SAMPLE_U24: priv->output_bits = 24; break;
309 case SAD_SAMPLE_S32:
310 case SAD_SAMPLE_U32: priv->output_bits = 32; break;
311 case SAD_SAMPLE_FLOAT: break;
312 default:
313 free(priv);
314 *error = SAD_ERROR_INCORRECT_OUTPUT_SAMPLEFORMAT;
315 return NULL;
316 }
317
318 switch(inbuf_format->sample_format){
319 case SAD_SAMPLE_S8:
320 case SAD_SAMPLE_U8: priv->input_bits = 8; break;
321 case SAD_SAMPLE_S16:
322 case SAD_SAMPLE_U16: priv->input_bits = 16; break;
323 case SAD_SAMPLE_S24:
324 case SAD_SAMPLE_U24: priv->input_bits = 24; break;
325 case SAD_SAMPLE_S32:
326 case SAD_SAMPLE_U32: priv->input_bits = 32; break;
327 case SAD_SAMPLE_FIXED32: priv->input_fracbits = inbuf_format->fracbits; break;
328 case SAD_SAMPLE_FLOAT: break;
329 default:
330 free(priv);
331 *error = SAD_ERROR_INCORRECT_INPUT_SAMPLEFORMAT;
332 return NULL;
333 }
334
335 *error = SAD_ERROR_OK;
336 return (SAD_dither_t*)priv;
337 }
338
339 int SAD_dither_free(SAD_dither_t* state) {
340 DEBUG_MSG("f: SAD_dither_free\n",0);
341 free(state);
342 return SAD_ERROR_OK;
343 }
344
345 /*
346 * Depend on format->channels_order inbuf and outbuf will be treated as
347 * smth* or smth** if channels_order = SAD_CHORDER_INTERLEAVED or SAD_CHORDER_SEPARATED
348 * accordingly
349 *
350 * frame is aggregate of format->channels samples
351 */
352
353 #define GET_FLOAT_SAMPLE(b,o,n,c,i) (o == SAD_CHORDER_INTERLEAVED ? (((float*)b)[i*n+c]) : (((float**)b)[c][i]))
354 #define PUT_FLOAT_SAMPLE(b,o,n,c,i,s) { \
355 if (o == SAD_CHORDER_INTERLEAVED) { \
356 ((float*)b)[i*n+c] = s; \
357 } else { \
358 ((float**)b)[c][i] = s; \
359 } \
360 }
361
362 int SAD_dither_process_buffer (SAD_dither_t *state, void *inbuf, void *outbuf, int frames)
363 {
364 SAD_state_priv *priv = (SAD_state_priv*) state;
365 int i, ch;
366 int channels = priv->channels;
367 int inbits = priv->input_bits;
368 int outbits = priv->output_bits;
369 int fracbits = priv->input_fracbits;
370 float scale = priv->scale * priv->rg_scale;
371 int dither = priv->dither;
372 int hardlimit = priv->hardlimit;
373 SAD_channels_order input_chorder = priv->input_chorder;
374 SAD_channels_order output_chorder = priv->output_chorder;
375
376 SAD_get_sample_proc get_sample = priv->get_sample;
377 SAD_put_sample_proc put_sample = priv->put_sample;
378
379 #ifdef DEEP_DEBUG
380 printf("f: SAD_process_buffer\n");
381 #endif
382
383 if (priv->input_sample_format == SAD_SAMPLE_FLOAT) {
384 if (priv->output_sample_format == SAD_SAMPLE_FLOAT) {
385 /* process buffer */
386 for(i=0; i<frames; i++) {
387 for(ch=0; ch<channels; ch++) {
388 float sample = GET_FLOAT_SAMPLE(inbuf, input_chorder, channels, ch ,i);
389 sample = __dither_sample_float_to_float(sample, scale, hardlimit);
390 PUT_FLOAT_SAMPLE(outbuf, output_chorder, channels, ch ,i, sample);
391 }
392 }
393 } else {
394 if (put_sample == NULL) return SAD_ERROR_CORRUPTED_PRIVATE_DATA;
395 /* process buffer */
396 for(i=0; i<frames; i++) {
397 for(ch=0; ch<channels; ch++) {
398 float sample = GET_FLOAT_SAMPLE(inbuf, input_chorder, channels, ch ,i);
399 sad_sint32 isample = __dither_sample_float_to_int(sample, outbits, scale, dither, hardlimit);
400 put_sample (outbuf, isample, channels, ch, i);
401 }
402 }
403 }
404 } else {
405 if (priv->output_sample_format == SAD_SAMPLE_FLOAT) {
406 if (get_sample == NULL) return SAD_ERROR_CORRUPTED_PRIVATE_DATA;
407 /* process buffer */
408 for(i=0; i<frames; i++) {
409 for(ch=0; ch<channels; ch++) {
410 sad_sint32 sample = get_sample (inbuf, channels, ch, i);
411 float fsample = __dither_sample_fixed_to_float (sample, inbits, fracbits, scale, hardlimit);
412 PUT_FLOAT_SAMPLE(outbuf, output_chorder, channels, ch ,i, fsample);
413 }
414 }
415 } else {
416 if (put_sample == NULL || get_sample == NULL) return SAD_ERROR_CORRUPTED_PRIVATE_DATA;
417 /* process buffer */
418 for(i=0; i<frames; i++) {
419 for(ch=0; ch<channels; ch++){
420 sad_sint32 sample = get_sample (inbuf, channels, ch, i);
421 sad_sint32 isample = __dither_sample_fixed_to_int (sample, inbits, fracbits, outbits, scale, dither, hardlimit);
422 put_sample (outbuf, isample, channels, ch, i);
423 }
424 }
425 }
426 }
427
428 return SAD_ERROR_OK;
429 }
430
431 int SAD_dither_apply_replaygain (SAD_dither_t *state, SAD_replaygain_info *rg_info, SAD_replaygain_mode *mode) {
432 SAD_state_priv *priv = (SAD_state_priv*) state;
433 float scale = -1.0, peak = 0.0;
434
435 DEBUG_MSG("f: SAD_dither_apply_replaygain\n",0);
436
437 if(!rg_info->present) {
438 priv->rg_scale = 1.0;
439 priv->hardlimit = FALSE;
440 return SAD_ERROR_OK;
441 }
442
443 switch(mode->mode) {
444 case SAD_RG_ALBUM:
445 scale = db2scale(rg_info->album_gain);
446 peak = rg_info->album_peak;
447 if (peak == 0.0) {
448 scale = db2scale(rg_info->track_gain); // fallback to per-track mode
449 peak = rg_info->track_peak;
450 DEBUG_MSG("f: SAD_dither_apply_replaygain: fallback to track mode\n",0);
451 }
452 break;
453 case SAD_RG_TRACK:
454 scale = db2scale(rg_info->track_gain);
455 peak = rg_info->track_peak;
456 break;
457 case SAD_RG_NONE:
458 scale = -1.0;
459 }
460
461 if (scale != -1.0 && peak != 0.0) {
462 DEBUG_MSG("f: SAD_dither_apply_replaygain: applying\n",0);
463 scale *= mode->preamp;
464 // Clipping prevention
465 if(mode->clipping_prevention) {
466 #ifdef DEBUG
467 if(scale * peak > 1.0) DEBUG_MSG("f: SAD_dither_apply_replaygain: clipping prevented\n",0);
468 #endif
469 scale = scale * peak > 1.0 ? 1.0 / peak : scale;
470 }
471 scale = scale > 15.0 ? 15.0 : scale; // safety
472 priv->rg_scale = scale;
473 priv->hardlimit = mode->hard_limit; // apply settings
474 } else {
475 priv->rg_scale = 1.0;
476 priv->hardlimit = FALSE;
477 }
478
479 return SAD_ERROR_OK;
480 }
481
482 int SAD_dither_set_scale (SAD_dither_t *state, float scale) {
483 SAD_state_priv *priv = (SAD_state_priv*) state;
484 priv->scale = scale;
485 return SAD_ERROR_OK;
486 }
487
488 int SAD_dither_set_dither (SAD_dither_t *state, int dither) {
489 SAD_state_priv *priv = (SAD_state_priv*) state;
490 priv->dither = dither;
491 return SAD_ERROR_OK;
492 }