// =============================================================== // // // // File : AP_Tree.cxx // // Purpose : // // // // Institute of Microbiology (Technical University Munich) // // http://www.arb-home.de/ // // // // =============================================================== // #include "AP_Tree.hxx" #include "AP_TreeShader.hxx" #include #include #include #include #include #include #include #include #include using namespace std; /*!*************************************************************************************** ************ Rates ********** *****************************************************************************************/ void AP_rates::print() { AP_FLOAT max; int i; max = 0.0; for (i=0; i max) max = rates[i]; } printf("rates:"); for (i=0; iget_timestamp() <= this->update) return NULp; rate_len = fil->get_filtered_length(); delete [] rates; rates = new AP_FLOAT[rate_len]; for (int i=0; iupdate = fil->get_timestamp(); return NULp; } char *AP_rates::init(AP_FLOAT * ra, AP_filter *fil) { if (fil->get_timestamp() <= this->update) return NULp; rate_len = fil->get_filtered_length(); delete [] rates; rates = new AP_FLOAT[rate_len]; int i, j; for (j=i=0; iuse_position(j)) { rates[i++] = ra[j]; } } this->update = fil->get_timestamp(); return NULp; } AP_rates::~AP_rates() { delete [] rates; } /*!*************************************************************************************** ************ AP_tree_root ********** *****************************************************************************************/ AP_tree_root::AP_tree_root(AliView *aliView, AP_sequence *seq_proto, bool add_delete_callbacks, const group_scaling *scaling) : ARB_seqtree_root(aliView, seq_proto, add_delete_callbacks), root_changed_cb(NULp), root_changed_cd(NULp), node_deleted_cb(NULp), node_deleted_cd(NULp), gScale(scaling), gb_tree_gone(NULp), gone_tree_name(NULp), tree_timer(0), species_timer(0), rates(NULp) { GBDATA *gb_main = get_gb_main(); GB_transaction ta(gb_main); gb_species_data = GBT_get_species_data(gb_main); } AP_tree_root::~AP_tree_root() { predelete(); free(gone_tree_name); ap_assert(!get_root_node()); } bool AP_tree_root::is_tree_updated() { GBDATA *gbtree = get_gb_tree(); if (gbtree) { GB_transaction ta(gbtree); return GB_read_clock(gbtree)>tree_timer; } return true; } bool AP_tree_root::is_species_updated() { if (gb_species_data) { GB_transaction ta(gb_species_data); return GB_read_clock(gb_species_data)>species_timer; } return true; } void AP_tree_root::update_timers() { if (gb_species_data) { GB_transaction ta(GB_get_root(gb_species_data)); GBDATA *gbtree = get_gb_tree(); if (gbtree) tree_timer = GB_read_clock(gbtree); species_timer = GB_read_clock(gb_species_data); } } /*!*************************************************************************************** ************ AP_tree ********** *****************************************************************************************/ static void ap_tree_node_deleted(GBDATA *, AP_tree *tree) { tree->gb_node = NULp; } AP_tree::~AP_tree() { if (gr.callback_exists && gb_node) { GB_remove_callback(gb_node, GB_CB_DELETE, makeDatabaseCallback(ap_tree_node_deleted, this)); } AP_tree_root *root = get_tree_root(); if (root) root->inform_about_delete(this); } void AP_tree::initial_insert(AP_tree *new_brother, AP_tree_root *troot) { ap_assert(troot); ap_assert(is_leaf()); ap_assert(new_brother->is_leaf()); ap_assert(!troot->get_root_node()); ASSERT_VALID_TREE(this); ASSERT_VALID_TREE(new_brother); AP_tree *new_root = DOWNCAST(AP_tree*, troot->makeNode()); new_root->leftson = this; new_root->rightson = new_brother; new_root->father = NULp; father = new_root; new_brother->father = new_root; new_root->leftlen = 0.5; new_root->rightlen = 0.5; troot->change_root(NULp, new_root); set_tree_root(troot); new_brother->set_tree_root(troot); } void AP_tree::insert(AP_tree *new_brother) { ASSERT_VALID_TREE(this); ASSERT_VALID_TREE(new_brother); ap_assert(new_brother->get_tree_root()->get_root_node()->has_valid_root_remarks()); AP_tree *new_tree = DOWNCAST(AP_tree*, new_brother->get_tree_root()->makeNode()); AP_FLOAT laenge; if (new_brother->is_son_of_root()) { new_brother->get_father()->remove_root_remark(); } new_tree->leftson = this; new_tree->rightson = new_brother; new_tree->father = new_brother->father; father = new_tree; if (new_brother->father) { if (new_brother->father->leftson == new_brother) { laenge = new_brother->father->leftlen = new_brother->father->leftlen*.5; new_brother->father->leftson = new_tree; } else { laenge = new_brother->father->rightlen = new_brother->father->rightlen*.5; new_brother->father->rightson = new_tree; } } else { laenge = 0.5; } new_tree->leftlen = laenge; new_tree->rightlen = laenge; new_brother->father = new_tree; AP_tree_root *troot = new_brother->get_tree_root(); ap_assert(troot); // Note: initial_insert() has to be used to build initial tree if (!new_tree->father) troot->change_root(new_brother, new_tree); new_tree->set_tree_root(troot); set_tree_root(troot); ASSERT_VALID_TREE(troot->get_root_node()); ap_assert(get_tree_root()->get_root_node()->has_valid_root_remarks()); } void AP_tree_root::change_root(TreeNode *oldroot, TreeNode *newroot) { if (root_changed_cb) { // @@@ better call after calling base::change_root? root_changed_cb(root_changed_cd, DOWNCAST(AP_tree*, oldroot), DOWNCAST(AP_tree*, newroot)); } if (!oldroot) { ap_assert(newroot); if (gb_tree_gone) { // when tree was temporarily deleted (e.g. by 'Remove & add all') set_gb_tree(gb_tree_gone); // re-use previous DB entry gb_tree_gone = NULp; } } if (!newroot) { // tree empty GBDATA *gbtree = get_gb_tree(); if (gbtree) { ap_assert(!gb_tree_gone); // no tree should be remembered yet gb_tree_gone = gbtree; // remember for deletion (done in AP_tree::save) } } ARB_seqtree_root::change_root(oldroot, newroot); } void AP_tree_root::inform_about_delete(AP_tree *del) { if (node_deleted_cb) node_deleted_cb(node_deleted_cd, del); } void AP_tree_root::set_root_changed_callback(AP_rootChangedCb cb, void *cd) { root_changed_cb = cb; root_changed_cd = cd; } void AP_tree_root::set_node_deleted_callback(AP_nodeDelCb cb, void *cd) { node_deleted_cb = cb; node_deleted_cd = cd; } AP_tree *AP_tree::REMOVE() { // Remove this + father from tree. Father node will be destroyed. // Caller has to destroy the removed node (if intended). // // Warning: when removing the 2nd to last node from the tree, // the whole tree will be removed. // In that case both leaf nodes remain undestroyed. ASSERT_VALID_TREE(this); if (!father) { get_tree_root()->change_root(this, NULp); // tell AP_tree_root that the root node has been removed forget_origin(); // removed nodes are rootless } else { AP_tree *brother = get_brother(); // brother remains in tree GBT_LEN brothersLen = brother->get_branchlength(); AP_tree *fath = get_father(); // fath of this is removed as well ARB_seqtree *grandfather = fath->get_father(); AP_tree_root *troot = get_tree_root(); if (fath->gb_node) { // move inner information to remaining subtree if (!brother->gb_node && !brother->is_leaf()) { brother->gb_node = fath->gb_node; fath->gb_node = NULp; } } remove_remarks_from_this_and_parent(); // remove remarks of this + father if (grandfather) { brother->unlink_from_father(); bool wasLeftSon = fath->is_leftson(); fath->unlink_from_father(); if (wasLeftSon) { ap_assert(!grandfather->leftson); grandfather->leftlen += brothersLen; grandfather->leftson = brother; } else { ap_assert(!grandfather->rightson); grandfather->rightlen += brothersLen; grandfather->rightson = brother; } brother->father = grandfather; if (!grandfather->father) { ap_assert(brother->is_son_of_root()); if (!brother->is_leaf()) brother->remove_remark(); } } else { // father is root, make brother the new root if (brother->is_leaf()) { troot->change_root(fath, NULp); // erase tree from root brother->unlink_from_father(); // do not automatically delete brother } else { brother->unlink_from_father(); troot->change_root(fath, brother); } } ap_assert(fath == father); ASSERT_VALID_TREE_OR_NULL(troot->get_root_node()); troot->inform_about_delete(fath); troot->inform_about_delete(this); fath->forget_origin(); ASSERT_VALID_TREE(fath); unlink_from_father(); destroy(fath, troot); ASSERT_VALID_TREE(this); } return this; } GB_ERROR AP_tree::cantMoveNextTo(AP_tree *new_brother) { GB_ERROR error = NULp; if (!father) error = "Can't move the root of the tree"; else if (!new_brother->father) error = "Can't move to the root of the tree"; else if (new_brother->father == father) error = "Already there"; else if (new_brother == father) error = "Already there"; else if (!father->father) error = "Can't move son of root"; else if (new_brother->is_inside(this)) error = "Can't move a subtree into itself"; return error; } void AP_tree::moveNextTo(AP_tree *new_brother, AP_FLOAT rel_pos) { // rel_pos: 0.0 -> branch at father; 1.0 -> branch at brother; 0.5 -> branch at half distance between father and brother // @@@ "move subtree" needs better correction for groups (esp. if moving from root-of-group or when keeled groups are involved; see #785) ap_assert(father); ap_assert(new_brother); ap_assert(new_brother->father); ap_assert(new_brother->father != father); // already there ap_assert(new_brother != father); // already there ap_assert(!new_brother->is_inside(this)); // can't move tree into itself ap_assert(get_tree_root()->get_root_node()->has_valid_root_remarks()); remove_remarks_from_this_and_parent(); get_brother()->smart_remove_remark(); new_brother->remove_remarks_from_this_and_parent(); if (father->leftson != this) get_father()->swap_sons(); AP_tree *new_root = NULp; if (!father->father) { // move son of root ap_assert(!father->has_group_info()); get_father()->remove_root_remark(); new_root = get_brother(); new_root->father = NULp; ap_assert(!new_root->is_leaf()); // a leaf is no valid tree (should be impossible, because new_brother!=brother) } else { AP_tree *grandfather = get_father()->get_father(); if (!grandfather->father) { // move grandchild of root grandfather->remove_root_remark(); } if (grandfather->leftson == father) { grandfather->leftlen += father->rightlen; grandfather->leftson = father->rightson; } else { grandfather->rightlen += father->rightlen; grandfather->rightson = father->rightson; } father->rightson->father = grandfather; } AP_tree *new_tree = get_father(); AP_tree *brother_father = new_brother->get_father(); AP_FLOAT laenge; if (brother_father->leftson == new_brother) { laenge = brother_father->leftlen; laenge -= brother_father->leftlen = laenge * rel_pos; brother_father->leftson = new_tree; } else { laenge = brother_father->rightlen; laenge -= brother_father->rightlen = laenge * rel_pos; brother_father->rightson = new_tree; } new_tree->rightlen = laenge; new_brother->father = new_tree; new_tree->rightson = new_brother; new_tree->father = brother_father; if (new_root) { new_tree->get_tree_root()->change_root(new_tree, new_root); new_root->remove_root_remark(); } ap_assert(new_tree->get_tree_root()->get_root_node()->has_valid_root_remarks()); } inline int tree_read_byte(GBDATA *tree, const char *key, int init) { if (tree) { GBDATA *gbd = GB_entry(tree, key); if (gbd) return GB_read_byte(gbd); } return init; } inline float tree_read_float(GBDATA *tree, const char *key, float init) { if (tree) { GBDATA *gbd = GB_entry(tree, key); if (gbd) return GB_read_float(gbd); } return init; } //! moves all node/leaf information from struct TreeNode to AP_tree void AP_tree::load_node_info() { gr.spread = tree_read_float(gb_node, "spread", 1.0); gr.left_angle = tree_read_float(gb_node, "left_angle", 0.0); gr.right_angle = tree_read_float(gb_node, "right_angle", 0.0); gr.left_linewidth = tree_read_byte (gb_node, "left_linewidth", 0); gr.right_linewidth = tree_read_byte (gb_node, "right_linewidth", 0); gr.grouped = tree_read_byte (gb_node, "grouped", 0); } void AP_tree::load_subtree_info() { load_node_info(); if (!is_leaf()) { get_leftson()->load_subtree_info(); get_rightson()->load_subtree_info(); } } #if defined(DEBUG) #if defined(DEVEL_RALF) #define DEBUG_tree_write_byte #endif // DEVEL_RALF #endif // DEBUG static GB_ERROR tree_write_byte(GBDATA *gb_tree, AP_tree *node, short i, const char *key, int init) { GBDATA *gbd; GB_ERROR error = NULp; if (i==init) { if (node->gb_node) { gbd = GB_entry(node->gb_node, key); if (gbd) { #if defined(DEBUG_tree_write_byte) printf("[tree_write_byte] deleting db entry %p\n", gbd); #endif // DEBUG_tree_write_byte GB_delete(gbd); } } } else { if (!node->gb_node) { node->gb_node = GB_create_container(gb_tree, "node"); #if defined(DEBUG_tree_write_byte) printf("[tree_write_byte] created node-container %p\n", node->gb_node); #endif // DEBUG_tree_write_byte } gbd = GB_entry(node->gb_node, key); if (!gbd) { gbd = GB_create(node->gb_node, key, GB_BYTE); #if defined(DEBUG_tree_write_byte) printf("[tree_write_byte] created db entry %p\n", gbd); #endif // DEBUG_tree_write_byte } error = GB_write_byte(gbd, i); } return error; } static GB_ERROR tree_write_float(GBDATA *gb_tree, AP_tree *node, float f, const char *key, float init) { GB_ERROR error = NULp; if (f==init) { if (node->gb_node) { GBDATA *gbd = GB_entry(node->gb_node, key); if (gbd) error = GB_delete(gbd); } } else { if (!node->gb_node) { node->gb_node = GB_create_container(gb_tree, "node"); if (!node->gb_node) error = GB_await_error(); } if (!error) error = GBT_write_float(node->gb_node, key, f); } return error; } GB_ERROR AP_tree::tree_write_tree_rek(GBDATA *gb_tree) { GB_ERROR error = NULp; if (!is_leaf()) { error = get_leftson()->tree_write_tree_rek(gb_tree); if (!error) error = get_rightson()->tree_write_tree_rek(gb_tree); if (!error) error = tree_write_float(gb_tree, this, gr.spread, "spread", 1.0); if (!error) error = tree_write_float(gb_tree, this, gr.left_angle, "left_angle", 0.0); if (!error) error = tree_write_float(gb_tree, this, gr.right_angle, "right_angle", 0.0); if (!error) error = tree_write_byte (gb_tree, this, gr.left_linewidth, "left_linewidth", 0); if (!error) error = tree_write_byte (gb_tree, this, gr.right_linewidth, "right_linewidth", 0); if (!error) error = tree_write_byte (gb_tree, this, gr.grouped, "grouped", 0); } return error; } GB_ERROR AP_tree_root::saveToDB() { GB_ERROR error = GB_push_transaction(get_gb_main()); if (get_gb_tree()) { error = get_root_node()->tree_write_tree_rek(get_gb_tree()); } else { ap_assert(!gb_tree_gone); // should have been handled by caller (e.g. in AWT_graphic_tree::save) } if (!error) { if (!get_gb_tree() && gone_tree_name) { // tree was deleted before GBDATA *gb_tree_exists = GBT_find_tree(get_gb_main(), gone_tree_name); if (gb_tree_exists) { error = "tree already exists"; } else { error = GBT_write_tree(get_gb_main(), gone_tree_name, get_root_node()); if (!error) { gb_tree_exists = GBT_find_tree(get_gb_main(), gone_tree_name); ap_assert(gb_tree_exists); if (gb_tree_exists) { set_gb_tree_and_name(GBT_find_tree(get_gb_main(), gone_tree_name), gone_tree_name); aw_message(GBS_global_string("Recreated previously deleted '%s'", gone_tree_name)); freenull(gone_tree_name); } } } if (error) aw_message(GBS_global_string("Failed to recreate '%s' (Reason: %s)", gone_tree_name, error)); } if (!error) error = ARB_seqtree_root::saveToDB(); } if (!error) update_timers(); return GB_end_transaction(get_gb_main(), error); } inline GBDATA *find_group_name_entry(TreeNode *node) { return node->has_group_info() ? GB_entry(node->gb_node, "group_name") : NULp; } inline void TreeNode::swap_node_info(TreeNode *other, bool ofKeeledGroups) { if (ofKeeledGroups) { // member 'inverseLeft' cannot be handled correctly w/o knowing son nodes get_father()->swap_node_info(other->get_father(), false); fixKeeledOrientation(); other->fixKeeledOrientation(); } else if (this == other) { gb_assert(keeledOver && other->keeledOver); inverseLeft = !inverseLeft; } else { std::swap(name, other->name); std::swap(gb_node, other->gb_node); std::swap(keeledOver, other->keeledOver); } } GB_ERROR AP_tree::swap_group_with(AP_tree *dest, bool swapKeeled) { GB_ERROR error = NULp; if (swapKeeled) { ap_assert(!is_leaf() && is_keeled_group()); AP_tree *parent = get_father(); AP_tree *dest_parent = dest->get_father(); if (parent != dest_parent && dest_parent->has_group_info() && dest_parent->keelTarget() != dest ) { error = GBS_global_string("cannot move group '%s' (would create partial overlap with '%s')", parent->name, dest_parent->name); } if (!error && dest_parent->is_root_node()) { error = "invalid move of keeled group to tree-root"; } if (!error) { swap_node_info(dest, true); parent->load_node_info(); dest_parent->load_node_info(); } } else { ap_assert(!is_leaf() && is_normal_group()); if (dest->has_group_info() && !dest->is_normal_group()) { error = GBS_global_string("cannot move group '%s' (would create partial overlap with '%s')", name, dest->name); } if (!error) { swap_node_info(dest, false); load_node_info(); dest->load_node_info(); } } return error; } GB_ERROR AP_tree::move_group_to(AP_tree *dest) { GB_ERROR error = NULp; bool src_normal = !is_leaf() && is_normal_group(); bool src_keeled = !is_leaf() && is_keeled_group(); if (!src_normal && !src_keeled) { error = "Please select a valid source group"; } else { if (dest->is_leaf()) { error = GBS_global_string("'%s' is no valid destination for a group", dest->name); } else { if (src_keeled) { error = swap_group_with(dest, true); if (error && src_normal) { GB_ERROR error1 = error; error = swap_group_with(dest, false); if (error) error = GBS_global_string("Neighter keeled nor normal group can be moved that way:\n%s\n%s", error1, error); } } else { error = swap_group_with(dest, false); } if (!error) { GBDATA *gb_retax = NULp; if (!gb_retax && src_normal) gb_retax = find_group_name_entry(dest); if (!gb_retax && src_keeled) gb_retax = find_group_name_entry(dest->get_father()); if (!gb_retax) gb_retax = find_group_name_entry(this); if (!gb_retax) gb_retax = find_group_name_entry(get_father()); if (gb_retax) GB_touch(gb_retax); // force taxonomy reload ap_assert(gb_retax); // should normally always find at least one name } } } return error; } #if defined(ASSERTION_USED) || defined(UNIT_TESTS) bool AP_tree::has_correct_mark_flags() const { if (is_leaf()) return true; if (!get_leftson() ->has_correct_mark_flags()) return false; if (!get_rightson()->has_correct_mark_flags()) return false; const AP_tree_members& left = get_leftson()->gr; const AP_tree_members& right = get_rightson()->gr; unsigned wanted_mark_sum = left.mark_sum + right.mark_sum; return gr.mark_sum == wanted_mark_sum; } #endif class AP_DefaultTreeShader: public AP_TreeShader { // default tree shader (as used by unit-tests and ARB_PARSIMONY) static void default_shader_never_shades() { ap_assert(0); } public: AP_DefaultTreeShader() {} void init() OVERRIDE {} void update_settings() OVERRIDE { colorize_marked = true; colorize_groups = AW_color_groups_active(); shade_species = false; } ShadedValue calc_shaded_leaf_GC(GBDATA */*gb_node*/) const OVERRIDE { default_shader_never_shades(); return NULp; } ShadedValue calc_shaded_inner_GC(const ShadedValue& /*left*/, float /*left_ratio*/, const ShadedValue& /*right*/) const OVERRIDE { default_shader_never_shades(); return NULp; } int to_GC(const ShadedValue& /*val*/) const OVERRIDE { default_shader_never_shades(); return -1; } }; AP_TreeShader *AP_tree::shader = new AP_DefaultTreeShader; void AP_tree::set_tree_shader(AP_TreeShader *new_shader) { delete shader; shader = new_shader; shader->init(); } inline void AP_tree::recalc_hidden() { // gr.hidden of father needs to be up-to-date! gr.hidden = get_father() && (get_father()->gr.hidden || get_father()->gr.grouped); } inline void AP_tree::recalc_view_sum(const group_scaling& gscaling) { // childs need to be up-to-date // gr.leaf_sum needs to be up-to-date ap_assert(!is_leaf()); if (is_folded_group()) { ap_assert(gscaling.has_been_set()); const unsigned MIN_GROUP_SIZE = 2U; unsigned squared_size = unsigned(pow(double(gr.leaf_sum), gscaling.pow) * gscaling.linear); gr.view_sum = std::max(squared_size, MIN_GROUP_SIZE); gr.view_sum = std::min(gr.leaf_sum, gr.view_sum); // folded group will never use more space than unfolded } else { gr.view_sum = get_leftson()->gr.view_sum + get_rightson()->gr.view_sum; } } GB_ERROR AP_tree::update_and_write_folding(GBDATA *gb_tree, const group_scaling& gscaling) { // Warning: only use if you know what you are doing! // // recalculates gr.hidden and gr.view_sum (after gr.grouped was modified) // writes changed gr.grouped to database GB_ERROR error = NULp; recalc_hidden(); if (!is_leaf()) { error = get_leftson()->update_and_write_folding(gb_tree, gscaling); if (!error) error = get_rightson()->update_and_write_folding(gb_tree, gscaling); if (!error) { recalc_view_sum(gscaling); error = tree_write_byte(gb_tree, this, gr.grouped, "grouped", 0); } } return error; } void AP_tree::recompute_and_write_folding() { // Warning: only use if you know what you are doing! // // Usable only if: folding changed (but nothing else!) // Call with root-node of subtree for which folding shall be recalculated // (needs to be the root of ALL folding changes!). // Need to do resize afterwards. // // Note: gr.hidden is assumed to be correct for parent nodes! AP_tree_root *troot = get_tree_root(); GBDATA *gb_tree = troot->get_gb_tree(); const group_scaling *gscaling = troot->get_group_scaling(); ap_assert(gb_tree); ap_assert(gscaling); GB_ERROR error = update_and_write_folding(gb_tree, *gscaling); if (error) aw_message(GBS_global_string("Error in folding-update: %s", error)); // update view_sum of parent-nodes { AP_tree *fath = get_father(); while (fath) { fath->recalc_view_sum(*gscaling); fath = fath->get_father(); } } } template<> ShadedValue AP_tree::update_subtree_information(const group_scaling& gscaling) { ShadedValue val; recalc_hidden(); if (is_leaf()) { gr.view_sum = 1; gr.leaf_sum = 1; gr.max_tree_depth = 0.0; gr.min_tree_depth = 0.0; bool is_marked = gb_node && GB_read_flag(gb_node); gr.mark_sum = int(is_marked); gr.gc = shader->calc_leaf_GC(gb_node, is_marked); if (shader->does_shade()) { val = shader->calc_shaded_leaf_GC(gb_node); if (gr.gc == AWT_GC_NONE_MARKED) { gr.gc = shader->to_GC(val); } } } else { ShadedValue leftVal = get_leftson()->update_subtree_information(gscaling); ShadedValue rightVal = get_rightson()->update_subtree_information(gscaling); { AP_tree_members& left = get_leftson()->gr; AP_tree_members& right = get_rightson()->gr; gr.leaf_sum = left.leaf_sum + right.leaf_sum; recalc_view_sum(gscaling); gr.min_tree_depth = std::min(leftlen+left.min_tree_depth, rightlen+right.min_tree_depth); gr.max_tree_depth = std::max(leftlen+left.max_tree_depth, rightlen+right.max_tree_depth); gr.mark_sum = left.mark_sum + right.mark_sum; gr.gc = shader->calc_inner_GC(left.gc, right.gc); if (shader->does_shade()) { float left_weight = left.leaf_sum / float(gr.leaf_sum); val = shader->calc_shaded_inner_GC(leftVal, left_weight, rightVal); if (gr.gc == AWT_GC_NONE_MARKED) { gr.gc = shader->to_GC(val); } } } } ap_assert(implicated(shader->does_shade(), val.isSet())); // expect we have shaded value (if shading is performed) return val; } unsigned AP_tree::count_leafs() const { return is_leaf() ? 1 : get_leftson()->count_leafs() + get_rightson()->count_leafs(); } int AP_tree::colorize(GB_HASH *hashptr) { // currently only used for multiprobes // colorizes the tree according to 'hashptr' // hashkey = species id // hashvalue = gc int res; if (is_leaf()) { if (gb_node) { res = GBS_read_hash(hashptr, name); if (!res && GB_read_flag(gb_node)) { // marked but not in hash -> black res = AWT_GC_BLACK; } } else { res = AWT_GC_ONLY_ZOMBIES; } } else { int l = get_leftson()->colorize(hashptr); int r = get_rightson()->colorize(hashptr); if (l == r) res = l; else if (l == AWT_GC_ONLY_ZOMBIES) res = r; else if (r == AWT_GC_ONLY_ZOMBIES) res = l; else res = AWT_GC_SOME_MARKED; } gr.gc = res; return res; } void AP_tree::compute_tree() { #if defined(DEVEL_RALF) && 0 fputs(" - AP_tree::compute_tree() called\n", stderr); #endif AP_tree_root *troot = get_tree_root(); const group_scaling *gscaling = troot->get_group_scaling(); ap_assert(gscaling && gscaling->has_been_set()); { GB_transaction ta(troot->get_gb_main()); shader->update_settings(); update_subtree_information(*gscaling); } } GB_ERROR AP_tree_root::loadFromDB(const char *name) { GB_ERROR error = ARB_seqtree_root::loadFromDB(name); if (!error) { get_root_node()->load_subtree_info(); } update_timers(); // maybe after link() ? // @@@ really do if error? return error; } GB_ERROR AP_tree::relink() { GB_transaction ta(get_tree_root()->get_gb_main()); // open close a transaction GB_ERROR error = GBT_link_tree(this, get_tree_root()->get_gb_main(), false, NULp, NULp); // no status get_tree_root()->update_timers(); return error; } AP_UPDATE_FLAGS AP_tree_root::check_update() { GBDATA *gb_main = get_gb_main(); if (!gb_main) { return AP_UPDATE_RELOADED; } else { GB_transaction ta(gb_main); if (is_tree_updated()) return AP_UPDATE_RELOADED; if (is_species_updated()) return AP_UPDATE_RELINKED; return AP_UPDATE_OK; } } void AP_tree::buildLeafList_rek(AP_tree **list, long& num) { // builds a list of all species if (!is_leaf()) { get_leftson()->buildLeafList_rek(list, num); get_rightson()->buildLeafList_rek(list, num); } else { list[num] = this; num++; } } void AP_tree::buildLeafList(AP_tree **&list, long &num) { num = count_leafs(); list = new AP_tree *[num+1]; list[num] = NULp; long count = 0; buildLeafList_rek(list, count); ap_assert(count == num); } long AP_tree_root::remove_leafs(AWT_RemoveType awt_remove_type) { // may remove the complete tree ASSERT_VALID_TREE(get_root_node()); AP_tree **list; long count; get_root_node()->buildLeafList(list, count); GB_transaction ta(get_gb_main()); long removed = 0; for (long i=0; igb_node) { if ((awt_remove_type & AWT_REMOVE_NO_SEQUENCE) && !leaf->get_seq()) { removeNode = true; } else if (awt_remove_type & (AWT_REMOVE_MARKED|AWT_REMOVE_UNMARKED)) { long flag = GB_read_flag(list[i]->gb_node); removeNode = (flag && (awt_remove_type&AWT_REMOVE_MARKED)) || (!flag && (awt_remove_type&AWT_REMOVE_UNMARKED)); } } else { if (awt_remove_type & AWT_REMOVE_ZOMBIES) { removeNode = true; } } if (removeNode) { destroyNode(list[i]->REMOVE()); removed++; if (!get_root_node()) { break; // tree has been deleted } } ASSERT_VALID_TREE(get_root_node()); } delete [] list; ASSERT_VALID_TREE_OR_NULL(get_root_node()); return removed; } // -------------------------------------------------------------------------------- template class ValueCounter { T min, max, sum; int count; char *mean_min_max_impl() const; char *mean_min_max_percent_impl() const; mutable char *buf; const char *set_buf(char *content) const { freeset(buf, content); return buf; } public: ValueCounter() : min(INT_MAX), max(INT_MIN), sum(0), count(0), buf(NULp) {} ValueCounter(const ValueCounter& other) : min(other.min), max(other.max), sum(other.sum), count(other.count), buf(NULp) {} ~ValueCounter() { free(buf); } DECLARE_ASSIGNMENT_OPERATOR(ValueCounter); void count_value(T val) { count++; min = std::min(min, val); max = std::max(max, val); sum += val; } int get_count() const { return count; } T get_min() const { return min; } T get_max() const { return max; } double get_mean() const { return sum/double(count); } const char *mean_min_max() const { return count ? set_buf(mean_min_max_impl()) : ""; } const char *mean_min_max_percent() const { return count ? set_buf(mean_min_max_percent_impl()) : ""; } ValueCounter& operator += (const T& inc) { min += inc; max += inc; sum += inc*count; return *this; } ValueCounter& operator += (const ValueCounter& other) { min = std::min(min, other.min); max = std::max(max, other.max); sum += other.sum; count += other.count; return *this; } }; template inline ValueCounter operator+(const ValueCounter& c1, const ValueCounter& c2) { return ValueCounter(c1) += c2; } template inline ValueCounter operator+(const ValueCounter& c, const T& inc) { return ValueCounter(c) += inc; } template<> char *ValueCounter::mean_min_max_impl() const { return GBS_global_string_copy("%.2f (range: %i .. %i)", get_mean(), get_min(), get_max()); } template<> char *ValueCounter::mean_min_max_impl() const { return GBS_global_string_copy("%.2f (range: %.2f .. %.2f)", get_mean(), get_min(), get_max()); } template<> char *ValueCounter::mean_min_max_percent_impl() const { return GBS_global_string_copy("%.2f%% (range: %.2f%% .. %.2f%%)", get_mean()*100.0, get_min()*100.0, get_max()*100.0); } class LongBranchMarker { double min_rel_diff; double min_abs_diff; int leafs; int nonzeroleafs; int multifurcs; ValueCounter absdiff; ValueCounter reldiff; ValueCounter absdiff_marked; ValueCounter reldiff_marked; double perform_marking(AP_tree *at, bool& marked) { marked = false; if (at->is_leaf()) { if (at->get_branchlength() != 0.0) { nonzeroleafs++; } leafs++; return 0.0; } if (!at->is_root_node()) { if (at->get_branchlength() == 0.0) { // is multifurcation if (!at->get_father()->is_root_node() || at->is_leftson()) { // do not count two multifurcs @ sons of root multifurcs++; } } } bool marked_left; bool marked_right; double max = perform_marking(at->get_leftson(), marked_left) + at->leftlen; double min = perform_marking(at->get_rightson(), marked_right) + at->rightlen; bool max_is_left = true; if (maxmin_abs_diff && rel_diff>min_rel_diff) { if (max_is_left) { if (!marked_left) { at->get_leftson()->mark_subtree(); marked = true; } } else { if (!marked_right) { at->get_rightson()->mark_subtree(); marked = true; } } } if (marked) { // just marked one of my subtrees absdiff_marked.count_value(abs_diff); reldiff_marked.count_value(rel_diff); } else { marked = marked_left||marked_right; } return min; // use minimal distance for whole subtree } static char *meanDiffs(const ValueCounter& abs, const ValueCounter& rel) { return GBS_global_string_copy( "Mean absolute diff: %s\n" "Mean relative diff: %s", abs.mean_min_max(), rel.mean_min_max_percent()); } public: LongBranchMarker(AP_tree *root, double min_rel_diff_, double min_abs_diff_) : min_rel_diff(min_rel_diff_), min_abs_diff(min_abs_diff_), leafs(0), nonzeroleafs(0), multifurcs(0) { bool dummy; perform_marking(root, dummy); } const char *get_report() const { char *diffs_all = meanDiffs(absdiff, reldiff); char *diffs_marked = meanDiffs(absdiff_marked, reldiff_marked); int nodes = leafs_2_nodes(leafs, UNROOTED); int edges = nodes_2_edges(nodes); int zeroleafs = leafs-nonzeroleafs; int zeroedges = multifurcs+zeroleafs; int realedges = edges-zeroedges; int furcs = nodes-leafs; // = inner nodes int realfurcs = furcs-multifurcs; int node_digits = calc_digits(nodes); ap_assert(zeroleafs<=leafs); ap_assert(zeroedges<=edges); ap_assert(realedges<=edges); ap_assert(multifurcs<=furcs); ap_assert(realfurcs<=furcs); const char *msg = GBS_global_string( "Unrooted tree contains %*i furcations,\n" " of which %*i are multifurcations,\n" " i.e. %*i are \"real\" furcations.\n" "\n" "Unrooted tree contains %*i edges,\n" " of which %*i have a length > zero.\n" "\n" "%s\n" "\n" "%i subtrees have been marked:\n" "%s\n" "\n", node_digits, furcs, node_digits, multifurcs, node_digits, realfurcs, node_digits, edges, node_digits, realedges, diffs_all, absdiff_marked.get_count(), diffs_marked); free(diffs_all); free(diffs_marked); return msg; } double get_max_abs_diff() const { return absdiff.get_max(); } }; struct DepthMarker { // limits (marked if depth and dist are above) int min_depth; double min_rootdist; // current values (for recursion) int depth; double dist; // results ValueCounter depths, depths_marked; ValueCounter distances, distances_marked; void perform_marking(AP_tree *at, AP_FLOAT atLen) { int depthInc = atLen == 0.0 ? 0 : 1; // do NOT increase depth at multifurcations depth += depthInc; dist += atLen; if (at->is_leaf()) { depths.count_value(depth); distances.count_value(dist); int mark = depth >= min_depth && dist >= min_rootdist; if (at->gb_node) { GB_write_flag(at->gb_node, mark); if (mark) { depths_marked.count_value(depth); distances_marked.count_value(dist); } } } else { perform_marking(at->get_leftson(), at->leftlen); perform_marking(at->get_rightson(), at->rightlen); } depth -= depthInc; dist -= atLen; } public: DepthMarker(AP_tree *root, int min_depth_, double min_rootdist_) : min_depth(min_depth_), min_rootdist(min_rootdist_), depth(0), dist(0.0) { perform_marking(root, 0.0); } const char *get_report() const { int leafs = depths.get_count(); int marked = depths_marked.get_count(); double balanced_depth = log10(leafs) / log10(2); const char *msg = GBS_global_string( "The optimal mean depth of a tree with %i leafs\n" " would be %.2f\n" "\n" "Your tree:\n" "mean depth: %s\n" "mean distance: %s\n" "\n" "%i species (%.2f%%) have been marked:\n" "mean depth: %s\n" "mean distance: %s\n" , leafs, balanced_depth, depths.mean_min_max(), distances.mean_min_max(), marked, marked/double(leafs)*100.0, depths_marked.mean_min_max(), distances_marked.mean_min_max() ); return msg; } int get_max_depth() const { return depths.get_max(); } double get_max_rootdist() const { return distances.get_max(); } }; const char *AP_tree::mark_long_branches(double min_rel_diff, double min_abs_diff, double& found_max_abs_diff) { // look for asymmetric parts of the tree and mark all species with long branches LongBranchMarker lmarker(this, min_rel_diff, min_abs_diff); found_max_abs_diff = lmarker.get_max_abs_diff(); return lmarker.get_report(); } const char *AP_tree::mark_deep_leafs(int min_depth, double min_rootdist, int& found_max_depth, double& found_max_rootdist) { // mark all leafs with min_depth and min_rootdist DepthMarker dmarker(this, min_depth, min_rootdist); found_max_depth = dmarker.get_max_depth(); found_max_rootdist = dmarker.get_max_rootdist(); return dmarker.get_report(); } // -------------------------------------------------------------------------------- typedef ValueCounter Distance; class DistanceCounter { Distance min, max, mean; public: void count_distance(const Distance& d) { mean.count_value(d.get_mean()); min.count_value(d.get_min()); max.count_value(d.get_max()); } int get_count() const { return mean.get_count(); } char *get_report() const { return GBS_global_string_copy( "Mean mean distance: %s\n" "Mean min. distance: %s\n" "Mean max. distance: %s", mean.mean_min_max(), min.mean_min_max(), max.mean_min_max() ); } }; class EdgeDistances { typedef map DistanceMap; DistanceMap downdist; // inclusive length of branch itself DistanceMap updist; // inclusive length of branch itself GBT_LEN distSum; // of all distances in tree arb_progress progress; const Distance& calc_downdist(AP_tree *at, AP_FLOAT len) { if (at->is_leaf()) { Distance d; d.count_value(len); downdist[at] = d; progress.inc(); } else { downdist[at] = calc_downdist(at->get_leftson(), at->leftlen) + calc_downdist(at->get_rightson(), at->rightlen) + len; } return downdist[at]; } const Distance& calc_updist(AP_tree *at, AP_FLOAT len) { ap_assert(at->father); // impossible - root has no updist! AP_tree *father = at->get_father(); AP_tree *brother = at->get_brother(); if (father->father) { ap_assert(updist.find(father) != updist.end()); ap_assert(downdist.find(brother) != downdist.end()); updist[at] = updist[father] + downdist[brother] + len; } else { ap_assert(downdist.find(brother) != downdist.end()); updist[at] = downdist[brother]+len; } if (!at->is_leaf()) { calc_updist(at->get_leftson(), at->leftlen); calc_updist(at->get_rightson(), at->rightlen); } else { progress.inc(); } return updist[at]; } DistanceCounter alldists, markeddists; void calc_distance_stats(AP_tree *at) { if (at->is_leaf()) { ap_assert(updist.find(at) != updist.end()); const Distance& upwards = updist[at]; alldists.count_distance(upwards); if (at->gb_node && GB_read_flag(at->gb_node)) { markeddists.count_distance(upwards); } progress.inc(); } else { calc_distance_stats(at->get_leftson()); calc_distance_stats(at->get_rightson()); } } public: EdgeDistances(AP_tree *root) : distSum(root->sum_child_lengths()), progress("Analysing distances", root->count_leafs()*3L) { calc_downdist(root->get_leftson(), root->leftlen); calc_downdist(root->get_rightson(), root->rightlen); calc_updist(root->get_leftson(), root->leftlen); calc_updist(root->get_rightson(), root->rightlen); calc_distance_stats(root); } const char *get_report() const { char *alldists_report = alldists.get_report(); char *markeddists_report = markeddists.get_report(); const char *msg = GBS_global_string( "Overall in-tree-distance (ITD): %.3f\n" " per-species-distance (PSD): %.6f\n" "\n" "Distance statistic for %i leafs:\n" "(each leaf to all other leafs)\n" "\n" "%s\n" "\n" "Distance statistic for %i marked leafs:\n" "\n" "%s\n", distSum, distSum / alldists.get_count(), alldists.get_count(), alldists_report, markeddists.get_count(), markeddists_report); free(markeddists_report); free(alldists_report); return msg; } }; const char *AP_tree::analyse_distances() { return EdgeDistances(this).get_report(); } // -------------------------------------------------------------------------------- static int ap_mark_degenerated(AP_tree *at, double degeneration_factor, double& max_degeneration) { // returns number of species in subtree if (at->is_leaf()) return 1; int lSons = ap_mark_degenerated(at->get_leftson(), degeneration_factor, max_degeneration); int rSons = ap_mark_degenerated(at->get_rightson(), degeneration_factor, max_degeneration); double this_degeneration = 0; if (lSons= degeneration_factor) { at->get_leftson()->mark_subtree(); } } else if (rSons= degeneration_factor) { at->get_rightson()->mark_subtree(); } } if (this_degeneration >= max_degeneration) { max_degeneration = this_degeneration; } return lSons+rSons; } double AP_tree::mark_degenerated_branches(double degeneration_factor) { // marks all species in degenerated branches. // For all nodes, where one branch contains 'degeneration_factor' more species than the // other branch, the smaller branch is considered degenerated. double max_degeneration = 0; ap_mark_degenerated(this, degeneration_factor, max_degeneration); return max_degeneration; } static int ap_mark_duplicates_rek(AP_tree *at, GB_HASH *seen_species) { if (at->is_leaf()) { if (at->name) { if (GBS_read_hash(seen_species, at->name)) { // already seen -> mark species if (at->gb_node) { GB_write_flag(at->gb_node, 1); } else { // duplicated zombie return 1; } } else { // first occurrence GBS_write_hash(seen_species, at->name, 1); } } } else { return ap_mark_duplicates_rek(at->get_leftson(), seen_species) + ap_mark_duplicates_rek(at->get_rightson(), seen_species); } return 0; } void AP_tree::mark_duplicates() { GB_HASH *seen_species = GBS_create_hash(gr.leaf_sum, GB_IGNORE_CASE); int dup_zombies = ap_mark_duplicates_rek(this, seen_species); if (dup_zombies) { aw_message(GBS_global_string("Warning: Detected %i duplicated zombies (can't mark them)", dup_zombies)); } GBS_free_hash(seen_species); } static double ap_just_tree_rek(AP_tree *at) { if (at->is_leaf()) { return 0.0; } else { double bl = ap_just_tree_rek(at->get_leftson()); double br = ap_just_tree_rek(at->get_rightson()); double l = at->leftlen + at->rightlen; double diff = fabs(bl - br); if (l < diff * 1.1) l = diff * 1.1; double go = (bl + br + l) * .5; at->leftlen = go - bl; at->rightlen = go - br; return go; } } void AP_tree::justify_branch_lenghs(GBDATA *gb_main) { // shift branches to create a symmetric looking tree GB_transaction ta(gb_main); ap_just_tree_rek(this); } static void relink_tree_rek(AP_tree *node, void (*relinker)(GBDATA *&ref_gb_node, char *&ref_name, GB_HASH *organism_hash), GB_HASH *organism_hash) { if (node->is_leaf()) { relinker(node->gb_node, node->name, organism_hash); } else { relink_tree_rek(node->get_leftson(), relinker, organism_hash); relink_tree_rek(node->get_rightson(), relinker, organism_hash); } } void AP_tree::relink_tree(GBDATA *gb_main, void (*relinker)(GBDATA *&ref_gb_node, char *&ref_name, GB_HASH *organism_hash), GB_HASH *organism_hash) { // relinks the tree using a relinker-function // every node in tree is passed to relinker, relinker might modify // these values (ref_gb_node and ref_name) and the modified values are written back into tree GB_transaction ta(gb_main); relink_tree_rek(this, relinker, organism_hash); } void AP_tree::reset_child_angles() { if (!is_leaf()) { gr.reset_both_child_angles(); get_leftson()->reset_child_angles(); get_rightson()->reset_child_angles(); } } void AP_tree::reset_child_linewidths() { if (!is_leaf()) { gr.reset_both_child_linewidths(); get_leftson()->reset_child_linewidths(); get_rightson()->reset_child_linewidths(); } } void AP_tree::set_linewidth_recursive(int width) { set_linewidth(width); if (!is_leaf()) { get_leftson()->set_linewidth_recursive(width); get_rightson()->set_linewidth_recursive(width); } } void AP_tree::reset_child_layout() { if (!is_leaf()) { gr.reset_child_spread(); gr.reset_both_child_angles(); gr.reset_both_child_linewidths(); get_leftson()->reset_child_layout(); get_rightson()->reset_child_layout(); } } void AP_tree::reset_subtree_spreads() { gr.reset_child_spread(); if (!is_leaf()) { get_leftson()->reset_subtree_spreads(); get_rightson()->reset_subtree_spreads(); } } void AP_tree::reset_subtree_angles() { reset_angle(); if (!is_leaf()) reset_child_angles(); } void AP_tree::reset_subtree_linewidths() { reset_linewidth(); if (!is_leaf()) reset_child_linewidths(); } void AP_tree::reset_subtree_layout() { reset_linewidth(); reset_angle(); if (!is_leaf()) reset_child_layout(); } bool AP_tree::is_inside_folded_group() const { if (!is_leaf() && is_folded_group()) return true; if (!father) return false; return get_father()->is_inside_folded_group(); }