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
143
144
145
146
147
148
149
150
151
152
|
/**
* Author: Mark Larkin
*
* Copyright (c) 2007 Des Higgins, Julie Thompson and Toby Gibson.
*/
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include "RootedTreeOutput.h"
namespace clustalw
{
RootedTreeOutput::RootedTreeOutput(SeqInfo* seqInfo)
{
firstSeq = seqInfo->firstSeq;
lastSeq = seqInfo->lastSeq;
numSeqs = seqInfo->numSeqs;
}
void RootedTreeOutput::printPhylipTree(RootedGuideTree* tree, ofstream* ptrToFile, Alignment *alignPtr,
DistMatrix* distMat)
{
if(!ptrToFile || !ptrToFile->is_open())
{
return;
}
//cerr << "\n******************* DEBUG AW RootedTreeOutput::printPhylipTree: we only get here when using UPGMA with >2 seqs. Why?\n";
// If we have only 2 sequences, use the distances in the distMat
if (lastSeq - firstSeq + 1 == 2)
{
(*ptrToFile) << "(" << alignPtr->getName(firstSeq) << ":" << fixed << setprecision(5)
<< (*distMat)(firstSeq, firstSeq + 1) << "," << alignPtr->getName(firstSeq + 1)
<< ":" << fixed << setprecision(5) << (*distMat)(firstSeq, firstSeq + 1)
;
// << ")";
// AW 2009-05-15: final ")" added here, but not tested
// when is this ever used? tried
// infile=test.aln -clustering=nj|upgma
// on pairwise input, but this is only called for multiple
// input. why?
} else {
// AW 2009-05-08: fixed bug 155 by removing unnecessary, extra
// outer parenthesis from output
phylipTraverse(ptrToFile, alignPtr, tree->getRoot());
}
(*ptrToFile) << ";\n";
}
void RootedTreeOutput::printNexusTree(RootedGuideTree* tree, ofstream* ptrToFile,
Alignment *alignPtr, DistMatrix* distMat)
{
if(!ptrToFile || !ptrToFile->is_open())
{
return;
}
(*ptrToFile) << "#NEXUS\n\n";
(*ptrToFile) << "BEGIN TREES;\n\n";
(*ptrToFile) << "\tTRANSLATE\n";
for(int j = 1; j < numSeqs; j++)
{
(*ptrToFile) << "\t\t" << j << "\t" << alignPtr->getName(j) <<",\n";
}
(*ptrToFile) << "\t\t" << numSeqs << "\t" << alignPtr->getName(numSeqs) << "\n";
(*ptrToFile) << "\t\t;\n";
(*ptrToFile) << "\tUTREE PAUP_1= ";
// IF we have only 2 seqs
if (lastSeq - firstSeq + 1 == 2)
{
(*ptrToFile) << "(" << alignPtr->getName(firstSeq) << ":" << fixed << setprecision(5)
<< (*distMat)(firstSeq, firstSeq + 1) << "," << alignPtr->getName(firstSeq + 1)
<< ":" << fixed << setprecision(5) << (*distMat)(firstSeq, firstSeq + 1);
}
else
{
(*ptrToFile) << "(";
nexusTraverse(ptrToFile, alignPtr, tree->getRoot());
}
(*ptrToFile) << ");\n";
(*ptrToFile) << "\nENDBLOCK;\n";
}
/**
* PRIVATE FUNCTIONS
*/
void RootedTreeOutput::phylipTraverse(ofstream* ptrToFile, Alignment *alignPtr, Node* t)
{
if(!ptrToFile)
{
return;
}
if(t != 0)
{
if(t->isLeafNode())
{
if(alignPtr)
{
(*ptrToFile) << alignPtr->getName(t->getSeqNum()) << ":" << t->getHeight();
}
else
{
(*ptrToFile) << t->getSeqNum() << ":" << t->getHeight();
}
}
else // Internal node
{
(*ptrToFile) << "(\n";
phylipTraverse(ptrToFile, alignPtr, t->getLeft());
(*ptrToFile) << ",\n";
phylipTraverse(ptrToFile, alignPtr, t->getRight());
(*ptrToFile) << "):" << t->getHeight();
}
}
}
void RootedTreeOutput::nexusTraverse(ofstream* ptrToFile, Alignment *alignPtr, Node* t)
{
if(t != 0)
{
if(!t->isLeafNode()) // Internal node
{
(*ptrToFile) << "(";
nexusTraverse(ptrToFile, alignPtr, t->getLeft());
(*ptrToFile) << ",";
nexusTraverse(ptrToFile, alignPtr, t->getRight());
(*ptrToFile) << "):" << t->getHeight();
}
else // Leaf node
{
if(alignPtr)
{
(*ptrToFile) << alignPtr->getName(t->getSeqNum()) << ":" << t->getHeight();
}
else
{
(*ptrToFile) << t->getSeqNum() << ":" << t->getHeight();
}
}
}
}
}
|
