stochasticencoder.h

Go to the documentation of this file.

#ifndef DYSCO_STOCHASTIC_ENCODER_H
#define DYSCO_STOCHASTIC_ENCODER_H

#include "uvector.h"

#include <algorithm>
#include <cmath>
#include <cstring>

namespace dyscostman {

template <typename ValueType = float>
class StochasticEncoder {
 public:
  StochasticEncoder(size_t quantCount, ValueType stddev,
                    bool gaussianMapping = true);

  static StochasticEncoder StudentTEncoder(size_t quantCount, double nu,
                                           double rms) {
    StochasticEncoder<ValueType> encoder(quantCount);
    encoder.initializeStudentT(nu, rms);
    return encoder;
  }

  static StochasticEncoder TruncatedGausEncoder(size_t quantCount, double trunc,
                                                double rms) {
    StochasticEncoder<ValueType> encoder(quantCount);
    encoder.initializeTruncatedGaussian(trunc, rms);
    return encoder;
  }

  typedef unsigned symbol_t;

  typedef ValueType value_t;

  symbol_t Encode(ValueType value) const {
    if (std::isfinite(value))
      return _encDictionary.symbol(_encDictionary.lower_bound(value));
    else
      return QuantizationCount() - 1;
  }

  static std::uniform_int_distribution<unsigned> GetDitherDistribution() {
    return std::uniform_int_distribution<unsigned>(0, ((1u << 31) - 1));
  }

  symbol_t EncodeWithDithering(ValueType value, unsigned ditherValue) const {
    if (std::isfinite(value)) {
      const typename Dictionary::const_iterator lowerBound =
          _decDictionary.lower_bound(value);
      if (lowerBound == _decDictionary.begin())
        return _decDictionary.symbol(lowerBound);
      if (lowerBound == _decDictionary.end())
        return _decDictionary.symbol(lowerBound - 1);
      const ValueType rightValue = _decDictionary.value(lowerBound);
      const ValueType leftValue = _decDictionary.value(lowerBound - 1);

      ValueType ditherMark =
          ValueType(1u << 31) * (value - leftValue) / (rightValue - leftValue);
      if (ditherMark > ditherValue)
        return _decDictionary.symbol(lowerBound);
      else
        return _decDictionary.symbol(lowerBound - 1);
    } else {
      return _encDictionary.size();
    }
  }

  value_t RightBoundary(symbol_t symbol) const {
    if (symbol != _encDictionary.size())
      return _encDictionary.value(symbol);
    else
      return 0.0;
  }

  ValueType Decode(symbol_t symbol) const {
    return _decDictionary.value(symbol);
  }

  size_t QuantizationCount() const { return _decDictionary.size() + 1; }

  ValueType MaxQuantity() const { return _decDictionary.largest_value(); }

  ValueType MinQuantity() const { return _decDictionary.smallest_value(); }

 private:
  explicit StochasticEncoder(size_t quantCount)
      : _encDictionary(quantCount - 1), _decDictionary(quantCount - 1) {}

  void initializeStudentT(double nu, double rms);

  void initializeTruncatedGaussian(double truncationValue, double rms);

  class Dictionary {
   public:
    typedef value_t *iterator;
    typedef const value_t *const_iterator;

    Dictionary() : _values() {}

    explicit Dictionary(size_t size) : _values(size) {}

    void reserve(size_t size) { _values.reserve(size); }

    void resize(size_t size) { _values.resize(size); }

    const_iterator lower_bound(value_t val) const {
      size_t p = 0, q = _values.size();
      size_t m = (p + q) / 2;
      if (_values[m] <= val)
        p = m;
      else
        q = m;
      while (p + 1 != q) {
        size_t m = (p + q) / 2;
        if (_values[m] <= val)
          p = m;
        else
          q = m;
      }
      return (_values[p] < val) ? (&_values[q]) : (&_values[p]);
    }

    const_iterator lower_bound_fast(value_t val) const {
      size_t p = 0, q = _values.size();
      while (p + 1 != q) {
        size_t m = (p + q) / 2;
        if (_values[m] <= val)
          p = m;
        else
          q = m;
      }
      return (_values[p] < val) ? (&_values[q]) : (&_values[p]);
    }

    const_iterator lower_bound_slow(value_t val) const {
      const value_t *p = &*_values.begin(), *q = p + _values.size();
      while (p + 1 != q) {
        // This is a bit inefficient, but (p + q)/2 was not allowed, because
        // operator+(ptr,ptr) is not allowed.
        const value_t *m = p + (q - p) / 2;
        if (*m <= val)
          p = m;
        else
          q = m;
      }
      return p;
    }

    iterator begin() { return &*_values.begin(); }
    const_iterator begin() const { return &*_values.begin(); }
    const_iterator end() const { return &*_values.end(); }
    symbol_t symbol(const_iterator iter) const { return (iter - begin()); }
    symbol_t largest_symbol() const { return _values.size() - 1; }
    value_t value(const_iterator iter) const { return *iter; }
    value_t value(symbol_t sym) const { return _values[sym]; }
    value_t largest_value() const { return _values.back(); }
    value_t smallest_value() const { return _values.front(); }
    size_t size() const { return _values.size(); }
    size_t capacity(size_t) const { return _values.capacity(); }

   private:
    ao::uvector<value_t> _values;
  };

  typedef long double num_t;

  static num_t cumulative(num_t x);
  static num_t invCumulative(num_t c, num_t err = num_t(1e-13));

  Dictionary _encDictionary;
  Dictionary _decDictionary;
};

}  // namespace dyscostman

#endif