libsequence  1.9.5
CoalescentSimTypes.cc
1 /*
2 
3  Copyright (C) 2003-2009 Kevin Thornton, krthornt[]@[]uci.edu
4 
5  Remove the brackets to email me.
6 
7  This file is part of libsequence.
8 
9  libsequence is free software: you can redistribute it and/or modify
10  it under the terms of the GNU General Public License as published by
11  the Free Software Foundation, either version 3 of the License, or
12  (at your option) any later version.
13 
14  libsequence is distributed in the hope that it will be useful,
15  but WITHOUT ANY WARRANTY; without even the implied warranty of
16  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17  GNU General Public License for more details.
18 
19  You should have received a copy of the GNU General Public License
20  long with libsequence. If not, see <http://www.gnu.org/licenses/>.
21 
22 */
23 
25 #include <cassert>
26 #include <iostream>
27 #include <cstdlib>
28 
29 namespace Sequence
30 {
31  namespace coalsim {
32  segment::segment() : beg(0),end(0),desc(0)
37  {
38  }
39 
40  segment::segment(const int & b,
41  const int & e,
42  const int & d = 0) : beg(b),end(e),desc(d)
49  {
50  assert ( end >= beg );
51  }
52 
53  chromosome::chromosome() : segs(NULL),
54  pop(0),
55  nsegs(0)
60  {
61  }
62 
63  chromosome::chromosome( const std::vector<segment> & initial_segs,
64  const int & population ) :
65  //segs((segment *)(malloc(initial_segs.size()*sizeof(segment)))),
66  segs( static_cast<segment*>(malloc(initial_segs.size()*sizeof(segment)))),
67  pop(population),
68  nsegs(unsigned(initial_segs.size()))
75  {
76  std::copy(initial_segs.begin(),initial_segs.end(),segs);
77  }
78 
83  {
84  if(nsegs>0) free(segs);
85  }
86 
88  //segs((segment *)malloc(ch.nsegs*sizeof(segment))),
89  segs(static_cast<segment *>(malloc(ch.nsegs*sizeof(segment)))),
90  pop(ch.pop),
91  nsegs(ch.nsegs)
95  {
96  std::copy(ch.segs,ch.segs+ch.nsegs,segs);
97  }
98 
103  {
104  if(this == &ch) return *this;
105  if(ch.nsegs>this->nsegs)
106  {
107  //this->segs = (segment *)realloc(this->segs,ch.nsegs*sizeof(segment));
108  this->segs = static_cast<segment *>(realloc(this->segs,ch.nsegs*sizeof(segment)));
109  }
110  std::copy(ch.segs,ch.segs+ch.nsegs,segs);
111  this->nsegs = ch.nsegs;
112  this->pop = ch.pop;
113  return *this;
114  }
115 
116 
127  {
128  std::swap(this->segs,ch.segs);
129  std::swap(this->nsegs,ch.nsegs);
130  std::swap(this->pop,ch.pop);
131  }
132 
134  const unsigned & new_nsegs )
140  {
141  free(segs);
142  segs=newsegs;
143  nsegs=new_nsegs;
144  }
145 
150  {
151  return segs;
152  }
153 
158  {
159  return segs+nsegs;
160  }
161 
166  {
167  return segs;
168  }
169 
174  {
175  return segs+nsegs;
176  }
177 
178  int chromosome::links() const
184  {
185  return( (segs+nsegs-1)->end - segs->beg );
186  }
187 
188  node::node(const double & t,
189  const int & a) : time(t),abv(a)
194  {
195  }
196 
197  marginal::marginal(const int & b, const int & ns,
198  const int & nn,
199  const std::vector<node> & t) :
200  beg(b),nsam(ns),nnodes(nn),tree(t)
207  {
208  }
209 
210  marginal::iterator marginal::begin()
214  {
215  return tree.begin();
216  }
217 
218  marginal::iterator marginal::end()
222  {
223  return tree.end();
224  }
225 
226  marginal::const_iterator marginal::begin() const
230  {
231  return tree.begin();
232  }
233 
234  marginal::const_iterator marginal::end() const
238  {
239  return tree.end();
240  }
241 
242  marginal::reference marginal::operator[](const std::vector<node>::size_type &i)
247  {
248  assert( int(i) <= 2*nsam-2 );
249  return tree[i];
250  }
251 
252  marginal::const_reference marginal::operator[](const std::vector<node>::size_type &i) const
257  {
258  assert( int(i) <= 2*nsam-2 );
259  return tree[i];
260  }
261 
262  bool marginal::operator<(const marginal & rhs) const
266  {
267  return this->beg < rhs.beg;
268  }
269 
270  std::ostream & operator<<(std::ostream & s,const chromosome & c)
276  {
277  s << '(';
279  for( ; seg < c.end() ;
280  ++seg)
281  {
282  s << seg->beg << ':' << seg->end << ',';
283  }
284  s << ' ' << c.nsegs << ' ' << c.links() << ')';
285  return s;
286  }
287 
288  std::ostream & operator<<(std::ostream & s, const marginal & m)
293  {
294  s << "marginal tree begins at: " << m.beg << '\n';
295  int seg=0;
296  for( ;seg < m.nnodes ; ++seg)
297  {
298  s << seg << ' '
299  << (m.begin()+seg)->time << ' '
300  << (m.begin()+seg)->abv << '\n';
301  }
302  s << seg << ' '
303  << (m.begin()+seg)->time << ' '
304  << (m.begin()+seg)->abv;
305  return s;
306  }
307 
312  {
313  public:
314  marginal::const_iterator mi;
315  const int nsam;
316  std::vector<int> left,right;
317  std::vector<node> tree;
318  void init();
319  std::ostream & parens( const int & noden,
320  std::ostream & o) const;
323  newick_stream_marginal_tree_impl(arg::const_iterator m);
324  newick_stream_marginal_tree_impl(arg::iterator m);
325  };
326 
327  newick_stream_marginal_tree_impl::newick_stream_marginal_tree_impl(const marginal & m) :
328  mi(m.begin()),nsam(m.nsam),
329  left(std::vector<int>::size_type(2*nsam-1),-1),right(left),
330  tree(std::vector<node>())
331  {
332  init();
333  }
334 
335  newick_stream_marginal_tree_impl::newick_stream_marginal_tree_impl(const marginal * m) :
336  mi(m->begin()),nsam(m->nsam),
337  left(std::vector<int>::size_type(2*nsam-1),-1),right(left),
338  tree(std::vector<node>())
339  {
340  init();
341  }
342 
343  newick_stream_marginal_tree_impl::newick_stream_marginal_tree_impl(arg::const_iterator m) :
344  mi(m->begin()),nsam(m->nsam),
345  left(std::vector<int>::size_type(2*nsam-1),-1),right(left),
346  tree(std::vector<node>())
347  {
348  init();
349  }
350 
351  newick_stream_marginal_tree_impl::newick_stream_marginal_tree_impl(arg::iterator m) :
352  mi(m->begin()),nsam(m->nsam),
353  left(std::vector<int>::size_type(2*nsam-1),-1),right(left),
354  tree(std::vector<node>())
355  {
356  init();
357  }
358 
359  std::ostream &
361  std::ostream & o) const
365  {
366  if( left[std::vector<int>::size_type(noden)] == -1 )
367  {
368  assert( (mi+((mi+noden)->abv))->time >= 0. );
369  o << noden+1 << ':' <<(mi+((mi+noden)->abv))->time;
370  }
371  else
372  {
373  o << '(';
374  parens( left[std::vector<int>::size_type(noden)], o ) ;
375  o << ',';
376  parens( right[std::vector<int>::size_type(noden)], o ) ;
377  if( (mi+noden)->abv == -1 )
378  {
379  o << ");";
380  }
381  else
382  {
383  double time = (mi + (mi+noden)->abv )->time - (mi+noden)->time ;
384  assert(time >= 0.);
385  o << "):"<<time;
386  }
387  }
388  return o;
389  }
390 
395  {
396  for(int i=0;i<2*nsam-2;++i)
397  {
398  if( left[ std::vector<int>::size_type( (mi+i)->abv ) ] == -1 )
399  {
400  left[std::vector<int>::size_type( (mi+i)->abv ) ] = i;
401  }
402  else
403  {
404  right[ std::vector<int>::size_type( (mi+i)->abv ) ] = i;
405  }
406  }
407 
408  }
409 
410  newick_stream_marginal_tree::newick_stream_marginal_tree( const marginal & m ) :
411  impl( new newick_stream_marginal_tree_impl(m) )
412  {
413  }
414 
415  newick_stream_marginal_tree::newick_stream_marginal_tree( const marginal * m ) :
416  impl( new newick_stream_marginal_tree_impl(m) )
417  {
418  }
419 
420  newick_stream_marginal_tree::newick_stream_marginal_tree( arg::const_iterator m ) :
421  impl( new newick_stream_marginal_tree_impl(m) )
422  {
423  }
424 
425  newick_stream_marginal_tree::newick_stream_marginal_tree( arg::iterator m ) :
426  impl( new newick_stream_marginal_tree_impl(m) )
427  {
428  }
429 
430  newick_stream_marginal_tree::~newick_stream_marginal_tree( )
431  {
432  }
433 
434  std::vector<node> newick_stream_marginal_tree::get_tree() const
439  {
440  return impl->tree;
441  }
442 
443  std::ostream & newick_stream_marginal_tree::print( std::ostream & o ) const
447  {
448  impl->parens( 2*(impl->nsam)-2,o);
449  return o;
450  }
451 
452  std::istream & newick_stream_marginal_tree::read( std::istream & i )
457  {
458  return i;
459  }
460 
461  std::ostream & operator<<(std::ostream & o, const newick_stream_marginal_tree & n)
465  {
466  return n.print(o);
467  }
468 
469  std::istream & operator>>(std::istream & i, newick_stream_marginal_tree & n)
473  {
474  return n.read(i);
475  }
476  }
477 }
std::ostream & parens(const int &noden, std::ostream &o) const
chromosome()
constructor sets segs to NULL, pop to 0, and nsegs to 0
void assign_allocated_segs(segment *newsegs, const unsigned &new_nsegs)
reference operator[](const std::vector< node >::size_type &i)
Class that provides a typecast-on-output of a marginal tree to a newick tree Example use: ...
Definition: SimTypes.hpp:221
std::ostream & operator<<(std::ostream &s, const chromosome &c)
output operator for chromosome types in coalescent simulation Outputs the segments contained by the c...
The namespace in which this library resides.
declaration of types for coalescent simulation
segment()
default constructor initializes all members to 0
std::istream & operator>>(std::istream &i, newick_stream_marginal_tree &n)
marginal(const int &b, const int &ns, const int &nn, const std::vector< node > &tree)
std::vector< node > tree
Definition: SimTypes.hpp:188
A chromosome is a container of segments.
Definition: SimTypes.hpp:92
chromosome & operator=(const chromosome &ch)
assignment operator
bool operator<(const marginal &rhs) const
node(const double &t=0., const int &a=-1)
A portion of a recombining chromosome.
Definition: SimTypes.hpp:83
std::ostream & print(std::ostream &o) const
The genealogy of a portion of a chromosome on which no recombination has occurred.
Definition: SimTypes.hpp:166