src/targets.h

#ifndef TRANSCRIPTS_H
#define TRANSCRIPTS_H

#include <boost/scoped_ptr.hpp>
#include "boost/shared_ptr.hpp"
#include <boost/thread.hpp>
#include <boost/unordered_map.hpp>
#include <iostream>
#include <fstream>
#include <map>
#include <string>
#include <vector>
#include "main.h"
#include "bundles.h"
#include "sequence.h"

class LengthDistribution;
class FragHit;
class BiasBoss;
class MismatchTable;
class Librarian;
class HaplotypeHandler;
class TargetTable;

struct RoundParams {
  double mass;
  double ambig_mass;
  double tot_ambig_mass;
  double mass_var;
  double var_sum;
  /***
   * A shared pointer to the target's HaplotypeHandler.  Null if target has no
   * haplotype partner.
   **/
  boost::shared_ptr<HaplotypeHandler> haplotype;
  RoundParams() : mass(LOG_0), ambig_mass(LOG_0), tot_ambig_mass(LOG_0),
                  mass_var(LOG_0), var_sum(LOG_0) {}
};

typedef size_t TargID;

class Target {
  friend class HaplotypeHandler;
  friend class TargetTable;
  const Librarian* _libs;
  TargID _id;
  std::string _name;
  SequenceFwd _seq_f;
  SequenceRev _seq_r;
  double _alpha;
  RoundParams _curr_params;
  RoundParams _last_params;
  RoundParams* _ret_params;
  size_t _uniq_counts;
  size_t _tot_counts;
  double _init_pseudo_mass;
  Bundle* _bundle;
  mutable boost::mutex _mutex;
  boost::scoped_ptr<std::vector<float> > _start_bias;
  boost::scoped_ptr<std::vector<float> > _start_bias_buffer;
  boost::scoped_ptr<std::vector<float> > _end_bias;
  boost::scoped_ptr<std::vector<float> > _end_bias_buffer;
  double _avg_bias;
  double _avg_bias_buffer;
  double _cached_eff_len;
  double _cached_eff_len_buffer;
  bool _solvable;

public:
  Target(TargID id, const std::string& name, const std::string& seq,
         bool prob_seq, double alpha, const Librarian* libs,
         const BiasBoss* known_bias_boss, const LengthDistribution* known_fld);
  void lock() const { _mutex.lock(); }
  void unlock() const { _mutex.unlock(); }
  const std::string& name() const { return _name; }
  TargID id() const { return _id; }
  const Sequence& seq(bool rev=false) const {
    if (rev) {
      return _seq_r;
    }
    return _seq_f;
  }
  Sequence& seq(bool rev) {
    if (rev) {
      return _seq_r;
    }
    return _seq_f;
  }
  void haplotype(boost::shared_ptr<HaplotypeHandler> hh) {
    _curr_params.haplotype = hh;
  }
  void alpha(double alpha) { _alpha = log(alpha); }
  size_t length() const { return _seq_f.length(); }
  double rho() const;
  double mass(bool with_pseudo=true) const;
  double mass_var() const;
  double var_sum() const { return _ret_params->var_sum; }
  double tot_ambig_mass() const { return _ret_params->tot_ambig_mass; }
  void round_reset();
  size_t tot_counts() const { return _tot_counts; }
  size_t uniq_counts() const { return _uniq_counts; }
  Bundle* bundle() const { return _bundle; }
  void bundle(Bundle* b) { _bundle = b; }
  void add_hit(const FragHit& h, double v, double mass);
  void incr_counts(bool uniq, size_t incr_amt = 1) {
    if (uniq) {
      _solvable = true;
    }
    _tot_counts += incr_amt;
    _uniq_counts += incr_amt * uniq;
  }
  double sample_likelihood(bool with_pseudo,
                           const std::vector<const Target*>* neighbors = NULL)
  const;
  double align_likelihood(const FragHit& frag) const;
  double est_effective_length(const LengthDistribution* fld = NULL,
                              bool with_bias=true) const;
  double cached_effective_length(bool with_bias=true) const;
  void update_target_bias_buffer(const BiasBoss* bias_table = NULL,
                                 const LengthDistribution* fld = NULL);
  void swap_bias_parameters();
  bool solvable() const { return _solvable; }
  void solvable(bool s) { _solvable = s; }
};

class HaplotypeHandler {
  std::vector<const Target*> _targets;
  FrequencyMatrix<double> _haplo_taus;
  std::string _frag_name_buff;
  std::vector<double> _align_likelihoods_buff;
  std::vector<double> _masses_buff;
  bool _committed;
  size_t find_target(const Target* targ);
  void commit_buffer();
 public:
  HaplotypeHandler(std::vector<Target*> targets, double alpha);
  double get_mass(const Target* targ, bool with_pseudo);
  void update_mass(const Target* targ, const std::string& frag_name,
                   double align_likelihood, double mass);
};


typedef std::vector<Target*> TransMap;
typedef boost::unordered_map<std::string, size_t> TransIndex;
typedef boost::unordered_map<size_t, float> CovarMap;
typedef boost::unordered_map<std::string, double> AlphaMap;
typedef boost::unordered_set<std::vector<Target*> > HaplotypeSet;

class TargetTable {
  const Librarian* _libs;
  TransMap _targ_map;
  BundleTable _bundle_table;
  CovarTable _covar_table;
  HaplotypeSet _haplotype_groups;
  double _total_fpb;
  mutable boost::mutex _fpb_mut;

  void add_targ(const std::string& name, const std::string& seq, bool prob_seqs,
                bool known_aux_params, double alpha,
                const TransIndex& targ_index, const TransIndex& targ_lengths);

public:
  TargetTable(std::string targ_fasta_file, std::string haplotype_file,
              bool prob_seqs, bool known_aux_params, double alpha,
              const AlphaMap* alpha_map, const Librarian* libs);
  ~TargetTable();
  Target* get_targ(TargID id);
  void round_reset();
  size_t size() const { return _targ_map.size(); }
  double total_fpb() const;
  void update_total_fpb(double incr_amt);
  void update_covar(TargID targ1, TargID targ2, double covar) {
    _covar_table.increment(targ1, targ2, covar);
  }
  double get_covar(TargID targ1, TargID targ2) {
    return _covar_table.get(targ1, targ2);
  }
  size_t covar_size() const { return _covar_table.size(); }
  Bundle* merge_bundles(Bundle* b1, Bundle* b2);
  size_t num_bundles() const { return _bundle_table.size(); }
  void masses_to_counts();
  void output_results(std::string output_dir, size_t tot_counts,
                      bool output_varcov=false, bool output_rdds=false);
  void asynch_bias_update(boost::mutex* mutex);
  void enable_bundle_threadsafety() { _bundle_table.threadsafe_mode(true); }
  void disable_bundle_threadsafety() { _bundle_table.threadsafe_mode(false); }
  void collapse_bundles() { _bundle_table.collapse(); }
};

#endif