source: branches/stable/PARSIMONY/AP_tree_nlen.cxx

Last change on this file was 17478, checked in by westram, 7 years ago
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 45.5 KB
Line 
1// =============================================================== //
2//                                                                 //
3//   File      : AP_tree_nlen.cxx                                  //
4//   Purpose   :                                                   //
5//                                                                 //
6//   Coded by Ralf Westram (coder@reallysoft.de) in Summer 1995    //
7//   Institute of Microbiology (Technical University Munich)       //
8//   http://www.arb-home.de/                                       //
9//                                                                 //
10// =============================================================== //
11
12#include "ap_tree_nlen.hxx"
13#include "ap_main.hxx"
14
15#include <AP_seq_dna.hxx>
16#include <aw_root.hxx>
17
18using namespace std;
19
20// ---------------------------------
21//      Section base operations:
22// ---------------------------------
23
24AP_UPDATE_FLAGS AP_pars_root::check_update() {
25    // disables load if tree changes in DB
26    // (ignore changes performed in arb_ntree while tree is loaded in arb_pars)
27    AP_UPDATE_FLAGS res = AP_tree_root::check_update();
28    return res == AP_UPDATE_RELOADED ? AP_UPDATE_OK : res;
29}
30
31ostream& operator<<(ostream& out, const AP_tree_nlen *node) {
32    out << ' ';
33
34    if (!node) {
35        out << "NULp";
36    }
37    if (node->is_leaf()) {
38        out << ((void *)node) << '(' << node->name << ')';
39    }
40    else {
41        static int notTooDeep;
42
43        if (notTooDeep) {
44            out << ((void *)node);
45            if (!node->father) out << " (ROOT)";
46        }
47        else {
48            notTooDeep = 1;
49
50            out << "NODE(" << ((void *)node);
51
52            if (!node->father) {
53                out << " (ROOT)";
54            }
55            else {
56                out << ", father=" << node->father;
57            }
58
59            out << ", leftson=" << node->leftson
60                << ", rightson=" << node->rightson
61                << ", edge[0]=" << node->edge[0]
62                << ", edge[1]=" << node->edge[1]
63                << ", edge[2]=" << node->edge[2]
64                << ")";
65
66            notTooDeep = 0;
67        }
68    }
69
70    return out << ' ';
71}
72
73int AP_tree_nlen::unusedEdgeIndex() const {
74    for (int e=0; e<3; e++) if (!edge[e]) return e;
75    return -1;
76}
77
78AP_tree_edge* AP_tree_nlen::edgeTo(const AP_tree_nlen *neighbour) const {
79    for (int e=0; e<3; e++) {
80        if (edge[e] && edge[e]->node[1-index[e]]==neighbour) {
81            return edge[e];
82        }
83    }
84    return NULp;
85}
86
87AP_tree_edge* AP_tree_nlen::nextEdge(const AP_tree_edge *afterThatEdge) const {
88    /*! @return one edge of 'this'
89     *
90     * @param afterThatEdge
91     * - if == NULp -> returns the "first" edge (edge[0])
92     * - otherwise -> returns the next edge following 'afterThatEdge' in the array edge[]
93     */
94    return edge[afterThatEdge ? ((indexOf(afterThatEdge)+1) % 3) : 0];
95}
96
97void AP_tree_nlen::unlinkAllEdges(AP_tree_edge **edgePtr1, AP_tree_edge **edgePtr2, AP_tree_edge **edgePtr3) {
98    ap_assert(edge[0]);
99    ap_assert(edge[1]);
100    ap_assert(edge[2]);
101
102    *edgePtr1 = edge[0]->unlink();
103    *edgePtr2 = edge[1]->unlink();
104    *edgePtr3 = edge[2]->unlink();
105}
106
107void AP_tree_nlen::linkAllEdges(AP_tree_edge *edge1, AP_tree_edge *edge2, AP_tree_edge *edge3) {
108    ap_assert(!edge[0]);
109    ap_assert(!edge[1]);
110    ap_assert(!edge[2]);
111
112    edge1->relink(this, get_father()->get_father() ? get_father() : get_brother());
113    edge2->relink(this, get_leftson());
114    edge3->relink(this, get_rightson());
115}
116
117// -----------------------------
118//      Check tree structure
119
120#if defined(PROVIDE_TREE_STRUCTURE_TESTS)
121
122#if defined(DEBUG)
123#define DUMP_INVALID_SUBTREES
124#endif
125
126#if defined(DEVEL_RALF)
127#define CHECK_CORRECT_INVALIDATION // recombines all up-to-date nodes to find missing invalidations (VERY slow)
128#endif
129
130
131#if defined(DUMP_INVALID_SUBTREES)
132inline void dumpSubtree(const char *msg, const AP_tree_nlen *node) {
133    fprintf(stderr, "%s:\n", msg);
134    char *printable = GBT_tree_2_newick(node, NewickFormat(nSIMPLE|nWRAP), true);
135    fputs(printable, stderr);
136    fputc('\n', stderr);
137    free(printable);
138}
139#endif
140
141inline const AP_tree_edge *edge_between(const AP_tree_nlen *node1, const AP_tree_nlen *node2) {
142    AP_tree_edge *edge_12 = node1->edgeTo(node2);
143
144#if defined(ASSERTION_USED)
145    AP_tree_edge *edge_21 = node2->edgeTo(node1);
146    ap_assert(edge_12 == edge_21); // nodes should agree about their edge
147#endif
148
149    return edge_12;
150}
151
152inline const char *no_valid_edge_between(const AP_tree_nlen *node1, const AP_tree_nlen *node2) {
153    AP_tree_edge *edge_12 = node1->edgeTo(node2);
154    AP_tree_edge *edge_21 = node2->edgeTo(node1);
155
156    if (edge_12 == edge_21) {
157        return edge_12 ? NULp : "edge missing";
158    }
159    return "edge inconsistent";
160}
161
162#if defined(DUMP_INVALID_SUBTREES)
163#define PRINT_BAD_EDGE(msg,node) dumpSubtree(msg,node)
164#else // !defined(DUMP_INVALID_SUBTREES)
165#define PRINT_BAD_EDGE(msg,node) fprintf(stderr, "Warning: %s (at node=%p)\n", (msg), (node))
166#endif
167
168#define SHOW_BAD_EDGE(format,str,node) do{              \
169        char *msg = GBS_global_string_copy(format,str); \
170        PRINT_BAD_EDGE(msg, node);                      \
171        free(msg);                                      \
172    }while(0)
173
174Validity AP_tree_nlen::has_valid_edges() const {
175    Validity valid;
176    if (get_father()) {                     // root has no edges
177        if (get_father()->is_root_node()) { // sons of root have one edge between them
178            if (is_leftson()) { // test root-edge only from one son
179                const char *invalid = no_valid_edge_between(this, get_brother());
180                if (invalid) {
181                    SHOW_BAD_EDGE("root-%s. root", invalid, get_father());
182                    valid = Validity(false, "no valid edge between sons of root");
183                }
184            }
185            const char *invalid = no_valid_edge_between(this, get_father());
186            if (!invalid || !strstr(invalid, "missing")) {
187                SHOW_BAD_EDGE("unexpected edge (%s) between root and son", invalid ? invalid : "valid", this);
188                valid = Validity(false, "unexpected edge between son-of-root and root");
189            }
190        }
191        else {
192            const char *invalid = no_valid_edge_between(this, get_father());
193            if (invalid) {
194                SHOW_BAD_EDGE("son-%s. father", invalid, get_father());
195                SHOW_BAD_EDGE("parent-%s. son", invalid, this);
196                valid = Validity(false, "invalid edge between son and father");
197            }
198        }
199    }
200
201    if (!is_leaf()) {
202        if (valid) valid = get_leftson()->has_valid_edges();
203        if (valid) valid = get_rightson()->has_valid_edges();
204    }
205    return valid;
206}
207
208Validity AP_tree_nlen::sequence_state_valid() const {
209    // if some node has a sequence, all son-nodes have to have sequences!
210
211    Validity valid;
212
213    const AP_combinableSeq *sequence = get_seq();
214    if (sequence) {
215        if (sequence->hasSequence()) {
216            if (!is_leaf()) {
217                bool leftson_hasSequence  = get_leftson()->hasSequence();
218                bool rightson_hasSequence = get_rightson()->hasSequence();
219
220#if defined(DUMP_INVALID_SUBTREES)
221                if (!leftson_hasSequence) dumpSubtree("left subtree has no sequence", get_leftson());
222                if (!rightson_hasSequence) dumpSubtree("right subtree has no sequence", get_rightson());
223                if (!(leftson_hasSequence && rightson_hasSequence)) {
224                    dumpSubtree("while father HAS sequence", this);
225                }
226#endif
227
228                valid = Validity(leftson_hasSequence && rightson_hasSequence, "node has sequence and son w/o sequence");
229
230#if defined(CHECK_CORRECT_INVALIDATION)
231                if (valid) {
232                    // check for missing invalidations
233                    // (if recalculating a node (via combine) does not reproduce the current sequence, it should have been invalidated)
234
235                    AP_combinableSeq *recombined             = sequence->dup();
236                    Mutations         mutations_from_combine = recombined->noncounting_combine_seq(get_leftson()->get_seq(), get_rightson()->get_seq());
237
238                    valid = Validity(recombined->combinedEquals(sequence), "recombining changed existing sequence (missing invalidation?)");
239                    if (valid) {
240                        Mutations expected_mutrate = mutations_from_combine + get_leftson()->mutations + get_rightson()->mutations;
241                        valid = Validity(expected_mutrate == mutations, "invalid mutation_rate");
242                    }
243
244                    delete recombined;
245
246#if defined(DUMP_INVALID_SUBTREES)
247                    if (!valid) {
248                        dumpSubtree(valid.why_not(), this);
249                    }
250#endif
251                }
252#endif
253            }
254        }
255#if defined(ASSERTION_USED)
256        else {
257            if (is_leaf()) ap_assert(sequence->is_bound_to_species()); // can do lazy load if needed
258        }
259#endif
260    }
261
262    if (!is_leaf()) {
263        if (valid) valid = get_leftson()->sequence_state_valid();
264        if (valid) valid = get_rightson()->sequence_state_valid();
265    }
266
267    return valid;
268}
269
270Validity AP_tree_nlen::is_valid() const {
271    ap_assert(knownNonNull(this));
272
273    Validity valid   = AP_tree::is_valid();
274    if (valid) valid = has_valid_edges();
275    if (valid) valid = sequence_state_valid();
276
277    return valid;
278}
279
280#endif // PROVIDE_TREE_STRUCTURE_TESTS
281
282// -------------------------
283//      Tree operations:
284//
285// insert
286// remove
287// swap
288// set_root
289// move
290// costs
291
292
293inline void push_all_upnode_sequences(AP_tree_nlen *nodeBelow) {
294    for  (AP_tree_nlen *upnode = nodeBelow->get_father();
295          upnode;
296          upnode = upnode->get_father())
297    {
298        ap_main->push_node(upnode, SEQUENCE);
299    }
300}
301
302inline void sortOldestFirst(AP_tree_edge **e1, AP_tree_edge **e2) {
303    if ((*e1)->Age() > (*e2)->Age()) {
304        swap(*e1, *e2);
305    }
306}
307
308inline void sortOldestFirst(AP_tree_edge **e1, AP_tree_edge **e2, AP_tree_edge **e3) {
309    sortOldestFirst(e1, e2);
310    sortOldestFirst(e2, e3);
311    sortOldestFirst(e1, e2);
312}
313
314void AP_tree_nlen::initial_insert(AP_tree_nlen *newBrother, AP_pars_root *troot) {
315    // construct initial tree from 'this' and 'newBrother'
316    // (both have to be leafs)
317
318    ap_assert(newBrother);
319    ap_assert(is_leaf());
320    ap_assert(newBrother->is_leaf());
321
322    AP_tree::initial_insert(newBrother, troot);
323    makeEdge(newBrother, this); // build the root edge
324
325    ASSERT_VALID_TREE(this->get_father());
326}
327
328void AP_tree_nlen::insert(AP_tree_nlen *newBrother) {
329    //  inserts 'this' (a new node) at the father-edge of 'newBrother'
330    ap_assert(newBrother);
331    ap_assert(rootNode()->has_valid_root_remarks());
332
333    ASSERT_VALID_TREE(this);
334    ASSERT_VALID_TREE(newBrother);
335
336    ap_main->push_node(this, STRUCTURE);
337    ap_main->push_node(newBrother, STRUCTURE);
338
339    AP_tree_nlen *brothersFather = newBrother->get_father();
340    if (brothersFather) {
341        ap_main->push_node(brothersFather, BOTH);
342        push_all_upnode_sequences(brothersFather);
343
344        if (brothersFather->get_father()) {
345            AP_tree_edge *oldEdge = newBrother->edgeTo(newBrother->get_father())->unlink();
346            AP_tree::insert(newBrother);
347            oldEdge->relink(get_father(), get_father()->get_father());
348        }
349        else { // insert to son of root
350            AP_tree_nlen *brothersOldBrother = newBrother->get_brother();
351            ap_main->push_node(brothersOldBrother, STRUCTURE);
352
353            AP_tree_edge *oldEdge = newBrother->edgeTo(brothersOldBrother)->unlink();
354            AP_tree::insert(newBrother);
355            oldEdge->relink(get_father(), get_father()->get_brother());
356        }
357
358        makeEdge(this, get_father());
359        makeEdge(get_father(), newBrother);
360
361        ASSERT_VALID_TREE(get_father()->get_father());
362    }
363    else { // insert at root
364        ap_assert(!newBrother->is_leaf()); // either swap 'this' and 'newBrother' or use initial_insert() to construct the initial tree
365
366        AP_tree_nlen *lson = newBrother->get_leftson();
367        AP_tree_nlen *rson = newBrother->get_rightson();
368
369        ap_main->push_node(lson, STRUCTURE);
370        ap_main->push_node(rson, STRUCTURE);
371
372        AP_tree_edge *oldEdge = lson->edgeTo(rson)->unlink();
373
374        AP_tree::insert(newBrother);
375
376        oldEdge->relink(this, newBrother);
377        makeEdge(newBrother, rson);
378        makeEdge(newBrother, lson);
379
380        ASSERT_VALID_TREE(get_father());
381    }
382    ap_assert(rootNode()->has_valid_root_remarks());
383}
384
385AP_tree_nlen *AP_tree_nlen::REMOVE() {
386    // Removes 'this' and its father from the tree:
387    //
388    //       grandpa                grandpa
389    //           /                    /
390    //          /                    /
391    //    father        =>        brother
392    //       /     \                                            .
393    //      /       \                                           .
394    //   this       brother
395    //
396    // One of the edges is relinked between brother and grandpa.
397    // 'father' is destroyed, 'this' is returned.
398
399    AP_tree_nlen *oldBrother = get_brother();
400
401    ASSERT_VALID_TREE(this);
402
403    ap_assert(father); // can't remove complete tree,
404
405    ap_main->push_node(this, STRUCTURE);
406    ap_main->push_node(oldBrother, STRUCTURE);
407    push_all_upnode_sequences(get_father());
408
409    AP_tree_edge *oldEdge;
410    AP_tree_nlen *grandPa = get_father()->get_father();
411    if (grandPa) {
412        ASSERT_VALID_TREE(grandPa);
413
414        ap_main->push_node(get_father(), BOTH);
415        ap_main->push_node(grandPa, STRUCTURE);
416
417        destroyEdge(edgeTo(get_father())->unlink());
418        destroyEdge(get_father()->edgeTo(oldBrother)->unlink());
419
420        if (grandPa->father) {
421            oldEdge = get_father()->edgeTo(grandPa)->unlink();
422            AP_tree::REMOVE();
423            oldEdge->relink(oldBrother, grandPa);
424        }
425        else { // remove grandson of root
426            AP_tree_nlen *uncle = get_father()->get_brother();
427            ap_main->push_node(uncle, STRUCTURE);
428
429            oldEdge = get_father()->edgeTo(uncle)->unlink();
430            AP_tree::REMOVE();
431            oldEdge->relink(oldBrother, uncle);
432        }
433        ASSERT_VALID_TREE(grandPa);
434    }
435    else {                                          // remove son of root
436        AP_tree_nlen *oldRoot = get_father();
437        ASSERT_VALID_TREE(oldRoot);
438
439        if (oldBrother->is_leaf()) {
440            //           root
441            //            oo
442            //           o  o
443            //          o    o
444            // oldBrother --- this         ----->   NULp
445            //
446            ap_main->push_node(oldRoot, ROOT);
447
448            destroyEdge(edgeTo(oldBrother)->unlink());
449
450#if defined(ASSERTION_USED)
451            AP_pars_root *troot = get_tree_root();
452#endif // ASSERTION_USED
453            AP_tree::REMOVE();
454            ap_assert(!troot->get_root_node()); // tree should have been removed
455        }
456        else {
457            //
458            //           root
459            //            oo                                                              .
460            //           o  o                                     root (=oldBrother)
461            //          o    o                                     oo                      .
462            // oldBrother --- this          ----->                o  o                     .
463            //       /\                                          o    o                    .
464            //      /  \                                     lson ----- rson
465            //     /    \                                                                .
466            //    lson  rson
467            //
468            AP_tree_nlen *lson = oldBrother->get_leftson();
469            AP_tree_nlen *rson = oldBrother->get_rightson();
470
471            ap_assert(lson && rson);
472
473            ap_main->push_node(lson, STRUCTURE);
474            ap_main->push_node(rson, STRUCTURE);
475            ap_main->push_node(oldRoot, ROOT);
476
477            destroyEdge(edgeTo(oldBrother)->unlink());
478            destroyEdge(oldBrother->edgeTo(lson)->unlink());
479
480            oldEdge = oldBrother->edgeTo(rson)->unlink();
481            AP_tree::REMOVE();
482            oldEdge->relink(lson, rson);
483
484            ap_assert(lson->get_tree_root()->get_root_node() == oldBrother);
485            ASSERT_VALID_TREE(oldBrother);
486        }
487    }
488
489    father = NULp;
490    set_tree_root(NULp);
491
492    ASSERT_VALID_TREE(this);
493    return this;
494}
495
496void AP_tree_nlen::swap_sons() {
497    ap_assert(!is_leaf()); // cannot swap sons at leafs
498
499    ap_main->push_node(this, STRUCTURE);
500    AP_tree::swap_sons();
501}
502
503void AP_tree_nlen::swap_assymetric(AP_TREE_SIDE mode) {
504    // mode AP_LEFT exchanges leftson with brother
505    // mode AP_RIGHT exchanges rightson with brother
506
507    // @@@ "NNI" works really bad for keeled groups (fixable?; #785)
508
509    ap_assert(!is_leaf());                          // cannot swap leafs
510    ap_assert(father);                              // cannot swap root (has no brother)
511    ap_assert(mode == AP_LEFT || mode == AP_RIGHT); // illegal mode
512
513    AP_tree_nlen *oldBrother = get_brother();
514    AP_tree_nlen *movedSon   = mode == AP_LEFT ? get_leftson() : get_rightson();
515
516    if (!father->father) {
517        // son of root case
518        // take leftson of brother to exchange with
519
520        if (!oldBrother->is_leaf()) { // swap needed ?
521            AP_tree_nlen *nephew = oldBrother->get_leftson();
522
523            ap_main->push_node(this, BOTH);
524            ap_main->push_node(movedSon, STRUCTURE);
525            ap_main->push_node(get_father(), SEQUENCE);
526            ap_main->push_node(nephew, STRUCTURE);
527            ap_main->push_node(oldBrother, BOTH);
528
529            AP_tree_edge *edge1 = edgeTo(movedSon)->unlink();
530            AP_tree_edge *edge2 = oldBrother->edgeTo(nephew)->unlink();
531
532            if (mode == AP_LEFT) {
533                swap(leftson->father, nephew->father);
534                swap(leftson, oldBrother->leftson);
535            }
536            else {
537                swap(rightson->father, nephew->father);
538                swap(rightson, oldBrother->leftson);
539            }
540
541            edge2->relink(this, nephew);
542            edge1->relink(oldBrother, movedSon);
543
544            if (nephew->gr.mark_sum != movedSon->gr.mark_sum) {
545                get_brother()->recalc_marked_from_sons();
546                this->recalc_marked_from_sons_and_forward_upwards();
547            }
548        }
549    }
550    else {
551        ap_main->push_node(this, BOTH);
552        ap_main->push_node(get_father(), BOTH);
553        ap_main->push_node(oldBrother, STRUCTURE);
554        ap_main->push_node(movedSon, STRUCTURE);
555
556        push_all_upnode_sequences(get_father());
557
558        AP_tree_edge *edge1 = edgeTo(movedSon)->unlink();
559        AP_tree_edge *edge2 = get_father()->edgeTo(oldBrother)->unlink();
560
561        if (mode == AP_LEFT) { // swap leftson with brother
562            swap(leftson->father, oldBrother->father);
563            if (father->leftson == this) {
564                swap(leftson, father->rightson);
565            }
566            else {
567                swap(leftson, father->leftson);
568            }
569        }
570        else { // swap rightson with brother
571            swap(rightson->father, oldBrother->father);
572            if (father->leftson == this) {
573                swap(rightson, father->rightson);
574            }
575            else {
576                swap(rightson, father->leftson);
577            }
578        }
579
580        edge2->relink(this, oldBrother);
581        edge1->relink(get_father(), movedSon);
582
583        if (oldBrother->gr.mark_sum != movedSon->gr.mark_sum) {
584            recalc_marked_from_sons_and_forward_upwards(); // father is done implicit
585        }
586    }
587}
588
589void AP_tree_nlen::set_root() {
590    if (at_root()) return; // already root
591
592    // from this to root buffer the nodes
593    ap_main->push_node(this,  STRUCTURE);
594    ap_assert(rootNode()->has_valid_root_remarks());
595
596    AP_tree_nlen *son_of_root = NULp; // in previous topology 'this' was contained inside 'son_of_root'
597    AP_tree_nlen *old_root    = NULp;
598    {
599        AP_tree_nlen *pntr;
600        for (pntr = get_father(); pntr->father; pntr = pntr->get_father()) {
601            ap_main->push_node(pntr, BOTH);
602            son_of_root = pntr;
603        }
604        old_root = pntr;
605    }
606
607    ap_assert(son_of_root); // always true
608
609    {
610        AP_tree_nlen *other_son_of_root = son_of_root->get_brother();
611        ap_main->push_node(other_son_of_root, STRUCTURE);
612    }
613
614    ap_main->push_node(old_root, ROOT);
615    ap_assert(old_root->has_valid_root_remarks());
616    AP_tree::set_root();
617
618    for (AP_tree_nlen *node = son_of_root; node ; node = node->get_father()) {
619        node->recalc_marked_from_sons();
620    }
621}
622
623void AP_tree_nlen::moveNextTo(AP_tree_nlen *newBrother, AP_FLOAT rel_pos) {
624    // Note: see http://bugs.arb-home.de/ticket/627#comment:8 for an experimental
625    // replacement of moveNextTo with REMOVE() + insert()
626
627    ap_assert(father);
628    ap_assert(newBrother);
629    ap_assert(newBrother->father);
630    ap_assert(newBrother->father != father); // already there
631    ap_assert(newBrother != father);         // already there
632
633    ASSERT_VALID_TREE(rootNode());
634
635    // push everything that will be modified onto stack
636    ap_main->push_node(this,  STRUCTURE);
637    ap_main->push_node(get_brother(), STRUCTURE);
638
639    if (father->father) {
640        AP_tree_nlen *grandpa = get_father()->get_father();
641
642        ap_main->push_node(get_father(), BOTH);
643
644        if (grandpa->father) {
645            ap_main->push_node(grandpa, BOTH);
646            push_all_upnode_sequences(grandpa);
647        }
648        else { // 'this' is grandson of root
649            ap_main->push_node(grandpa, ROOT);
650            ap_main->push_node(get_father()->get_brother(), STRUCTURE);
651        }
652    }
653    else { // 'this' is son of root
654        ap_main->push_node(get_father(), ROOT);
655
656        if (!get_brother()->is_leaf()) {
657            ap_main->push_node(get_brother()->get_leftson(), STRUCTURE);
658            ap_main->push_node(get_brother()->get_rightson(), STRUCTURE);
659        }
660    }
661
662    ap_main->push_node(newBrother,  STRUCTURE);
663    if (newBrother->father) {
664        AP_tree_nlen *newBrothersFather = newBrother->get_father();
665        ap_main->push_node(newBrothersFather, BOTH);
666        if (newBrothersFather->is_son_of_root()) {
667            ap_main->push_node(newBrothersFather->get_brother(), STRUCTURE);
668        }
669
670        if (newBrother->is_son_of_root()) { // move to son of root
671            ap_main->push_node(newBrother->get_brother(), STRUCTURE);
672        }
673        push_all_upnode_sequences(newBrothersFather);
674    }
675
676    AP_tree_nlen *thisFather        = get_father();
677    AP_tree_nlen *grandFather       = thisFather->get_father();
678    AP_tree_nlen *oldBrother        = get_brother();
679    AP_tree_nlen *newBrothersFather = newBrother->get_father();
680    AP_tree_edge *e1, *e2, *e3;
681
682    if (thisFather==newBrothersFather->get_father()) { // son -> son of brother
683        if (grandFather) {
684            if (grandFather->get_father()) {
685                // covered by test at PARS_main.cxx@COVER3
686                thisFather->unlinkAllEdges(&e1, &e2, &e3);
687                AP_tree_edge *e4 = newBrother->edgeTo(oldBrother)->unlink();
688
689                AP_tree::moveNextTo(newBrother, rel_pos);
690
691                sortOldestFirst(&e1, &e2, &e3);
692                e1->relink(oldBrother, grandFather); // use oldest edge at remove position
693                thisFather->linkAllEdges(e2, e3, e4);
694            }
695            else { // grandson of root -> son of brother
696                // covered by test at PARS_main.cxx@COVER2
697                AP_tree_nlen *uncle = thisFather->get_brother();
698
699                thisFather->unlinkAllEdges(&e1, &e2, &e3);
700                AP_tree_edge *e4 = newBrother->edgeTo(oldBrother)->unlink();
701
702                AP_tree::moveNextTo(newBrother, rel_pos);
703
704                sortOldestFirst(&e1, &e2, &e3);
705                e1->relink(oldBrother, uncle);
706                thisFather->linkAllEdges(e2, e3, e4);
707            }
708        }
709        else { // son of root -> grandson of root
710            // covered by test at PARS_main.cxx@COVER1
711            oldBrother->unlinkAllEdges(&e1, &e2, &e3);
712            AP_tree::moveNextTo(newBrother, rel_pos);
713            thisFather->linkAllEdges(e1, e2, e3);
714        }
715    }
716    else if (grandFather==newBrothersFather) { // son -> brother of father
717        if (grandFather->father) {
718            // covered by test at PARS_main.cxx@COVER4
719            thisFather->unlinkAllEdges(&e1, &e2, &e3);
720            AP_tree_edge *e4 = grandFather->edgeTo(newBrother)->unlink();
721
722            AP_tree::moveNextTo(newBrother, rel_pos);
723
724            sortOldestFirst(&e1, &e2, &e3);
725            e1->relink(oldBrother, grandFather);
726            thisFather->linkAllEdges(e2, e3, e4);
727        }
728        else { // no edges change if we move grandson of root -> son of root
729            AP_tree::moveNextTo(newBrother, rel_pos);
730        }
731    }
732    else {
733        //  now we are sure, the minimal distance
734        //  between 'this' and 'newBrother' is 4 edges
735        //  or if the root-edge is between them, the
736        //  minimal distance is 3 edges
737
738        if (!grandFather) { // son of root
739            oldBrother->unlinkAllEdges(&e1, &e2, &e3);
740            AP_tree_edge *e4 = newBrother->edgeTo(newBrothersFather)->unlink();
741
742            AP_tree::moveNextTo(newBrother, rel_pos);
743
744            sortOldestFirst(&e1, &e2, &e3);
745            e1->relink(oldBrother->get_leftson(), oldBrother->get_rightson()); // new root-edge
746            thisFather->linkAllEdges(e2, e3, e4);   // old root
747        }
748        else if (!grandFather->get_father()) { // grandson of root
749            if (newBrothersFather->is_son_of_root()) { // grandson of root -> grandson of root
750                thisFather->unlinkAllEdges(&e1, &e2, &e3);
751                AP_tree_edge *e4 = newBrother->edgeTo(newBrothersFather)->unlink();
752
753                AP_tree::moveNextTo(newBrother, rel_pos);
754
755                sortOldestFirst(&e1, &e2, &e3);
756                e1->relink(oldBrother, newBrothersFather);  // new root-edge
757                thisFather->linkAllEdges(e2, e3, e4);
758            }
759            else {
760                AP_tree_nlen *uncle = thisFather->get_brother();
761
762                thisFather->unlinkAllEdges(&e1, &e2, &e3);
763                AP_tree_edge *e4 = newBrother->edgeTo(newBrothersFather)->unlink();
764
765                AP_tree::moveNextTo(newBrother, rel_pos);
766
767                sortOldestFirst(&e1, &e2, &e3);
768                e1->relink(oldBrother, uncle);
769                thisFather->linkAllEdges(e2, e3, e4);
770            }
771        }
772        else {
773            if (!newBrothersFather->get_father()) { // move to son of root
774                AP_tree_nlen *newBrothersBrother = newBrother->get_brother();
775
776                thisFather->unlinkAllEdges(&e1, &e2, &e3);
777                AP_tree_edge *e4 = newBrother->edgeTo(newBrothersBrother)->unlink();
778
779                AP_tree::moveNextTo(newBrother, rel_pos);
780
781                sortOldestFirst(&e1, &e2, &e3);
782                e1->relink(oldBrother, grandFather);
783                thisFather->linkAllEdges(e2, e3, e4);
784            }
785            else { // simple independent move
786                thisFather->unlinkAllEdges(&e1, &e2, &e3);
787                AP_tree_edge *e4 = newBrother->edgeTo(newBrothersFather)->unlink();
788
789                AP_tree::moveNextTo(newBrother, rel_pos);
790
791                sortOldestFirst(&e1, &e2, &e3);
792                e1->relink(oldBrother, grandFather);
793                thisFather->linkAllEdges(e2, e3, e4);
794            }
795        }
796    }
797
798    ASSERT_VALID_TREE(this);
799    ASSERT_VALID_TREE(rootNode());
800
801    ap_assert(is_leftson());
802    ap_assert(get_brother() == newBrother);
803}
804
805void AP_tree_nlen::unhash_sequence() {
806    /*! removes the parsimony sequence from an inner node
807     * (has no effect for leafs)
808     */
809
810    AP_sequence *sequence = get_seq();
811    if (sequence && !is_leaf()) sequence->forget_sequence();
812}
813
814bool AP_tree_nlen::acceptCurrentState(Level frame_level) {
815    // returns
816    // - true           if the top state has been removed
817    // - false          if the top state was kept/extended for possible revert at lower frame_level
818
819    if (remembered_for_frame != frame_level) {
820        ap_assert(0); // internal control number check failed
821        return false;
822    }
823
824    NodeState *state   = states.pop();
825    bool       removed = true;
826
827    Level next_frame_level   = frame_level-1;
828    Level stored_frame_level = state->frameNr;
829
830    if (!next_frame_level) { // accept() called at top-level
831        delete state;
832    }
833    else if (stored_frame_level == next_frame_level) {
834        // node already is buffered for next_frame_level
835
836        // if the currently accepted state->mode is not completely covered by previous state->mode
837        // => a future revert() would only restore partially
838        // To avoid that, move missing state information to previous NodeState
839        {
840            NodeState     *prev_state = states.top();
841            AP_STACK_MODE  prev_mode  = prev_state->mode;
842            AP_STACK_MODE  common     = AP_STACK_MODE(prev_mode & state->mode);
843
844            if (common != state->mode) {
845                AP_STACK_MODE missing = AP_STACK_MODE(state->mode & ~common); // previous is missing this state information
846
847                ap_assert((prev_mode&missing) == NOTHING);
848                state->move_info_to(*prev_state, missing);
849            }
850        }
851
852        delete state;
853    }
854    else {
855        // keep state for future revert
856        states.push(state);
857        removed = false;
858    }
859    remembered_for_frame = next_frame_level;
860
861    return removed;
862}
863
864
865bool AP_tree_nlen::rememberState(AP_STACK_MODE mode, Level frame_level) {
866    // according to mode
867    // tree_structure or sequence is buffered in the node
868
869    NodeState *store;
870    bool       ret;
871
872    if (is_leaf() && !(STRUCTURE & mode)) return false;    // tips push only structure
873
874    if (remembered_for_frame == frame_level) { // node already has a push (at current frame_level)
875        NodeState *is_stored = states.top();
876
877        if (0 == (mode & ~is_stored->mode)) { // already buffered
878            AP_sequence *sequence = get_seq();
879            if (sequence && (mode & SEQUENCE)) sequence->forget_sequence();
880            return false;
881        }
882        store = is_stored;
883        ret = false;
884    }
885    else { // first push for this node (in current stack frame)
886        store = new NodeState(remembered_for_frame);
887        states.push(store);
888
889        remembered_for_frame = frame_level;
890        ret                  = true;
891    }
892
893    if ((mode & (STRUCTURE|SEQUENCE)) && !(store->mode & (STRUCTURE|SEQUENCE))) {
894        store->mark_sum = gr.mark_sum;
895    }
896    if ((mode & STRUCTURE) && !(store->mode & STRUCTURE)) {
897        store->father    = get_father();
898        store->leftson   = get_leftson();
899        store->rightson  = get_rightson();
900        store->leftlen   = leftlen;
901        store->rightlen  = rightlen;
902        store->root      = get_tree_root();
903        store->gb_node   = gb_node;
904        store->keelState = keeledStateInfo();
905
906        if (!is_leaf()) store->remark = get_remark_ptr();
907
908        for (int e=0; e<3; e++) {
909            store->edge[e]      = edge[e];
910            store->edgeIndex[e] = index[e];
911        }
912    }
913
914    if (mode & SEQUENCE) {
915        ap_assert(!is_leaf()); // only allowed to push SEQUENCE for inner nodes
916        if (!(store->mode & SEQUENCE)) {
917            AP_sequence *sequence = take_seq();
918            store->sequence       = sequence;
919            store->mutations      = mutations;
920            mutations             = 0;
921        }
922        else {
923            AP_sequence *sequence = get_seq();
924            if (sequence) sequence->forget_sequence();
925        }
926    }
927
928    store->mode = (AP_STACK_MODE)(store->mode|mode);
929
930    return ret;
931}
932
933void NodeState::move_info_to(NodeState& target, AP_STACK_MODE what) {
934    // rescue partial NodeState information
935
936    ap_assert((mode&what)        == what); // this has to contain 'what' is moved
937    ap_assert((target.mode&what) == NOTHING); // target shall not already contain 'what' is moved
938
939    if ((what & (STRUCTURE|SEQUENCE)) && !(target.mode & (STRUCTURE|SEQUENCE))) {
940        target.mark_sum = mark_sum;
941    }
942    if (what & STRUCTURE) {
943        target.father    = father;
944        target.leftson   = leftson;
945        target.rightson  = rightson;
946        target.leftlen   = leftlen;
947        target.rightlen  = rightlen;
948        target.root      = root;
949        target.gb_node   = gb_node;
950        target.keelState = keelState;
951        target.remark    = remark;
952
953        for (int e=0; e<3; e++) {
954            target.edge[e]      = edge[e];
955            target.edgeIndex[e] = edgeIndex[e];
956        }
957    }
958    if (what & SEQUENCE) {
959        target.sequence  = sequence;
960        target.mutations = mutations;
961        sequence         = NULp;
962
963    }
964    // nothing needs to be done for ROOT
965    target.mode = AP_STACK_MODE(target.mode|what);
966}
967
968void AP_tree_nlen::restore_structure(const NodeState& state) {
969    father   = state.father;
970    leftson  = state.leftson;
971    rightson = state.rightson;
972    leftlen  = state.leftlen;
973    rightlen = state.rightlen;
974    set_tree_root(state.root);
975    gb_node  = state.gb_node;
976    setKeeledState(state.keelState);
977
978    if (!is_leaf()) use_as_remark(state.remark);
979
980    gr.mark_sum = state.mark_sum;
981
982    for (int e=0; e<3; e++) {
983        edge[e]  = state.edge[e];
984        index[e] = state.edgeIndex[e];
985        if (edge[e]) {
986            edge[e]->index[index[e]] = e;
987            edge[e]->node[index[e]]  = this;
988        }
989    }
990}
991void AP_tree_nlen::restore_sequence(NodeState& state) {
992    replace_seq(state.sequence);
993    state.sequence = NULp;
994    mutations      = state.mutations;
995    gr.mark_sum    = state.mark_sum;
996}
997void AP_tree_nlen::restore_sequence_nondestructive(const NodeState& state) {
998    replace_seq(state.sequence ? state.sequence->dup() : NULp);
999    mutations = state.mutations;
1000}
1001void AP_tree_nlen::restore_root(const NodeState& state) {
1002    state.root->change_root(state.root->get_root_node(), this);
1003}
1004
1005void AP_tree_nlen::restore(NodeState& state) {
1006    //! restore 'this' from NodeState (cheap; only call once for each 'state')
1007    AP_STACK_MODE mode = state.mode;
1008    if (mode&STRUCTURE) restore_structure(state);
1009    if (mode&SEQUENCE) restore_sequence(state);
1010    if (ROOT==mode) restore_root(state);
1011}
1012void AP_tree_nlen::restore_nondestructive(const NodeState& state) {
1013    //! restore 'this' from NodeState (expensive; may be called multiple times for each 'state')
1014    AP_STACK_MODE mode = state.mode;
1015    if (mode&STRUCTURE) restore_structure(state);
1016    if (mode&SEQUENCE) restore_sequence_nondestructive(state);
1017    if (ROOT==mode) restore_root(state);
1018}
1019
1020void AP_tree_nlen::revertToPreviousState(Level IF_ASSERTION_USED(curr_frameLevel), bool& IF_ASSERTION_USED(rootPopped)) { // pop old tree costs
1021    ap_assert(remembered_for_frame == curr_frameLevel); // error in node stack (node wasnt remembered in current frame!)
1022
1023    NodeState *previous = states.pop();
1024#if defined(ASSERTION_USED)
1025    if (previous->mode == ROOT) { // @@@ remove test code later
1026        ap_assert(!rootPopped); // only allowed once
1027        rootPopped = true;
1028    }
1029#endif
1030    restore(*previous);
1031
1032    remembered_for_frame = previous->frameNr;
1033    delete previous;
1034}
1035
1036void AP_tree_nlen::parsimony_rec(char *mutPerSite) {
1037    AP_combinableSeq *sequence = get_seq();
1038
1039    if (is_leaf()) {
1040        ap_assert(sequence); // tree w/o aliview?
1041        sequence->ensure_sequence_loaded();
1042    }
1043    else {
1044        if (!sequence) {
1045            sequence = set_seq(get_tree_root()->get_seqTemplate()->dup());
1046            ap_assert(sequence);
1047        }
1048
1049        if (!sequence->hasSequence()) {
1050            AP_tree_nlen *lson = get_leftson();
1051            AP_tree_nlen *rson = get_rightson();
1052
1053            ap_assert(lson);
1054            ap_assert(rson);
1055
1056            lson->parsimony_rec(mutPerSite);
1057            rson->parsimony_rec(mutPerSite);
1058
1059            AP_combinableSeq *lseq = lson->get_seq();
1060            AP_combinableSeq *rseq = rson->get_seq();
1061
1062            ap_assert(lseq);
1063            ap_assert(rseq);
1064
1065            Mutations mutations_for_combine = sequence->combine_seq(lseq, rseq, mutPerSite);
1066            mutations                   = lson->mutations + rson->mutations + mutations_for_combine;
1067        }
1068    }
1069}
1070
1071Mutations AP_tree_nlen::costs(char *mutPerSite) {
1072    // returns costs of a tree ( = number of mutations)
1073
1074    ap_assert(get_tree_root()->get_seqTemplate());  // forgot to set_seqTemplate() ?  (previously returned 0.0 in this case)
1075    ap_assert(sequence_state_valid());
1076
1077    parsimony_rec(mutPerSite);
1078    return mutations;
1079}
1080
1081Mutations AP_tree_nlen::nn_interchange_rec(EdgeSpec whichEdges, AP_BL_MODE mode) {
1082    if (!father) {
1083        return rootEdge()->nni_rec(whichEdges, mode, NULp, true);
1084    }
1085    if (!father->father) {
1086        AP_tree_edge *e = rootEdge();
1087        return e->nni_rec(whichEdges, mode, e->otherNode(this), false);
1088    }
1089    return edgeTo(get_father())->nni_rec(whichEdges, mode, get_father(), false);
1090}
1091
1092CONSTEXPR_INLINE AP_TREE_SIDE idx2side(const int idx) {
1093    return idx&1 ? AP_RIGHT : AP_LEFT;
1094}
1095
1096bool AP_tree_edge::kl_rec(const KL_params& KL, const int rec_depth, Mutations pars_best) {
1097    /*! does K.L. recursion
1098     * @param KL                parameters defining how recursion is done
1099     * @param rec_depth         current recursion depth (starts with 0)
1100     * @param pars_best         current parsimony value of topology
1101     */
1102
1103    ap_assert(!is_leaf_edge());
1104    if (rec_depth >= KL.max_rec_depth) return false;
1105
1106    ap_assert(implicated(rec_depth>0, kl_visited));
1107
1108    int           order[8];
1109    AP_tree_edge *descend[8];
1110
1111    {
1112        if (rec_depth == 0) {
1113            descend[0] = this;
1114            descend[2] = NULp;
1115            descend[4] = NULp;
1116            descend[6] = NULp;
1117        }
1118        else {
1119            AP_tree_nlen *son    = sonNode();
1120            AP_tree_nlen *notSon = otherNode(son); // brother or father
1121
1122            descend[0] = notSon->nextEdge(this);
1123            descend[2] = notSon->nextEdge(descend[0]);
1124
1125            ap_assert(descend[2] != this);
1126
1127            descend[4] = son->nextEdge(this);
1128            descend[6] = son->nextEdge(descend[4]);
1129
1130            ap_assert(descend[6] != this);
1131        }
1132
1133        descend[1] = descend[0];
1134        descend[3] = descend[2];
1135        descend[5] = descend[4];
1136        descend[7] = descend[6];
1137    }
1138
1139    // ---------------------------------
1140    //      detect parsimony values
1141
1142    ap_main->remember(); // @@@ i think this is unneeded. better reset root after all done in caller
1143    set_root();
1144    rootNode()->costs();
1145
1146    int rec_width_dynamic = 0;
1147    int visited_subtrees  = 0;
1148    int better_subtrees   = 0;
1149
1150    Mutations pars[8]; // eight parsimony values (produced by 2*swap_assymetric at each adjacent edge)
1151
1152#if defined(ASSERTION_USED)
1153    int forbidden_descends = 0;
1154#endif
1155    {
1156        AP_FLOAT schwellwert = KL.thresFunctor.calculate(rec_depth); // @@@ skip if not needed
1157        for (int i = 0; i < 8; i++) {
1158            order[i] = i;
1159            AP_tree_edge * const subedge = descend[i];
1160
1161            if (subedge                  &&
1162                !subedge->is_leaf_edge() &&
1163                !subedge->kl_visited     &&
1164                (!KL.stopAtFoldedGroups || !subedge->next_to_folded_group())
1165                )
1166            {
1167                ap_main->remember();
1168                subedge->sonNode()->swap_assymetric(idx2side(i));
1169                pars[i] = rootNode()->costs();
1170                if (pars[i] < pars_best) {
1171                    better_subtrees++;
1172                    pars_best = pars[i]; // @@@ do not overwrite yet; store and overwrite when done with this loop
1173                }
1174                if (pars[i] < schwellwert) {
1175                    rec_width_dynamic++;
1176                }
1177                ap_main->revert();
1178                visited_subtrees++;
1179            }
1180            else {
1181                pars[i] = -1;
1182#if defined(ASSERTION_USED)
1183                if (subedge && subedge->kl_visited) {
1184                    forbidden_descends++;
1185                }
1186#endif
1187            }
1188        }
1189    }
1190
1191    // bubblesort pars[]+order[], such that pars[0] contains best (=smallest) parsimony value
1192    {
1193        for (int i=7, t=0; t<i; t++) { // move negative (=unused) parsimony values to the end
1194            if (pars[t] <0) {
1195                pars[t]  = pars[i];
1196                order[t] = i;
1197                t--;
1198                i--;
1199            }
1200        }
1201
1202        for (int t = visited_subtrees - 1; t > 0; t--) {
1203            bool bubbled = false;
1204            for (int i = 0; i < t; i++) {
1205                if (pars[i] > pars[i+1]) {
1206                    std::swap(order[i], order[i+1]);
1207                    std::swap(pars[i],  pars[i+1]);
1208                    bubbled = true;
1209                }
1210            }
1211            if (!bubbled) break;
1212        }
1213    }
1214
1215#if defined(ASSERTION_USED)
1216    // rec_depth == 0 (called with start-node)
1217    // rec_depth == 1 (called twice with start-node (swap_assymetric AP_LEFT + AP_RIGHT))
1218    // rec_depth == 2 (called twice with each adjacent node -> 8 calls)
1219    // rec_depth == 3 (called twice with each adjacent node, but not with those were recursion came from -> 6 calls)
1220
1221    if (!is_root_edge()) {
1222        switch (rec_depth) {
1223            case 0:
1224                ap_assert(visited_subtrees == 2);
1225                ap_assert(forbidden_descends == 0);
1226                break;
1227            case 1:
1228                ap_assert(visited_subtrees <= 8);
1229                ap_assert(forbidden_descends == 0);
1230                break;
1231            default:
1232                ap_assert(visited_subtrees <= 6);
1233                ap_assert(forbidden_descends == 2);
1234                break;
1235        }
1236    }
1237    else { // at root
1238        switch (rec_depth) {
1239            case 0:
1240                ap_assert(visited_subtrees <= 2);
1241                break;
1242            case 1:
1243                ap_assert(visited_subtrees <= 8);
1244                ap_assert(forbidden_descends <= 2); // in case of subtree-optimization, 2 descends may be forbidden
1245                break;
1246            default:
1247                ap_assert(visited_subtrees <= 8);
1248                break;
1249        }
1250    }
1251#endif
1252
1253    int rec_width;
1254    if (better_subtrees) {
1255        rec_width = better_subtrees; // @@@ wrong if static/dynamic reduction would allow more
1256
1257        // @@@ IMO the whole concept of incrementing depth when a better topology was found has no positive effect
1258        // (the better topology is kept anyway and next recursive KL will do a full optimization starting from that edge as well)
1259
1260    }
1261    else {
1262        rec_width = visited_subtrees;
1263        if (KL.rec_type & AP_STATIC) {
1264            int rec_width_static = (rec_depth < CUSTOM_DEPTHS) ? KL.rec_width[rec_depth] : 1;
1265            rec_width            = std::min(rec_width, rec_width_static);
1266        }
1267        if (KL.rec_type & AP_DYNAMIK) {
1268            rec_width = std::min(rec_width, rec_width_dynamic);
1269        }
1270    }
1271    ap_assert(rec_width<=visited_subtrees);
1272
1273    bool found_better = false;
1274    for (int i=0; i<rec_width && !found_better; i++) {
1275        AP_tree_edge * const subedge = descend[order[i]];
1276
1277        ap_main->remember();
1278        subedge->kl_visited = true; // mark
1279        subedge->sonNode()->swap_assymetric(idx2side(order[i])); // swap
1280        rootNode()->parsimony_rec();
1281
1282        if (better_subtrees) {
1283            KL_params modified      = KL;
1284            modified.rec_type       = AP_STATIC;
1285            modified.max_rec_depth += KL.inc_rec_depth;
1286
1287            subedge->kl_rec(modified, rec_depth+1, pars_best);
1288            found_better = true;
1289        }
1290        else {
1291            found_better = subedge->kl_rec(KL, rec_depth+1, pars_best);
1292        }
1293
1294        subedge->kl_visited = false;      // unmark
1295        ap_main->accept_if(found_better); // revert
1296    }
1297
1298    ap_main->accept_if(found_better); // undo set_root otherwise
1299    return found_better;
1300}
1301
1302void AP_tree_nlen::exchange(AP_tree_nlen *other) {
1303    // exchange 'this' with other
1304    // 'this' has to be in tree; other has to be a "single node"
1305    //
1306    // Used by quick-add.
1307
1308    AP_tree_root *root = get_tree_root();
1309    ap_assert(root);
1310    ap_assert(!other->get_tree_root());
1311
1312    ASSERT_VALID_TREE(rootNode());
1313    ASSERT_VALID_TREE(other);
1314
1315    ap_main->push_node(this, STRUCTURE);
1316    ap_main->push_node(other, STRUCTURE);
1317    ap_main->push_node(get_father(), BOTH);
1318    push_all_upnode_sequences(get_father());
1319
1320    AP_tree_edge *myEdge    = nextEdge(NULp);
1321    AP_tree_nlen *connected = myEdge->otherNode(this);
1322    myEdge->unlink();
1323
1324    if (is_leftson()) {
1325        father->leftson = other;
1326    }
1327    else {
1328        father->rightson = other;
1329    }
1330    other->father = father;
1331    father        = NULp;
1332
1333    other->set_tree_root(root);
1334    set_tree_root(NULp);
1335
1336    myEdge->relink(other, connected);
1337
1338    ASSERT_VALID_TREE(rootNode());
1339    ASSERT_VALID_TREE(this);
1340}
1341
1342const char* AP_tree_nlen::sortByName() {
1343    if (name) return name;  // leafs
1344
1345    const char *n1 = get_leftson()->sortByName();
1346    const char *n2 = get_rightson()->sortByName();
1347
1348    if (strcmp(n1, n2)<0) return n1;
1349
1350    AP_tree::swap_sons();
1351
1352    return n2;
1353}
1354
1355const char *AP_tree_nlen::fullname() const {
1356    if (!name) {
1357        static char *buffer;
1358        char        *lName = ARB_strdup(get_leftson()->fullname());
1359        char        *rName = ARB_strdup(get_rightson()->fullname());
1360        int          len   = strlen(lName)+strlen(rName)+4;
1361
1362        if (buffer) free(buffer);
1363
1364        ARB_alloc(buffer, len);
1365
1366        strcpy(buffer, "[");
1367        strcat(buffer, lName);
1368        strcat(buffer, ",");
1369        strcat(buffer, rName);
1370        strcat(buffer, "]");
1371
1372        free(lName);
1373        free(rName);
1374
1375        return buffer;
1376    }
1377
1378    return name;
1379}
1380
1381
1382char* AP_tree_nlen::getSequenceCopy() {
1383    costs();
1384
1385    AP_sequence_parsimony *pseq = DOWNCAST(AP_sequence_parsimony*, get_seq());
1386    ap_assert(pseq->hasSequence());
1387
1388    size_t  len = pseq->get_sequence_length();
1389    char   *s   = new char[len];
1390    memcpy(s, pseq->get_sequence(), len);
1391
1392    return s;
1393}
1394
1395
1396GB_ERROR AP_pars_root::saveToDB() {
1397    has_been_saved = true;
1398    return AP_tree_root::saveToDB();
1399}
1400
Note: See TracBrowser for help on using the repository browser.