RNAMake
basepair_state.h
1 //
2 // BasepairState.h
3 // REDESIGNC
4 //
5 // Created by Joseph Yesselman on 9/28/14.
6 // Copyright (c) 2014 Joseph Yesselman. All rights reserved.
7 //
8 
9 #ifndef __REDESIGNC__BasepairState__
10 #define __REDESIGNC__BasepairState__
11 
12 #include <iostream>
13 #include <vector>
14 
15 //RNAMake Headers
16 #include "math/xyz_vector.h"
17 #include "math/xyz_matrix.h"
18 #include "math/transform.h"
19 #include "math/numerical.h"
20 #include "structure/basepair_state.fwd.h"
21 
22 class BasepairState {
23 
24 public:
25 
26  inline
27  BasepairState():
28  d_( Point(0.0) ),
29  r_( Matrix(0.0) ),
30  r_T_( Matrix(0.0) ),
31  sugars_ ( Points(2) ),
32  diff_ ( Vector(0.0) ),
33  diff_sugars_( Vectors(2) )
34  {}
35 
36  inline
37  BasepairState(
38  Point const & d,
39  Matrix const & r,
40  Points const & sugars):
41  d_( d ),
42  r_( r ),
43  r_T_( Matrix(0) ),
44  sugars_( sugars ),
45  diff_ ( Vector(0.0) ),
46  diff_sugars_( Vectors(2) )
47  { transpose(r_, r_T_); }
48 
49  inline
50  BasepairState(BasepairState const & b):
51  d_( b.d_ ),
52  r_( b.r_ ),
53  r_T_( Matrix(0) ),
54  sugars_ (b.sugars_),
55  diff_ ( Vector(0.0) ),
56  diff_sugars_( Vectors(2) )
57  { transpose(r_, r_T_); }
58 
59 
60  ~BasepairState()
61  {}
62 
63  BasepairState
64  copy() const {
65  //Point d = d_;
66  //Matrix r = r_;
67  //Points sugars = sugars_;
68  return BasepairState(d_, r_, sugars_);
69  }
70 
71 
72 public:
73  inline
74  void
75  calculate_r_T() { transpose(r_, r_T_); }
76 
77  inline
78  void
79  get_transforming_r_and_t (
80  BasepairState const & o_state, //state with desired rotation and translation
81  BasepairState & r_state) {
82 
83  //TODO figure out why I have to do this each time???
84  calculate_r_T();
85 
86  //calculate transforming rotation matrix and store it in r_state (result state)
87  dot(r_T_,o_state.r_,r_state.r_);
88  r_state.calculate_r_T();
89 
90  //rotate sugars to new position and store them in r_state
91  dot_vectors(r_state.r_T_,o_state.sugars_,r_state.sugars_);
92 
93  diff_ = -o_state.d() + d_;
94  int i;
95 
96  for(i = 0; i < 2; i+=1) {
97  r_state.sugars_[i] += diff_;
98  }
99 
100  for(i = 0; i < 2; i+=1) {
101  diff_sugars_[i] = sugars_[i] - r_state.sugars_[i];
102  }
103 
104  diff_ = (diff_sugars_[0] + diff_sugars_[1]) / 2.0f;
105  r_state.d_ = -o_state.d() + diff_ + d_;
106 
107  }
108 
109  inline
110  void
111  get_transformed_state(
112  BasepairState const & o_state,
113  BasepairState & r_state) {
114 
115  dot (r_, o_state.r_T_, r_state.r_);
116  dot_vector (o_state.r_T_, d_, r_state.d_);
117  dot_vectors(o_state.r_T_, sugars_, r_state.sugars_);
118 
119  int i;
120  for(i = 0; i < 2; i++) {
121  r_state.sugars_[i] += o_state.d_;
122  }
123  r_state.d_ += o_state.d_;
124 
125  }
126 
127  inline
128  void
129  flip() {
130  r_ = transform_1(r_);
131  calculate_r_T();
132  }
133 
134  inline
135  float
136  diff(
137  BasepairStateOP const & state) {
138  float diff = d_.distance(state->d());
139  diff += _rot_diff(state)*2;
140  return diff;
141  }
142 
143  inline
144  float
145  _rot_diff(
146  BasepairStateOP const & state) {
147  float r_diff = r_.difference(state->r());
148  state->flip();
149  float r_diff_2 = r_.difference(state->r());
150  state->flip();
151 
152  if( r_diff > r_diff_2) { r_diff = r_diff_2;}
153 
154  return r_diff;
155  }
156 
157 public:
158 
159  String const
160  to_str() const{
161  String s = vector_to_str(d_) + ";" + matrix_to_str(r_) + ";" + vector_to_str(sugars_[0]) + vector_to_str(sugars_[1]);
162  return s;
163  }
164 
165 
166 public: //getters
167 
168  inline
169  const
170  Point &
171  d() {
172  return d_;
173  }
174 
175  inline
176  const
177  Point &
178  d() const {
179  return d_;
180  }
181 
182  inline
183  const
184  Matrix &
185  r() {
186  return r_;
187  }
188 
189  inline
190  const
191  Matrix &
192  r() const {
193  return r_;
194  }
195 
196  inline
197  const
198  Matrix &
199  r_T() {
200  return r_T_;
201  }
202 
203  inline
204  const
205  Matrix &
206  r_T() const {
207  return r_T_;
208  }
209 
210  inline
211  const
212  Points &
213  sugars() {
214  return sugars_;
215  }
216 
217  inline
218  const
219  Points &
220  sugars() const {
221  return sugars_;
222  }
223 
224 
225 public: //setters
226 
227  inline
228  void
229  d(Point const & newd) { d_ = newd; }
230 
231  inline
232  void
233  r( Matrix const & newr) { r_ = newr; }
234 
235  inline
236  void
237  sugars( Points const & newsug) { sugars_ = newsug; }
238 
239  inline
240  void
241  set (BasepairState const & nstate_) {
242  d_ = nstate_.d_;
243  r_ = nstate_.r_;
244  sugars_ = nstate_.sugars_;
245  calculate_r_T();
246  }
247 
248 
249 private:
250  Point d_;
251  Matrix r_;
252  Matrix r_T_; //holds the transpose of r_ for speed
253  Points sugars_;
254 
255  /*
256  Inclusion of these two tempory variables increases the speed of
257  get_transforming_r_and_t() by 240 percent seems worth it for just
258  a few more bytes of memory since this a bulk of the calculations
259  done in the algorithm -JDY 2014.10.4
260  */
261 
262  Vector diff_; //stores partial products to speed up algorithm
263  Vectors diff_sugars_; //stores partial products to speed up algorithm
264 
265 
266 };
267 
268 BasepairState
269 str_to_basepairstate(
270  String const &);
271 
272 BasepairState
273 get_ref_bp_state();
274 
275 float
276 get_bpstate_rotation_diff(
277  BasepairState const &,
278  BasepairState const &);
279 
280 
281 std::ostream&
282 operator <<(
283  std::ostream&,
284  const BasepairState&);
285 
286 
287 inline
288 const
289 float
290 frame_distance(
291  BasepairStateOP const & current,
292  BasepairStateOP const & end,
293  BasepairStateOP const & endflip) {
294 
295  float score = current->d().distance(end->d());
296 
297  float r_diff = end->r().difference(current->r());
298  float r_diff_flip = endflip->r().difference(current->r());
299 
300  if(r_diff > r_diff_flip) {
301  r_diff = r_diff_flip;
302  }
303 
304  score += r_diff;
305 
306  return score;
307 
308 }
309 
310 inline
311 const
312 float
313 new_score_function(
314  BasepairStateOP const & current,
315  BasepairStateOP const & end,
316  BasepairStateOP const & endflip) {
317 
318  float d_diff = current->d().distance(end->d());
319 
320  if(d_diff > 25) { return d_diff; }
321 
322  float r_diff = current->r().difference(end->r());
323  float r_diff_flip = current->r().difference(endflip->r());
324 
325  if(r_diff > r_diff_flip) {
326  r_diff = r_diff_flip;
327  }
328 
329  float scale = (log(150/d_diff) - 1);
330  if (scale > 2) { scale = 2; }
331  if (scale < 0) { scale = 0; }
332 
333  return d_diff + scale*r_diff;
334 }
335 
336 inline
337 const
338 float
339 new_score_function_new(
340  BasepairStateOP const & current,
341  BasepairStateOP const & end,
342  BasepairStateOP const & endflip) {
343 
344  float d_diff = (current->sugars()[0].distance(end->sugars()[1]) +
345  current->sugars()[1].distance(end->sugars()[0]))*0.50;
346 
347  if(d_diff > 25) { return d_diff; }
348 
349  float r_diff = current->r().difference(end->r());
350  float r_diff_flip = current->r().difference(endflip->r());
351 
352  if(r_diff > r_diff_flip) {
353  r_diff = r_diff_flip;
354  }
355 
356  float scale = (log(150/d_diff) - 1);
357  if (scale > 2) { scale = 2; }
358  if (scale < 0) { scale = 0; }
359 
360  return d_diff + scale*r_diff;
361 }
362 
363 
364 int
365 are_BasepairStates_equal(
366  BasepairState const,
367  BasepairState const);
368 
369 
370 
371 #endif /* defined(__REDESIGNC__BasepairState__) */
BasepairState
Definition: basepair_state.h:22
Generated by   doxygen 1.8.10