File:Delta sigma dithering 1st order results.svg

Quantizing performance of 1st order sigma-delta dithering.

Own work made using the following scilab script (also public domain)

//--------------------------------------------------------------------------------------
// Investigate equivalent bit resolution of first-order delta-sigma dithered
// quantized samples for various oversampling and quantizing amounts on random data.
// The delta-sigma dithering used is:
//
//   out[n] = quantize(in[n]+errsum[n])
//   errsum[n+1] = errsum[n] + in[n] - out[n]
//
// The output is then "perfectly" filtered, decimated back to the original sample rate
// and the standard deviation of the error (difference between input & output) is found.
//
// The number of bits resolution at the original sample rate giving the same error is
// then calculated.
//
// e.g. if we oversample 4 bits by 8 times using delta-sigma dithering and then "perfectly"
// filter the output, what is the equivalent bit resolution at the original sample rate?
//
// Notes:
//   * "perfectly filtering" means using the FFT and zeroing unwanted bins
//
//   * Input data is not quite uniform distribution because of oversampling and then
//     normalizing to ensure no out of range data
//
//   * Calaculation of effective bit resolution assumes evenly distributed input data
//     (which isn't exactly the case as stated in the above note)
//
// This script is public domain
//
// version 1: June 2008, darrell.barrell
//

//--------------------------------------------------------------------------------------
function out = real_fft(x)
// real to complex fft
  out = fft(x(:));
  out = [ out(1); out(2:$/2)*2; out($/2+1)];
endfunction

//--------------------------------------------------------------------------------------
function out = real_ifft(x)
// complex to real fft
  t = [x(1); x(2:$-1)*.5; x($); conj(x($-1:-1:2))*.5];
  out = real(ifft(t));
endfunction

//--------------------------------------------------------------------------------------
function out = clip(x, lower, upper)
// clip "x" to be within bounds "lower" and "upper"
  out = max(min(x, upper), lower);
endfunction


//--------------------------------------------------------------------------------------
function n_bits = equiv_bits(stdev_quan_err)
//
// find the equivalent number of sampling bits resolution that gives a quantizing
// error std deviation of "stdev_quan_err"
//
// Assume
//    linear quantizing
//    rounding to closest quantized value
//    uniform input distribution
//    full scale conversion (i.e. quantized output is from 0 to 1)
//
// stdev_quan_err = sqrt(1/12) ./ (2.^quan_bits-1);
//
// Note, sqrt(1/12) is the stdev of uniform distribution with
// an output range of 1 unit.
//
  n_bits = log(1 ./ (sqrt(12)*stdev_quan_err) + 1) / log(2);
endfunction


// resample factor
oversample_v = [2 4 8 16 32 64];
//oversample_v = 8;

// data length
len_d = 2048;

// dithered bit resolutions to test
dith_bits_v = [1 4 8 12];

// ordinary quantized bit resolutions to plot
quan_bits_v = 1:16;

// create some random data
d = rand(len_d, 1);

// maximum oversample
max_over = max(oversample_v);

// use FFT resample
fft_d = real_fft(d);

// pad spectrum with 0
fft_d(length(d)/2*max_over+1) = 0;

// resampled data at maximum resample rate
d_max_resamp = real_ifft(fft_d)*max_over;

// normalize resampled data between 0 & 1
// note that the distribution is no longer uniform
min_r = min(d_max_resamp);
max_r = max(d_max_resamp);
d_max_resamp = (d_max_resamp-min_r)/(max_r-min_r);
d = (d-min_r)/(max_r-min_r);

//printf("resample error=%f\n", sum((d_max_resamp(1:max_over:$)-d).^2));

// do ordinary quantizing to test
quan_err = [];
for n = 2.^quan_bits-1
  d_quan = round(d*n)/n;
  quan_err($+1) = stdev(d_quan - d);
end
equiv_bits(quan_err)


title("Quantizing performance of 1st order sigma-delta dithering");
xlabel("Times oversample");
ylabel("Equivalent bit resolution");
xgrid;

// dithered bit resolution
for bits = dith_bits_v

  // number of quantizing levels
  levels = 2^bits-1;

  // oversampling
  over = 2;

  // std deviation of dithering error after "perfect" resample
  dith_err = [];

  for over = oversample_v
    printf("levels=%d, over=%d\n", levels, over);

    // decimate resampled data to get resampling that we want
    decimate_stp = max_over/over;
    d_resamp = d_max_resamp(1:decimate_stp:$);

    // current error sum
    err = 0;

    // do the dithering to
    dith = zeros(len_d * over, 1);
    for idx = 1:length(d_resamp)
      // input sample
      in = d_resamp(idx);

      // output sample gets first order sigma-delta dithered & rounded
      out = round((in+err)*levels)/levels;
      out = clip(out, 0, 1);

      // update error
      err = err + in - out;

      // output
      dith(idx) = out;
    end

    // "perfect" filter & decimate using FFT
    dith_fft = real_fft(dith)/over;
    dith_dec = real_ifft(dith_fft(1:len_d/2+1));

    // "noise" power due to quantizing errors (in dB) using the
    // standard deviation (which is sqrt(power))
    dith_err($+1) = stdev(dith_dec - d);
  end

  dith_equiv_bits = equiv_bits(dith_err);
  plot(dith_equiv_bits);
  plot(dith_equiv_bits, 'x');
  if length(dith_equiv_bits) > 1
    printf("%f ", diff(dith_equiv_bits));
    printf("\n");
  end

end

printf("Done\n");

Darrell.barrell

See below.