~ubuntu-branches/ubuntu/trusty/clustalx/trusty

« back to all changes in this revision

Viewing changes to clustalW/tree/UPGMA/RootedClusterTree.cpp

  • Committer: Bazaar Package Importer
  • Author(s): Charles Plessy, Steffen Moeller, Charles Plessy
  • Date: 2009-10-21 13:25:44 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20091021132544-r4hbcnjxp354wxh0
Tags: 2.0.12-1
* New upstream release (LP: #423648, #393769):
  - Uses Qt instead of lesstif.
  - Includes new code for UPGMA guide trees.
  - Includes iterative alignment facility.

[ Steffen Moeller ]
* New upstream release.
* Updated watch file (Closes: #550893).
* Removed LICENSE from debian/clustalx.docs
* rename to clustalx seems no longer required in debian/rules
* moved clustalx.1 into debian folder (eases working with svn-buildpackage)
* added German translation to desktop file

[ Charles Plessy ]
* Updated my email address.
* debian/copyright made machine-readable.
* Added various informations in debian/upstream-metadata.yaml.
* Switched to Debhelper 7.
  (debian/rules, debian/control, debian/patches, debian/compat)
* Removed useless Debhelper file debian/clustalx.dirs.
* Updated package description.
* Hardcoded the localisation of accessory files in /usr/share/clustalx.
  (debian/patches/hardcode-accessory-file-locations.patch)
* Documented in debian/README.source that the documentation for quilt
  is in /usr/share/doc/quilt.
* Added upstream changelog downloaded from upstream website
  (debian/rules, debian/CHANGELOG.upstream).
* Incremented Standards-Version to reflect conformance with Policy 3.8.3
  (debian/control, no other changes needed).
* Updated homepage in debian/clustalw.1.

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/**
 
2
 * Author: Mark Larkin
 
3
 * 
 
4
 * Copyright (c) 2007 Des Higgins, Julie Thompson and Toby Gibson.  
 
5
 */
 
6
#ifdef HAVE_CONFIG_H
 
7
    #include "config.h"
 
8
#endif
 
9
#include "RootedClusterTree.h"
 
10
#include "../../general/OutputFile.h"
 
11
#include "UPGMAAlgorithm.h"
 
12
#include "RootedTreeOutput.h"
 
13
//#include "../RandomGenerator.h"
 
14
namespace clustalw
 
15
{
 
16
 
 
17
auto_ptr<AlignmentSteps> RootedClusterTree::treeFromDistMatrix(RootedGuideTree* phyloTree,DistMatrix* distMat, Alignment *alignPtr, 
 
18
                                           int seq1, int nSeqs, string& phylipName)
 
19
{
 
20
    OutputFile phylipPhyTreeFile;
 
21
    auto_ptr<AlignmentSteps> progSteps;
 
22
    try
 
23
    {
 
24
        // Test to see if the inputs are valid
 
25
        if(seq1 < 1 || nSeqs < 1)
 
26
        {
 
27
            cerr << "Invalid inputs into treeFromDistMatrix \n"
 
28
                 << "seq1 = " << seq1 << " nSeqs = " << nSeqs << "\n"
 
29
                 << "Need to end program!\n";
 
30
            exit(1); 
 
31
            return progSteps;
 
32
        }
 
33
 
 
34
        float dist;
 
35
        string path;
 
36
        verbose = false;
 
37
        firstSeq = seq1;
 
38
        lastSeq = firstSeq + nSeqs - 1;
 
39
    
 
40
        SeqInfo info;
 
41
        info.firstSeq = firstSeq;
 
42
        info.lastSeq = lastSeq;
 
43
        info.numSeqs = nSeqs;
 
44
 
 
45
        utilityObject->getPath(userParameters->getSeqName(), &path);
 
46
        
 
47
        if(nSeqs >= 2)
 
48
        {
 
49
            string name = phylipName;
 
50
            if(!phylipPhyTreeFile.openFile(&name, 
 
51
                             "\nEnter name for new GUIDE TREE           file  ", &path, "dnd",
 
52
                             "Guide tree"))
 
53
            {
 
54
                return progSteps;
 
55
            }
 
56
            phylipName = name;                    
 
57
        }
 
58
        else
 
59
        {
 
60
            return progSteps;
 
61
        }
 
62
                
 
63
        RootedTreeOutput outputTree(&info);
 
64
        
 
65
        ofstream* ptrToFile = phylipPhyTreeFile.getPtrToFile();
 
66
        
 
67
        if (nSeqs == 2)
 
68
        {
 
69
            dist = (*distMat)(firstSeq, firstSeq + 1) / 2.0;
 
70
            if(ptrToFile->is_open())
 
71
            {
 
72
                (*ptrToFile) <<  "(" << alignPtr->getName(firstSeq) << ":" 
 
73
                             << setprecision(5)
 
74
                             << dist << "," << alignPtr->getName(firstSeq + 1) << ":" 
 
75
                             << setprecision(5) << dist <<");\n";
 
76
            }
 
77
            progSteps.reset(new AlignmentSteps);
 
78
            vector<int> groups;
 
79
            groups.resize(nSeqs + 1, 0);
 
80
            groups[1] = 1;
 
81
            groups[2] = 2;
 
82
        }
 
83
        else
 
84
        {
 
85
            UPGMAAlgorithm clusAlgorithm;
 
86
            progSteps = clusAlgorithm.generateTree(phyloTree, distMat, &info, false);
 
87
            outputTree.printPhylipTree(phyloTree, ptrToFile, alignPtr, distMat);
 
88
        }
 
89
        return progSteps;
 
90
    }
 
91
    catch(const exception &ex)
 
92
    {
 
93
        cerr << "ERROR: Error has occured in treeFromDistMatrix. " 
 
94
             << "Need to terminate program.\n"
 
95
             << ex.what();
 
96
        exit(1);   
 
97
    }
 
98
    catch(...)
 
99
    {
 
100
        cerr << "ERROR: Error has occured in treeFromDistMatrix. " 
 
101
             << "Need to terminate program.\n";      
 
102
        exit(1);  
 
103
    }
 
104
}
 
105
 
 
106
void RootedClusterTree::treeFromAlignment(TreeNames* treeNames, Alignment *alignPtr)
 
107
{
 
108
    try
 
109
    {
 
110
        OutputFile phylipPhyTreeFile;   
 
111
        OutputFile clustalPhyTreeFile;
 
112
        OutputFile distancesPhyTreeFile;
 
113
        OutputFile nexusPhyTreeFile;
 
114
        OutputFile pimFile;
 
115
        
 
116
        RootedGuideTree phyloTree;
 
117
        
 
118
        string path;
 
119
        int j;
 
120
        int overspill = 0;
 
121
        int totalDists;
 
122
        numSeqs = alignPtr->getNumSeqs(); // NOTE class variable
 
123
        
 
124
        /**
 
125
         * Check if numSeqs is ok
 
126
         */
 
127
        if(!checkIfConditionsMet(numSeqs, 2))
 
128
        {
 
129
            return;
 
130
        }
 
131
              
 
132
        firstSeq = 1;
 
133
        lastSeq = numSeqs;
 
134
        
 
135
        // The SeqInfo struct is passed to reduce the number of parameters passed!
 
136
        SeqInfo info;
 
137
        info.firstSeq = firstSeq;
 
138
        info.lastSeq = lastSeq;
 
139
        info.numSeqs = numSeqs;
 
140
    
 
141
        RootedTreeOutput outputTree(&info); // No bootstrap!
 
142
 
 
143
        utilityObject->getPath(userParameters->getSeqName(), &path);
 
144
 
 
145
        /**
 
146
         * Open the required output files.
 
147
         */
 
148
        if(!openFilesForTreeFromAlignment(&clustalPhyTreeFile, &phylipPhyTreeFile, 
 
149
                        &distancesPhyTreeFile, &nexusPhyTreeFile, &pimFile, treeNames, &path))
 
150
        {
 
151
            return; // Problem opeing one of the files, cannot continue!
 
152
        } 
 
153
    
 
154
        int _lenFirstSeq = alignPtr->getSeqLength(firstSeq);
 
155
    
 
156
        bootPositions.clear();
 
157
        bootPositions.resize(_lenFirstSeq + 2);
 
158
 
 
159
        for (j = 1; j <= _lenFirstSeq; ++j)
 
160
        {
 
161
            bootPositions[j] = j;
 
162
        }
 
163
 
 
164
        /**
 
165
         * Calculate quickDist and overspill
 
166
         */
 
167
        overspill = calcQuickDistMatForAll(clustalPhyTreeFile.getPtrToFile(),
 
168
                       phylipPhyTreeFile.getPtrToFile(), nexusPhyTreeFile.getPtrToFile(),
 
169
                       pimFile.getPtrToFile(), distancesPhyTreeFile.getPtrToFile(), alignPtr);
 
170
 
 
171
        // check if any distances overflowed the distance corrections 
 
172
        if (overspill > 0)
 
173
        {
 
174
            totalDists = (numSeqs *(numSeqs - 1)) / 2;
 
175
            overspillMessage(overspill, totalDists);
 
176
        }
 
177
 
 
178
        if (userParameters->getOutputTreeClustal())
 
179
        {
 
180
            verbose = true;
 
181
        }
 
182
        // Turn on file output 
 
183
    
 
184
 
 
185
        if (userParameters->getOutputTreeClustal() ||
 
186
            userParameters->getOutputTreePhylip() 
 
187
            || userParameters->getOutputTreeNexus())
 
188
        {
 
189
            UPGMAAlgorithm clusAlgorithm;
 
190
            clusAlgorithm.setVerbose(true);
 
191
            clusAlgorithm.generateTree(&phyloTree, quickDistMat.get(), &info, false,
 
192
                                       clustalPhyTreeFile.getPtrToFile());
 
193
            clusAlgorithm.setVerbose(false);
 
194
        }
 
195
 
 
196
        if (userParameters->getOutputTreePhylip())
 
197
        {
 
198
            outputTree.printPhylipTree(&phyloTree, phylipPhyTreeFile.getPtrToFile(), alignPtr,
 
199
                                        quickDistMat.get());
 
200
        }
 
201
 
 
202
        if (userParameters->getOutputTreeNexus())
 
203
        {
 
204
            outputTree.printNexusTree(&phyloTree, nexusPhyTreeFile.getPtrToFile(), alignPtr, 
 
205
                                       quickDistMat.get());
 
206
        }
 
207
    
 
208
        /** Free up resources!!!!! */
 
209
    
 
210
        treeGaps.clear();
 
211
        bootPositions.clear();
 
212
          
 
213
    }
 
214
    catch(const exception& ex)
 
215
    {
 
216
        cerr << ex.what() << endl;
 
217
        utilityObject->error("Terminating program. Cannot continue\n");
 
218
        exit(1);    
 
219
    }
 
220
}
 
221
                                
 
222
}