A toolkit for working with phylogenetic data.
v0.18.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
circular_layout.cpp
Go to the documentation of this file.
1 /*
2  Genesis - A toolkit for working with phylogenetic data.
3  Copyright (C) 2014-2017 Lucas Czech
4 
5  This program is free software: you can redistribute it and/or modify
6  it under the terms of the GNU General Public License as published by
7  the Free Software Foundation, either version 3 of the License, or
8  (at your option) any later version.
9 
10  This program is distributed in the hope that it will be useful,
11  but WITHOUT ANY WARRANTY; without even the implied warranty of
12  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13  GNU General Public License for more details.
14 
15  You should have received a copy of the GNU General Public License
16  along with this program. If not, see <http://www.gnu.org/licenses/>.
17 
18  Contact:
19  Lucas Czech <lucas.czech@h-its.org>
20  Exelixis Lab, Heidelberg Institute for Theoretical Studies
21  Schloss-Wolfsbrunnenweg 35, D-69118 Heidelberg, Germany
22 */
23 
32 
41 
44 
45 #include <algorithm>
46 #include <assert.h>
47 #include <cmath>
48 #include <stdexcept>
49 
50 namespace genesis {
51 namespace tree {
52 
53 // =================================================================================================
54 // Constructors
55 // =================================================================================================
56 
57 CircularLayout::CircularLayout( Tree const& orig_tree, Type const drawing_type )
58 {
59  // Need to set drawing type first, so that init_tree_() can use it.
60  type( drawing_type );
61  tree( orig_tree );
62 }
63 
64 // =================================================================================================
65 // Options
66 // =================================================================================================
67 
69 {
70  scaler_r_ = value;
71  return *this;
72 }
73 
75 {
76  return scaler_r_;
77 }
78 
79 // =================================================================================================
80 // Virtual Functions
81 // =================================================================================================
82 
83 void CircularLayout::init_tree_( Tree const& orig_tree )
84 {
85  // Init with proper data types.
86  init_nodes_( orig_tree );
87  init_edges_( orig_tree );
88 
89  // Set radiuses of nodes.
90  if( type() == Type::kCladogram ) {
91  set_node_r_cladogram_();
92  } else if( type() == Type::kPhylogram ) {
93  set_node_r_phylogram_();
94  } else {
95  assert( false );
96  }
97 
98  // Set angles of nodes.
99  set_node_a_();
100 }
101 
102 std::string svg_arc(
103  double center_x, double center_y, double radius, double start_angle, double end_angle
104 ) {
105  std::string large_arc;
106  if( start_angle > end_angle ) {
107  large_arc = ( end_angle - start_angle <= utils::PI ? "1" : "0" );
108  } else {
109  large_arc = ( end_angle - start_angle <= utils::PI ? "0" : "1" );
110  }
111 
112  double start_x = center_x + ( radius * cos( end_angle ));
113  double start_y = center_y + ( radius * sin( end_angle ));
114  double end_x = center_x + ( radius * cos( start_angle ));
115  double end_y = center_y + ( radius * sin( start_angle ));
116 
117  std::ostringstream os;
118  os << "M " << start_x << " " << start_y << " ";
119  os << "A " << radius << " " << radius << " " << 0 << " " << large_arc << " " << 0 << " ";
120  os << end_x << " " << end_y;
121  return os.str();
122 }
123 
124 utils::SvgDocument CircularLayout::to_svg_document_() const
125 {
126  using namespace utils;
127  SvgDocument doc;
128  SvgGroup tree_lines;
129  SvgGroup taxa_names;
130  SvgGroup node_shapes;
131 
132  for( auto const& node_it : tree().nodes() ) {
133  auto const& node = *node_it;
134 
135  auto const& node_data = node.data<CircularNodeData>();
136  auto const& parent_data = tree().node_at( node_data.parent_index ).data<CircularNodeData>();
137 
138  auto const node_x = node_data.r * cos( node_data.a );
139  auto const node_y = node_data.r * sin( node_data.a );
140 
141  // Get the edge between the node and its parent.
142  auto edge_data_ptr = edge_between( node, tree().node_at( node_data.parent_index ) );
143 
144  // If there is an edge (i.e., we are not at the root), draw lines between the nodes.
145  if( edge_data_ptr ) {
146  auto stroke = edge_data_ptr->data<CircularEdgeData>().stroke;
147  stroke.line_cap = utils::SvgStroke::LineCap::kRound;
148 
149  auto start_a = parent_data.a;
150  auto end_a = node_data.a;
151  if( parent_data.a > node_data.a ) {
152  std::swap( start_a, end_a );
153  }
154 
155  tree_lines << SvgPath(
156  { svg_arc( 0, 0, parent_data.r, start_a, end_a ) },
157  stroke,
158  SvgFill( SvgFill::Type::kNone )
159  );
160 
161  tree_lines << SvgLine(
162  parent_data.r * cos( node_data.a ), parent_data.r * sin( node_data.a ),
163  node_x, node_y,
164  stroke
165  );
166 
167  } else {
168 
169  // If there is no edge, it must be the root.
170  assert( node.is_root() );
171  }
172 
173  // If the node has a name, print it.
174  if( node_data.name != "" ) {
175  // auto label = SvgText( node_data.name, SvgPoint( node_data.x + 5, node_data.y ) );
176  // label.dy = "0.4em";
177 
178  auto label = text_template();
179  label.text = node_data.name;
180  label.transform.append( SvgTransform::Translate(
181  ( node_data.r + 10 ) * cos( node_data.a ),
182  ( node_data.r + 10 ) * sin( node_data.a )
183  ));
184  label.transform.append( SvgTransform::Rotate(
185  360 * node_data.a / ( 2.0 * utils::PI )
186  ));
187  taxa_names << std::move( label );
188  }
189 
190  // If there is a node shape, draw it.
191  if( ! node_data.shape.empty() ) {
192  auto ns = node_data.shape;
193  ns.transform.append( SvgTransform::Translate( node_x, node_y ));
194  node_shapes << std::move( ns );
195  }
196  }
197 
198  // We are sure that we won't use the groups again, so let's move them!
199  doc << std::move( tree_lines );
200  if( ! taxa_names.empty() ) {
201  doc << std::move( taxa_names );
202  }
203  if( ! node_shapes.empty() ) {
204  doc << std::move( node_shapes );
205  }
206  return doc;
207 }
208 
209 // =================================================================================================
210 // Internal Functions
211 // =================================================================================================
212 
213 void CircularLayout::init_nodes_( Tree const& orig_tree )
214 {
215  // Init nodes.
216  for( size_t i = 0; i < tree().node_count(); ++i ) {
217  // Safety: correct indices.
218  assert( tree().node_at(i).index() == i && orig_tree.node_at(i).index() == i );
219 
220  // Set the tree node data.
222  auto& node_data = tree().node_at(i).data<CircularNodeData>();
223 
224  // If the original tree has node names, use them.
225  auto orig_node_data_ptr = orig_tree.node_at(i).data_cast<DefaultNodeData>();
226  if( orig_node_data_ptr ) {
227  node_data.name = orig_node_data_ptr->name;
228  }
229  }
230 }
231 
232 void CircularLayout::init_edges_( Tree const& orig_tree )
233 {
234  // Init edges.
235  for( size_t i = 0; i < tree().edge_count(); ++i ) {
236  // Safety: correct indices.
237  assert( tree().edge_at(i).index() == i && orig_tree.edge_at(i).index() == i );
238 
239  // Set the tree edge data.
241  auto& edge_data = tree().edge_at(i).data<CircularEdgeData>();
242 
243  // If the original tree has edge branch lengths, use them.
244  auto orig_edge_data_ptr = orig_tree.edge_at(i).data_cast<DefaultEdgeData>();
245  if( orig_edge_data_ptr ) {
246  edge_data.branch_length = orig_edge_data_ptr->branch_length;
247  }
248  }
249 }
250 
251 void CircularLayout::set_node_a_()
252 {
253  auto const num_leaves = static_cast<double>( leaf_node_count( tree() ));
254 
255  // Set node parents and angles of leaves.
256  size_t leaf_count = 0;
257  size_t parent = 0;
258  for( auto it : eulertour( tree() )) {
259  auto& node_data = tree().node_at( it.node().index() ).data<CircularNodeData>();
260 
261  if( node_data.parent_index == -1 ) {
262  node_data.parent_index = parent;
263  }
264  if( it.node().is_leaf() ) {
265  node_data.a = 2.0 * utils::PI * static_cast<double>(leaf_count) / num_leaves;
266  ++leaf_count;
267  }
268 
269  parent = it.node().index();
270  }
271 
272  // Set remaining angles of inner nodes to mid-points of their children.
273  for( auto it : postorder( tree() )) {
274  auto& node_data = tree().node_at( it.node().index() ).data<CircularNodeData>();
275  auto& parent_data = tree().node_at( node_data.parent_index ).data<CircularNodeData>();
276 
277  if( node_data.a < 0.0 ) {
278  auto min_max_diff = node_data.children_max_a - node_data.children_min_a;
279  node_data.a = node_data.children_min_a + min_max_diff / 2.0;
280  }
281 
282  if( parent_data.children_min_a < 0.0 || parent_data.children_min_a > node_data.a ) {
283  parent_data.children_min_a = node_data.a;
284  }
285  if( parent_data.children_max_a < 0.0 || parent_data.children_max_a < node_data.a ) {
286  parent_data.children_max_a = node_data.a;
287  }
288  }
289 }
290 
291 void CircularLayout::set_node_r_phylogram_()
292 {
293  auto node_dists = node_branch_length_distance_vector( tree() );
294 
295  for( size_t i = 0; i < node_dists.size(); ++i ) {
296  tree().node_at(i).data<CircularNodeData>().r = node_dists[i] * scaler_r_ * 20.0;
297  }
298 }
299 
300 void CircularLayout::set_node_r_cladogram_()
301 {
302  // Set root r to 0.
303  tree().root_node().data<CircularNodeData>().r = 0.0;
304 
305  // Get the heights of all subtrees starting from the root.
306  auto heights = subtree_max_path_heights( tree() );
307 
308  // Get the height of the tree, i.e. longest path from root to any leaf.
309  auto root_height = heights[ tree().root_node().index() ];
310 
311  for( auto it : preorder( tree() )) {
312  // The subtree height calculation does not work for the root, so skip it.
313  // Also, we already set the leaf.
314  if( it.is_first_iteration() ) {
315  continue;
316  }
317 
318  // Get the height of the subtree starting at the current node.
319  auto height = heights[ it.node().index() ];
320  assert( height <= root_height );
321 
322  // Set the radius.
323  auto r = ( root_height - height ) * scaler_r_;
324  tree().node_at( it.node().index() ).data<CircularNodeData>().r = r;
325  }
326 }
327 
328 } // namespace tree
329 } // namespace genesis
static std::unique_ptr< CircularEdgeData > create()
size_t index() const
Return the index of this Edge.
Definition: edge.cpp:45
double height(Tree const &tree)
Get the height of the tree, i.e., the longest distance from the root to a leaf, measured using the br...
void swap(SequenceSet &lhs, SequenceSet &rhs)
TreeNode & node_at(size_t index)
Return the TreeNode at a certain index.
Definition: tree/tree.cpp:304
Tree operator functions.
static std::unique_ptr< CircularNodeData > create()
size_t leaf_node_count(Tree const &tree)
Count the number of leaf Nodes of a Tree.
TreeEdge & reset_data(std::unique_ptr< BaseEdgeData > data)
Reset the data pointer of this TreeEdge.
Definition: edge.cpp:200
TreeEdge * edge_between(TreeNode &lhs, TreeNode &rhs)
Return the TreeEdge between two TreeNode&s, if they are neighbours, or nullptr otherwise.
size_t index() const
Return the index of this Node.
Definition: node.cpp:48
utils::Range< IteratorPreorder< TreeLink const, TreeNode const, TreeEdge const > > preorder(ElementType const &element)
TreeNode & root_node()
Return the TreeNode at the current root of the Tree.
Definition: tree/tree.cpp:258
size_t node_count() const
Return the number of TreeNodes of the Tree.
Definition: tree/tree.cpp:350
constexpr double PI
Make the world go round.
Definition: common.hpp:49
Tree const & tree() const
Definition: layout_base.cpp:51
utils::Range< IteratorEulertour< TreeLink const, TreeNode const, TreeEdge const > > eulertour(ElementType const &element)
Definition: eulertour.hpp:190
Class for representing phylogenetic trees.
Definition: tree/tree.hpp:95
std::string svg_arc(double center_x, double center_y, double radius, double start_angle, double end_angle)
std::vector< size_t > subtree_max_path_heights(Tree const &tree, TreeNode const &node)
Header of Default Tree distance methods.
Provides easy and fast logging functionality.
size_t edge_count() const
Return the number of TreeEdges of the Tree.
Definition: tree/tree.cpp:358
TreeNode & reset_data(std::unique_ptr< BaseNodeData > data)
Reset the data pointer of this TreeNode.
Definition: node.cpp:163
TreeEdge & edge_at(size_t index)
Return the TreeEdge at a certain index.
Definition: tree/tree.cpp:324
NodeDataType & data()
Definition: node.hpp:108
std::vector< double > node_branch_length_distance_vector(Tree const &tree, TreeNode const *node)
Return a vector containing the distance of all nodes with respect to the given start node...
Header of Tree distance methods.
utils::Range< IteratorPostorder< TreeLink const, TreeNode const, TreeEdge const > > postorder(ElementType const &element)
EdgeDataType & data()
Definition: edge.hpp:118
utils::SvgText & text_template()
Definition: layout_base.cpp:85