A quantiser is a non-smooth function and no inverse function exist that can be applied to correct for the error it introduces. Applying a dither signal and averaging to a quantiser produces a smooth image for which an inverse function does exist. This article describes methods for minimising the error after inverse compensation. We show that there is an optimal dither variance that minimises the error after inversion. Simple rules for choosing the optimal dither variance are presented. The error after inversion can be made arbitrarily small by increasing the averaging length. This can be done by oversampling the signal by the same factor as the number of averages. Quantisation of a dither signal with a continuous probability distribution, produces a discrete probability mass function. We discuss a method for recovering an unknown continuous probability distribution from the empirical discrete probability mass function of the quantised dither signal. This enables inverse compensation in systems where exact control of the dither signal is not possible, as inverse compensation requires information about the continuous probability distribution of the dither signal and the step-size of the quantiser. Keywords: Analog digital conversion, Quantisation error, Nonlinear distortion, Dither, Averaging, Gaussian noise 1. Introduction Quantisation is the process of mapping a large set of values to a smaller set of values; thereby discarding some values and introducing a signal dependent error, as illus- trated in Fig. 1. However, by applying a dithering signal as illustrated in Fig. 2, a bijective correspondence between the expected value of the quantised signal and the refer- ence signal can be obtained; that is, Efyg=N(x), where y=y(t) is the quantized version of x=x(t), as illustrated in Fig. 4 (see also Fig. 5) By applying averaging, a resolution below the step-size of the dithered quantiser can then be achieved [1]. Im- proving the resolution in this manner is particularly useful when a signal is sampled using coarse quantisation or in applications where a large dynamic range is required. In general, the function N(x) is non-linear, and de- termined by the distribution of the dither signal and the step-size of the quantiser. If this functional relationship between the reference and the expected value of the quan- tised signal is known, the error in the expected value in- troduced by N(x) compared to xcan be compensated by using the inverse of the function [2]. When the probability distribution for a dither signal ful ls certain criteria, the moments of the quantised signal Corresponding author. Email addresses: eielsen@ux.uis.no (Arn nn A. Eielsen), jjl@es.aau.dk (John Leth), andrew.fleming@newcastle.edu.au (Andrew J. Fleming)can be made independent of the reference signal [3, 4, 5, 6, 7, 8]. In terms of the rst moment, the expected value, the quantiser can then be perfectly linearised by applying a uniformly distributed dither signal with a range that is exactly=2, where  denotes the step-size of the quan- tiser [7]. However, dither distributions that render the moments of the quantised signal independent of the reference, e.g. a uniform distribution, are infeasible to produce in a physical system as they require a high degree of precision in the dither signal reproduction. This is because any dynamics in the signal path, by design or due to parasitic e ects, will cause the reproduced dither signal to tend to a normal distribution. This e ect is explained by the central limit theorem [9]. If the realisation of the distribution is not exact, the linearising e ect is signi cantly reduced [10]. In contrast, since introducing dynamics in a signal path can be done using e.g. linear lters, a normally distributed dither signal is almost trivial to produce. In the normally distributed case, the quantiser can only be perfectly linearised in terms of the expected value when the standard deviation, or variance, of the dither is in nite, as shown in Sec. 3.4. When applying averaging, consider- ing the combined error introduced by N(x) (a determinis- tic e ect) as well as the variance of the dither (a stochastic e ect), the minimum error of the dithered and averaged quantiser output signal is obtained when the standard de- viation is approximately  =2 [10, 11]. This error is il- lustrated in Fig. 5 using a non-optimal dither variance to Preprint submitted to Signal Processing December 28, 2019 Figure 1: A quantised signal will have a signal dependent error, referred to as the quantisation error: y(t)x(t). Q + G(z)x(t)d(t) u(t)y(t)y[n] y[n] Figure 2: Quantiser dithering and averaging using an LTI- lter. emphasise the e ect of N(x). As noted above, if the dither distribution parameters and quantiser step-size are known, then applying the in- verse of the non-linear function N1(^y) to the output of the dithered and averaged quantiser reduces the overall error [2]. This is illustrated in Fig. 6. By using the inverse, it may be intuitively expected that the standard deviation of the dither can be reduced signi cantly below  =2, asN(x) becomes bijective for any standard deviation greater than zero. As the results in this paper show, this is not the case: The derivative of N(x) approaches zero in the midpoints between the steps when the dither variance is reduced, implying that the derivative ofN1(^y) at these points approaches in nity. A large derivative, or sensitivity, for N1(^y) ampli es noise, and this noise contributes to the overall error more than the removal of the systematic error due to N(x) can reduce it. An optimal dither variance that minimises the overall error in the inverted signal therefore exists. 1.1. Contributions The main contributions are an analysis of the overall error in the signal when applying N1(^y) and a method for identifying the parameters of the probability distribution using the quantised dither signal. The error analysis makes it possible to determine the performance limits when ap- plying the inverse, and the identi cation method enables the application of the inverse to systems where a dither with a known type of probability distribution can be in- troduced, but where the exact distribution parameters are unknown. The functional relationship between the reference sig- nal and the expected value of the quantised signal as well as the expression for the variance of the output signal is developed directly using convolution integrals rather thanvia characteristic functions, such as in [10, 2, 11]. This ap- proach makes it simpler to prove convergence to a linear function, as compared to the approaches in [3, 4, 5, 7, 8]. It is demonstrated that when using N1(^y) in combi- nation with sucient averaging, it is possible to reduce the error compared to the case when not using an inverse. It is also demonstrated that there exists an optimal dither vari- ance, and that the optimal dither depends on assumptions made about the reference signal. 1.2. Notation The set of natural numbers f0;1;2;:::gis denoted N. A de nition is denoted by ,andindicate the convo- lution product. The Laplace operator is denoted Land the Z-transformZ. Functions of time t(or discrete-time n) are usually denoted by lower case, e.g. g(t) (org[n]), and the Laplace transform (or Z-transform) by upper case, e.g.G(s) =Lfgg(s) (orG[z] =Zfgg[z]). The standard notationLp(R) andlpindicate the Lp-space and lp-space p= 1;2;:::;1andkgkpthep-norm ofg(t) org[n]. For a stochastic process d(t;!) the dependency on the sam- ple variable !is omitted, = Efd(t)gdenote the mean value,Cd(t;s),Efd(t)d(s)gEfd(t)gEfd(s)gdenote the auto-covariance and Cd(t),Cd(t;t) the variance. The indicator function for the set Sis denoted S(z) and the oor operator by bc. 2. The uniform quantiser A uniform quantiser can represent values separated by the quantisation step-size  >0 and can be de ned using the truncation operator trunc( u) as trunc(u),u +1 2 2Z: (1) The quantised output y=y(t) given input u=u(t) is then y=n(u), trunc(u): (2) The uniform quantiser is demonstrated in Fig. 1. 2.1. Alternative formulation The uniform quantiser can equivalently be expressed as y=n(u) =X k2Nnk(u) (3) where each (odd) step-function, as shown in Fig.3, nk(u),n+ k(u)n k(u) (4a) in the quantiser, is described by the step functions n+ k(u),H(uTk); (4b) n k(u),H(uTk); (4c) 2 Tk Tk  unk(u) Figure 3: Step-function. Figure 4: A quantiser is not bijective as it maps the continuum to a discrete set, but the smoothed quantiser (17) is; thus giving a one- to-one correspondence between reference and expected value. Here a zero-mean, normally distributed dither (with d= =6) has been used. with the threshold of the step Tkgiven by Tk, (k+ 1=2) (5) andH(u), [0;1)(u) the Heaviside step-function. Note that the sum in (3) converges for every tsince the sum- mands become zero when TkuTkor equivalently whenkjuj=1=2. 3. Smoothing the quantiser The input signal uwhen a stochastic dithering signal is present is u(t),x(t) +d(t); (6) wherex(t) is the deterministic reference signal that should be recovered from a quantised and averaged time-series. It is assumed that the stochastic dither signal d(t) is an iden- tically distributed stochastic process having the strictly white noise property, i.e. d(ti) andd(tj) are independent for each pair tiandtj[12]. In the sequel the input uis viewed as a function of xandtviz., u=u(x;t) =x+d(t); when the time dependency of xis irrelevant. If this is not the case we write u(t) =u(x(t);t). The dither signal d(t) is either already present natu- rally due to sampling a noisy measurement, or it is added arti cially. In either case, the e ect is that by averag- ing the quantised signal y(x;t) =n(u(x;t)), the e ect ofthe discontinuous quantiser can be smoothed, imbuing the quantiser with continuity, as illustrated in Fig. 5. How- ever, apart from the smoothing e ect, the dither signal is otherwise unwanted in the recovered signal. LetFd(z) be the cumulative distribution function (CDF) andfd(z) the probability density function (PDF) corre- sponding to d(t), and note that neither of these functions depend on time tby the i.d., property of d(t). Moreover, let() denotes the Dirac delta function and x2R. Then fx(z;x) =(zx) is the PDF corresponding to xwhich depends on time tifx=x(t) does;fx(z;t),fx(z;x(t)). The signal u(x;t) will then have the PDF fu(z;x) =Z Rfx(zw;x)fd(w) dw=fd(zx):(7) withfu(z;t),fu(z;x(t)) whenx=x(t). The application of the dither signal makes it possible to de ne an averaged and smoothed step-function Nk(z),Z Rnk(z+w)fd(w) dw; (8) related tonk(u) and the stochastic dither signal d(t). Note that in general, when nk(z) is any arbitrary func- tion of bounded variation, kNkk1kfdk1knkk1=knkk1 (9) sincenk(u)2L1(R), and that if LD,supz2Rjfd(z)j< 1then by [13, Lemma A.1] LNkLDTV(nk) (10) withLNka Lipshitz constant of Nk(z) andTV(nk) the total variation of nk(u). In the case considered in this arti- cleTV(nk) = 2, hence the dither PDF and quantisation step-size determines an upper bound on the smoothing ef- fect on the quantiser. Moreover, using (4) Nk(z) =  Z1 Tkzfd(w)dwZTkz 1fd(w)dw! (11a) =  (1Fd(Tkz)Fd(Tkz)): (11b) It will be assumed throughout that the CDF Fdis such that N(z),X k2NNk(z) (12a) =X n2Z1 2Fd(Tnz) (12b) is nite for each z. The summand kin (12a) correspond to the two summands n=kandn=(k+ 1) in (12b). In the special case where fd(z) is symmetric around zero1(11) can be written as Nk(z) = ZTk+z Tkzfd(w)dw (13a) = (Fd(Tk+z)Fd(Tkz)): (13b) 1fd(+x) =fd(x) with=Efd(t)g= 0 3 Figure 5: Applying dither makes is possible to reduce the quanti- sation error if the output is averaged. Here a zero-mean, normally distributed dither (with d= =6) has been used. The averaged signal was obtained with M= 1000 samples. Figure 6: The inverse of the non-linear function relating the reference and expected value can be used to reduce the quantisation error of the averaged signal. Using a zero-mean, normally distributed dither with (d= =6), the averaged signal ^ y=hy(t)iMwas obtained with M= 1000 samples. Hence sign ( Nk(z)) = sign(z). 3.1. First and Second Order Moments Consider the strictly white noise output for each step yk(x;t),nk(u(x;t)); (14) and the (strictly white noise) quantised output signal y(x;t),X k2Nyk(x;t); (15) withyk(t),yk(x(t);t) andy(t),y(x(t);t) ifx=x(t). Using (8), the expected value for each step is Efyk(x;t)g=Z Rnk(u)fd(ux) du (16a) =Z Rnk(x+w)fd(w) dw (16b) =Nk(x): (16c) The expected value of the quantiser output is therefore Efy(x;t)g=X k2NNk(x) =N(x): (17)With (17) and x=x(t), this leads to an output auto- covariance given by Cy(t;s) = Efy(t)y(s)gEfy(t)gEfy(s)g (18a) =( 0 t6=s Efy2(t)gN2(x(t))t=s(18b) ,v(t) ftg(s) (18c) where the rst case in (18b) follows by the strictly white noise property of y(t) and the second case is by (17). It is remarked that the right-hand side of (16) and (17) depend on time only when x=x(t), implying the same for the left-hand sides. That is, the time evolution of the stochastic dither does not a ect the expected value of the output. This is in general not the case for the auto- covariance (18). However, as shown below the variance v=v(t) can be viewed as a function of xwhen time de- pendency of x(t) is irrelevant. Consider the zero-mean strictly white noise error signal "(x;t) =y(x;t)N(x); (19) having, with x=x(t), auto-covariance C"(t;s) =Cy(t;s). In order to nd the time-varying variance C"(t) =Cy(t) =v(t) (20) of the error signal "(t),"(x(t);t), the output of the terms (4) in the summation (3) is considered individually as correlated variables. The variance and covariance of the terms can then be found individually and then summed using the formula v(t) =X k2Zvar(yk(t)) + 2X k2ZX l>kcov(yk(t);yl(t));(21) to nd the variance of the error. For this note rst that the Heaviside step-function has the properties H2(u) =H(u); (22) and (more generally) H(ua)H(ub) =H(ub) (23) ifab. The variance for yk(t) can then be found as var(yk(t)) ,Efy2 k(t)gNk(x(t))2 = 2 1Fd(Tkx(t)) +Fd(Tkx(t)) Nk(x(t))2(24a) = 22Fd(Tkx(t)) + Nk(x(t)) Nk(x(t))2(24b) using (22) to obtain (24a), and (11b) to obtain (24b). For k<l , the covariance between terms can be found as cov(yk(t);yl(t)),Efyk(t)yl(t)gEfyk(t)gEfyl(t)g = 2 1Fd(Tlx(t)) +Fd(Tlx(t)) Nk(x(t))Nl(x(t)) (25) 4 Figure 7: The normalised standard deviation for the normally dis- tributed case with d= = 1=4. using (23). Note that when the time dependency of x(t) is irrelevant then, from (24) and (25), the variance vmay be viewed as a function of xviz.,v=v(x) with the special casev(t) =v(x(t)). For this reason (20) is generalised to C"(x) =Cy(x) =v(x) (26) with (20) obtained from (26) by letting x=x(t). An example of the variance v(x) and its dependence on xis shown in Fig. 7. 3.2. Convergence to a linear function We show thatd dzN(z) is constant, and therefore that N(z) is ane linear, when the characteristic function 'd(v) of the PDF fd(z) has a speci c property (see (31)). From (12b) we get d dzN(z) =X n2Zfd(Tnz) (27a) =X n2Zfd(=2z+n); (27b) which by the Poisson sum formula [9] yields: d dzN(z) =X n2Z^fd 2n  ei2n (z=2)(28) =X n2Z'd2n  ei2n (z=2)(29) = 1 +X n2Znf0g'd2n  ei2n (z=2)(30) It can then be seen that if 'd2n  = 08n2Znf0g; (31) then the derivative of N(z) is constant, and hence N(z) is ane linear, that is, when (31) is ful lled, then the quantisation error ex(x(t)),N(x(t))x(t) =N(0); (32)is constant. Note that the quantisation error (32) is in terms of the reference x(t), whereas the output error (19), withx=x(t) is in terms of the output y(t). IfN(x(t)) intersects the origin, N(0) = 0. The above derivation leads to the same well-known result as found in [3, 4, 5, 7, 8]. Below are two examples of distributions that can ful ll the criterion (31) and in principle linearise the quantiser and remove the quantisation error. As will be seen, neither have practical use in measurement systems. 3.3. Application of a uniformly distributed dither The simplest distribution that can ful l (31) is the uni- form distribution [7, 8]. In this case fd(z) =( 1 ;jzj 2 0;jzj> 2(33) which has characteristic function 'd(v) =eiv 2eiv 2 iv=sin 2v  2v; hence 'd2n  =sin(n) n= sinc (n) = 08n2Znf0g: If the dither signal has a uniform distribution, the dither signal reconstruction must be exact in order to linearise the quantiser. Engineering a device that can reproduce a physical white noise signal with a uniform distribution is in general infeasible due to the precision requirements. Hence, the uniformly distributed case has little practical use in most physical measurement systems. 3.4. Application of a normally distributed dither Producing a normally distributed signal is compara- tively trivial compared to a uniformly distributed white noise signal. If the process dis normally distributed with meandand standard deviation d, then fd(z) =1p 2de(zd)2 22 d (34) 'd(v) = eidve1 22 dv2; (35) hence 'd2n  = eid2n e2 d22n2 2; (36) andN(z) will therefore approach an ane linear function asd! 1 , i.e. this distribution only ful lls (31) for an in nite standard deviation. An in nite variance is of course impractical, and in general, dithering with a large variance has limited practical application. 5 4. Averaging 4.1. Sample-averaging Considering (18b) and (19), the output y(x;t) can be seen to be a linear combination of the e ective non-linearity N(x), and the white noise error term "(x;t), that is y(x;t) =N(x) +"(x;t) = Efy(x;t)g+"(x;t):(37) The expression Efy(x;t)g=N(x) represents the expected value in terms of an average over the sample space. We can in principle nd this average as follows. Fori= 1;2;:::;M letdi(t) denote a stochastic dither signal and assume that d1(t);d2(t);:::;dM(t) are indepen- dent and identically distributed (i.i.d.) for each xed time t. As above we construct, for each i, the (strictly white noise) quantised output signal yi(x;t),X k2Nyi k(x;t),X k2Nnk(x+di(t)): (38) Hencey(x;t);y1(x;t);:::;yM(x;t) are i.i.d., and the sam- ple average hy(x;t)iM=1 MMX i=1yi(x;t); (39) converges almost surely to E fy(x;t)g=N(x) by the law of large numbers. That is, hy(x)i=N(x) withhy(x)i,limM!1hy(x;t)iM. Note that by direct calculations the variance of hy(x;t)iMis ChyiM(x) =1 MCy(x) =1 Mv(x); (40) which converge to zero, implying that hy(x;t)iMconverge toEfy(x;t)g=N(x) in the mean square sense since Ef(hy(x;t)iMEfy(x;t)g)2g =Ef(hy(x;t)iMEfhy(x;t)iMg)2g=ChyiM(x):(41) Moreover, for h"(x;t)iM,1 MMX i=1"i(x;t) =hy(x;t)iMEfy(x;t)g;(42) with"i(x;t),yi(x;t)Efyi(x;t)g=yi(x;t)Efy(x;t)g, we get h"i,lim M!1h"(x;t)iM= 0; and Ch"iM(x) =1 MC"(x) =1 MCy(x); (43) with the time dependent case x=x(t) as a special case.The sample averaging is often not implementable in applications since, in general, a large value of Mis re- quired (that is, a large number of physical channels is required). A common con guration in measurement sys- tems is to have a single channel. Hence, time-averaging may therefore be the only option for producing an average value. 4.2. Time-averaging Time-averaging the output of the quantiser is in prac- tice done using discrete time-samples of y(t), at the in- stancestn= nwhere is the sampling-time. The sam- ples will then be averaged using a linear time-invariant (LTI) low-pass lter. The discrete-time output is denoted y[n],y(tn), and the LTI- lter impulse response is de- notedg[n]. Denoting the time-averaged output  y[n], (gy)[n], we have y[n] = (gN(x))[n] + (g")[n] (44a) =N(x[n]) + (hN(x))[n] + (g")[n] (44b) whereh[n] =g[n][n]. Thez-domain representation is then y(z) =G(z)y(z) (45a) =ZfN(x)g(z) +H(z)ZfN(x)g(z) +G(z)E(z); (45b) whereH(z) =G(z)1. The lter H(z) is used to describe the error introduced by the time-averaging operation. Since the error signal (19) is not stationary, common systems norms [14] do not apply. However, an LTI- lter will still reduce the variance of the error signal "(t) as will be shown next. First, the ltered error signal is de ned as "[n],(g")[n]; (46) whereg[n] is the output lter and "[n] is the error signal. The variance of the ltered error signal is then C"[n;n] = Ef"2[n]g (47a) = (vg2)[n]; (47b) where (47b) follows by the discrete-time version of [12, Theorem 9-3]. Hence if v2l1andg22l2then kC"k2=kvg2k2 (48a) kvk1kg2k2 (48b) where (48b) follows by Young's convolution theorem for sequences. By (48) the variance of the ltered error signal will be attenuated by the energy of the impulse response squared g2 2. For strictly proper, unity gain low-pass l- ters, the energy of the impulse response is directly related to the cut-o frequency: A lower cut-o frequency yields lower energy, hence smaller variance for the ltered error signal ". 6 4.2.1. Moving average lter Applying a moving average lter g[n] =1 MM1X i=0[ni]; gives y[n] = (gy)[n] =1 MM1X i=0y[ni]: The variance of the ltered error signal will then be Ef"2[n]g=1 M2M1X i=0v[ni] (49) according to (47), and it will be attenuated by kg2k2=1 M; (50) according to (48). Note that a moving average lter in the z-domain is a uniformly weighted sum of delay elements, which can simpli ed using the geometric series sum formula to yield G(z) =1 MM1X i=0zi=1 M1zM 1z1: (51) By setting z= ei2fthe frequency response of G(z) can be found as G(f) =sinc(Mf) sinc(f)ei2(M 21)f: (52) The magnitude response is shown in Fig. 8. The frequency response can be used to compensate for some of the e ects the moving average lter has on the reference signal x[n], as explained in Sec. 5.2. 4.2.2. Accumulate-and-dump lter A large number of samples is desireable in order to reduce the variace due to the dither signal. However, in terms of implementability, it may be necceary to decimate, or downsample, the output of the moving average lter in order to reduce the sampling rate. If the moving lter with Msamples is subsequently downsampled by the same fac- torM, it is known as an accumulate-and-dump lter. This operation is attractive due to its simplicity. The e ect of the decimation can be analysed by considering the inter- polation formula y0(t) =X n2Zy[n]sinc(tn): Decimating  y[n] is achieved by evaluating  y0[Mn],y0(Mn).4.3. Variance for Small  If it can be assumed that  is small in the sense that x(t) varies approximately linearly over the quantisation step-size, the variance v(x) can be approximated by an average [2, 11]; the limiting case being that xvaries linearly across one step. In the sample-averaged case the average variance is de- ned as v,1 Z=2 =2v(x) dx (53) and the variance Ch"iM(x) from (43) can be approximated by Ch"iM(x)1 Mv: (54) using (40). In the time-averaged case the average variance is sim- ilarly de ned as vd=1 MM1X k=0v  2+ Mk (55) and, assuming a moving average lter, the variance of  "[n] will then be Ef"2[n]g1 Mvd; (56) by (49). Note that if  =Min (55) is interpreted as the sampling- time then vvdfor suciently small sampling-time (or equivalently for suciently large M). An example of the average variance is shown in Fig. 7, and compared to the case when it is dependent on x. By assuming that the vari- ance is constant, the process describing the error becomes stationary and information about the worst-case error vari- ance is neglected. 5. Inverting the E ective Non-linearity Consider again the output of the quantiser (37) with x=x(t). If there was no error signal, "(t) = 0, then the reference signal x(t) could directly be recovered by application of the inverse of N(z), that is x(t) =N1(y(t)): (57) Clearly (57) can not be true since we always have "(t)6= 0 due to the presence of the dither signal. However, we may obtain an estimate of x(t) as follows. Let ^ ydenote an aver- age of the measured output (e.g. hy(t)iM,hy(x(t);t)iM or y[n]). If the dither dis given, then N(z) is known analytically and the estimate ^ x,N1(^y) can then be computed numerically by solving F(^x) = 0 where F(z),N(z)^y: (58) One such method is the bisection method [15], which is guaranteed to converge if F(z) is a continuous function 7 in a domain [ zl;zh] andF(zl) andF(zh) have opposite signs. SincejN(z)^yj, it is always possible to choose a constant z0such thatwl= ^yz0andzh= ^y+z0 will ful l this encompassing condition. See Fig. 6 for an example of N1(y(t)), solved using the bisection method. It is clearly of interest to obtain an estimate of the variance of ^ x. This can be done using a rst-order Taylor series ofN1(w), atN(x(t)), evaluated along the quan- tised output signal y(t). That is, N1(y(t))x(t) + (x(t))"(t); (59) where we have de ned the sensitivity as (x),dN1 dw w=N(x)=dN dz z=x1 : (60) Using (59) directly, or the delta-method [16], an approxi- mation to the variance of N1(y(t)) can be found as CN1(y)(x) 2(x)v(x); (61) withx(t) replaced by xsince time dependence is irrelevant. 5.1. Inverting the sample-average Sample-averaging reduces the variance of the quantised signal, and hence also the variance of the signal obtained by application of the inverse N1(z). Using (59) and (42) we get N1(hy(t)iM)x(t) +"e;M(t); (62) where "e;M(t) = (x(t))h"(t)iM (63) denotes the stochastic error due to the dithering signal, see Fig. 6 for an example of N1(hy(t)iM), found solving (58). Using (43) and (40), the variance of (62) is CN1(hyiM)(x)C"e;M(x) =1 M 2(x)v(x)1 MCN1(y)(x);(64) withx(t) replaced by xsince time dependence is irrele- vant. Hence the sample-average reduces the variance (61) by approximately a factor of 1 =M. 5.2. Inverting the time-average Time-averaging using an LTI low-pass lter G(z) such as (51) reduces the variance of the stochastic error signal "[n] but introduces a frequency dependent error for the deterministic reference signal x[n]. At low frequencies the error forx[n] is small but increases toward the Nyquist- frequency. This can be seen from the high-pass character- istic of the lter H(z) =G(z)1 in (45b). From (59) and (19), the inverse can be written as N1(y[n])x[n] +ef[n] +"t[n]; (65) where ef[n] = [n](hN(x))[n]; (66) Figure 8: The frequency response of the moving average lter G(z) in the passband can be compensated for bye using a pre- lter W(z). By using the pre- lter, the amplitude response can be made to be unity in the passband, compensating for the amplitude attenuation by the moving average lter. The lter W(z) has been synthesised using the least-squares method [17]. denotes the deterministic error due to ltering, and "t[n] = [n](g")[n]: (67) denotes the stochastic error due to the dithering signal. Under the assumptions of Sec. 4.3, the variance of "t[n] will be reduced by a factor of 1 =M. Note that if x(t) is bandwidth limited, it is straight- forward to compensate for the error ef[n] by designing a - nite impulse-response (FIR) compensating pre- lter W(z) which mininises k1G(f)W(f)k2(68) in the passband for x(t), and minimises k1W(f)k2(69) elsewhere, as demonstrated in Fig. 8 for the moving aver- age lter. Given (68), the error H(f) =G(f)W(f)10 in the passband. If the passband is suciently small, the energy due to the di erence kG(z)G(z)W(z)k2will be small. Hence, the e ect of H(z) can be ignored (e.g. it was neglected for the results obtained in Fig. 11). 5.3. Optimal dither variance when inverting The variance CN1(hyiM)(x) ofN1(hy(t)iM) from (64) dependes on the sensitivity (x) from (60) and the vari- ance of the error v(x) from (26). Consider (27) for a nor- mal distribution: The variance v(x) will then have maxima at the thresholds (5), that is at x=Tk, as seen in Fig. 7. Moreover, the sensitivity (x) will have maxima at k, as demonstrated in Fig. 10, with the maximum sensitivity tending to in nity as dgoes to zero (and (x) tending to a constant as dgoes to1). It then follows that the varianceCN1(hyiM)(x) will be dominated by (x) when d= 1 and byv(x) whend= grows larger. As shown below it turns out that there is an optimal d, balancing the sensitivity and the magnitude of d. Two 8 Figure 9: The variance 2(x)v(x) from (61). As Fig. 7, the plot shows the normalised standard deviation for the normally distributed case withd= = 1=4. Figure 10: The sensitivity (x) from (60) and its dependency on the variancedfor a normally distributed dither. metrics for CN1(hyiM)(x) can be devised, suitable for dif- ferent conditions on the reference signal: An averaged case where it is assumed that x=x(t) varies linearly over the quantisation step-size (small , as suggested in [2, 11]), and a worst case where the maximum variance is consid- ered. 5.3.1. Small  As in Sec. 4.3, it is assumed that  is small in the sense thatx(t) varies approximately linearly over the quantisa- tion step-size. Then the variance CN1(hyiM)(x) can be approximated by an average variance  v with CN1(hyiM)(x)v ,1 Z=2 =21 M 2(x)v(x) dx(70) v MZ=2 =2 2(x) dx,aC;(71) where the Cauchy-Schwarz inequality and (53) has been used to obtain (71). Note that aC=aC(M;d), so for each Man approximative smallest upper bound aC(M; d) on the variance CN1(hyiM)(x) can be obtained through (71) by varying das in Fig. 10. This is shown in Fig. 11 for a scaled version of aC(M;d=).It is also informative to look at a mean square error (MSE) estimate of the quantisation error eq(x(t)),y(t) x(t), which is the case without inversion [11]. Consider the sample average heq(x(t))iM,1 MX yi(t) +x(t) =hy(t)iMx(t) =ex(x(t)) +h"(t)iM (72) where (72) follows from (42) and (32). The second moment of the sample average is therefore Efheq(x(t))i2 Mg=e2 x(x(t)) +CheqiM(x(t)) =e2 x(x(t)) +Ch"iM(x(t)) =e2 x(x(t)) +1 Mv(x(t)) (73) Settingx=x(t) and averaging yields an MSE mse =1 Z=2 =2Efheq(x)i2 Mgdx= e2 x+1 Mv (74) where the MSE of the deterministic component is e2 x,1 Z=2 =2e2 x(x) dx; (75) and is in addition to the average variance  vrepresenting the stochastic part of eq. The MSE can be computed nu- merically, as in Fig. 11, or using the expressions prestented in [11]. 5.3.2. Worst case If no assumptions can be made about x(t), the worst- case variance mC= max xCN1(hyiM)(x)1 Mmax x 2(x)v(x):(76) will be a more reasonable estimate of CN1(hyiM)(x). Note thatmC=mC(M;d), similar to aC, which is shown in Fig. 11 for a scaled version. 5.3.3. Optimal variances Upper bounds for the standard deviationp CN1(hyiM)(x) from (71) and (76), and the root mean squared error (RMSE) using (74) are shown in Fig. 11 for a normally distributed dither signal for di erent standard deviations dand three values ofM. Both axes have been normalised relative to the step-size . The numerical values of the optimal cases are presented in Tab. 1. Considering (74), it can be seen that increasing Mde- creases the error due to stocastic dither,  v. The contri- bution to the overall error due to the deterministic com- ponent extherefore becomes more signi cant relative to the stochastic component when using a large number of averages. This is the reason the minimal RMSE depends onM, in addition to the dither variance d. When apply- ing the inverse N1(z), the deterministic component  exis 9 Figure 11: Upper bounds for the standard deviationp CN1(hyiM)(x), using normally distributed dither. Dia- mond, square and circle indicate optimal d= for the three cases. Table 1: Optimal d= values for the three cases: The averaged standard deviation aC, the worst-case standard deviation mCand the root mean squared error (mse) without using inverse compensa- tion. M= 1 (number of averages)paC= 0.390 for d= = 0:230pmC= 0.466 for d= = 0:310pmse= 0.289 for d= = 0:00 M= 100paC= 0.0390 for d= = 0:230pmC= 0.0466 for d= = 0:310pmse= 0.0493 for d= = 0:370 M= 10000paC= 0.00390 for d= = 0:230pmC= 0.00466 for d= = 0:310pmse= 0.00600 for d= = 0:500 removed from the error, and the variances aCandmCare then only scaled by the factor 1 =M. The minimal error variance is determined by the relative signi cance of the size ofdand the sensitivity which depends on d. By numerical evaluation, it was found that aC(M;d) gave a smaller minimizer than mse( M;d) whenM 4, and that mC(M;d) gave a smaller minimizer than mse(M;d) whenM38. Hence, for M4 inverting the average output provides improved performance over only averaging an optimally dithered quantiser. 6. Probability mass function of the quantised dither Consider the case when the reference signal is a con- stant,x(t) =x. The PDF for the input will then be fu(z) =fd(zx): (77)In a similar manner as described in [7], the probability that the input signal u(t) =x+d(t) will truncate to a given integerk2Zcan be found by considering the probability ofu(t) in the domain [ k1=2;k+ 1=2] f(k),Z(k+1=2) (k1=2)fd(zx) dz (78) =Fd((k+ 1=2)x)Fd((k1=2)x) (79) withFdthe CDF of d. Ifd(t) is normally distributed, then the input distribution is fu(z;u;d) =1p 2de(zu)2 2d2(80) withu=x+d, andffrom (78) can be written as f(k;u;d) = erf  k+1 2 up 2d! erf  k1 2 up 2d! :(81) 6.1. Identi cation of the dither probability density function Assuming that the reference signal can be held con- stant and that the dither signal has a PDF of a known type, the parameters 2Rpof the input PDF fu(k;) can be determined from experiment and used to construct an accurate inverse N1(z). By sampling a quantised sig- nal, sayfs1;s2;:::;sng, when the reference is constant, an empirical probabilty mass function (PMF) ^fn(k),1 nnX i=1(sik) can be produced. This can then be used to nd the param- eters, by solving the parameter identi cation problem min  f(k;)^fn(k) 2; (82) using the analytic expression (78). In the case of a normal distribution, robust estimates for ^ uand ^dare found by solving a problem of the form min u;d p f(k;u;d)q ^fn(k) 2; (83) where the problem has been scaled by the square root, as the masses associated with values away from ktend to be very small. In this case it is only possible to identify u, which means that the bias xcan not be compensated for usingN1(z). The latter causes a constant bias to be introduced when applying N1(z). 6.1.1. Numerical example Using (83) to nd the parameters for a normal dis- tribution is illustrated in Fig. 12. In this case  = 1, u==4 andd= =3, andn= 106samples have been used to produce the estimate ^fn(k). Solving (83) with the data used in Fig. 12 yields ^ u=0:2502 and ^d= 0:3330, which are very accurate estimates. 10 -2 -1 0 1 200.20.40.60.81Figure 12: Histogram plot of f(k) and ^fn(k) forn= 106samples of a quantiser with  = 1 driven with a normally distributed signal withu==4 andd= =3. 7. Conclusions When dithering a quantiser, the expected value of the output has a one-to-one correspondence with the reference value, and this function can be determined if the proba- bility density function of the dither and the step-size of the quantiser are known. The deterministic error intro- duced by the quantiser can then be compensated by using the inverse of this function, and the stochastic error due to the dither can be reduced by averaging. Sample and time-averaging achieves approximately the same perfor- mance gain. Time-averaging is more feasible in a physical system, but may require a compensation lter to reduce amplitude and phase distortion of the reference signal due to the averaging operation. It was demonstrated that when applying inversion there is an optimal dither variance that minimises the error af- ter averaging. The optimal variance depends on assump- tions about the reference signal: Either the step-size  is assumed to be small compared to the reference and an average error variance is considered, or the worst-case is considerd by nding the maximum error variance. The variance of the error when using inversion scales with a fac- tor 1=M, whereMis the number of averages. For M4, using an optimal dither, inverting the average output pro- vides improved performance over only averaging an opti- mally dithered quantiser. With a very small amount of averaging, applying inversion is better than mere averag- ing. A method for identifying the probability density func- tion for the dither using empirical data was presented and demonstrated numerically. This enables compensation by inversion in systems where it is dicult to have exact con- trol over the generated dither signal.