Mercurial > audlegacy

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/README.new_fmts	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,23 @@
+Since 38449056c293 input plugins can open audio output in their's native format,
+without worrying of output plugin's capabilities. Audio stream will be converted to format specified
+in prefs win (currently FMT_S16_NE and FMT_S24_NE). Dithering will be applied automatically if
+precision loss expected (except FMT_*32_*, on which dithering can cause clipping).
+
+For full list of available formats see plugin.h
+
+Known issues:
+
+1. postproc_flow can deal only with S16_NE streams, so it will be bypassed in 24bit mode, i.e.: no vis, no equalizer, no softvolume.
+2. Currently none of output plugins supports 24bit output :)
+
+Possible solution:
+
+1. It's easy to produce S16_NE stream directly from decoder's output special for vis_flow, as far as the latter has no output.
+2. Switch equalizer and softvolume (and, in perspective, DSP PAPI) to floating-point. Another solution can cause precision
+   loss and avoid dithering. Why floating-point, not fixed? Because:
+     1. libsamplerate requires it.
+     2. Floating-point is standard for inter-component and inter-process audio exchange (LADSPA, JACK etc.).
+
+Have a look at proposed_pipeline.dia
+
+-- asphyx
--- a/configure.ac	Wed Jan 30 20:32:50 2008 +0300
+++ b/configure.ac	Mon Feb 04 05:55:43 2008 +0300
@@ -282,6 +282,54 @@
    [AC_MSG_ERROR([Cannot find libmcs])]
 )

+dnl AltiVec support
+
+AC_ARG_ENABLE(altivec,
+    [  --disable-altivec       disable AltiVec support. (default=enabled) ],
+    [enable_altivec=$enableval],
+    [enable_altivec=yes]
+)
+if test "x$enable_altivec" = "xyes"; then
+  AC_CHECK_HEADERS(altivec.h,
+      [AC_DEFINE(HAVE_ALTIVEC, 1, [Define to 1 if your system has AltiVec.])
+       AC_DEFINE(HAVE_ALTIVEC_H, 1, [Define to 1 if your system has an altivec.h file.])
+       AC_DEFINE(ARCH_POWERPC, 1, [Define to 1 if your system is a PowerPC.])
+       SIMD_CFLAGS=-maltivec
+       AC_SUBST(SIMD_CFLAGS)],
+      [enable_altivec="no"]
+  )
+fi
+
+dnl SSE2 support
+
+AC_ARG_ENABLE(sse2,
+    [  --disable-sse2          disable SSE2 support. (default=enabled) ],
+    [enable_sse2=$enableval],
+    [enable_sse2=yes]
+)
+if test "x$enable_sse2" = "xyes"; then
+  AC_MSG_CHECKING(SSE2)
+  ac_save_CFLAGS="$CFLAGS"
+  CFLAGS="-msse2"
+  AC_TRY_RUN([
+  #include <emmintrin.h>
+  int main()
+  {
+    _mm_setzero_pd();
+    return 0;
+  }
+  ],
+  [AC_MSG_RESULT(SSE2 yes)
+  AC_DEFINE(HAVE_SSE2, 1, [Define to 1 if your system has SSE2])
+  SIMD_CFLAGS=-msse2],
+  [AC_MSG_RESULT(SSE2 no)
+  enable_sse2="no"],
+  )
+fi
+AC_SUBST(SIMD_CFLAGS)
+CFLAGS="$ac_save_CFLAGS"
+
+
 AC_CHECK_HEADERS(unistd.h)
 AC_CHECK_HEADERS(fcntl.h)
 AC_CHECK_HEADERS(sys/ioctl.h)
@@ -464,6 +512,9 @@
 echo "  Automatic character code detection:     $enable_chardet"
 echo "  Sample rate upconversion:               $enable_samplerate"
 echo "  D-Bus support:                          $enable_dbus"
+echo
+echo "  SSE2:                                   $enable_sse2"
+echo "  AltiVec:                                $enable_altivec"
 echo

 if test "$beep_cv_lib_xlibs_threadsafe" = "no"; then
--- a/extra.mk.in	Wed Jan 30 20:32:50 2008 +0300
+++ b/extra.mk.in	Mon Feb 04 05:55:43 2008 +0300
@@ -338,3 +338,4 @@
 TRANSPORT_PLUGINS ?= @TRANSPORT_PLUGINS@
 PCH ?= @PCH@
 IMPLIB_LIBS = @IMPLIB_LIBS@
+SIMD_CFLAGS ?= @SIMD_CFLAGS@
Binary file proposed_pipeline.dia has changed
--- a/src/Makefile	Wed Jan 30 20:32:50 2008 +0300
+++ b/src/Makefile	Mon Feb 04 05:55:43 2008 +0300
@@ -1,4 +1,4 @@
-SUBDIRS = libguess audacious libid3tag
+SUBDIRS = libguess libSAD audacious libid3tag

 include ../extra.mk
--- a/src/audacious/Makefile	Wed Jan 30 20:32:50 2008 +0300
+++ b/src/audacious/Makefile	Mon Feb 04 05:55:43 2008 +0300
@@ -134,13 +134,15 @@
        images/playback.png					\
        images/playlist.png					\
        images/plugins.png					\
+       images/replay_gain.png					\
        ui/equalizer.ui						\
        ui/mainwin.ui						\
        ui/playlist.ui						\
        ui/carbon-menubar.ui

 CLEAN = build_stamp.c libaudacious.exe.a
-EXT_DEPS = ../libguess/libguess.a
+EXT_DEPS = ../libguess/libguess.a \
+           ../libSAD/libSAD.a

 include ../../buildsys.mk

@@ -161,7 +163,9 @@
 	    ${samplerate_CFLAGS}	\
 	    ${REGEX_CFLAGS}		\
 	    ${LIBMCS_CFLAGS}		\
+	    ${SIMD_CFLAGS}		\
 	    -D_AUDACIOUS_CORE		\
+	    -I../libSAD			\
 	    -I.. -I../..		\
 	    -I./intl

@@ -169,6 +173,7 @@
 	${LIBINTL}		\
 	${samplerate_LIBS}	\
 	../libguess/libguess.a	\
+	../libSAD/libSAD.a	\
 	${GTK_LIBS}		\
 	${DBUS_LIBS}		\
 	${MOWGLI_LIBS}		\
Binary file src/audacious/images/replay_gain.png has changed
--- a/src/audacious/input.c	Wed Jan 30 20:32:50 2008 +0300
+++ b/src/audacious/input.c	Mon Feb 04 05:55:43 2008 +0300
@@ -250,6 +250,8 @@
                 right[i] = (*ptr16++);
             }
         break;
+    default:
+        AUDDBG("Incorrect sample format!\n");
     }
 }
--- a/src/audacious/main.c	Wed Jan 30 20:32:50 2008 +0300
+++ b/src/audacious/main.c	Mon Feb 04 05:55:43 2008 +0300
@@ -79,6 +79,7 @@
 #include "ui_fileinfo.h"
 #include "signals.h"
 #include "ui_skinned_window.h"
+#include "libSAD.h"

 #include "icons-stock.h"
 #include "images/audacious_player.xpm"
@@ -230,6 +231,14 @@
     .warn_about_broken_gtk_engines = TRUE,           /* warn about broken gtk themes */
     FALSE,          /* disable inline themes */
     TRUE,           /* remember jtf text entry */
+    16,             /* output bit depth */
+    TRUE,           /* enable replay gain */
+    TRUE,           /* enable clipping prevention */
+    TRUE,           /* track mode */
+    FALSE,          /* album mode */
+    FALSE,          /* enable hard limiter */
+    0.0,            /* preamp */
+    -9.0,           /* default gain */
 };

 typedef struct bmp_cfg_boolent_t {
@@ -341,6 +350,11 @@
     {"warn_about_broken_gtk_engines", &cfg.warn_about_broken_gtk_engines, TRUE},
     {"disable_inline_gtk", &cfg.disable_inline_gtk, TRUE},
     {"remember_jtf_entry", &cfg.remember_jtf_entry, TRUE},
+    {"enable_replay_gain",         &cfg.enable_replay_gain, TRUE},
+    {"enable_clipping_prevention", &cfg.enable_clipping_prevention, TRUE},
+    {"replay_gain_track",          &cfg.replay_gain_track, TRUE},
+    {"replay_gain_album",          &cfg.replay_gain_album, TRUE},
+    {"enable_hard_limiter",        &cfg.enable_hard_limiter, TRUE},
 };

 static gint ncfgbent = G_N_ELEMENTS(bmp_boolents);
@@ -378,6 +392,7 @@
     {"colorize_r", &cfg.colorize_r, TRUE},
     {"colorize_g", &cfg.colorize_g, TRUE},
     {"colorize_b", &cfg.colorize_b, TRUE},
+    {"output_bit_depth", &cfg.output_bit_depth, TRUE},
 };

 static gint ncfgient = G_N_ELEMENTS(bmp_numents);
@@ -628,6 +643,10 @@
         }
     }

+    /* RG settings */
+    cfg_db_get_float(db, NULL, "replay_gain_preamp", &cfg.replay_gain_preamp);
+    cfg_db_get_float(db, NULL, "default_gain", &cfg.default_gain);
+
     cfg_db_close(db);


@@ -880,6 +899,10 @@
     }

     cfg_db_set_float(db, NULL, "equalizer_preamp", cfg.equalizer_preamp);
+
+    /* RG settings */
+    cfg_db_set_float(db, NULL, "replay_gain_preamp", cfg.replay_gain_preamp);
+    cfg_db_set_float(db, NULL, "default_gain",       cfg.default_gain);

     for (i = 0; i < 10; i++) {
         str = g_strdup_printf("equalizer_band%d", i);
@@ -1427,6 +1450,7 @@
     signal_handlers_init();

     g_random_set_seed(time(NULL));
+    SAD_dither_init_rand((gint32)time(NULL));

     bmp_config_load();
--- a/src/audacious/main.h	Wed Jan 30 20:32:50 2008 +0300
+++ b/src/audacious/main.h	Mon Feb 04 05:55:43 2008 +0300
@@ -142,6 +142,14 @@
     gboolean warn_about_broken_gtk_engines;
     gboolean disable_inline_gtk;
     gboolean remember_jtf_entry;
+    gint output_bit_depth;
+    gboolean enable_replay_gain;
+    gboolean enable_clipping_prevention;
+    gboolean replay_gain_track;
+    gboolean replay_gain_album;
+    gboolean enable_hard_limiter;
+    gfloat replay_gain_preamp;
+    gfloat default_gain;
 };

 typedef struct _BmpConfig BmpConfig;
--- a/src/audacious/output.c	Wed Jan 30 20:32:50 2008 +0300
+++ b/src/audacious/output.c	Mon Feb 04 05:55:43 2008 +0300
@@ -1,5 +1,5 @@
 /*  Audacious - Cross-platform multimedia player
- *  Copyright (C) 2005-2007  Audacious team
+ *  Copyright (C) 2005-2008  Audacious team
  *
  *  Based on BMP:
  *  Copyright (C) 2003-2004  BMP development team.
@@ -23,6 +23,8 @@
  *  Audacious or using our public API to be a derived work.
  */

+/* #define AUD_DEBUG */
+
 #ifdef HAVE_CONFIG_H
 #  include "config.h"
 #endif
@@ -42,12 +44,16 @@
 #include "volumecontrol.h"
 #include "visualization.h"

+#include "libSAD.h"
+
 #include <math.h>

 #ifdef USE_SRC
 #include <samplerate.h>
 #endif

+#define FMT_FRACBITS(a) ( (a) == FMT_FIXED32 ? __AUDACIOUS_ASSUMED_MAD_F_FRACBITS__ : 0 )
+
 OutputPluginData op_data = {
     NULL,
     NULL
@@ -77,6 +83,82 @@
     .written_time = get_written_time,
 };

+static const struct {
+    AFormat afmt;
+    SAD_sample_format sadfmt;
+} format_table[] = {
+    {FMT_U8, SAD_SAMPLE_U8},
+    {FMT_S8, SAD_SAMPLE_S8},
+
+    {FMT_S16_LE, SAD_SAMPLE_S16_LE},
+    {FMT_S16_BE, SAD_SAMPLE_S16_BE},
+    {FMT_S16_NE, SAD_SAMPLE_S16},
+
+    {FMT_U16_LE, SAD_SAMPLE_U16_LE},
+    {FMT_U16_BE, SAD_SAMPLE_U16_BE},
+    {FMT_U16_NE, SAD_SAMPLE_U16},
+
+    {FMT_S24_LE, SAD_SAMPLE_S24_LE},
+    {FMT_S24_BE, SAD_SAMPLE_S24_BE},
+    {FMT_S24_NE, SAD_SAMPLE_S24},
+
+    {FMT_U24_LE, SAD_SAMPLE_U24_LE},
+    {FMT_U24_BE, SAD_SAMPLE_U24_BE},
+    {FMT_U24_NE, SAD_SAMPLE_U24},
+
+    {FMT_S32_LE, SAD_SAMPLE_S32_LE},
+    {FMT_S32_BE, SAD_SAMPLE_S32_BE},
+    {FMT_S32_NE, SAD_SAMPLE_S32},
+
+    {FMT_U32_LE, SAD_SAMPLE_U32_LE},
+    {FMT_U32_BE, SAD_SAMPLE_U32_BE},
+    {FMT_U32_NE, SAD_SAMPLE_U32},
+
+    {FMT_FLOAT, SAD_SAMPLE_FLOAT},
+    {FMT_FIXED32, SAD_SAMPLE_FIXED32},
+};
+
+static void apply_replaygain_info (ReplayGainInfo *rg_info);
+
+static inline unsigned sample_size(AFormat fmt) {
+  switch(fmt) {
+    case FMT_S8:
+    case FMT_U8: return sizeof(gint8);
+    case FMT_S16_NE:
+    case FMT_S16_LE:
+    case FMT_S16_BE:
+    case FMT_U16_NE:
+    case FMT_U16_LE:
+    case FMT_U16_BE: return sizeof(gint16);
+    case FMT_S24_NE:
+    case FMT_S24_LE:
+    case FMT_S24_BE:
+    case FMT_U24_NE:
+    case FMT_U24_LE:
+    case FMT_U24_BE:
+    case FMT_S32_NE:
+    case FMT_S32_LE:
+    case FMT_S32_BE:
+    case FMT_U32_NE:
+    case FMT_U32_LE:
+    case FMT_U32_BE:
+    case FMT_FIXED32: return sizeof(gint32);
+    case FMT_FLOAT: return sizeof(float);
+    default: return 0;
+  }
+}
+
+static SAD_sample_format
+sadfmt_from_afmt(AFormat fmt)
+{
+    int i;
+    for (i = 0; i < sizeof(format_table) / sizeof(format_table[0]); i++) {
+        if (format_table[i].afmt == fmt) return format_table[i].sadfmt;
+    }
+
+    return -1;
+}
+
 OutputPlugin *
 get_current_output_plugin(void)
 {
@@ -251,16 +333,44 @@

 #endif

+static SAD_dither_t *sad_state = NULL;
+static SAD_dither_t *sad_state_to_float = NULL;
+static SAD_dither_t *sad_state_from_float = NULL;
+static void *sad_out_buf = NULL;
+static int sad_out_buf_length = 0;
+static ReplayGainInfo replay_gain_info = {
+    .track_gain = 0.0,
+    .track_peak = 0.0,
+    .album_gain = 0.0,
+    .album_peak = 0.0,
+};
+
+static void
+freeSAD()
+{
+    if (sad_state != NULL)            {SAD_dither_free(sad_state);            sad_state = NULL;}
+    if (sad_state_from_float != NULL) {SAD_dither_free(sad_state_from_float); sad_state_from_float = NULL;}
+    if (sad_state_to_float != NULL)   {SAD_dither_free(sad_state_to_float);   sad_state_to_float = NULL;}
+    if (sad_out_buf != NULL)          {free(sad_out_buf);                     sad_out_buf = NULL; sad_out_buf_length = 0;}
+}
+
 gint
 output_open_audio(AFormat fmt, gint rate, gint nch)
 {
     gint ret;
     OutputPlugin *op;
-
+    AUDDBG("\n");
+
+    AFormat output_fmt;
+    int bit_depth;
+    SAD_buffer_format input_sad_fmt;
+    SAD_buffer_format output_sad_fmt;
+
 #ifdef USE_SRC
-    ConfigDb *db;
     gboolean src_enabled;
     gint src_rate, src_type;
+    ConfigDb *db;
+
     db = cfg_db_open();

     if (cfg_db_get_bool(db, NULL, "enable_src", &src_enabled) == FALSE)
@@ -280,13 +390,9 @@
     else
       converter_type = src_type;

-    cfg_db_close(db);
-
     freeSRC();

-    if(src_enabled&&
-       (fmt == FMT_S16_NE||(fmt == FMT_S16_LE && G_BYTE_ORDER == G_LITTLE_ENDIAN)||
-	(fmt == FMT_S16_BE && G_BYTE_ORDER == G_BIG_ENDIAN)))
+    if(src_enabled)
       {
 	src_state = src_new(converter_type, nch, &srcError);
 	if (src_state != NULL)
@@ -295,22 +401,119 @@
 	    rate = overSamplingFs;
 	  }
 	else
+        {
 	  fprintf(stderr, "src_new(): %s\n\n", src_strerror(srcError));
+          src_enabled = FALSE;
+        }
       }
+
+    cfg_db_close(db);
 #endif

+    /*if (cfg_db_get_int(db, NULL, "output_bit_depth", &bit_depth) == FALSE) bit_depth = 16;*/
+    bit_depth = cfg.output_bit_depth;
+
+    AUDDBG("bit depth: %d\n", bit_depth);
+    output_fmt = (bit_depth == 24) ? FMT_S24_NE : FMT_S16_NE;
+    /*output_fmt = (bit_depth == 24) ? FMT_S24_LE : FMT_S16_LE;*/ /* no reason to support other output formats --asphyx */
+
+    freeSAD();
+
+#ifdef USE_SRC
+    if (src_enabled) {
+        AUDDBG("initializing dithering engine for 2 stage conversion\n");
+        input_sad_fmt.sample_format = sadfmt_from_afmt(fmt);
+        if (input_sad_fmt.sample_format < 0) return FALSE;
+        input_sad_fmt.fracbits = FMT_FRACBITS(fmt);
+        input_sad_fmt.channels = nch;
+        input_sad_fmt.channels_order = SAD_CHORDER_INTERLEAVED;
+        input_sad_fmt.samplerate = 0;
+
+        output_sad_fmt.sample_format = SAD_SAMPLE_FLOAT;
+        output_sad_fmt.fracbits = 0;
+        output_sad_fmt.channels = nch;
+        output_sad_fmt.channels_order = SAD_CHORDER_INTERLEAVED;
+        output_sad_fmt.samplerate = 0;
+
+        sad_state_to_float = SAD_dither_init(&input_sad_fmt, &output_sad_fmt, &ret);
+        if (sad_state_to_float == NULL) {
+            AUDDBG("ditherer init failed (decoder's native --> float)\n");
+            return FALSE;
+        }
+        SAD_dither_set_dither (sad_state_to_float, FALSE);
+
+        input_sad_fmt.sample_format = SAD_SAMPLE_FLOAT;
+        input_sad_fmt.fracbits = 0;
+        input_sad_fmt.channels = nch;
+        input_sad_fmt.channels_order = SAD_CHORDER_INTERLEAVED;
+        input_sad_fmt.samplerate = 0;
+
+        output_sad_fmt.sample_format = sadfmt_from_afmt(output_fmt);
+        if (output_sad_fmt.sample_format < 0) return FALSE;
+        output_sad_fmt.fracbits = FMT_FRACBITS(output_fmt);
+        output_sad_fmt.channels = nch;
+        output_sad_fmt.channels_order = SAD_CHORDER_INTERLEAVED;
+        output_sad_fmt.samplerate = 0;
+
+        sad_state_from_float = SAD_dither_init(&input_sad_fmt, &output_sad_fmt, &ret);
+        if (sad_state_from_float == NULL) {
+            SAD_dither_free(sad_state_to_float);
+            AUDDBG("ditherer init failed (float --> output)\n");
+            return FALSE;
+        }
+        SAD_dither_set_dither (sad_state_from_float, TRUE);
+
+        fmt = output_fmt;
+    } else
+#endif /* USE_SRC */
+    {   /* needed for RG processing !*/
+        AUDDBG("initializing dithering engine for direct conversion\n");
+
+        input_sad_fmt.sample_format = sadfmt_from_afmt(fmt);
+        if (input_sad_fmt.sample_format < 0) return FALSE;
+        input_sad_fmt.fracbits = FMT_FRACBITS(fmt);
+        input_sad_fmt.channels = nch;
+        input_sad_fmt.channels_order = SAD_CHORDER_INTERLEAVED;
+        input_sad_fmt.samplerate = 0; /* resampling not implemented yet in libSAD */
+
+        output_sad_fmt.sample_format = sadfmt_from_afmt(output_fmt);
+        output_sad_fmt.fracbits = FMT_FRACBITS(output_fmt);
+        output_sad_fmt.channels = nch;
+        output_sad_fmt.channels_order = SAD_CHORDER_INTERLEAVED;
+        output_sad_fmt.samplerate = 0;
+
+        sad_state = SAD_dither_init(&input_sad_fmt, &output_sad_fmt, &ret);
+        if (sad_state == NULL) {
+            AUDDBG("ditherer init failed\n");
+            return FALSE;
+        }
+        SAD_dither_set_dither (sad_state, TRUE);
+
+        fmt = output_fmt;
+    }
+
+    if(replay_gain_info.album_peak == 0.0 && replay_gain_info.track_peak == 0.0) {
+        AUDDBG("RG info isn't set yet. Filling replay_gain_info with default values.\n");
+        replay_gain_info.track_gain = cfg.default_gain;
+        replay_gain_info.track_peak = 1.0;
+        replay_gain_info.album_gain = cfg.default_gain;
+        replay_gain_info.album_peak = 1.0;
+    }
+    apply_replaygain_info(&replay_gain_info);
+
     op = get_current_output_plugin();

     if (op == NULL)
-        return -1;
+        return FALSE;

     /* Is our output port already open? */
     if ((op_state.rate != 0 && op_state.nch != 0) &&
         (op_state.rate == rate && op_state.nch == nch && op_state.fmt == fmt))
     {
         /* Yes, and it's the correct sampling rate. Reset the counter and go. */
+        AUDDBG("flushing output instead of reopening\n");
         op->flush(0);
-        return 1;
+        return TRUE;
     }
     else if (op_state.rate != 0 && op_state.nch != 0)
         op->close_audio();
@@ -319,6 +522,7 @@

     if (ret == 1)            /* Success? */
     {
+        AUDDBG("opened audio: fmt=%d, rate=%d, nch=%d\n", fmt, rate, nch);
         op_state.fmt = fmt;
         op_state.rate = rate;
         op_state.nch = nch;
@@ -344,9 +548,18 @@
 {
     OutputPlugin *op = get_current_output_plugin();

+    AUDDBG("\n");
+
 #ifdef USE_SRC
     freeSRC();
 #endif
+    freeSAD();
+
+    AUDDBG("clearing RG settings\n");
+    replay_gain_info.track_gain = 0.0;
+    replay_gain_info.track_peak = 0.0;
+    replay_gain_info.album_gain = 0.0;
+    replay_gain_info.album_peak = 0.0;

     /* Do not close if there are still songs to play and the user has
      * not requested a stop.  --nenolod
@@ -423,6 +636,8 @@
     static Flow *postproc_flow = NULL;
     OutputPlugin *op = playback->output;
     gint writeoffs;
+
+    if (length <= 0) return;
     gint time = playback->output->written_time();

     if (postproc_flow == NULL)
@@ -435,9 +650,10 @@
     }

 #ifdef USE_SRC
-    if(src_state != NULL&&length > 0)
+    if(src_state != NULL)
       {
-        int lrLength = length / nch;
+        /*int lrLength = length / nch;*/
+        int lrLength = length / sample_size(fmt);
         int overLrLength = (int)floor(lrLength*(src_data.src_ratio+1));
 	if(lengthOfSrcIn < lrLength)
 	  {
@@ -451,9 +667,10 @@
 	    free(srcOut);
 	    free(wOut);
 	    srcOut = (float*)malloc(sizeof(float)*overLrLength);
-	    wOut = (short int*)malloc(sizeof(short int)*overLrLength);
+	    wOut = (short int*)malloc(sample_size(op_state.fmt) * overLrLength);
 	  }
-        src_short_to_float_array((short int*)ptr, srcIn, lrLength);
+        /*src_short_to_float_array((short int*)ptr, srcIn, lrLength);*/
+        SAD_dither_process_buffer(sad_state_to_float, ptr, srcIn, lrLength / nch);
         src_data.data_in = srcIn;
         src_data.data_out = srcOut;
         src_data.end_of_input = 0;
@@ -465,14 +682,32 @@
           }
         else
           {
-            src_float_to_short_array(srcOut, wOut, src_data.output_frames_gen*2);
+            /*src_float_to_short_array(srcOut, wOut, src_data.output_frames_gen*2);*/
+            SAD_dither_process_buffer(sad_state_from_float, srcOut, wOut, src_data.output_frames_gen);
             ptr = wOut;
-            length = src_data.output_frames_gen * (nch * 2);
+            length = src_data.output_frames_gen *  op_state.nch * sample_size(op_state.fmt);
           }
-      }
+      } else
 #endif
-
-    length = flow_execute(postproc_flow, time, &ptr, length, op_state.fmt, op_state.rate, op_state.nch);
+    if(sad_state != NULL) {
+        int frames  = length / nch / sample_size(fmt);
+        int len =  frames * op_state.nch * sample_size(op_state.fmt);
+        if(sad_out_buf == NULL || sad_out_buf_length < len ) {
+            if(sad_out_buf != NULL) free (sad_out_buf);
+            sad_out_buf = malloc(len);
+            sad_out_buf_length = len;
+        }
+        SAD_dither_process_buffer(sad_state, ptr, sad_out_buf, frames);
+        ptr = sad_out_buf;
+        length = len;
+    }
+
+    if (op_state.fmt == FMT_S16_NE || (op_state.fmt == FMT_S16_LE && G_BYTE_ORDER == G_LITTLE_ENDIAN) ||
+                                      (op_state.fmt == FMT_S16_BE && G_BYTE_ORDER == G_BIG_ENDIAN)) {
+        length = flow_execute(postproc_flow, time, &ptr, length, op_state.fmt, op_state.rate, op_state.nch);
+    } else {
+        AUDDBG("postproc_flow can deal only with S16_NE streams\n"); /*FIXME*/
+    }

     writeoffs = 0;
     while (writeoffs < length)
@@ -521,3 +756,58 @@
         writeoffs += writable;
     }
 }
+
+/* called by input plugin when RG info available --asphyx */
+void
+output_set_replaygain_info (InputPlayback *pb, ReplayGainInfo *rg_info)
+{
+    replay_gain_info = *rg_info;
+    apply_replaygain_info(rg_info);
+}
+
+static void
+apply_replaygain_info (ReplayGainInfo *rg_info)
+{
+    SAD_replaygain_mode mode;
+    SAD_replaygain_info info;
+    /*ConfigDb *db;*/
+    gboolean rg_enabled;
+    gboolean album_mode;
+    SAD_dither_t *active_state;
+
+
+    if(sad_state == NULL && sad_state_from_float == NULL) {
+        AUDDBG("SAD not initialized!\n");
+        return;
+    }
+
+    rg_enabled = cfg.enable_replay_gain;
+    album_mode = cfg.replay_gain_album;
+    mode.clipping_prevention = cfg.enable_clipping_prevention;
+    mode.hard_limit = cfg.enable_hard_limiter;
+
+    if(!rg_enabled) return;
+
+    mode.mode = album_mode ? SAD_RG_ALBUM : SAD_RG_TRACK;
+    mode.preamp = cfg.replay_gain_preamp;
+
+    info.present = TRUE;
+    info.track_gain = rg_info->track_gain;
+    info.track_peak = rg_info->track_peak;
+    info.album_gain = rg_info->album_gain;
+    info.album_peak = rg_info->album_peak;
+
+    AUDDBG("Applying Replay Gain settings:\n");
+    AUDDBG("* mode:                %s\n",     mode.mode == SAD_RG_ALBUM ? "album" : "track");
+    AUDDBG("* clipping prevention: %s\n",     mode.clipping_prevention ? "yes" : "no");
+    AUDDBG("* hard limit:          %s\n",     mode.hard_limit ? "yes" : "no");
+    AUDDBG("* preamp:              %+f dB\n", mode.preamp);
+    AUDDBG("Replay Gain info for current track:\n");
+    AUDDBG("* track gain:          %+f dB\n", info.track_gain);
+    AUDDBG("* track peak:          %f\n",     info.track_peak);
+    AUDDBG("* album gain:          %+f dB\n", info.album_gain);
+    AUDDBG("* album peak:          %f\n",     info.album_peak);
+
+    active_state = sad_state != NULL ? sad_state : sad_state_from_float;
+    SAD_dither_apply_replaygain(active_state, &info, &mode);
+}
--- a/src/audacious/output.h	Wed Jan 30 20:32:50 2008 +0300
+++ b/src/audacious/output.h	Mon Feb 04 05:55:43 2008 +0300
@@ -63,6 +63,7 @@
 gint output_buffer_playing(void);

 void output_pass_audio(InputPlayback *, AFormat, gint, gint, gpointer, int *);
+void output_set_replaygain_info(InputPlayback *, ReplayGainInfo *);

 gint get_written_time(void);
 gint get_output_time(void);
--- a/src/audacious/playback.c	Wed Jan 30 20:32:50 2008 +0300
+++ b/src/audacious/playback.c	Mon Feb 04 05:55:43 2008 +0300
@@ -386,6 +386,7 @@
     playback->set_params = playback_set_pb_params;
     playback->set_title = playback_set_pb_title;
     playback->pass_audio = output_pass_audio;
+    playback->set_replaygain_info = output_set_replaygain_info;

     return playback;
 }
--- a/src/audacious/plugin.h	Wed Jan 30 20:32:50 2008 +0300
+++ b/src/audacious/plugin.h	Mon Feb 04 05:55:43 2008 +0300
@@ -62,14 +62,56 @@
 typedef enum {
     FMT_U8,
     FMT_S8,
+
     FMT_U16_LE,
     FMT_U16_BE,
     FMT_U16_NE,
     FMT_S16_LE,
     FMT_S16_BE,
-    FMT_S16_NE
+    FMT_S16_NE,
+
+    /* added in Audacious 1.5 */
+    FMT_U24_LE, /* stored in lower 3 bytes of 32-bit value, highest byte must be always set to 0 --asphyx */
+    FMT_U24_BE,
+    FMT_U24_NE,
+    FMT_S24_LE,
+    FMT_S24_BE,
+    FMT_S24_NE,
+
+    FMT_U32_LE,
+    FMT_U32_BE,
+    FMT_U32_NE,
+    FMT_S32_LE,
+    FMT_S32_BE,
+    FMT_S32_NE,
+
+    FMT_FLOAT,
+    FMT_FIXED32, /* equivalent of libmad's mad_fixed_t explained below */
 } AFormat;

+/* From mad.h:
+ *
+ * Fixed-point format: 0xABBBBBBB
+ * A == whole part      (sign + 3 bits)
+ * B == fractional part (28 bits)
+ *
+ * Values are signed two's complement, so the effective range is:
+ * 0x80000000 to 0x7fffffff
+ *       -8.0 to +7.9999999962747097015380859375
+ *
+ * The smallest representable value is:
+ * 0x00000001 == 0.0000000037252902984619140625 (i.e. about 3.725e-9)
+ *
+ * 28 bits of fractional accuracy represent about
+ * 8.6 digits of decimal accuracy.
+ *
+ * Fixed-point numbers can be added or subtracted as normal
+ * integers, but multiplication requires shifting the 64-bit result
+ * from 56 fractional bits back to 28 (and rounding.)
+ */
+
+#define __AUDACIOUS_ASSUMED_MAD_F_FRACBITS__ 28 /* useful for build time check for plugins linked against libmad, i.e. madplug */
+
 typedef enum {
     INPUT_VIS_ANALYZER,
     INPUT_VIS_SCOPE,
@@ -94,6 +136,13 @@
     InputPlugin *ip;
 } ProbeResult;

+typedef struct {
+    gfloat track_gain; /* in dB !!! --asphyx */
+    gfloat track_peak;
+    gfloat album_gain;
+    gfloat album_peak;
+} ReplayGainInfo;
+
 typedef GHashTable INIFile;

 #include "audacious/playlist.h"
@@ -1053,6 +1102,10 @@
     void (*set_title) (InputPlayback *, gchar * text);

     void (*pass_audio) (InputPlayback *, AFormat, gint, gint, gpointer, gint *);
+
+    /* added in Audacious 1.5 */
+    /* called by input plugin when RG info available --asphyx */
+    void (*set_replaygain_info) (InputPlayback *, ReplayGainInfo *);
 };

 struct _InputPlugin {
--- a/src/audacious/ui_preferences.c	Wed Jan 30 20:32:50 2008 +0300
+++ b/src/audacious/ui_preferences.c	Mon Feb 04 05:55:43 2008 +0300
@@ -130,6 +130,7 @@
 static Category categories[] = {
     {DATA_DIR "/images/appearance.png",   N_("Appearance")},
     {DATA_DIR "/images/audio.png",        N_("Audio")},
+    {DATA_DIR "/images/replay_gain.png",  N_("Replay Gain")},
     {DATA_DIR "/images/connectivity.png", N_("Connectivity")},
     {DATA_DIR "/images/mouse.png",        N_("Mouse")},
     {DATA_DIR "/images/playback.png",     N_("Playback")},
@@ -191,7 +192,8 @@
 static void mainwin_font_set_cb();
 static void playlist_font_set_cb();
 GtkWidget *ui_preferences_chardet_table_populate(void);
-static gint titlestring_timeout_counter = 0;
+static GtkWidget *ui_preferences_bit_depth(void);
+static GtkWidget *ui_preferences_rg_params(void);

 static PreferencesWidget appearance_misc_widgets[] = {
     {WIDGET_LABEL, N_("<b>_Fonts</b>"), NULL, NULL, NULL, FALSE},
@@ -219,6 +221,20 @@
         N_("When checked, Audacious will detect file formats on demand. This can result in a messier playlist, but delivers a major speed benefit."), FALSE},
     {WIDGET_CHK_BTN, N_("Detect file formats by extension."), &cfg.use_extension_probing, NULL,
         N_("When checked, Audacious will detect file formats based by extension. Only files with extensions of supported formats will be loaded."), FALSE},
+    {WIDGET_LABEL, N_("<b>Bit Depth</b>"), NULL, NULL, NULL, FALSE},
+    {WIDGET_CUSTOM, NULL, NULL, NULL, NULL, FALSE, ui_preferences_bit_depth},
+};
+
+static PreferencesWidget replay_gain_page_widgets[] = {
+    {WIDGET_LABEL, N_("<b>Replay Gain configuration</b>"), NULL, NULL, NULL, FALSE},
+    {WIDGET_CHK_BTN, N_("Enable Replay Gain"), &cfg.enable_replay_gain, NULL, NULL, FALSE},
+    {WIDGET_LABEL, N_("<b>Replay Gain mode</b>"), NULL, NULL, NULL, TRUE},
+    {WIDGET_RADIO_BTN, N_("Track gain/peak"), &cfg.replay_gain_track, NULL, NULL, TRUE},
+    {WIDGET_RADIO_BTN, N_("Album gain/peak"), &cfg.replay_gain_album, NULL, NULL, TRUE},
+    {WIDGET_LABEL, N_("<b>Miscellaneous</b>"), NULL, NULL, NULL, TRUE},
+    {WIDGET_CHK_BTN, N_("Enable clipping prevention"), &cfg.enable_clipping_prevention, NULL, NULL, TRUE},
+    {WIDGET_CHK_BTN, N_("Enable 6 dB hard limiter"), &cfg.enable_hard_limiter, NULL, NULL, TRUE},
+    {WIDGET_CUSTOM, NULL, NULL, NULL, NULL, TRUE, ui_preferences_rg_params},
 };

 static PreferencesWidget playback_page_widgets[] = {
@@ -1657,6 +1673,97 @@
     return widget;
 }

+static void
+on_bit_depth_cbox_changed(GtkWidget *cbox, gpointer data)
+{
+    gint active = gtk_combo_box_get_active(GTK_COMBO_BOX(cbox));
+    cfg.output_bit_depth = (active == 1) ? 24 : 16;
+}
+
+GtkWidget *
+ui_preferences_bit_depth(void)
+{
+    GtkWidget *box = gtk_hbox_new(FALSE, 10);
+    GtkWidget *label = gtk_label_new(_("Output bit depth:"));
+    gtk_misc_set_alignment(GTK_MISC(label), 1, 0.5);
+    gtk_box_pack_start(GTK_BOX(box), label, FALSE, FALSE, 0);
+
+    GtkWidget *combo = gtk_combo_box_new_text ();
+    gtk_combo_box_append_text(GTK_COMBO_BOX (combo), "16");
+    gtk_combo_box_append_text(GTK_COMBO_BOX (combo), "24");
+    gtk_box_pack_start(GTK_BOX(box), combo, FALSE, FALSE, 0);
+
+    gint active = (cfg.output_bit_depth == 24) ? 1 : 0;
+    gtk_combo_box_set_active(GTK_COMBO_BOX(combo), active);
+    g_signal_connect(combo, "changed", G_CALLBACK(on_bit_depth_cbox_changed), NULL);
+
+    gtk_tooltips_set_tip (tooltips, box,
+                          _("All streams will be converted to this bit depth.\n"
+                          "This should be the max supported bit depth of\nthe sound card or output plugin."),
+                          NULL);
+
+    return box;
+}
+
+static void
+on_rg_spin_changed(GtkSpinButton *spinbutton, gpointer user_data)
+{
+    *((gfloat*) user_data) = gtk_spin_button_get_value(spinbutton);
+}
+
+static GtkWidget *
+ui_preferences_rg_params(void)
+{
+    GtkWidget *alignment = gtk_alignment_new(0.5, 0.5, 1, 1);
+    GtkWidget *table = gtk_table_new(2, 3, FALSE);
+    gtk_table_set_row_spacings (GTK_TABLE (table), 6);
+    gtk_table_set_col_spacings (GTK_TABLE (table), 6);
+
+    GtkWidget *label = gtk_label_new(_("Preamp:"));
+    gtk_table_attach(GTK_TABLE(table), label, 0, 1, 0, 1,
+                     (GtkAttachOptions) (GTK_FILL),
+                     (GtkAttachOptions) (0), 0, 0);
+    gtk_misc_set_alignment(GTK_MISC(label), 0, 0.5);
+
+    GtkWidget *spin = gtk_spin_button_new_with_range(-15, 15, 0.01);
+    gtk_table_attach(GTK_TABLE(table), spin, 1, 2, 0, 1,
+                     (GtkAttachOptions) (0),
+                     (GtkAttachOptions) (0), 0, 0);
+    gtk_spin_button_set_value(GTK_SPIN_BUTTON(spin), cfg.replay_gain_preamp);
+    g_signal_connect(G_OBJECT(spin), "value_changed", G_CALLBACK(on_rg_spin_changed), &cfg.replay_gain_preamp);
+
+    label = gtk_label_new(_("dB"));
+    gtk_table_attach(GTK_TABLE(table), label, 2, 3, 0, 1,
+                     (GtkAttachOptions) (GTK_FILL),
+                     (GtkAttachOptions) (0), 0, 0);
+    gtk_misc_set_alignment(GTK_MISC(label), 0, 0.5);
+
+    label = gtk_label_new(_("Default gain:"));
+    gtk_table_attach(GTK_TABLE(table), label, 0, 1, 1, 2,
+                     (GtkAttachOptions) (GTK_FILL),
+                     (GtkAttachOptions) (0), 0, 0);
+    gtk_misc_set_alignment(GTK_MISC(label), 0, 0.5);
+
+    spin = gtk_spin_button_new_with_range(-15, 15, 0.01);
+    gtk_table_attach(GTK_TABLE(table), spin, 1, 2, 1, 2,
+                     (GtkAttachOptions) (0),
+                     (GtkAttachOptions) (0), 0, 0);
+    gtk_spin_button_set_value(GTK_SPIN_BUTTON(spin), cfg.default_gain);
+    g_signal_connect(G_OBJECT(spin), "value_changed", G_CALLBACK(on_rg_spin_changed), &cfg.default_gain);
+    gtk_tooltips_set_tip (tooltips, spin, _("This gain will be used if file doesn't contain Replay Gain metadata."), NULL);
+
+    label = gtk_label_new(_("dB"));
+    gtk_table_attach(GTK_TABLE(table), label, 2, 3, 1, 2,
+                     (GtkAttachOptions) (GTK_FILL),
+                     (GtkAttachOptions) (0), 0, 0);
+    gtk_misc_set_alignment(GTK_MISC(label), 0, 0.5);
+
+    gtk_container_add(GTK_CONTAINER(alignment), table);
+    gtk_alignment_set_padding(GTK_ALIGNMENT(alignment), 0, 0, 12, 0);
+
+    return alignment;
+}
+
 /* it's at early stage */
 static void
 create_widgets(GtkBox *box, PreferencesWidget *widgets, gint amt)
@@ -1705,6 +1812,7 @@
                 gtk_misc_set_alignment(GTK_MISC(widget), 0, 0.5);
                 break;
             case WIDGET_RADIO_BTN:
+                gtk_alignment_set_padding(GTK_ALIGNMENT(alignment), 0, 0, 12, 0);
                 widget = gtk_radio_button_new_with_mnemonic(radio_btn_group, _(widgets[x].label));
                 radio_btn_group = gtk_radio_button_get_group (GTK_RADIO_BUTTON (widget));
                 g_signal_connect(G_OBJECT(widget), "toggled",
@@ -1966,6 +2074,20 @@
 }

 static void
+create_replay_gain_category(void)
+{
+    GtkWidget *rg_page_vbox;
+    GtkWidget *widgets_vbox;
+
+    rg_page_vbox = gtk_vbox_new (FALSE, 0);
+    gtk_container_add (GTK_CONTAINER (category_notebook), rg_page_vbox);
+
+    widgets_vbox = gtk_vbox_new (FALSE, 0);
+    create_widgets(GTK_BOX(widgets_vbox), replay_gain_page_widgets, G_N_ELEMENTS(replay_gain_page_widgets));
+    gtk_box_pack_start (GTK_BOX (rg_page_vbox), widgets_vbox, TRUE, TRUE, 0);
+}
+
+static void
 create_playlist_category(void)
 {
     GtkWidget *playlist_page_vbox;
@@ -2981,6 +3103,7 @@

     create_appearence_category();
     create_audio_category();
+    create_replay_gain_category();
     create_connectivity_category();
     create_mouse_category();
     create_playback_category();
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/Makefile	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,8 @@
+STATIC_LIB_NOINST = libSAD.a
+SRCS = dither.c \
+       dither_ops.c \
+       noicegen.c
+
+include ../../buildsys.mk
+
+CPPFLAGS += ${SIMD_CFLAGS}
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/SFMT-alti.h	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,156 @@
+/**
+ * @file SFMT-alti.h
+ *
+ * @brief SIMD oriented Fast Mersenne Twister(SFMT)
+ * pseudorandom number generator
+ *
+ * @author Mutsuo Saito (Hiroshima University)
+ * @author Makoto Matsumoto (Hiroshima University)
+ *
+ * Copyright (C) 2007 Mutsuo Saito, Makoto Matsumoto and Hiroshima
+ * University. All rights reserved.
+ *
+ * The new BSD License is applied to this software.
+ * see LICENSE.txt
+ */
+
+#ifndef SFMT_ALTI_H
+#define SFMT_ALTI_H
+
+inline static vector unsigned int vec_recursion(vector unsigned int a,
+						vector unsigned int b,
+						vector unsigned int c,
+						vector unsigned int d)
+    ALWAYSINLINE;
+
+/**
+ * This function represents the recursion formula in AltiVec and BIG ENDIAN.
+ * @param a a 128-bit part of the interal state array
+ * @param b a 128-bit part of the interal state array
+ * @param c a 128-bit part of the interal state array
+ * @param d a 128-bit part of the interal state array
+ * @return output
+ */
+inline static vector unsigned int vec_recursion(vector unsigned int a,
+						vector unsigned int b,
+						vector unsigned int c,
+						vector unsigned int d) {
+
+    const vector unsigned int sl1 = ALTI_SL1;
+    const vector unsigned int sr1 = ALTI_SR1;
+#ifdef ONLY64
+    const vector unsigned int mask = ALTI_MSK64;
+    const vector unsigned char perm_sl = ALTI_SL2_PERM64;
+    const vector unsigned char perm_sr = ALTI_SR2_PERM64;
+#else
+    const vector unsigned int mask = ALTI_MSK;
+    const vector unsigned char perm_sl = ALTI_SL2_PERM;
+    const vector unsigned char perm_sr = ALTI_SR2_PERM;
+#endif
+    vector unsigned int v, w, x, y, z;
+    x = vec_perm(a, (vector unsigned int)perm_sl, perm_sl);
+    v = a;
+    y = vec_sr(b, sr1);
+    z = vec_perm(c, (vector unsigned int)perm_sr, perm_sr);
+    w = vec_sl(d, sl1);
+    z = vec_xor(z, w);
+    y = vec_and(y, mask);
+    v = vec_xor(v, x);
+    z = vec_xor(z, y);
+    z = vec_xor(z, v);
+    return z;
+}
+
+/**
+ * This function fills the internal state array with pseudorandom
+ * integers.
+ */
+inline static void gen_rand_all(void) {
+    int i;
+    vector unsigned int r, r1, r2;
+
+    r1 = sfmt[N - 2].s;
+    r2 = sfmt[N - 1].s;
+    for (i = 0; i < N - POS1; i++) {
+	r = vec_recursion(sfmt[i].s, sfmt[i + POS1].s, r1, r2);
+	sfmt[i].s = r;
+	r1 = r2;
+	r2 = r;
+    }
+    for (; i < N; i++) {
+	r = vec_recursion(sfmt[i].s, sfmt[i + POS1 - N].s, r1, r2);
+	sfmt[i].s = r;
+	r1 = r2;
+	r2 = r;
+    }
+}
+
+/**
+ * This function fills the user-specified array with pseudorandom
+ * integers.
+ *
+ * @param array an 128-bit array to be filled by pseudorandom numbers.
+ * @param size number of 128-bit pesudorandom numbers to be generated.
+ */
+inline static void gen_rand_array(w128_t *array, int size) {
+    int i, j;
+    vector unsigned int r, r1, r2;
+
+    r1 = sfmt[N - 2].s;
+    r2 = sfmt[N - 1].s;
+    for (i = 0; i < N - POS1; i++) {
+	r = vec_recursion(sfmt[i].s, sfmt[i + POS1].s, r1, r2);
+	array[i].s = r;
+	r1 = r2;
+	r2 = r;
+    }
+    for (; i < N; i++) {
+	r = vec_recursion(sfmt[i].s, array[i + POS1 - N].s, r1, r2);
+	array[i].s = r;
+	r1 = r2;
+	r2 = r;
+    }
+    /* main loop */
+    for (; i < size - N; i++) {
+	r = vec_recursion(array[i - N].s, array[i + POS1 - N].s, r1, r2);
+	array[i].s = r;
+	r1 = r2;
+	r2 = r;
+    }
+    for (j = 0; j < 2 * N - size; j++) {
+	sfmt[j].s = array[j + size - N].s;
+    }
+    for (; i < size; i++) {
+	r = vec_recursion(array[i - N].s, array[i + POS1 - N].s, r1, r2);
+	array[i].s = r;
+	sfmt[j++].s = r;
+	r1 = r2;
+	r2 = r;
+    }
+}
+
+#ifndef ONLY64
+#if defined(__APPLE__)
+#define ALTI_SWAP (vector unsigned char) \
+	(4, 5, 6, 7, 0, 1, 2, 3, 12, 13, 14, 15, 8, 9, 10, 11)
+#else
+#define ALTI_SWAP {4, 5, 6, 7, 0, 1, 2, 3, 12, 13, 14, 15, 8, 9, 10, 11}
+#endif
+/**
+ * This function swaps high and low 32-bit of 64-bit integers in user
+ * specified array.
+ *
+ * @param array an 128-bit array to be swaped.
+ * @param size size of 128-bit array.
+ */
+inline static void swap(w128_t *array, int size) {
+    int i;
+    const vector unsigned char perm = ALTI_SWAP;
+
+    for (i = 0; i < size; i++) {
+	array[i].s = vec_perm(array[i].s, (vector unsigned int)perm, perm);
+    }
+}
+#endif
+
+#endif
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/SFMT-params.h	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,97 @@
+#ifndef SFMT_PARAMS_H
+#define SFMT_PARAMS_H
+
+#if !defined(MEXP)
+#ifdef __GNUC__
+  #warning "MEXP is not defined. I assume MEXP is 19937."
+#endif
+  #define MEXP 19937
+#endif
+/*-----------------
+  BASIC DEFINITIONS
+  -----------------*/
+/** Mersenne Exponent. The period of the sequence
+ *  is a multiple of 2^MEXP-1.
+ * #define MEXP 19937 */
+/** SFMT generator has an internal state array of 128-bit integers,
+ * and N is its size. */
+#define N (MEXP / 128 + 1)
+/** N32 is the size of internal state array when regarded as an array
+ * of 32-bit integers.*/
+#define N32 (N * 4)
+/** N64 is the size of internal state array when regarded as an array
+ * of 64-bit integers.*/
+#define N64 (N * 2)
+
+/*----------------------
+  the parameters of SFMT
+  following definitions are in paramsXXXX.h file.
+  ----------------------*/
+/** the pick up position of the array.
+#define POS1 122
+*/
+
+/** the parameter of shift left as four 32-bit registers.
+#define SL1 18
+ */
+
+/** the parameter of shift left as one 128-bit register.
+ * The 128-bit integer is shifted by (SL2 * 8) bits.
+#define SL2 1
+*/
+
+/** the parameter of shift right as four 32-bit registers.
+#define SR1 11
+*/
+
+/** the parameter of shift right as one 128-bit register.
+ * The 128-bit integer is shifted by (SL2 * 8) bits.
+#define SR2 1
+*/
+
+/** A bitmask, used in the recursion.  These parameters are introduced
+ * to break symmetry of SIMD.
+#define MSK1 0xdfffffefU
+#define MSK2 0xddfecb7fU
+#define MSK3 0xbffaffffU
+#define MSK4 0xbffffff6U
+*/
+
+/** These definitions are part of a 128-bit period certification vector.
+#define PARITY1	0x00000001U
+#define PARITY2	0x00000000U
+#define PARITY3	0x00000000U
+#define PARITY4	0xc98e126aU
+*/
+
+#if MEXP == 607
+  #include "SFMT-params607.h"
+#elif MEXP == 1279
+  #include "SFMT-params1279.h"
+#elif MEXP == 2281
+  #include "SFMT-params2281.h"
+#elif MEXP == 4253
+  #include "SFMT-params4253.h"
+#elif MEXP == 11213
+  #include "SFMT-params11213.h"
+#elif MEXP == 19937
+  #include "SFMT-params19937.h"
+#elif MEXP == 44497
+  #include "SFMT-params44497.h"
+#elif MEXP == 86243
+  #include "SFMT-params86243.h"
+#elif MEXP == 132049
+  #include "SFMT-params132049.h"
+#elif MEXP == 216091
+  #include "SFMT-params216091.h"
+#else
+#ifdef __GNUC__
+  #error "MEXP is not valid."
+  #undef MEXP
+#else
+  #undef MEXP
+#endif
+
+#endif
+
+#endif /* SFMT_PARAMS_H */
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/SFMT-params19937.h	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,46 @@
+#ifndef SFMT_PARAMS19937_H
+#define SFMT_PARAMS19937_H
+
+#define POS1	122
+#define SL1	18
+#define SL2	1
+#define SR1	11
+#define SR2	1
+#define MSK1	0xdfffffefU
+#define MSK2	0xddfecb7fU
+#define MSK3	0xbffaffffU
+#define MSK4	0xbffffff6U
+#define PARITY1	0x00000001U
+#define PARITY2	0x00000000U
+#define PARITY3	0x00000000U
+#define PARITY4	0x13c9e684U
+
+
+/* PARAMETERS FOR ALTIVEC */
+#if defined(__APPLE__)	/* For OSX */
+    #define ALTI_SL1	(vector unsigned int)(SL1, SL1, SL1, SL1)
+    #define ALTI_SR1	(vector unsigned int)(SR1, SR1, SR1, SR1)
+    #define ALTI_MSK	(vector unsigned int)(MSK1, MSK2, MSK3, MSK4)
+    #define ALTI_MSK64 \
+	(vector unsigned int)(MSK2, MSK1, MSK4, MSK3)
+    #define ALTI_SL2_PERM \
+	(vector unsigned char)(1,2,3,23,5,6,7,0,9,10,11,4,13,14,15,8)
+    #define ALTI_SL2_PERM64 \
+	(vector unsigned char)(1,2,3,4,5,6,7,31,9,10,11,12,13,14,15,0)
+    #define ALTI_SR2_PERM \
+	(vector unsigned char)(7,0,1,2,11,4,5,6,15,8,9,10,17,12,13,14)
+    #define ALTI_SR2_PERM64 \
+	(vector unsigned char)(15,0,1,2,3,4,5,6,17,8,9,10,11,12,13,14)
+#else	/* For OTHER OSs(Linux?) */
+    #define ALTI_SL1	{SL1, SL1, SL1, SL1}
+    #define ALTI_SR1	{SR1, SR1, SR1, SR1}
+    #define ALTI_MSK	{MSK1, MSK2, MSK3, MSK4}
+    #define ALTI_MSK64	{MSK2, MSK1, MSK4, MSK3}
+    #define ALTI_SL2_PERM	{1,2,3,23,5,6,7,0,9,10,11,4,13,14,15,8}
+    #define ALTI_SL2_PERM64	{1,2,3,4,5,6,7,31,9,10,11,12,13,14,15,0}
+    #define ALTI_SR2_PERM	{7,0,1,2,11,4,5,6,15,8,9,10,17,12,13,14}
+    #define ALTI_SR2_PERM64	{15,0,1,2,3,4,5,6,17,8,9,10,11,12,13,14}
+#endif	/* For OSX */
+#define IDSTR	"SFMT-19937:122-18-1-11-1:dfffffef-ddfecb7f-bffaffff-bffffff6"
+
+#endif /* SFMT_PARAMS19937_H */
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/SFMT-sse2.h	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,121 @@
+/**
+ * @file  SFMT-sse2.h
+ * @brief SIMD oriented Fast Mersenne Twister(SFMT) for Intel SSE2
+ *
+ * @author Mutsuo Saito (Hiroshima University)
+ * @author Makoto Matsumoto (Hiroshima University)
+ *
+ * @note We assume LITTLE ENDIAN in this file
+ *
+ * Copyright (C) 2006, 2007 Mutsuo Saito, Makoto Matsumoto and Hiroshima
+ * University. All rights reserved.
+ *
+ * The new BSD License is applied to this software, see LICENSE.txt
+ */
+
+#ifndef SFMT_SSE2_H
+#define SFMT_SSE2_H
+
+inline static __m128i mm_recursion(__m128i *a, __m128i *b, __m128i c,
+				   __m128i d, __m128i mask) ALWAYSINLINE;
+
+/**
+ * This function represents the recursion formula.
+ * @param a a 128-bit part of the interal state array
+ * @param b a 128-bit part of the interal state array
+ * @param c a 128-bit part of the interal state array
+ * @param d a 128-bit part of the interal state array
+ * @param mask 128-bit mask
+ * @return output
+ */
+inline static __m128i mm_recursion(__m128i *a, __m128i *b,
+				   __m128i c, __m128i d, __m128i mask) {
+    __m128i v, x, y, z;
+
+    x = _mm_load_si128(a);
+    y = _mm_srli_epi32(*b, SR1);
+    z = _mm_srli_si128(c, SR2);
+    v = _mm_slli_epi32(d, SL1);
+    z = _mm_xor_si128(z, x);
+    z = _mm_xor_si128(z, v);
+    x = _mm_slli_si128(x, SL2);
+    y = _mm_and_si128(y, mask);
+    z = _mm_xor_si128(z, x);
+    z = _mm_xor_si128(z, y);
+    return z;
+}
+
+/**
+ * This function fills the internal state array with pseudorandom
+ * integers.
+ */
+inline static void gen_rand_all(void) {
+    int i;
+    __m128i r, r1, r2, mask;
+    mask = _mm_set_epi32(MSK4, MSK3, MSK2, MSK1);
+
+    r1 = _mm_load_si128(&sfmt[N - 2].si);
+    r2 = _mm_load_si128(&sfmt[N - 1].si);
+    for (i = 0; i < N - POS1; i++) {
+	r = mm_recursion(&sfmt[i].si, &sfmt[i + POS1].si, r1, r2, mask);
+	_mm_store_si128(&sfmt[i].si, r);
+	r1 = r2;
+	r2 = r;
+    }
+    for (; i < N; i++) {
+	r = mm_recursion(&sfmt[i].si, &sfmt[i + POS1 - N].si, r1, r2, mask);
+	_mm_store_si128(&sfmt[i].si, r);
+	r1 = r2;
+	r2 = r;
+    }
+}
+
+/**
+ * This function fills the user-specified array with pseudorandom
+ * integers.
+ *
+ * @param array an 128-bit array to be filled by pseudorandom numbers.
+ * @param size number of 128-bit pesudorandom numbers to be generated.
+ */
+inline static void gen_rand_array(w128_t *array, int size) {
+    int i, j;
+    __m128i r, r1, r2, mask;
+    mask = _mm_set_epi32(MSK4, MSK3, MSK2, MSK1);
+
+    r1 = _mm_load_si128(&sfmt[N - 2].si);
+    r2 = _mm_load_si128(&sfmt[N - 1].si);
+    for (i = 0; i < N - POS1; i++) {
+	r = mm_recursion(&sfmt[i].si, &sfmt[i + POS1].si, r1, r2, mask);
+	_mm_store_si128(&array[i].si, r);
+	r1 = r2;
+	r2 = r;
+    }
+    for (; i < N; i++) {
+	r = mm_recursion(&sfmt[i].si, &array[i + POS1 - N].si, r1, r2, mask);
+	_mm_store_si128(&array[i].si, r);
+	r1 = r2;
+	r2 = r;
+    }
+    /* main loop */
+    for (; i < size - N; i++) {
+	r = mm_recursion(&array[i - N].si, &array[i + POS1 - N].si, r1, r2,
+			 mask);
+	_mm_store_si128(&array[i].si, r);
+	r1 = r2;
+	r2 = r;
+    }
+    for (j = 0; j < 2 * N - size; j++) {
+	r = _mm_load_si128(&array[j + size - N].si);
+	_mm_store_si128(&sfmt[j].si, r);
+    }
+    for (; i < size; i++) {
+	r = mm_recursion(&array[i - N].si, &array[i + POS1 - N].si, r1, r2,
+			 mask);
+	_mm_store_si128(&array[i].si, r);
+	_mm_store_si128(&sfmt[j++].si, r);
+	r1 = r2;
+	r2 = r;
+    }
+}
+
+#endif
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/SFMT.c	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,618 @@
+/**
+ * @file  SFMT.c
+ * @brief SIMD oriented Fast Mersenne Twister(SFMT)
+ *
+ * @author Mutsuo Saito (Hiroshima University)
+ * @author Makoto Matsumoto (Hiroshima University)
+ *
+ * Copyright (C) 2006,2007 Mutsuo Saito, Makoto Matsumoto and Hiroshima
+ * University. All rights reserved.
+ *
+ * The new BSD License is applied to this software, see LICENSE.txt
+ */
+#include <string.h>
+#include <assert.h>
+#include "SFMT.h"
+#include "SFMT-params.h"
+
+#if defined(__BIG_ENDIAN__) && !defined(__amd64) && !defined(BIG_ENDIAN64)
+#define BIG_ENDIAN64 1
+#endif
+#if defined(HAVE_ALTIVEC) && !defined(BIG_ENDIAN64)
+#define BIG_ENDIAN64 1
+#endif
+#if defined(ONLY64) && !defined(BIG_ENDIAN64)
+  #if defined(__GNUC__)
+    #error "-DONLY64 must be specified with -DBIG_ENDIAN64"
+  #endif
+#undef ONLY64
+#endif
+/*------------------------------------------------------
+  128-bit SIMD data type for Altivec, SSE2 or standard C
+  ------------------------------------------------------*/
+#if defined(HAVE_ALTIVEC)
+  #if !defined(__APPLE__)
+    #include <altivec.h>
+  #endif
+/** 128-bit data structure */
+union W128_T {
+    vector unsigned int s;
+    uint32_t u[4];
+};
+/** 128-bit data type */
+typedef union W128_T w128_t;
+
+#elif defined(HAVE_SSE2)
+  #include <emmintrin.h>
+
+/** 128-bit data structure */
+union W128_T {
+    __m128i si;
+    uint32_t u[4];
+};
+/** 128-bit data type */
+typedef union W128_T w128_t;
+
+#else
+
+/** 128-bit data structure */
+struct W128_T {
+    uint32_t u[4];
+};
+/** 128-bit data type */
+typedef struct W128_T w128_t;
+
+#endif
+
+/*--------------------------------------
+  FILE GLOBAL VARIABLES
+  internal state, index counter and flag
+  --------------------------------------*/
+/** the 128-bit internal state array */
+static w128_t sfmt[N];
+/** the 32bit integer pointer to the 128-bit internal state array */
+static uint32_t *psfmt32 = &sfmt[0].u[0];
+#if !defined(BIG_ENDIAN64) || defined(ONLY64)
+/** the 64bit integer pointer to the 128-bit internal state array */
+static uint64_t *psfmt64 = (uint64_t *)&sfmt[0].u[0];
+#endif
+/** index counter to the 32-bit internal state array */
+static int idx;
+/** a flag: it is 0 if and only if the internal state is not yet
+ * initialized. */
+static int initialized = 0;
+/** a parity check vector which certificate the period of 2^{MEXP} */
+static uint32_t parity[4] = {PARITY1, PARITY2, PARITY3, PARITY4};
+
+/*----------------
+  STATIC FUNCTIONS
+  ----------------*/
+inline static int idxof(int i);
+inline static void rshift128(w128_t *out,  w128_t const *in, int shift);
+inline static void lshift128(w128_t *out,  w128_t const *in, int shift);
+inline static void gen_rand_all(void);
+inline static void gen_rand_array(w128_t *array, int size);
+inline static uint32_t func1(uint32_t x);
+inline static uint32_t func2(uint32_t x);
+static void period_certification(void);
+#if defined(BIG_ENDIAN64) && !defined(ONLY64)
+inline static void swap(w128_t *array, int size);
+#endif
+
+#if defined(HAVE_ALTIVEC)
+  #include "SFMT-alti.h"
+#elif defined(HAVE_SSE2)
+  #include "SFMT-sse2.h"
+#endif
+
+/**
+ * This function simulate a 64-bit index of LITTLE ENDIAN
+ * in BIG ENDIAN machine.
+ */
+#ifdef ONLY64
+inline static int idxof(int i) {
+    return i ^ 1;
+}
+#else
+inline static int idxof(int i) {
+    return i;
+}
+#endif
+/**
+ * This function simulates SIMD 128-bit right shift by the standard C.
+ * The 128-bit integer given in in is shifted by (shift * 8) bits.
+ * This function simulates the LITTLE ENDIAN SIMD.
+ * @param out the output of this function
+ * @param in the 128-bit data to be shifted
+ * @param shift the shift value
+ */
+#ifdef ONLY64
+inline static void rshift128(w128_t *out, w128_t const *in, int shift) {
+    uint64_t th, tl, oh, ol;
+
+    th = ((uint64_t)in->u[2] << 32) | ((uint64_t)in->u[3]);
+    tl = ((uint64_t)in->u[0] << 32) | ((uint64_t)in->u[1]);
+
+    oh = th >> (shift * 8);
+    ol = tl >> (shift * 8);
+    ol |= th << (64 - shift * 8);
+    out->u[0] = (uint32_t)(ol >> 32);
+    out->u[1] = (uint32_t)ol;
+    out->u[2] = (uint32_t)(oh >> 32);
+    out->u[3] = (uint32_t)oh;
+}
+#else
+inline static void rshift128(w128_t *out, w128_t const *in, int shift) {
+    uint64_t th, tl, oh, ol;
+
+    th = ((uint64_t)in->u[3] << 32) | ((uint64_t)in->u[2]);
+    tl = ((uint64_t)in->u[1] << 32) | ((uint64_t)in->u[0]);
+
+    oh = th >> (shift * 8);
+    ol = tl >> (shift * 8);
+    ol |= th << (64 - shift * 8);
+    out->u[1] = (uint32_t)(ol >> 32);
+    out->u[0] = (uint32_t)ol;
+    out->u[3] = (uint32_t)(oh >> 32);
+    out->u[2] = (uint32_t)oh;
+}
+#endif
+/**
+ * This function simulates SIMD 128-bit left shift by the standard C.
+ * The 128-bit integer given in in is shifted by (shift * 8) bits.
+ * This function simulates the LITTLE ENDIAN SIMD.
+ * @param out the output of this function
+ * @param in the 128-bit data to be shifted
+ * @param shift the shift value
+ */
+#ifdef ONLY64
+inline static void lshift128(w128_t *out, w128_t const *in, int shift) {
+    uint64_t th, tl, oh, ol;
+
+    th = ((uint64_t)in->u[2] << 32) | ((uint64_t)in->u[3]);
+    tl = ((uint64_t)in->u[0] << 32) | ((uint64_t)in->u[1]);
+
+    oh = th << (shift * 8);
+    ol = tl << (shift * 8);
+    oh |= tl >> (64 - shift * 8);
+    out->u[0] = (uint32_t)(ol >> 32);
+    out->u[1] = (uint32_t)ol;
+    out->u[2] = (uint32_t)(oh >> 32);
+    out->u[3] = (uint32_t)oh;
+}
+#else
+inline static void lshift128(w128_t *out, w128_t const *in, int shift) {
+    uint64_t th, tl, oh, ol;
+
+    th = ((uint64_t)in->u[3] << 32) | ((uint64_t)in->u[2]);
+    tl = ((uint64_t)in->u[1] << 32) | ((uint64_t)in->u[0]);
+
+    oh = th << (shift * 8);
+    ol = tl << (shift * 8);
+    oh |= tl >> (64 - shift * 8);
+    out->u[1] = (uint32_t)(ol >> 32);
+    out->u[0] = (uint32_t)ol;
+    out->u[3] = (uint32_t)(oh >> 32);
+    out->u[2] = (uint32_t)oh;
+}
+#endif
+
+/**
+ * This function represents the recursion formula.
+ * @param r output
+ * @param a a 128-bit part of the internal state array
+ * @param b a 128-bit part of the internal state array
+ * @param c a 128-bit part of the internal state array
+ * @param d a 128-bit part of the internal state array
+ */
+#ifdef ONLY64
+inline static void do_recursion(w128_t *r, w128_t *a, w128_t *b, w128_t *c,
+				w128_t *d) {
+    w128_t x;
+    w128_t y;
+
+    lshift128(&x, a, SL2);
+    rshift128(&y, c, SR2);
+    r->u[0] = a->u[0] ^ x.u[0] ^ ((b->u[0] >> SR1) & MSK2) ^ y.u[0]
+	^ (d->u[0] << SL1);
+    r->u[1] = a->u[1] ^ x.u[1] ^ ((b->u[1] >> SR1) & MSK1) ^ y.u[1]
+	^ (d->u[1] << SL1);
+    r->u[2] = a->u[2] ^ x.u[2] ^ ((b->u[2] >> SR1) & MSK4) ^ y.u[2]
+	^ (d->u[2] << SL1);
+    r->u[3] = a->u[3] ^ x.u[3] ^ ((b->u[3] >> SR1) & MSK3) ^ y.u[3]
+	^ (d->u[3] << SL1);
+}
+#else
+inline static void do_recursion(w128_t *r, w128_t *a, w128_t *b, w128_t *c,
+				w128_t *d) {
+    w128_t x;
+    w128_t y;
+
+    lshift128(&x, a, SL2);
+    rshift128(&y, c, SR2);
+    r->u[0] = a->u[0] ^ x.u[0] ^ ((b->u[0] >> SR1) & MSK1) ^ y.u[0]
+	^ (d->u[0] << SL1);
+    r->u[1] = a->u[1] ^ x.u[1] ^ ((b->u[1] >> SR1) & MSK2) ^ y.u[1]
+	^ (d->u[1] << SL1);
+    r->u[2] = a->u[2] ^ x.u[2] ^ ((b->u[2] >> SR1) & MSK3) ^ y.u[2]
+	^ (d->u[2] << SL1);
+    r->u[3] = a->u[3] ^ x.u[3] ^ ((b->u[3] >> SR1) & MSK4) ^ y.u[3]
+	^ (d->u[3] << SL1);
+}
+#endif
+
+#if (!defined(HAVE_ALTIVEC)) && (!defined(HAVE_SSE2))
+/**
+ * This function fills the internal state array with pseudorandom
+ * integers.
+ */
+inline static void gen_rand_all(void) {
+    int i;
+    w128_t *r1, *r2;
+
+    r1 = &sfmt[N - 2];
+    r2 = &sfmt[N - 1];
+    for (i = 0; i < N - POS1; i++) {
+	do_recursion(&sfmt[i], &sfmt[i], &sfmt[i + POS1], r1, r2);
+	r1 = r2;
+	r2 = &sfmt[i];
+    }
+    for (; i < N; i++) {
+	do_recursion(&sfmt[i], &sfmt[i], &sfmt[i + POS1 - N], r1, r2);
+	r1 = r2;
+	r2 = &sfmt[i];
+    }
+}
+
+/**
+ * This function fills the user-specified array with pseudorandom
+ * integers.
+ *
+ * @param array an 128-bit array to be filled by pseudorandom numbers.
+ * @param size number of 128-bit pseudorandom numbers to be generated.
+ */
+inline static void gen_rand_array(w128_t *array, int size) {
+    int i, j;
+    w128_t *r1, *r2;
+
+    r1 = &sfmt[N - 2];
+    r2 = &sfmt[N - 1];
+    for (i = 0; i < N - POS1; i++) {
+	do_recursion(&array[i], &sfmt[i], &sfmt[i + POS1], r1, r2);
+	r1 = r2;
+	r2 = &array[i];
+    }
+    for (; i < N; i++) {
+	do_recursion(&array[i], &sfmt[i], &array[i + POS1 - N], r1, r2);
+	r1 = r2;
+	r2 = &array[i];
+    }
+    for (; i < size - N; i++) {
+	do_recursion(&array[i], &array[i - N], &array[i + POS1 - N], r1, r2);
+	r1 = r2;
+	r2 = &array[i];
+    }
+    for (j = 0; j < 2 * N - size; j++) {
+	sfmt[j] = array[j + size - N];
+    }
+    for (; i < size; i++, j++) {
+	do_recursion(&array[i], &array[i - N], &array[i + POS1 - N], r1, r2);
+	r1 = r2;
+	r2 = &array[i];
+	sfmt[j] = array[i];
+    }
+}
+#endif
+
+#if defined(BIG_ENDIAN64) && !defined(ONLY64) && !defined(HAVE_ALTIVEC)
+inline static void swap(w128_t *array, int size) {
+    int i;
+    uint32_t x, y;
+
+    for (i = 0; i < size; i++) {
+	x = array[i].u[0];
+	y = array[i].u[2];
+	array[i].u[0] = array[i].u[1];
+	array[i].u[2] = array[i].u[3];
+	array[i].u[1] = x;
+	array[i].u[3] = y;
+    }
+}
+#endif
+/**
+ * This function represents a function used in the initialization
+ * by init_by_array
+ * @param x 32-bit integer
+ * @return 32-bit integer
+ */
+static uint32_t func1(uint32_t x) {
+    return (x ^ (x >> 27)) * (uint32_t)1664525UL;
+}
+
+/**
+ * This function represents a function used in the initialization
+ * by init_by_array
+ * @param x 32-bit integer
+ * @return 32-bit integer
+ */
+static uint32_t func2(uint32_t x) {
+    return (x ^ (x >> 27)) * (uint32_t)1566083941UL;
+}
+
+/**
+ * This function certificate the period of 2^{MEXP}
+ */
+static void period_certification(void) {
+    int inner = 0;
+    int i, j;
+    uint32_t work;
+
+    for (i = 0; i < 4; i++)
+	inner ^= psfmt32[idxof(i)] & parity[i];
+    for (i = 16; i > 0; i >>= 1)
+	inner ^= inner >> i;
+    inner &= 1;
+    /* check OK */
+    if (inner == 1) {
+	return;
+    }
+    /* check NG, and modification */
+    for (i = 0; i < 4; i++) {
+	work = 1;
+	for (j = 0; j < 32; j++) {
+	    if ((work & parity[i]) != 0) {
+		psfmt32[idxof(i)] ^= work;
+		return;
+	    }
+	    work = work << 1;
+	}
+    }
+}
+
+/*----------------
+  PUBLIC FUNCTIONS
+  ----------------*/
+/**
+ * This function returns the identification string.
+ * The string shows the word size, the Mersenne exponent,
+ * and all parameters of this generator.
+ */
+const char *get_idstring(void) {
+    return IDSTR;
+}
+
+/**
+ * This function returns the minimum size of array used for \b
+ * fill_array32() function.
+ * @return minimum size of array used for fill_array32() function.
+ */
+int get_min_array_size32(void) {
+    return N32;
+}
+
+/**
+ * This function returns the minimum size of array used for \b
+ * fill_array64() function.
+ * @return minimum size of array used for fill_array64() function.
+ */
+int get_min_array_size64(void) {
+    return N64;
+}
+
+#ifndef ONLY64
+/**
+ * This function generates and returns 32-bit pseudorandom number.
+ * init_gen_rand or init_by_array must be called before this function.
+ * @return 32-bit pseudorandom number
+ */
+uint32_t gen_rand32(void) {
+    uint32_t r;
+
+    assert(initialized);
+    if (idx >= N32) {
+	gen_rand_all();
+	idx = 0;
+    }
+    r = psfmt32[idx++];
+    return r;
+}
+#endif
+/**
+ * This function generates and returns 64-bit pseudorandom number.
+ * init_gen_rand or init_by_array must be called before this function.
+ * The function gen_rand64 should not be called after gen_rand32,
+ * unless an initialization is again executed.
+ * @return 64-bit pseudorandom number
+ */
+uint64_t gen_rand64(void) {
+#if defined(BIG_ENDIAN64) && !defined(ONLY64)
+    uint32_t r1, r2;
+#else
+    uint64_t r;
+#endif
+
+    assert(initialized);
+    assert(idx % 2 == 0);
+
+    if (idx >= N32) {
+	gen_rand_all();
+	idx = 0;
+    }
+#if defined(BIG_ENDIAN64) && !defined(ONLY64)
+    r1 = psfmt32[idx];
+    r2 = psfmt32[idx + 1];
+    idx += 2;
+    return ((uint64_t)r2 << 32) | r1;
+#else
+    r = psfmt64[idx / 2];
+    idx += 2;
+    return r;
+#endif
+}
+
+#ifndef ONLY64
+/**
+ * This function generates pseudorandom 32-bit integers in the
+ * specified array[] by one call. The number of pseudorandom integers
+ * is specified by the argument size, which must be at least 624 and a
+ * multiple of four.  The generation by this function is much faster
+ * than the following gen_rand function.
+ *
+ * For initialization, init_gen_rand or init_by_array must be called
+ * before the first call of this function. This function can not be
+ * used after calling gen_rand function, without initialization.
+ *
+ * @param array an array where pseudorandom 32-bit integers are filled
+ * by this function.  The pointer to the array must be \b "aligned"
+ * (namely, must be a multiple of 16) in the SIMD version, since it
+ * refers to the address of a 128-bit integer.  In the standard C
+ * version, the pointer is arbitrary.
+ *
+ * @param size the number of 32-bit pseudorandom integers to be
+ * generated.  size must be a multiple of 4, and greater than or equal
+ * to (MEXP / 128 + 1) * 4.
+ *
+ * @note \b memalign or \b posix_memalign is available to get aligned
+ * memory. Mac OSX doesn't have these functions, but \b malloc of OSX
+ * returns the pointer to the aligned memory block.
+ */
+void fill_array32(uint32_t *array, int size) {
+    assert(initialized);
+    assert(idx == N32);
+    assert(size % 4 == 0);
+    assert(size >= N32);
+
+    gen_rand_array((w128_t *)array, size / 4);
+    idx = N32;
+}
+#endif
+
+/**
+ * This function generates pseudorandom 64-bit integers in the
+ * specified array[] by one call. The number of pseudorandom integers
+ * is specified by the argument size, which must be at least 312 and a
+ * multiple of two.  The generation by this function is much faster
+ * than the following gen_rand function.
+ *
+ * For initialization, init_gen_rand or init_by_array must be called
+ * before the first call of this function. This function can not be
+ * used after calling gen_rand function, without initialization.
+ *
+ * @param array an array where pseudorandom 64-bit integers are filled
+ * by this function.  The pointer to the array must be "aligned"
+ * (namely, must be a multiple of 16) in the SIMD version, since it
+ * refers to the address of a 128-bit integer.  In the standard C
+ * version, the pointer is arbitrary.
+ *
+ * @param size the number of 64-bit pseudorandom integers to be
+ * generated.  size must be a multiple of 2, and greater than or equal
+ * to (MEXP / 128 + 1) * 2
+ *
+ * @note \b memalign or \b posix_memalign is available to get aligned
+ * memory. Mac OSX doesn't have these functions, but \b malloc of OSX
+ * returns the pointer to the aligned memory block.
+ */
+void fill_array64(uint64_t *array, int size) {
+    assert(initialized);
+    assert(idx == N32);
+    assert(size % 2 == 0);
+    assert(size >= N64);
+
+    gen_rand_array((w128_t *)array, size / 2);
+    idx = N32;
+
+#if defined(BIG_ENDIAN64) && !defined(ONLY64)
+    swap((w128_t *)array, size /2);
+#endif
+}
+
+/**
+ * This function initializes the internal state array with a 32-bit
+ * integer seed.
+ *
+ * @param seed a 32-bit integer used as the seed.
+ */
+void init_gen_rand(uint32_t seed) {
+    int i;
+
+    psfmt32[idxof(0)] = seed;
+    for (i = 1; i < N32; i++) {
+	psfmt32[idxof(i)] = 1812433253UL * (psfmt32[idxof(i - 1)]
+					    ^ (psfmt32[idxof(i - 1)] >> 30))
+	    + i;
+    }
+    idx = N32;
+    period_certification();
+    initialized = 1;
+}
+
+/**
+ * This function initializes the internal state array,
+ * with an array of 32-bit integers used as the seeds
+ * @param init_key the array of 32-bit integers, used as a seed.
+ * @param key_length the length of init_key.
+ */
+void init_by_array(uint32_t *init_key, int key_length) {
+    int i, j, count;
+    uint32_t r;
+    int lag;
+    int mid;
+    int size = N * 4;
+
+    if (size >= 623) {
+	lag = 11;
+    } else if (size >= 68) {
+	lag = 7;
+    } else if (size >= 39) {
+	lag = 5;
+    } else {
+	lag = 3;
+    }
+    mid = (size - lag) / 2;
+
+    memset(sfmt, 0x8b, sizeof(sfmt));
+    if (key_length + 1 > N32) {
+	count = key_length + 1;
+    } else {
+	count = N32;
+    }
+    r = func1(psfmt32[idxof(0)] ^ psfmt32[idxof(mid)]
+	      ^ psfmt32[idxof(N32 - 1)]);
+    psfmt32[idxof(mid)] += r;
+    r += key_length;
+    psfmt32[idxof(mid + lag)] += r;
+    psfmt32[idxof(0)] = r;
+
+    count--;
+    for (i = 1, j = 0; (j < count) && (j < key_length); j++) {
+	r = func1(psfmt32[idxof(i)] ^ psfmt32[idxof((i + mid) % N32)]
+		  ^ psfmt32[idxof((i + N32 - 1) % N32)]);
+	psfmt32[idxof((i + mid) % N32)] += r;
+	r += init_key[j] + i;
+	psfmt32[idxof((i + mid + lag) % N32)] += r;
+	psfmt32[idxof(i)] = r;
+	i = (i + 1) % N32;
+    }
+    for (; j < count; j++) {
+	r = func1(psfmt32[idxof(i)] ^ psfmt32[idxof((i + mid) % N32)]
+		  ^ psfmt32[idxof((i + N32 - 1) % N32)]);
+	psfmt32[idxof((i + mid) % N32)] += r;
+	r += i;
+	psfmt32[idxof((i + mid + lag) % N32)] += r;
+	psfmt32[idxof(i)] = r;
+	i = (i + 1) % N32;
+    }
+    for (j = 0; j < N32; j++) {
+	r = func2(psfmt32[idxof(i)] + psfmt32[idxof((i + mid) % N32)]
+		  + psfmt32[idxof((i + N32 - 1) % N32)]);
+	psfmt32[idxof((i + mid) % N32)] ^= r;
+	r -= i;
+	psfmt32[idxof((i + mid + lag) % N32)] ^= r;
+	psfmt32[idxof(i)] = r;
+	i = (i + 1) % N32;
+    }
+
+    idx = N32;
+    period_certification();
+    initialized = 1;
+}
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/SFMT.h	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,125 @@
+/**
+ * @file SFMT.h
+ *
+ * @brief SIMD oriented Fast Mersenne Twister(SFMT) pseudorandom
+ * number generator
+ *
+ * @author Mutsuo Saito (Hiroshima University)
+ * @author Makoto Matsumoto (Hiroshima University)
+ *
+ * Copyright (C) 2006, 2007 Mutsuo Saito, Makoto Matsumoto and Hiroshima
+ * University. All rights reserved.
+ *
+ * The new BSD License is applied to this software.
+ * see LICENSE.txt
+ *
+ * @note We assume that your system has inttypes.h.  If your system
+ * doesn't have inttypes.h, you have to typedef uint32_t and uint64_t,
+ * and you have to define PRIu64 and PRIx64 in this file as follows:
+ * @verbatim
+ typedef unsigned int uint32_t
+ typedef unsigned long long uint64_t
+ #define PRIu64 "llu"
+ #define PRIx64 "llx"
+@endverbatim
+ * uint32_t must be exactly 32-bit unsigned integer type (no more, no
+ * less), and uint64_t must be exactly 64-bit unsigned integer type.
+ * PRIu64 and PRIx64 are used for printf function to print 64-bit
+ * unsigned int and 64-bit unsigned int in hexadecimal format.
+ */
+
+#ifndef SFMT_H
+#define SFMT_H
+
+#include <stdio.h>
+
+#if defined(__STDC_VERSION__) && (__STDC_VERSION__ >= 199901L)
+  #include <inttypes.h>
+#elif defined(_MSC_VER) || defined(__BORLANDC__)
+  typedef unsigned int uint32_t;
+  typedef unsigned __int64 uint64_t;
+  #define inline __inline
+#else
+  #include <inttypes.h>
+  #if defined(__GNUC__)
+    #define inline __inline__
+  #endif
+#endif
+
+#ifndef PRIu64
+  #if defined(_MSC_VER) || defined(__BORLANDC__)
+    #define PRIu64 "I64u"
+    #define PRIx64 "I64x"
+  #else
+    #define PRIu64 "llu"
+    #define PRIx64 "llx"
+  #endif
+#endif
+
+#if defined(__GNUC__)
+#define ALWAYSINLINE __attribute__((always_inline))
+#endif
+
+uint32_t gen_rand32(void);
+uint64_t gen_rand64(void);
+void fill_array32(uint32_t *array, int size);
+void fill_array64(uint64_t *array, int size);
+void init_gen_rand(uint32_t seed);
+void init_by_array(uint32_t *init_key, int key_length);
+const char *get_idstring(void);
+int get_min_array_size32(void);
+int get_min_array_size64(void);
+
+/* These real versions are due to Isaku Wada */
+/** generates a random number on [0,1]-real-interval */
+inline static double to_real1(uint32_t v)
+{
+    return v * (1.0/4294967295.0);
+    /* divided by 2^32-1 */
+}
+
+/** generates a random number on [0,1]-real-interval */
+inline static double genrand_real1(void)
+{
+    return to_real1(gen_rand32());
+}
+
+/** generates a random number on [0,1)-real-interval */
+inline static double to_real2(uint32_t v)
+{
+    return v * (1.0/4294967296.0);
+    /* divided by 2^32 */
+}
+
+/** generates a random number on [0,1)-real-interval */
+inline static double genrand_real2(void)
+{
+    return to_real2(gen_rand32());
+}
+
+/** generates a random number on (0,1)-real-interval */
+inline static double to_real3(uint32_t v)
+{
+    return (((double)v) + 0.5)*(1.0/4294967296.0);
+    /* divided by 2^32 */
+}
+
+/** generates a random number on (0,1)-real-interval */
+inline static double genrand_real3(void)
+{
+    return to_real3(gen_rand32());
+}
+/** These real versions are due to Isaku Wada */
+
+/** generates a random number on [0,1) with 53-bit resolution*/
+inline static double to_res53(uint64_t v)
+{
+    return v * (1.0/18446744073709551616.0L);
+}
+
+/** generates a random number on [0,1) with 53-bit resolution*/
+inline static double genrand_res53(void)
+{
+    return to_res53(gen_rand64());
+}
+#endif
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/dither.c	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,527 @@
+/* Scale & Dither library (libSAD)
+ * High-precision bit depth converter with ReplayGain support
+ *
+ * Copyright (c) 2007-2008 Eugene Zagidullin (e.asphyx@gmail.com)
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation; either version 2 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+ */
+
+/* #define CLIPPING_DEBUG */
+/* #define DITHER_DEBUG */
+
+#include "common.h"
+#include "dither_ops.h"
+#include "noicegen.h"
+
+#include <assert.h>
+#include <math.h>
+
+/*
+ * Supported conversions:
+ *
+ *					O U T P U T
+ *   ,------------------+-----------------------------------------------.
+ *   |			|S8 U8 S16 U16 S24 U24 S32 U32 FLOAT FIXED-POINT|
+ *   +------------------+-----------------------------------------------+
+ *   | S8		|X  X  X   X   X   X   X   X   -     -          |
+ *   | U8		|X  X  X   X   X   X   X   X   -     -          |
+ * I | S16		|X  X  X   X   X   X   X   X   -     -          |
+ * N | U16		|X  X  X   X   X   X   X   X   -     -          |
+ * P | S24		|X  X  X   X   X   X   X   X   -     -          |
+ * U | U24		|X  X  X   X   X   X   X   X   -     -          |
+ * T | S32		|X  X  X   X   X   X   X   X   -     -          |
+ *   | U32		|X  X  X   X   X   X   X   X   -     -          |
+ *   | FLOAT		|X  X  X   X   X   X   X   X   X     -          |
+ *   | FIXED-POINT	|X  X  X   X   X   X   X   X   X     -          |
+ *   `------------------+-----------------------------------------------'
+ */
+
+#define SCALE(x,s) (s != 1.0 ? x * s : x)
+#define MAXINT(a) (1L << ((a)-1))
+#define CLIP(x,m) (x > m-1 ? m-1 : (x < -m ? -m : x))
+
+/* private object */
+typedef struct {
+  SAD_sample_format input_sample_format;
+  SAD_sample_format output_sample_format;
+  int input_bits;
+  int input_fracbits;
+  int output_bits;
+  int output_fracbits;
+  int channels;
+  SAD_channels_order input_chorder;
+  SAD_channels_order output_chorder;
+  SAD_get_sample_proc get_sample;
+  SAD_put_sample_proc put_sample;
+  int dither;
+  int hardlimit;
+  float scale;
+  float rg_scale;
+} SAD_state_priv;
+
+/* error code */
+
+//static SAD_error SAD_last_error = SAD_ERROR_OK;
+
+static inline double compute_hardlimit (double sample, float scale) {
+  sample *= scale;
+  const double k = 0.5;    /* -6dBFS */
+  if (sample > k) {
+    return tanh((sample - k) / (1 - k)) * (1 - k) + k;
+  }
+  else if (sample < -k) {
+    return tanh((sample + k) / (1 - k)) * (1 - k) - k;
+  }
+  return sample;
+}
+
+/*
+ * Dither fixed-point normalized or integer sample to n-bits integer
+ * samples < -1 and > 1 will be clipped
+ */
+
+static inline int32_t __dither_sample_fixed_to_int (int32_t sample, int inbits, int fracbits, int outbits, float scale, int dither,
+							int hardlimit)
+{
+  int n_bits_to_loose, bitwidth, precision_loss;
+  int32_t maxint = MAXINT(outbits);
+
+  n_bits_to_loose = 0;
+  bitwidth = inbits;
+  precision_loss = FALSE;
+
+/*#ifdef DEEP_DEBUG
+  printf("f: __dither_sample_fixed_to_int\n");
+#endif*/
+
+  if (fracbits == 0) {
+    if (inbits<29) {
+      /* convert to 4.28 fixed-point */
+      n_bits_to_loose = 29 - inbits;
+      sample <<= n_bits_to_loose;
+      bitwidth += n_bits_to_loose;
+    }
+
+    n_bits_to_loose += inbits - outbits;
+
+    if (inbits > outbits) {
+      precision_loss = TRUE;
+#ifdef PRECISION_DEBUG
+      printf("Precision loss, reason: bitwidth loss %d --> %d\n", inbits, outbits);
+#endif
+    }
+  } else {
+    n_bits_to_loose = fracbits + 1 - outbits;
+    bitwidth = fracbits;
+    precision_loss = TRUE;
+#ifdef PRECISION_DEBUG
+    printf("Precision loss, reason: fixed-point input\n", inbits, outbits);
+#endif
+  }
+
+  assert(n_bits_to_loose >=0 );
+
+  if (hardlimit) {
+    sample = (int32_t)(compute_hardlimit((double)sample/(double)MAXINT(bitwidth), scale) * (double)MAXINT(bitwidth));
+#ifdef PRECISION_DEBUG
+    printf("Precision loss, reason: hard limiter\n", inbits, outbits);
+#endif
+    precision_loss = TRUE;
+  } else {
+    sample = SCALE(sample, scale);
+  }
+
+  if (scale != 1.0){
+    precision_loss = TRUE;
+#ifdef PRECISION_DEBUG
+    printf("Precision loss, reason: scale\n", inbits, outbits);
+#endif
+  }
+
+  if (precision_loss && (n_bits_to_loose >= 1) && (inbits < 32 || fracbits != 0)) sample += (1L << (n_bits_to_loose - 1));
+
+#ifdef DITHER_DEBUG
+  int32_t val_wo_dither = sample >> n_bits_to_loose;
+  val_wo_dither = CLIP(val_wo_dither, maxint);
+#endif
+  if (dither && precision_loss && (n_bits_to_loose >= 1) && (inbits < 32 || fracbits != 0)) {
+    int32_t dither_num = triangular_dither_noise(n_bits_to_loose + 1);
+    sample += dither_num;
+  }
+
+  sample >>= n_bits_to_loose;
+
+  /* Clipping */
+#ifdef CLIPPING_DEBUG
+  int32_t val_wo_clip = sample;
+#endif
+  sample = CLIP(sample, maxint);
+#ifdef CLIPPING_DEBUG
+  if (val_wo_clip != sample) {
+    printf("Clipping: %d --> %d\n", val_wo_clip, sample);
+  }
+#endif
+#ifdef DITHER_DEBUG
+  if (dither && precision_loss && (n_bits_to_loose >= 1)) printf("%d --> %d, noise: %d\n", val_wo_dither, sample, sample - val_wo_dither);
+#endif
+  return sample;
+}
+
+/*
+ * Dither floating-point normalized sample to n-bits integer
+ * samples < -1 and > 1 will be clipped
+ */
+static inline int32_t __dither_sample_float_to_int (float sample, int nbits, float scale, int dither, int hardlimit) {
+
+#ifdef DEEP_DEBUG
+  printf("f: __dither_sample_float_to_int\n");
+#endif
+
+  int32_t maxint = MAXINT(nbits);
+
+  if (hardlimit) {
+    sample = compute_hardlimit((double)sample, scale);
+  } else {
+    sample = SCALE(sample, scale);
+  }
+
+  sample *= maxint;
+  /* we want to round precisely */
+  sample = (sample < 0 ? sample - 0.5 : sample + 0.5);
+
+#ifdef DITHER_DEBUG
+  int32_t val_wo_dither = (int32_t) sample;
+  val_wo_dither = CLIP(val_wo_dither, maxint);
+#endif
+  if (dither) {
+    float dither_num = triangular_dither_noise_f();
+    sample += dither_num;
+  }
+
+  /* Round and clipping */
+  int32_t value = (int32_t) sample;
+#ifdef CLIPPING_DEBUG
+  int32_t val_wo_clip = value;
+#endif
+  value = CLIP(value, maxint);
+#ifdef CLIPPING_DEBUG
+  if (val_wo_clip != value) {
+    printf("Clipping: %d --> %d\n", val_wo_clip, value);
+  }
+#endif
+
+#ifdef DITHER_DEBUG
+  printf("%d --> %d, noise: %d\n", val_wo_dither, value, value - val_wo_dither);
+#endif
+  return value;
+}
+
+static inline float __dither_sample_float_to_float (float sample, float scale, int hardlimit) {
+#ifdef DEEP_DEBUG
+  printf("f: __dither_sample_float_to_float\n");
+#endif
+  if (hardlimit) {
+    sample = compute_hardlimit((double)sample, scale);
+  } else {
+    sample = SCALE(sample, scale);
+  }
+  return sample;
+}
+
+static inline float __dither_sample_fixed_to_float (int32_t sample, int inbits, int fracbits, float scale, int hardlimit) {
+  float fsample;
+
+#ifdef DEEP_DEBUG
+  printf("f: __dither_sample_fixed_to_float\n");
+#endif
+  if (fracbits == 0) {
+     fsample = (float)sample / (float)MAXINT(inbits);
+  } else {
+     fsample = (float)sample / (float)MAXINT(fracbits+1);
+  }
+  return __dither_sample_float_to_float (fsample, scale, hardlimit);
+}
+
+
+
+
+
+SAD_dither_t* SAD_dither_init(SAD_buffer_format *inbuf_format, SAD_buffer_format *outbuf_format, int *error) {
+  SAD_state_priv *priv;
+
+  DEBUG_MSG("f: SAD_dither_init\n",0);
+
+  priv = calloc(sizeof(SAD_state_priv), 1);
+
+  /* Check buffer formats and assign buffer ops */
+  SAD_buffer_ops* inops = SAD_assign_buf_ops(inbuf_format);
+
+  if (inbuf_format->sample_format != SAD_SAMPLE_FLOAT) {
+    if (inops != NULL) {
+      priv->get_sample = inops->get_sample;
+    } else {
+      free(priv);
+      *error = SAD_ERROR_INCORRECT_INPUT_SAMPLEFORMAT;
+      return NULL;
+    }
+  }
+
+  SAD_buffer_ops* outops = SAD_assign_buf_ops(outbuf_format);
+
+  if (outbuf_format->sample_format != SAD_SAMPLE_FLOAT) {
+    if (outops != NULL) {
+      priv->put_sample = outops->put_sample;
+    } else {
+      free(priv);
+      *error = SAD_ERROR_INCORRECT_OUTPUT_SAMPLEFORMAT;
+      return NULL;
+    }
+  }
+
+  priv->input_fracbits = 0;
+  priv->output_fracbits = 0;
+  priv->input_sample_format = inbuf_format->sample_format;
+  priv->output_sample_format = outbuf_format->sample_format;
+  priv->input_chorder = inbuf_format->channels_order;
+  priv->output_chorder = outbuf_format->channels_order;
+  priv->channels = inbuf_format->channels;
+  priv->scale = 1.0;
+  priv->rg_scale = 1.0;
+  priv->dither = TRUE;
+  priv->hardlimit = FALSE;
+
+  switch(outbuf_format->sample_format){
+    case SAD_SAMPLE_S8:
+    case SAD_SAMPLE_U8: priv->output_bits = 8; break;
+    case SAD_SAMPLE_S16:
+    case SAD_SAMPLE_S16_LE:
+    case SAD_SAMPLE_S16_BE:
+    case SAD_SAMPLE_U16:
+    case SAD_SAMPLE_U16_LE:
+    case SAD_SAMPLE_U16_BE: priv->output_bits = 16; break;
+    case SAD_SAMPLE_S24:
+    case SAD_SAMPLE_S24_LE:
+    case SAD_SAMPLE_S24_BE:
+    case SAD_SAMPLE_U24:
+    case SAD_SAMPLE_U24_LE:
+    case SAD_SAMPLE_U24_BE: priv->output_bits = 24; break;
+    case SAD_SAMPLE_S32:
+    case SAD_SAMPLE_S32_LE:
+    case SAD_SAMPLE_S32_BE:
+    case SAD_SAMPLE_U32:
+    case SAD_SAMPLE_U32_LE:
+    case SAD_SAMPLE_U32_BE: priv->output_bits = 32; break;
+    case SAD_SAMPLE_FLOAT: break;
+    default:
+      free(priv);
+      *error = SAD_ERROR_INCORRECT_OUTPUT_SAMPLEFORMAT;
+      return NULL;
+  }
+
+  switch(inbuf_format->sample_format){
+    case SAD_SAMPLE_S8:
+    case SAD_SAMPLE_U8: priv->input_bits = 8; break;
+    case SAD_SAMPLE_S16:
+    case SAD_SAMPLE_S16_LE:
+    case SAD_SAMPLE_S16_BE:
+    case SAD_SAMPLE_U16:
+    case SAD_SAMPLE_U16_LE:
+    case SAD_SAMPLE_U16_BE: priv->input_bits = 16; break;
+    case SAD_SAMPLE_S24:
+    case SAD_SAMPLE_S24_LE:
+    case SAD_SAMPLE_S24_BE:
+    case SAD_SAMPLE_U24:
+    case SAD_SAMPLE_U24_LE:
+    case SAD_SAMPLE_U24_BE: priv->input_bits = 24; break;
+    case SAD_SAMPLE_S32:
+    case SAD_SAMPLE_S32_LE:
+    case SAD_SAMPLE_S32_BE:
+    case SAD_SAMPLE_U32:
+    case SAD_SAMPLE_U32_LE:
+    case SAD_SAMPLE_U32_BE: priv->input_bits = 32; break;
+    case SAD_SAMPLE_FIXED32: priv->input_fracbits = inbuf_format->fracbits; break;
+    case SAD_SAMPLE_FLOAT: break;
+    default:
+      free(priv);
+      *error = SAD_ERROR_INCORRECT_INPUT_SAMPLEFORMAT;
+      return NULL;
+  }
+
+  *error = SAD_ERROR_OK;
+  return (SAD_dither_t*)priv;
+}
+
+int SAD_dither_free(SAD_dither_t* state) {
+  DEBUG_MSG("f: SAD_dither_free\n",0);
+  free(state);
+  return SAD_ERROR_OK;
+}
+
+/*
+ * Depend on format->channels_order inbuf and outbuf will be treated as
+ * smth* or smth** if channels_order = SAD_CHORDER_INTERLEAVED or SAD_CHORDER_SEPARATED
+ * accordingly
+ *
+ * frame is aggregate of format->channels samples
+ */
+
+#define GET_FLOAT_SAMPLE(b,o,n,c,i) (o == SAD_CHORDER_INTERLEAVED ? (((float*)b)[i*n+c]) : (((float**)b)[c][i]))
+#define PUT_FLOAT_SAMPLE(b,o,n,c,i,s) { \
+    if (o == SAD_CHORDER_INTERLEAVED) { \
+      ((float*)b)[i*n+c] = s;		\
+    } else {				\
+      ((float**)b)[c][i] = s;		\
+    }					\
+  }
+
+int SAD_dither_process_buffer (SAD_dither_t *state, void *inbuf, void *outbuf, int frames)
+{
+  SAD_state_priv *priv = (SAD_state_priv*) state;
+  int i, ch;
+  int channels = priv->channels;
+  int inbits = priv->input_bits;
+  int outbits = priv->output_bits;
+  int fracbits = priv->input_fracbits;
+  float scale = priv->scale * priv->rg_scale;
+  int dither = priv->dither;
+  int hardlimit = priv->hardlimit;
+  SAD_channels_order input_chorder = priv->input_chorder;
+  SAD_channels_order output_chorder = priv->output_chorder;
+
+  SAD_get_sample_proc get_sample = priv->get_sample;
+  SAD_put_sample_proc put_sample = priv->put_sample;
+
+#ifdef DEEP_DEBUG
+  printf("f: SAD_process_buffer\n");
+#endif
+
+  if (priv->input_sample_format == SAD_SAMPLE_FLOAT) {
+      if (priv->output_sample_format == SAD_SAMPLE_FLOAT) {
+          /* process buffer */
+          for(i=0; i<frames; i++) {
+	      for(ch=0; ch<channels; ch++) {
+	          float sample = GET_FLOAT_SAMPLE(inbuf, input_chorder, channels, ch ,i);
+	          sample = __dither_sample_float_to_float(sample, scale, hardlimit);
+                  PUT_FLOAT_SAMPLE(outbuf, output_chorder, channels, ch ,i, sample);
+	      }
+	  }
+      } else {
+          if (put_sample == NULL) return SAD_ERROR_CORRUPTED_PRIVATE_DATA;
+          /* process buffer */
+          for(i=0; i<frames; i++) {
+	      for(ch=0; ch<channels; ch++) {
+	          float sample = GET_FLOAT_SAMPLE(inbuf, input_chorder, channels, ch ,i);
+	          int32_t isample = __dither_sample_float_to_int(sample, outbits, scale, dither, hardlimit);
+                  put_sample (outbuf, isample, channels, ch, i);
+	      }
+	  }
+      }
+  } else {
+      if (priv->output_sample_format == SAD_SAMPLE_FLOAT) {
+          if (get_sample == NULL) return SAD_ERROR_CORRUPTED_PRIVATE_DATA;
+          /* process buffer */
+          for(i=0; i<frames; i++) {
+	      for(ch=0; ch<channels; ch++) {
+  	          int32_t sample = get_sample (inbuf, channels, ch, i);
+                  float fsample = __dither_sample_fixed_to_float (sample, inbits, fracbits, scale, hardlimit);
+                  PUT_FLOAT_SAMPLE(outbuf, output_chorder, channels, ch ,i, fsample);
+	      }
+	  }
+      } else {
+          if (put_sample == NULL || get_sample == NULL) return SAD_ERROR_CORRUPTED_PRIVATE_DATA;
+          /* process buffer */
+          for(i=0; i<frames; i++) {
+	      for(ch=0; ch<channels; ch++){
+  	          int32_t sample = get_sample (inbuf, channels, ch, i);
+	          int32_t isample = __dither_sample_fixed_to_int (sample, inbits, fracbits, outbits, scale, dither, hardlimit);
+                  put_sample (outbuf, isample, channels, ch, i);
+	      }
+	  }
+      }
+  }
+
+  return SAD_ERROR_OK;
+}
+
+int SAD_dither_apply_replaygain (SAD_dither_t *state, SAD_replaygain_info *rg_info, SAD_replaygain_mode *mode) {
+  SAD_state_priv *priv = (SAD_state_priv*) state;
+  float scale = -1.0, peak = 0.0;
+
+  DEBUG_MSG("f: SAD_dither_apply_replaygain\n",0);
+
+  if(!rg_info->present) {
+    priv->rg_scale = 1.0;
+    priv->hardlimit = FALSE;
+    return SAD_ERROR_OK;
+  }
+
+  switch(mode->mode) {
+    case SAD_RG_ALBUM:
+      scale = db2scale(rg_info->album_gain);
+      peak = rg_info->album_peak;
+      if (peak == 0.0) {
+        scale = db2scale(rg_info->track_gain); // fallback to per-track mode
+	peak = rg_info->track_peak;
+	DEBUG_MSG("f: SAD_dither_apply_replaygain: fallback to track mode\n",0);
+      }
+      break;
+    case SAD_RG_TRACK:
+      scale = db2scale(rg_info->track_gain);
+      peak = rg_info->track_peak;
+      if (peak == 0.0) {
+        scale = db2scale(rg_info->album_gain); // fallback to per-album mode
+	peak = rg_info->album_peak;
+	DEBUG_MSG("f: SAD_dither_apply_replaygain: fallback to album mode\n",0);
+      }
+      break;
+    case SAD_RG_NONE:
+      scale = -1.0;
+  }
+
+  if (scale != -1.0 && peak != 0.0) {
+    DEBUG_MSG("f: SAD_dither_apply_replaygain: applying\n",0);
+    scale *= db2scale(mode->preamp);
+    // Clipping prevention
+    if(mode->clipping_prevention) {
+#ifdef DEBUG
+      if(scale * peak > 1.0) DEBUG_MSG("f: SAD_dither_apply_replaygain: clipping prevented\n",0);
+#endif
+      scale = scale * peak > 1.0 ? 1.0 / peak : scale;
+    }
+    scale = scale > 15.0 ? 15.0 : scale; // safety
+    priv->rg_scale = scale;
+    priv->hardlimit = mode->hard_limit; // apply settings
+  } else {
+    priv->rg_scale = 1.0;
+    priv->hardlimit = FALSE;
+  }
+
+  return SAD_ERROR_OK;
+}
+
+int SAD_dither_set_scale (SAD_dither_t *state, float scale) {
+  SAD_state_priv *priv = (SAD_state_priv*) state;
+  priv->scale = scale;
+  return SAD_ERROR_OK;
+}
+
+int SAD_dither_set_dither (SAD_dither_t *state, int dither) {
+  SAD_state_priv *priv = (SAD_state_priv*) state;
+  priv->dither = dither;
+  return SAD_ERROR_OK;
+}
+
+void SAD_dither_init_rand(uint32_t seed) {
+  noicegen_init_rand(seed);
+}
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/dither.h	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,49 @@
+/* Scale & Dither library (libSAD)
+ * High-precision bit depth converter with ReplayGain support
+ *
+ * Copyright (c) 2007-2008 Eugene Zagidullin (e.asphyx@gmail.com)
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation; either version 2 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+ */
+
+#ifndef DITHER_H
+#define DITHER_H
+
+#include "common.h"
+
+#define SAD_ERROR_INCORRECT_INPUT_SAMPLEFORMAT -2
+#define SAD_ERROR_INCORRECT_OUTPUT_SAMPLEFORMAT -3
+#define SAD_ERROR_CORRUPTED_PRIVATE_DATA -4
+
+typedef int32_t (*SAD_get_sample_proc) (void *buf, int nch, int ch, int i);
+typedef void (*SAD_put_sample_proc) (void *buf, int32_t sample, int nch, int ch, int i);
+
+typedef struct {
+  SAD_get_sample_proc get_sample;
+  SAD_put_sample_proc put_sample;
+} SAD_buffer_ops;
+
+/* private data */
+typedef struct {} SAD_dither_t;
+
+void SAD_dither_init_rand(uint32_t seed);
+
+SAD_dither_t* SAD_dither_init(SAD_buffer_format *inbuf_format, SAD_buffer_format *outbuf_format, int *error);
+int SAD_dither_free(SAD_dither_t* state);
+int SAD_dither_process_buffer (SAD_dither_t *state, void *inbuf, void *outbuf, int frames);
+int SAD_dither_apply_replaygain (SAD_dither_t *state, SAD_replaygain_info *rg_info, SAD_replaygain_mode *mode);
+int SAD_dither_set_scale (SAD_dither_t *state, float scale);
+int SAD_dither_set_dither (SAD_dither_t *state, int dither);
+
+#endif
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/dither_ops.c	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,716 @@
+/* Scale & Dither library (libSAD)
+ * High-precision bit depth converter with ReplayGain support
+ *
+ * Copyright (c) 2007-2008 Eugene Zagidullin (e.asphyx@gmail.com)
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation; either version 2 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+ */
+
+/* WARNING: reading this can damage your brain */
+
+#include "common.h"
+#include "dither_ops.h"
+#include "dither.h"
+
+#define SAD_GET_LE16(a) ( (uint16_t)(((uint8_t*)(a))[0])      | (uint16_t)(((uint8_t*)(a))[1]) << 8 )
+#define SAD_GET_BE16(a) ( (uint16_t)(((uint8_t*)(a))[0]) << 8 | (uint16_t)(((uint8_t*)(a))[1]) )
+
+#define SAD_GET_LE32(a) ( (uint32_t)(((uint8_t*)(a))[0])       | (uint32_t)(((uint8_t*)(a))[1]) << 8 | \
+                          (uint32_t)(((uint8_t*)(a))[2]) << 16 | (uint32_t)(((uint8_t*)(a))[3]) << 24 )
+#define SAD_GET_BE32(a) ( (uint32_t)(((uint8_t*)(a))[0]) << 24 | (uint32_t)(((uint8_t*)(a))[1]) << 16 | \
+                          (uint32_t)(((uint8_t*)(a))[2]) << 8 | (uint32_t)(((uint8_t*)(a))[3]) )
+
+#define SAD_PUT_LE16(a,b) { \
+          ((uint8_t*)(a))[0] = (uint8_t)((uint32_t)(b) &  0x000000ff);        \
+          ((uint8_t*)(a))[1] = (uint8_t)(((uint32_t)(b) & 0x0000ff00) >> 8); \
+        }
+
+#define SAD_PUT_BE16(a,b) { \
+          ((uint8_t*)(a))[0] = (uint8_t)(((uint32_t)(b) & 0x0000ff00) >> 8); \
+          ((uint8_t*)(a))[1] = (uint8_t)((uint32_t)(b) &  0x000000ff);        \
+        }
+
+#define SAD_PUT_LE32(a,b) { \
+          ((uint8_t*)(a))[0] = (uint8_t)((uint32_t)(b) &  0x000000ff);        \
+          ((uint8_t*)(a))[1] = (uint8_t)(((uint32_t)(b) & 0x0000ff00) >> 8); \
+          ((uint8_t*)(a))[2] = (uint8_t)(((uint32_t)(b) & 0x00ff0000) >> 16); \
+          ((uint8_t*)(a))[3] = (uint8_t)(((uint32_t)(b) & 0xff000000) >> 24); \
+        }
+
+#define SAD_PUT_BE32(a,b) { \
+          ((uint8_t*)(a))[0] = (uint8_t)(((uint32_t)(b) & 0xff000000) >> 24); \
+          ((uint8_t*)(a))[1] = (uint8_t)(((uint32_t)(b) & 0x00ff0000) >> 16); \
+          ((uint8_t*)(a))[2] = (uint8_t)(((uint32_t)(b) & 0x0000ff00) >> 8); \
+          ((uint8_t*)(a))[3] = (uint8_t)((uint32_t)(b) &  0x000000ff);       \
+        }
+
+
+/* buffer ops: */
+/****************************************************************************************************************
+ * 8-bit                                                                                                        *
+ ****************************************************************************************************************/
+
+/* signed */
+static int32_t get_s8_i_sample (void *buf, int nch, int ch, int i) {
+  return ((int8_t*)buf)[i*nch+ch];
+}
+
+static int32_t get_s8_s_sample (void *buf, int nch, int ch, int i) {
+  return ((int8_t**)buf)[ch][i];
+}
+
+static void put_s8_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((int8_t*)buf)[i*nch+ch] = (int8_t)sample;
+}
+
+static void put_s8_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((int8_t**)buf)[ch][i] = (int8_t)sample;
+}
+
+/* unsigned */
+static int32_t get_u8_i_sample (void *buf, int nch, int ch, int i) {
+  return (int32_t)(((uint8_t*)buf)[i*nch+ch]) - 128;
+}
+
+static int32_t get_u8_s_sample (void *buf, int nch, int ch, int i) {
+  return (int32_t)(((uint8_t**)buf)[ch][i]) - 128;
+}
+
+static void put_u8_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((uint8_t*)buf)[i*nch+ch] = (uint8_t)sample + 128;
+}
+
+static void put_u8_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((uint8_t**)buf)[ch][i] = (uint8_t)sample + 128;
+}
+
+static SAD_buffer_ops buf_s8_i_ops = {
+  &get_s8_i_sample,
+  &put_s8_i_sample
+};
+
+static SAD_buffer_ops buf_s8_s_ops = {
+  &get_s8_s_sample,
+  &put_s8_s_sample
+};
+
+static SAD_buffer_ops buf_u8_i_ops = {
+  &get_u8_i_sample,
+  &put_u8_i_sample
+};
+
+static SAD_buffer_ops buf_u8_s_ops = {
+  &get_u8_s_sample,
+  &put_u8_s_sample
+};
+
+/****************************************************************************************************************
+ * 16-bit                                                                                                       *
+ ****************************************************************************************************************/
+
+/* signed */
+static int32_t get_s16_i_sample (void *buf, int nch, int ch, int i) {
+  return (int32_t)(((int16_t*)buf)[i*nch+ch]);
+}
+
+static int32_t get_s16_s_sample (void *buf, int nch, int ch, int i) {
+  return (int32_t)(((int16_t**)buf)[ch][i]);
+}
+
+static void put_s16_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((int16_t*)buf)[i*nch+ch] = (int16_t)sample;
+}
+
+static void put_s16_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((int16_t**)buf)[ch][i] = (int16_t)sample;
+}
+
+/* unsigned */
+static int32_t get_u16_i_sample (void *buf, int nch, int ch, int i) {
+  return ((int32_t)(((uint16_t*)buf)[i*nch+ch])) - 32768;
+}
+
+static int32_t get_u16_s_sample (void *buf, int nch, int ch, int i) {
+  return ((int32_t)(((uint16_t**)buf)[ch][i])) - 32768;
+}
+
+static void put_u16_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((uint16_t*)buf)[i*nch+ch] = (uint16_t)(sample + 32768);
+}
+
+static void put_u16_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((uint16_t**)buf)[ch][i] = (uint16_t)(sample + 32768);
+}
+
+/* LE: signed */
+static int32_t get_s16_le_i_sample (void *buf, int nch, int ch, int i) {
+  int16_t *tmp = (int16_t*)buf+i*nch+ch;
+  return (int16_t)SAD_GET_LE16(tmp);
+}
+
+static int32_t get_s16_le_s_sample (void *buf, int nch, int ch, int i) {
+  int16_t *tmp = ((int16_t**)buf)[ch]+i;
+  return (int16_t)SAD_GET_LE16(tmp);
+}
+
+static void put_s16_le_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int16_t *tmp = (int16_t*)buf+i*nch+ch;
+  SAD_PUT_LE16(tmp, sample);
+}
+
+static void put_s16_le_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int16_t *tmp = ((int16_t**)buf)[ch]+i;
+  SAD_PUT_LE16(tmp, sample);
+}
+
+/* BE: signed */
+static int32_t get_s16_be_i_sample (void *buf, int nch, int ch, int i) {
+  int16_t *tmp = (int16_t*)buf+i*nch+ch;
+  return (int16_t)SAD_GET_BE16(tmp);
+}
+
+static int32_t get_s16_be_s_sample (void *buf, int nch, int ch, int i) {
+  int16_t *tmp = ((int16_t**)buf)[ch]+i;
+  return (int16_t)SAD_GET_BE16(tmp);
+}
+
+static void put_s16_be_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int16_t *tmp = (int16_t*)buf+i*nch+ch;
+  SAD_PUT_BE16(tmp, sample);
+}
+
+static void put_s16_be_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int16_t *tmp = ((int16_t**)buf)[ch]+i;
+  SAD_PUT_BE16(tmp, sample);
+}
+
+/* LE: unsigned */
+static int32_t get_u16_le_i_sample (void *buf, int nch, int ch, int i) {
+  int16_t *tmp = (int16_t*)buf+i*nch+ch;
+  return (int16_t)SAD_GET_LE16(tmp) - 32768;
+}
+
+static int32_t get_u16_le_s_sample (void *buf, int nch, int ch, int i) {
+  int16_t *tmp = ((int16_t**)buf)[ch]+i;
+  return (int16_t)SAD_GET_LE16(tmp) - 32768;
+}
+
+static void put_u16_le_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int16_t *tmp = (int16_t*)buf+i*nch+ch;
+  SAD_PUT_LE16(tmp, sample + 32768);
+}
+
+static void put_u16_le_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int16_t *tmp = ((int16_t**)buf)[ch]+i;
+  SAD_PUT_LE16(tmp, sample + 32768);
+}
+
+/* BE: unsigned */
+static int32_t get_u16_be_i_sample (void *buf, int nch, int ch, int i) {
+  int16_t *tmp = (int16_t*)buf+i*nch+ch;
+  return (int16_t)SAD_GET_BE16(tmp) - 32768;
+}
+
+static int32_t get_u16_be_s_sample (void *buf, int nch, int ch, int i) {
+  int16_t *tmp = ((int16_t**)buf)[ch]+i;
+  return (int16_t)SAD_GET_BE16(tmp) - 32768;
+}
+
+static void put_u16_be_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int16_t *tmp = (int16_t*)buf+i*nch+ch;
+  SAD_PUT_BE16(tmp, sample + 32768);
+}
+
+static void put_u16_be_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int16_t *tmp = ((int16_t**)buf)[ch]+i;
+  SAD_PUT_BE16(tmp, sample + 32768);
+}
+
+
+static SAD_buffer_ops buf_s16_i_ops = {
+  &get_s16_i_sample,
+  &put_s16_i_sample
+};
+
+static SAD_buffer_ops buf_s16_s_ops = {
+  &get_s16_s_sample,
+  &put_s16_s_sample
+};
+
+static SAD_buffer_ops buf_s16_le_i_ops = {
+  &get_s16_le_i_sample,
+  &put_s16_le_i_sample
+};
+
+static SAD_buffer_ops buf_s16_le_s_ops = {
+  &get_s16_le_s_sample,
+  &put_s16_le_s_sample
+};
+
+static SAD_buffer_ops buf_s16_be_i_ops = {
+  &get_s16_be_i_sample,
+  &put_s16_be_i_sample
+};
+
+static SAD_buffer_ops buf_s16_be_s_ops = {
+  &get_s16_be_s_sample,
+  &put_s16_be_s_sample
+};
+
+/* unsigned */
+
+static SAD_buffer_ops buf_u16_i_ops = {
+  &get_u16_i_sample,
+  &put_u16_i_sample
+};
+
+static SAD_buffer_ops buf_u16_s_ops = {
+  &get_u16_s_sample,
+  &put_u16_s_sample
+};
+
+static SAD_buffer_ops buf_u16_le_i_ops = {
+  &get_u16_le_i_sample,
+  &put_u16_le_i_sample
+};
+
+static SAD_buffer_ops buf_u16_le_s_ops = {
+  &get_u16_le_s_sample,
+  &put_u16_le_s_sample
+};
+
+static SAD_buffer_ops buf_u16_be_i_ops = {
+  &get_u16_be_i_sample,
+  &put_u16_be_i_sample
+};
+
+static SAD_buffer_ops buf_u16_be_s_ops = {
+  &get_u16_be_s_sample,
+  &put_u16_be_s_sample
+};
+
+/****************************************************************************************************************
+ * 24-bit                                                                                                       *
+ ****************************************************************************************************************/
+
+/*expand 24-bit signed value to 32-bit*/
+#define EXPAND_S24_TO_32(x) (((int32_t)(((x) & 0x00ffffff) << 8)) >> 8)
+#define EXPAND_U24_TO_32(x) ((int32_t)(x) & 0x00ffffff)
+
+/* signed */
+static int32_t get_s24_i_sample (void *buf, int nch, int ch, int i) {
+  return (int32_t)EXPAND_S24_TO_32(((int32_t*)buf)[i*nch+ch]);
+}
+
+static int32_t get_s24_s_sample (void *buf, int nch, int ch, int i) {
+  return (int32_t)EXPAND_S24_TO_32(((int32_t**)buf)[ch][i]);
+}
+
+static void put_s24_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((int32_t*)buf)[i*nch+ch] = (int32_t)sample & 0x00ffffff;
+}
+
+static void put_s24_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((int32_t**)buf)[ch][i] = (int32_t)sample & 0x00ffffff;
+}
+
+/* LE signed */
+
+static int32_t get_s24_le_i_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  return (int32_t)EXPAND_S24_TO_32(SAD_GET_LE32(tmp));
+}
+
+static int32_t get_s24_le_s_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  return (int32_t)EXPAND_S24_TO_32(SAD_GET_LE32(tmp));
+}
+
+static void put_s24_le_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  SAD_PUT_LE32(tmp, sample & 0x00ffffff);
+}
+
+static void put_s24_le_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  SAD_PUT_LE32(tmp, sample & 0x00ffffff);
+}
+
+/* BE signed */
+
+static int32_t get_s24_be_i_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  return (int32_t)EXPAND_S24_TO_32(SAD_GET_BE32(tmp));
+}
+
+static int32_t get_s24_be_s_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  return (int32_t)EXPAND_S24_TO_32(SAD_GET_BE32(tmp));
+}
+
+static void put_s24_be_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  SAD_PUT_BE32(tmp, sample & 0x00ffffff);
+}
+
+static void put_s24_be_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  SAD_PUT_BE32(tmp, sample & 0x00ffffff);
+}
+
+/* unsigned */
+static int32_t get_u24_i_sample (void *buf, int nch, int ch, int i) {
+  return (int32_t)EXPAND_U24_TO_32(((uint32_t*)buf)[i*nch+ch]) - 8388608;
+}
+
+static int32_t get_u24_s_sample (void *buf, int nch, int ch, int i) {
+  return (int32_t)EXPAND_U24_TO_32(((uint32_t**)buf)[ch][i]) - 8388608;
+}
+
+static void put_u24_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((uint32_t*)buf)[i*nch+ch] = ((uint32_t)sample + 8388608) & 0x00ffffff;
+}
+
+static void put_u24_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((uint32_t**)buf)[ch][i] = ((uint32_t)sample + 8388608) & 0x00ffffff;
+}
+
+/* LE unsigned */
+
+static int32_t get_u24_le_i_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  /*fprintf(stderr, "%d\n", (int32_t)EXPAND_U24_TO_32(SAD_GET_LE32(tmp)) - 8388608);*/
+  return (int32_t)EXPAND_U24_TO_32(SAD_GET_LE32(tmp)) - 8388608;
+}
+
+static int32_t get_u24_le_s_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  return (int32_t)EXPAND_U24_TO_32(SAD_GET_LE32(tmp)) - 8388608;
+}
+
+static void put_u24_le_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  SAD_PUT_LE32(tmp, (uint32_t)(sample + 8388608) & 0x00ffffff);
+}
+
+static void put_u24_le_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  SAD_PUT_LE32(tmp, (uint32_t)(sample + 8388608) & 0x00ffffff);
+}
+
+/* BE unsigned */
+
+static int32_t get_u24_be_i_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  return (int32_t)EXPAND_U24_TO_32(SAD_GET_BE32(tmp)) - 8388608;
+}
+
+static int32_t get_u24_be_s_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  return (int32_t)EXPAND_U24_TO_32(SAD_GET_BE32(tmp)) - 8388608;
+}
+
+static void put_u24_be_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  SAD_PUT_BE32(tmp, (uint32_t)(sample + 8388608) & 0x00ffffff);
+}
+
+static void put_u24_be_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  SAD_PUT_BE32(tmp, (uint32_t)(sample + 8388608) & 0x00ffffff);
+}
+
+static SAD_buffer_ops buf_s24_i_ops = {
+  &get_s24_i_sample,
+  &put_s24_i_sample
+};
+
+static SAD_buffer_ops buf_s24_s_ops = {
+  &get_s24_s_sample,
+  &put_s24_s_sample
+};
+
+static SAD_buffer_ops buf_s24_le_i_ops = {
+  &get_s24_le_i_sample,
+  &put_s24_le_i_sample
+};
+
+static SAD_buffer_ops buf_s24_le_s_ops = {
+  &get_s24_le_s_sample,
+  &put_s24_le_s_sample
+};
+
+static SAD_buffer_ops buf_s24_be_i_ops = {
+  &get_s24_be_i_sample,
+  &put_s24_be_i_sample
+};
+
+static SAD_buffer_ops buf_s24_be_s_ops = {
+  &get_s24_be_s_sample,
+  &put_s24_be_s_sample
+};
+
+static SAD_buffer_ops buf_u24_i_ops = {
+  &get_u24_i_sample,
+  &put_u24_i_sample
+};
+
+static SAD_buffer_ops buf_u24_s_ops = {
+  &get_u24_s_sample,
+  &put_u24_s_sample
+};
+
+static SAD_buffer_ops buf_u24_le_i_ops = {
+  &get_u24_le_i_sample,
+  &put_u24_le_i_sample
+};
+
+static SAD_buffer_ops buf_u24_le_s_ops = {
+  &get_u24_le_s_sample,
+  &put_u24_le_s_sample
+};
+
+static SAD_buffer_ops buf_u24_be_i_ops = {
+  &get_u24_be_i_sample,
+  &put_u24_be_i_sample
+};
+
+static SAD_buffer_ops buf_u24_be_s_ops = {
+  &get_u24_be_s_sample,
+  &put_u24_be_s_sample
+};
+
+/****************************************************************************************************************
+ * 32-bit                                                                                                       *
+ ****************************************************************************************************************/
+
+/* signed */
+static int32_t get_s32_i_sample (void *buf, int nch, int ch, int i) {
+  return ((int32_t*)buf)[i*nch+ch];
+}
+
+static int32_t get_s32_s_sample (void *buf, int nch, int ch, int i) {
+  return ((int32_t**)buf)[ch][i];
+}
+
+static void put_s32_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((int32_t*)buf)[i*nch+ch] = (int32_t)sample;
+}
+
+static void put_s32_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((int32_t**)buf)[ch][i] = (int32_t)sample;
+}
+
+/* LE: signed */
+static int32_t get_s32_le_i_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  return (int32_t)SAD_GET_LE32(tmp);
+}
+
+static int32_t get_s32_le_s_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  return (int32_t)SAD_GET_LE32(tmp);
+}
+
+static void put_s32_le_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  SAD_PUT_LE32(tmp, sample);
+}
+
+static void put_s32_le_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  SAD_PUT_LE32(tmp, sample);
+}
+
+/* BE: signed */
+static int32_t get_s32_be_i_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  return (int32_t)SAD_GET_BE32(tmp);
+}
+
+static int32_t get_s32_be_s_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  return (int32_t)SAD_GET_BE32(tmp);
+}
+
+static void put_s32_be_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  SAD_PUT_BE32(tmp, sample);
+}
+
+static void put_s32_be_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  SAD_PUT_BE32(tmp, sample);
+}
+
+/* unsigned */
+static int32_t get_u32_i_sample (void *buf, int nch, int ch, int i) {
+  return ((int32_t*)buf)[i*nch+ch] - (int32_t)(1L<<31);
+}
+
+static int32_t get_u32_s_sample (void *buf, int nch, int ch, int i) {
+  return ((int32_t**)buf)[ch][i] - (int32_t)(1L<<31);
+}
+
+static void put_u32_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((uint32_t*)buf)[i*nch+ch] = (uint32_t)(sample + (int32_t)(1L<<31));
+}
+
+static void put_u32_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  ((uint32_t**)buf)[ch][i] = (uint32_t)(sample + (int32_t)(1L<<31));
+}
+
+/* LE: unsigned */
+static int32_t get_u32_le_i_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  return (int32_t)SAD_GET_LE32(tmp) - (int32_t)(1L<<31);
+}
+
+static int32_t get_u32_le_s_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  return (int32_t)SAD_GET_LE32(tmp) - (int32_t)(1L<<31);
+}
+
+static void put_u32_le_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  SAD_PUT_LE32(tmp, sample + (int32_t)(1L<<31));
+}
+
+static void put_u32_le_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  SAD_PUT_LE32(tmp, sample + (int32_t)(1L<<31));
+}
+
+/* BE: unsigned */
+static int32_t get_u32_be_i_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  return (int32_t)SAD_GET_BE32(tmp) - (int32_t)(1L<<31);
+}
+
+static int32_t get_u32_be_s_sample (void *buf, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  return (int32_t)SAD_GET_BE32(tmp) - (int32_t)(1L<<31);
+}
+
+static void put_u32_be_i_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = (int32_t*)buf+i*nch+ch;
+  SAD_PUT_BE32(tmp, sample + (int32_t)(1L<<31));
+}
+
+static void put_u32_be_s_sample (void *buf, int32_t sample, int nch, int ch, int i) {
+  int32_t *tmp = ((int32_t**)buf)[ch]+i;
+  SAD_PUT_BE32(tmp, sample + (int32_t)(1L<<31));
+}
+
+static SAD_buffer_ops buf_s32_i_ops = {
+  &get_s32_i_sample,
+  &put_s32_i_sample
+};
+
+static SAD_buffer_ops buf_s32_s_ops = {
+  &get_s32_s_sample,
+  &put_s32_s_sample
+};
+
+static SAD_buffer_ops buf_s32_le_i_ops = {
+  &get_s32_le_i_sample,
+  &put_s32_le_i_sample
+};
+
+static SAD_buffer_ops buf_s32_le_s_ops = {
+  &get_s32_le_s_sample,
+  &put_s32_le_s_sample
+};
+
+static SAD_buffer_ops buf_s32_be_i_ops = {
+  &get_s32_be_i_sample,
+  &put_s32_be_i_sample
+};
+
+static SAD_buffer_ops buf_s32_be_s_ops = {
+  &get_s32_be_s_sample,
+  &put_s32_be_s_sample
+};
+
+static SAD_buffer_ops buf_u32_i_ops = {
+  &get_u32_i_sample,
+  &put_u32_i_sample
+};
+
+static SAD_buffer_ops buf_u32_s_ops = {
+  &get_u32_s_sample,
+  &put_u32_s_sample
+};
+
+static SAD_buffer_ops buf_u32_le_i_ops = {
+  &get_u32_le_i_sample,
+  &put_u32_le_i_sample
+};
+
+static SAD_buffer_ops buf_u32_le_s_ops = {
+  &get_u32_le_s_sample,
+  &put_u32_le_s_sample
+};
+
+static SAD_buffer_ops buf_u32_be_i_ops = {
+  &get_u32_be_i_sample,
+  &put_u32_be_i_sample
+};
+
+static SAD_buffer_ops buf_u32_be_s_ops = {
+  &get_u32_be_s_sample,
+  &put_u32_be_s_sample
+};
+
+static SAD_buffer_ops *SAD_buffer_optable[SAD_SAMPLE_MAX][SAD_CHORDER_MAX] = {
+  {&buf_s8_i_ops,     &buf_s8_s_ops},	  /* SAD_SAMPLE_S8     */
+  {&buf_u8_i_ops,     &buf_u8_s_ops},	  /* SAD_SAMPLE_U8     */
+
+  {&buf_s16_i_ops,    &buf_s16_s_ops},	  /* SAD_SAMPLE_S16    */
+  {&buf_s16_le_i_ops, &buf_s16_le_s_ops}, /* SAD_SAMPLE_S16_LE */
+  {&buf_s16_be_i_ops, &buf_s16_be_s_ops}, /* SAD_SAMPLE_S16_BE */
+  {&buf_u16_i_ops,    &buf_u16_s_ops},	  /* SAD_SAMPLE_U16    */
+  {&buf_u16_le_i_ops, &buf_u16_le_s_ops}, /* SAD_SAMPLE_U16_LE */
+  {&buf_u16_be_i_ops, &buf_u16_be_s_ops}, /* SAD_SAMPLE_U16_BE */
+
+  {&buf_s24_i_ops,    &buf_s24_s_ops},	  /* SAD_SAMPLE_S24    */
+  {&buf_s24_le_i_ops, &buf_s24_le_s_ops}, /* SAD_SAMPLE_S24_LE */
+  {&buf_s24_be_i_ops, &buf_s24_be_s_ops}, /* SAD_SAMPLE_S24_BE */
+  {&buf_u24_i_ops,    &buf_u24_s_ops},	  /* SAD_SAMPLE_U24    */
+  {&buf_u24_le_i_ops, &buf_u24_le_s_ops}, /* SAD_SAMPLE_U24_LE */
+  {&buf_u24_be_i_ops, &buf_u24_be_s_ops}, /* SAD_SAMPLE_U24_BE */
+
+  {&buf_s32_i_ops,    &buf_s32_s_ops}, 	  /* SAD_SAMPLE_S32    */
+  {&buf_s32_le_i_ops, &buf_s32_le_s_ops}, /* SAD_SAMPLE_S32_LE */
+  {&buf_s32_be_i_ops, &buf_s32_be_s_ops}, /* SAD_SAMPLE_S32_BE */
+  {&buf_u32_i_ops,    &buf_u32_s_ops}, 	  /* SAD_SAMPLE_U32    */
+  {&buf_u32_le_i_ops, &buf_u32_le_s_ops}, /* SAD_SAMPLE_U32_LE */
+  {&buf_u32_be_i_ops, &buf_u32_be_s_ops}, /* SAD_SAMPLE_U32_BE */
+
+  {&buf_s32_i_ops,    &buf_s32_s_ops}, 	  /* SAD_SAMPLE_FIXED32*/
+
+  {NULL,              NULL}		  /* SAD_SAMPLE_FLOAT  */
+};
+
+SAD_buffer_ops* SAD_assign_buf_ops (SAD_buffer_format *format) {
+#ifdef DEBUG
+  printf("f: SAD_assign_buf_ops\n");
+#endif
+  if (format->sample_format < SAD_SAMPLE_MAX && format->channels_order < SAD_CHORDER_MAX) {
+    return SAD_buffer_optable[format->sample_format][format->channels_order];
+  } else {
+    return NULL;
+  }
+}
+
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/dither_ops.h	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,28 @@
+/* Scale & Dither library (libSAD)
+ * High-precision bit depth converter with ReplayGain support
+ *
+ * Copyright (c) 2007-2008 Eugene Zagidullin (e.asphyx@gmail.com)
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation; either version 2 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+ */
+
+#ifndef DITHER_OPS_H
+#define DITHER_OPS_H
+
+#include "common.h"
+#include "dither.h"
+
+SAD_buffer_ops* SAD_assign_buf_ops (SAD_buffer_format *format);
+
+#endif /*DITHER_OPS_H*/
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/libSAD.h	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,26 @@
+/* Scale & Dither library (libSAD)
+ * High-precision bit depth converter with ReplayGain support
+ *
+ * Copyright (c) 2007-2008 Eugene Zagidullin (e.asphyx@gmail.com)
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation; either version 2 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+ */
+
+#ifndef LIBSAD_H
+#define LIBSAD_H
+
+#include "common.h"
+#include "dither.h"
+
+#endif /* LIBSAD_H */
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/noicegen.c	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,46 @@
+#include <stdio.h>
+#include <assert.h>
+#include "../../config.h"
+
+#ifdef HAVE_SSE2
+#  define SSE2 1
+#endif
+
+#ifdef HAVE_ALTIVEC
+#  define ALTIVEC 1
+#endif
+
+#define MEXP 19937
+
+/*#include "SFMT.h"*/
+#include "SFMT.c"
+
+#include "noicegen.h"
+
+int triangular_dither_noise(int nbits)
+{
+    // parameter nbits : the peak-to-peak amplitude desired (in bits)
+    //  use with nbits set to    2 + nber of bits to be trimmed.
+    // (because triangular is made from two uniformly distributed processes,
+    // it starts at 2 bits peak-to-peak amplitude)
+    // see The Theory of Dithered Quantization by Robert Alexander Wannamaker
+    // for complete proof of why that's optimal
+
+    int v = (gen_rand32() / 2 - gen_rand32() / 2);   // in ]-2^31, 2^31[
+    //int signe = (v>0) ? 1 : -1;
+    int P = 1 << (32 - nbits); // the power of 2
+    v /= P;
+    // now v in ]-2^(nbits-1), 2^(nbits-1) [
+
+    return v;
+}
+
+float triangular_dither_noise_f() {
+  // Сonditionally assume we have 16 bits in fractional part
+  // Please, check it thoroughly: is this assumption correct in floatin-point arithmetic?
+  return (float) triangular_dither_noise(17) / 65536.0;
+}
+
+void noicegen_init_rand(uint32_t seed) {
+  init_gen_rand(seed);
+}
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/libSAD/noicegen.h	Mon Feb 04 05:55:43 2008 +0300
@@ -0,0 +1,10 @@
+#ifndef NOICEGEN_H
+#define NOICEGEN_H
+
+#include <inttypes.h>
+
+int triangular_dither_noise(int nbits);
+float triangular_dither_noise_f(void);
+void noicegen_init_rand(uint32_t seed);
+
+#endif /*NOICEGEN_H*/