1 /*
2 * multireads.cpp
3 * cufflinks
4 *
5 * Created by Adam Roberts on 3/6/11.
6 * Copyright 2011 Adam Roberts. All rights reserved.
7 *
8 */
9
10
11 #include "hits.h"
12 #include "multireads.h"
13
add_hit(RefID r_id,int left,int right)14 void MultiRead::add_hit(RefID r_id, int left, int right)
15 {
16 _hits.push_back(MultiHit(r_id, left, right));
17 }
18
get_hit(RefID r_id,int left,int right)19 MultiHit* MultiRead::get_hit(RefID r_id, int left, int right)
20 {
21 for (size_t i = 0; i < num_hits(); ++i)
22 {
23 MultiHit& hit = _hits[_curr_index];
24 if (hit.r_id == r_id && hit.left == left && hit.right == right)
25 {
26 return &hit;
27 }
28 _curr_index = (_curr_index + 1) % num_hits();
29 }
30 fprintf(stderr, "\nWARNING: Multi-Hit not found (%d,%d).\n", left, right);
31 return NULL;
32 }
33
add_expr(RefID r_id,int left,int right,double expr)34 void MultiRead::add_expr(RefID r_id, int left, int right, double expr)
35 {
36 MultiHit* hit = get_hit(r_id, left, right);
37 if (hit)
38 {
39 hit->expr += expr;
40 _tot_expr += expr;
41 }
42 }
43
get_mass(RefID r_id,int left,int right,bool valid_mass)44 double MultiRead::get_mass(RefID r_id, int left, int right, bool valid_mass)
45 {
46 if (!valid_mass)
47 {
48 return 1.0/num_hits();
49 }
50
51 if (_tot_expr == 0.0)
52 return 0.0;
53
54 MultiHit* hit = get_hit(r_id, left, right);
55 if (hit)
56 return hit->expr/_tot_expr;
57 else
58 return 1.0/num_hits();
59 }
60
get_read(InsertID mr_id)61 MultiRead* MultiReadTable::get_read(InsertID mr_id)
62 {
63 MultiReadMap::iterator it;
64 it = _read_map.find(mr_id);
65 if (it == _read_map.end())
66 {
67 return NULL;
68 }
69 else
70 {
71 return &(it->second);
72 }
73 }
74
add_hit(const MateHit & hit)75 void MultiReadTable::add_hit(const MateHit& hit)
76 {
77 add_hit(hit.ref_id(), hit.left(), hit.right(), hit.insert_id(), hit.num_hits());
78 }
79
add_hit(RefID r_id,int left,int right,InsertID mr_id,int exp_num_hits)80 void MultiReadTable::add_hit(RefID r_id, int left, int right, InsertID mr_id, int exp_num_hits)
81 {
82 #if ENABLE_THREADS
83 boost::mutex::scoped_lock lock(_lock);
84 #endif
85 MultiRead* mr = get_read(mr_id);
86 if (!mr)
87 {
88 mr = &((_read_map.insert(std::make_pair(mr_id, MultiRead(mr_id, exp_num_hits)))).first->second);
89 }
90 mr->add_hit(r_id, left, right);
91 }
92
add_expr(const MateHit & hit,double expr)93 void MultiReadTable::add_expr(const MateHit& hit, double expr)
94 {
95 add_expr(hit.ref_id(), hit.left(), hit.right(), hit.insert_id(), expr);
96 }
97
add_expr(RefID r_id,int left,int right,InsertID mr_id,double expr)98 void MultiReadTable::add_expr(RefID r_id, int left, int right, InsertID mr_id, double expr)
99 {
100 #if ENABLE_THREADS
101 boost::mutex::scoped_lock lock(_lock);
102 #endif
103 MultiRead* mr = get_read(mr_id);
104 if (mr)
105 mr->add_expr(r_id, left, right, expr);
106 }
107
get_mass(const MateHit & hit)108 double MultiReadTable::get_mass(const MateHit& hit)
109 {
110 #if ENABLE_THREADS
111 boost::mutex::scoped_lock lock(_lock);
112 #endif
113 MultiRead* mr = get_read(hit.insert_id());
114 if(!mr)
115 return 1.0;
116 return mr->get_mass(hit.ref_id(), hit.left(), hit.right(), _valid_mass);
117 }
118
num_multireads()119 size_t MultiReadTable::num_multireads()
120 {
121 return (int)_read_map.size();
122 }
123
num_multihits()124 size_t MultiReadTable::num_multihits()
125 {
126 size_t count = 0;
127 for (MultiReadMap::iterator it=_read_map.begin() ; it != _read_map.end(); it++ )
128 {
129 count += it->second.num_hits();
130 }
131 return count;
132 }
133