1 
/* 
2 
* $Id: pa_dither.c 1418 20091012 21:00:53Z philburk $ 
3 
* Portable Audio I/O Library triangular dither generator 
4 
* 
5 
* Based on the Open Source API proposed by Ross Bencina 
6 
* Copyright (c) 19992002 Phil Burk, Ross Bencina 
7 
* 
8 
* Permission is hereby granted, free of charge, to any person obtaining 
9 
* a copy of this software and associated documentation files 
10 
* (the "Software"), to deal in the Software without restriction, 
11 
* including without limitation the rights to use, copy, modify, merge, 
12 
* publish, distribute, sublicense, and/or sell copies of the Software, 
13 
* and to permit persons to whom the Software is furnished to do so, 
14 
* subject to the following conditions: 
15 
* 
16 
* The above copyright notice and this permission notice shall be 
17 
* included in all copies or substantial portions of the Software. 
18 
* 
19 
* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, 
20 
* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF 
21 
* MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. 
22 
* IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR 
23 
* ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF 
24 
* CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION 
25 
* WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. 
26 
*/ 
27 

28 
/* 
29 
* The text above constitutes the entire PortAudio license; however, 
30 
* the PortAudio community also makes the following nonbinding requests: 
31 
* 
32 
* Any person wishing to distribute modifications to the Software is 
33 
* requested to send the modifications to the original developer so that 
34 
* they can be incorporated into the canonical version. It is also 
35 
* requested that these nonbinding requests be included along with the 
36 
* license above. 
37 
*/ 
38 

39 
/** @file 
40 
@ingroup common_src 
41 

42 
@brief Functions for generating dither noise 
43 
*/ 
44 

45 
#include "pa_types.h" 
46 
#include "pa_dither.h" 
47 

48 

49 
/* Note that the linear congruential algorithm requires 32 bit integers 
50 
* because it uses arithmetic overflow. So use PaUint32 instead of 
51 
* unsigned long so it will work on 64 bit systems. 
52 
*/ 
53 

54 
#define PA_DITHER_BITS_ (15) 
55 

56 

57 
void PaUtil_InitializeTriangularDitherState( PaUtilTriangularDitherGenerator *state ) 
58 
{ 
59 
state>previous = 0; 
60 
state>randSeed1 = 22222; 
61 
state>randSeed2 = 5555555; 
62 
} 
63 

64 

65 
PaInt32 PaUtil_Generate16BitTriangularDither( PaUtilTriangularDitherGenerator *state ) 
66 
{ 
67 
PaInt32 current, highPass; 
68 

69 
/* Generate two random numbers. */ 
70 
state>randSeed1 = (state>randSeed1 * 196314165) + 907633515; 
71 
state>randSeed2 = (state>randSeed2 * 196314165) + 907633515; 
72 

73 
/* Generate triangular distribution about 0. 
74 
* Shift before adding to prevent overflow which would skew the distribution. 
75 
* Also shift an extra bit for the high pass filter. 
76 
*/ 
77 
#define DITHER_SHIFT_ ((sizeof(PaInt32)*8  PA_DITHER_BITS_) + 1) 
78 

79 
current = (((PaInt32)state>randSeed1)>>DITHER_SHIFT_) + 
80 
(((PaInt32)state>randSeed2)>>DITHER_SHIFT_); 
81 

82 
/* High pass filter to reduce audibility. */ 
83 
highPass = current  state>previous; 
84 
state>previous = current; 
85 
return highPass; 
86 
} 
87 

88 

89 
/* Multiply by PA_FLOAT_DITHER_SCALE_ to get a float between 2.0 and +1.99999 */ 
90 
#define PA_FLOAT_DITHER_SCALE_ (1.0f / ((1<<PA_DITHER_BITS_)1)) 
91 
static const float const_float_dither_scale_ = PA_FLOAT_DITHER_SCALE_; 
92 

93 
float PaUtil_GenerateFloatTriangularDither( PaUtilTriangularDitherGenerator *state ) 
94 
{ 
95 
PaInt32 current, highPass; 
96 

97 
/* Generate two random numbers. */ 
98 
state>randSeed1 = (state>randSeed1 * 196314165) + 907633515; 
99 
state>randSeed2 = (state>randSeed2 * 196314165) + 907633515; 
100 

101 
/* Generate triangular distribution about 0. 
102 
* Shift before adding to prevent overflow which would skew the distribution. 
103 
* Also shift an extra bit for the high pass filter. 
104 
*/ 
105 
current = (((PaInt32)state>randSeed1)>>DITHER_SHIFT_) + 
106 
(((PaInt32)state>randSeed2)>>DITHER_SHIFT_); 
107 

108 
/* High pass filter to reduce audibility. */ 
109 
highPass = current  state>previous; 
110 
state>previous = current; 
111 
return ((float)highPass) * const_float_dither_scale_; 
112 
} 
113 

114 

115 
/* 
116 
The following alternate dither algorithms (from musicdsp.org) could be 
117 
considered 
118 
*/ 
119 

120 
/*Noise shaped dither (March 2000) 
121 
 
122 

123 
This is a simple implementation of highpass triangularPDF dither with 
124 
2ndorder noise shaping, for use when truncating floating point audio 
125 
data to fixed point. 
126 

127 
The noise shaping lowers the noise floor by 11dB below 5kHz (@ 44100Hz 
128 
sample rate) compared to triangularPDF dither. The code below assumes 
129 
input data is in the range +1 to 1 and doesn't check for overloads! 
130 

131 
To save time when generating dither for multiple channels you can do 
132 
things like this: r3=(r1 & 0x7F)<<8; instead of calling rand() again. 
133 

134 

135 

136 
int r1, r2; //rectangularPDF random numbers 
137 
float s1, s2; //error feedback buffers 
138 
float s = 0.5f; //set to 0.0f for no noise shaping 
139 
float w = pow(2.0,bits1); //word length (usually bits=16) 
140 
float wi= 1.0f/w; 
141 
float d = wi / RAND_MAX; //dither amplitude (2 lsb) 
142 
float o = wi * 0.5f; //remove dc offset 
143 
float in, tmp; 
144 
int out; 
145 

146 

147 
//for each sample... 
148 

149 
r2=r1; //can make HPTRI dither by 
150 
r1=rand(); //subtracting previous rand() 
151 

152 
in += s * (s1 + s1  s2); //error feedback 
153 
tmp = in + o + d * (float)(r1  r2); //dc offset and dither 
154 

155 
out = (int)(w * tmp); //truncate downwards 
156 
if(tmp<0.0f) out; //this is faster than floor() 
157 

158 
s2 = s1; 
159 
s1 = in  wi * (float)out; //error 
160 

161 

162 

163 
 
164 
paul.kellett@maxim.abel.co.uk 
165 
http://www.maxim.abel.co.uk 
166 
*/ 
167 

168 

169 
/* 
170 
16to8bit firstorder dither 
171 

172 
Type : First order error feedforward dithering code 
173 
References : Posted by Jon Watte 
174 

175 
Notes : 
176 
This is about as simple a dithering algorithm as you can implement, but it's 
177 
likely to sound better than just truncating to N bits. 
178 

179 
Note that you might not want to carry forward the full difference for infinity. 
180 
It's probably likely that the worst performance hit comes from the saturation 
181 
conditionals, which can be avoided with appropriate instructions on many DSPs 
182 
and integer SIMD type instructions, or CMOV. 
183 

184 
Last, if sound quality is paramount (such as when going from > 16 bits to 16 
185 
bits) you probably want to use a higherorder dither function found elsewhere 
186 
on this site. 
187 

188 

189 
Code : 
190 
// This code will downconvert and dither a 16bit signed short 
191 
// mono signal into an 8bit unsigned char signal, using a first 
192 
// order forwardfeeding error term dither. 
193 

194 
#define uchar unsigned char 
195 

196 
void dither_one_channel_16_to_8( short * input, uchar * output, int count, int * memory ) 
197 
{ 
198 
int m = *memory; 
199 
while( count > 0 ) { 
200 
int i = *input++; 
201 
i += m; 
202 
int j = i + 32768  128; 
203 
uchar o; 
204 
if( j < 0 ) { 
205 
o = 0; 
206 
} 
207 
else if( j > 65535 ) { 
208 
o = 255; 
209 
} 
210 
else { 
211 
o = (uchar)((j>>8)&0xff); 
212 
} 
213 
m = ((j32768+128)i); 
214 
*output++ = o; 
215 
} 
216 
*memory = m; 
217 
} 
218 
*/ 