src/main.h

#ifndef MAIN_H
#define MAIN_H

#include "config.h"
#include "logger.h"
#include <algorithm>
#include <limits>
#include <boost/foreach.hpp>
#include <cmath>
#include <cassert>
#include <limits>

#define foreach BOOST_FOREACH

class MapParser;
class TargetTable;
class BiasBoss;
class MismatchTable;
class LengthDistribution;

extern Logger logger;
extern bool running;
extern bool burned_out;
extern bool edit_detect;
extern size_t max_indel_size;
enum Direction { FR, RF, R, F, BOTH };
extern Direction direction;
extern size_t max_read_len;
const size_t NUM_NUCS = 4;
const char NUCS[] = {'A','C','G','T'};
const double LOG_0 = HUGE_VAL;
const double LOG_1 = 0;
const double LOG_QUARTER = log(0.25);
const double EPSILON = 0.000001;
const double LOG_EPSILON = log(EPSILON);
const double LOG_MAX = log(std::numeric_limits<double>::max());

inline bool approx_eq(double a, double b, double eps=EPSILON) {
  return fabs(a-b) <= eps;
}

inline double log_add(double x, double y) {
  if (fabs(x) == LOG_0) {
    return y;
  }
  if (fabs(y) == LOG_0) {
    return x;
  }

  if (y > x) {
    std::swap(x,y);
  }

  double sum = x+log(1+exp(y-x));
  return sum;
}
inline double log_sub(double x, double y) {
  if (fabs(y) == LOG_0) {
    return x;
  }
  
  // Have to be careful of numerical issues, so we allow y to be slightly
  // greater than x.
  if (x <= y) {
    assert(approx_eq(x, y));
    return LOG_0;
  }
  

  double diff = x+log(1-exp(y-x));
  return diff;
}
inline double islzero(double x)
{
  return (fabs(x) == LOG_0);
}

inline double sexp(double x)
{
  if (islzero(x)) {
    return 0.0;
  }
  return exp(x);
}

#endif