1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
|
#include "muscle.h"
#include "msa.h"
#include "tree.h"
#include "profile.h"
#include "pwpath.h"
#include "seqvect.h"
#include "estring.h"
#define TRACE 0
void DeleteProgNode(ProgNode &Node)
{
delete[] Node.m_Prof;
delete[] Node.m_EstringL;
delete[] Node.m_EstringR;
Node.m_Prof = 0;
Node.m_EstringL = 0;
Node.m_EstringR = 0;
}
static void MakeNode(ProgNode &OldNode, ProgNode &NewNode, bool bSwapLR)
{
if (bSwapLR)
{
NewNode.m_EstringL = OldNode.m_EstringR;
NewNode.m_EstringR = OldNode.m_EstringL;
}
else
{
NewNode.m_EstringL = OldNode.m_EstringL;
NewNode.m_EstringR = OldNode.m_EstringR;
}
NewNode.m_Prof = OldNode.m_Prof;
NewNode.m_uLength = OldNode.m_uLength;
NewNode.m_Weight = OldNode.m_Weight;
OldNode.m_Prof = 0;
OldNode.m_EstringL = 0;
OldNode.m_EstringR = 0;
}
void RealignDiffsE(const MSA &msaIn, const SeqVect &v,
const Tree &NewTree, const Tree &OldTree,
const unsigned uNewNodeIndexToOldNodeIndex[],
MSA &msaOut, ProgNode *OldProgNodes)
{
assert(OldProgNodes != 0);
const unsigned uNodeCount = NewTree.GetNodeCount();
if (uNodeCount%2 == 0)
Quit("RealignDiffs: Expected odd number of nodes");
const unsigned uMergeCount = (uNodeCount - 1)/2;
ProgNode *NewProgNodes = new ProgNode[uNodeCount];
for (unsigned uNewNodeIndex = 0; uNewNodeIndex < uNodeCount; ++uNewNodeIndex)
{
if (NODE_CHANGED == uNewNodeIndexToOldNodeIndex[uNewNodeIndex])
continue;
unsigned uOldNodeIndex = uNewNodeIndexToOldNodeIndex[uNewNodeIndex];
assert(uNewNodeIndex < uNodeCount);
assert(uOldNodeIndex < uNodeCount);
ProgNode &NewNode = NewProgNodes[uNewNodeIndex];
ProgNode &OldNode = OldProgNodes[uOldNodeIndex];
bool bSwapLR = false;
if (!NewTree.IsLeaf(uNewNodeIndex))
{
unsigned uNewLeft = NewTree.GetLeft(uNewNodeIndex);
unsigned uNewRight = NewTree.GetRight(uNewNodeIndex);
unsigned uOld = uNewNodeIndexToOldNodeIndex[uNewNodeIndex];
unsigned uOldLeft = OldTree.GetLeft(uOld);
unsigned uOldRight = OldTree.GetRight(uOld);
assert(uOldLeft < uNodeCount && uOldRight < uNodeCount);
if (uOldLeft != uNewNodeIndexToOldNodeIndex[uNewLeft])
{
assert(uOldLeft == uNewNodeIndexToOldNodeIndex[uNewRight]);
bSwapLR = true;
}
}
MakeNode(OldNode, NewNode, bSwapLR);
#if TRACE
Log("MakeNode old=%u new=%u swap=%d length=%u weight=%.3g\n",
uOldNodeIndex, uNewNodeIndex, bSwapLR, NewNode.m_uLength, NewNode.m_Weight);
#endif
}
unsigned uJoin = 0;
SetProgressDesc("Refine tree");
for (unsigned uNewNodeIndex = NewTree.FirstDepthFirstNode();
NULL_NEIGHBOR != uNewNodeIndex;
uNewNodeIndex = NewTree.NextDepthFirstNode(uNewNodeIndex))
{
if (NODE_CHANGED != uNewNodeIndexToOldNodeIndex[uNewNodeIndex])
continue;
Progress(uJoin, uMergeCount - 1);
++uJoin;
const unsigned uMergeNodeIndex = uNewNodeIndex;
ProgNode &Parent = NewProgNodes[uMergeNodeIndex];
const unsigned uLeft = NewTree.GetLeft(uNewNodeIndex);
const unsigned uRight = NewTree.GetRight(uNewNodeIndex);
ProgNode &Node1 = NewProgNodes[uLeft];
ProgNode &Node2 = NewProgNodes[uRight];
AlignTwoProfs(
Node1.m_Prof, Node1.m_uLength, Node1.m_Weight,
Node2.m_Prof, Node2.m_uLength, Node2.m_Weight,
Parent.m_Path,
&Parent.m_Prof, &Parent.m_uLength);
PathToEstrings(Parent.m_Path, &Parent.m_EstringL, &Parent.m_EstringR);
Parent.m_Weight = Node1.m_Weight + Node2.m_Weight;
delete[] Node1.m_Prof;
delete[] Node2.m_Prof;
Node1.m_Prof = 0;
Node2.m_Prof = 0;
}
ProgressStepsDone();
if (g_bBrenner)
MakeRootMSABrenner((SeqVect &) v, NewTree, NewProgNodes, msaOut);
else
MakeRootMSA(v, NewTree, NewProgNodes, msaOut);
#if DEBUG
AssertMSAEqIgnoreCaseAndGaps(msaIn, msaOut);
#endif
for (unsigned uNodeIndex = 0; uNodeIndex < uNodeCount; ++uNodeIndex)
DeleteProgNode(NewProgNodes[uNodeIndex]);
delete[] NewProgNodes;
}
|
