+++ /dev/null
-#include "muscle.h"\r
-#include "tree.h"\r
-#include "seqvect.h"\r
-#include "profile.h"\r
-#include "msa.h"\r
-#include "pwpath.h"\r
-#include "distfunc.h"\r
-#include "textfile.h"\r
-#include "estring.h"\r
-\r
-#define TRACE 0\r
-#define VALIDATE 0\r
-#define TRACE_LENGTH_DELTA 0\r
-\r
-static void LogLeafNames(const Tree &tree, unsigned uNodeIndex)\r
- {\r
- const unsigned uNodeCount = tree.GetNodeCount();\r
- unsigned *Leaves = new unsigned[uNodeCount];\r
- unsigned uLeafCount;\r
- GetLeaves(tree, uNodeIndex, Leaves, &uLeafCount);\r
- for (unsigned i = 0; i < uLeafCount; ++i)\r
- {\r
- if (i > 0)\r
- Log(",");\r
- Log("%s", tree.GetLeafName(Leaves[i]));\r
- }\r
- delete[] Leaves;\r
- }\r
-\r
-ProgNode *ProgressiveAlignE(const SeqVect &v, const Tree &GuideTree, MSA &a)\r
- {\r
- assert(GuideTree.IsRooted());\r
-\r
-#if TRACE\r
- Log("GuideTree:\n");\r
- GuideTree.LogMe();\r
-#endif\r
-\r
- const unsigned uSeqCount = v.Length();\r
- const unsigned uNodeCount = 2*uSeqCount - 1;\r
- const unsigned uIterCount = uSeqCount - 1;\r
-\r
- WEIGHT *Weights = new WEIGHT[uSeqCount];\r
- CalcClustalWWeights(GuideTree, Weights);\r
-\r
- ProgNode *ProgNodes = new ProgNode[uNodeCount];\r
-\r
- unsigned uJoin = 0;\r
- unsigned uTreeNodeIndex = GuideTree.FirstDepthFirstNode();\r
- SetProgressDesc("Align node");\r
- do\r
- {\r
- if (GuideTree.IsLeaf(uTreeNodeIndex))\r
- {\r
- if (uTreeNodeIndex >= uNodeCount)\r
- Quit("TreeNodeIndex=%u NodeCount=%u\n", uTreeNodeIndex, uNodeCount);\r
- ProgNode &Node = ProgNodes[uTreeNodeIndex];\r
- unsigned uId = GuideTree.GetLeafId(uTreeNodeIndex);\r
- if (uId >= uSeqCount)\r
- Quit("Seq index out of range");\r
- const Seq &s = *(v[uId]);\r
- Node.m_MSA.FromSeq(s);\r
- Node.m_MSA.SetSeqId(0, uId);\r
- Node.m_uLength = Node.m_MSA.GetColCount();\r
- Node.m_Weight = Weights[uId];\r
- // TODO: Term gaps settable\r
- Node.m_Prof = ProfileFromMSA(Node.m_MSA);\r
- Node.m_EstringL = 0;\r
- Node.m_EstringR = 0;\r
-#if TRACE\r
- Log("Leaf id=%u\n", uId);\r
- Log("MSA=\n");\r
- Node.m_MSA.LogMe();\r
- Log("Profile (from MSA)=\n");\r
- ListProfile(Node.m_Prof, Node.m_uLength, &Node.m_MSA);\r
-#endif\r
- }\r
- else\r
- {\r
- Progress(uJoin, uSeqCount - 1);\r
- ++uJoin;\r
-\r
- const unsigned uMergeNodeIndex = uTreeNodeIndex;\r
- ProgNode &Parent = ProgNodes[uMergeNodeIndex];\r
-\r
- const unsigned uLeft = GuideTree.GetLeft(uTreeNodeIndex);\r
- const unsigned uRight = GuideTree.GetRight(uTreeNodeIndex);\r
-\r
- if (g_bVerbose)\r
- {\r
- Log("Align: (");\r
- LogLeafNames(GuideTree, uLeft);\r
- Log(") (");\r
- LogLeafNames(GuideTree, uRight);\r
- Log(")\n");\r
- }\r
-\r
- ProgNode &Node1 = ProgNodes[uLeft];\r
- ProgNode &Node2 = ProgNodes[uRight];\r
-\r
-#if TRACE\r
- Log("AlignTwoMSAs:\n");\r
-#endif\r
- AlignTwoProfs(\r
- Node1.m_Prof, Node1.m_uLength, Node1.m_Weight,\r
- Node2.m_Prof, Node2.m_uLength, Node2.m_Weight,\r
- Parent.m_Path,\r
- &Parent.m_Prof, &Parent.m_uLength);\r
-#if TRACE_LENGTH_DELTA\r
- {\r
- unsigned L = Node1.m_uLength;\r
- unsigned R = Node2.m_uLength;\r
- unsigned P = Parent.m_Path.GetEdgeCount();\r
- unsigned Max = L > R ? L : R;\r
- unsigned d = P - Max;\r
- Log("LD%u;%u;%u;%u\n", L, R, P, d);\r
- }\r
-#endif\r
- PathToEstrings(Parent.m_Path, &Parent.m_EstringL, &Parent.m_EstringR);\r
-\r
- Parent.m_Weight = Node1.m_Weight + Node2.m_Weight;\r
-\r
-#if VALIDATE\r
- {\r
-#if TRACE\r
- Log("AlignTwoMSAs:\n");\r
-#endif\r
- PWPath TmpPath;\r
- AlignTwoMSAs(Node1.m_MSA, Node2.m_MSA, Parent.m_MSA, TmpPath);\r
- ProfPos *P1 = ProfileFromMSA(Node1.m_MSA, true);\r
- ProfPos *P2 = ProfileFromMSA(Node2.m_MSA, true);\r
- unsigned uLength = Parent.m_MSA.GetColCount();\r
- ProfPos *TmpProf = ProfileFromMSA(Parent.m_MSA, true);\r
-\r
-#if TRACE\r
- Log("Node1 MSA=\n");\r
- Node1.m_MSA.LogMe();\r
-\r
- Log("Node1 prof=\n");\r
- ListProfile(Node1.m_Prof, Node1.m_MSA.GetColCount(), &Node1.m_MSA);\r
- Log("Node1 prof (from MSA)=\n");\r
- ListProfile(P1, Node1.m_MSA.GetColCount(), &Node1.m_MSA);\r
-\r
- AssertProfsEq(Node1.m_Prof, Node1.m_uLength, P1, Node1.m_MSA.GetColCount());\r
-\r
- Log("Node2 prof=\n");\r
- ListProfile(Node2.m_Prof, Node2.m_MSA.GetColCount(), &Node2.m_MSA);\r
-\r
- Log("Node2 MSA=\n");\r
- Node2.m_MSA.LogMe();\r
-\r
- Log("Node2 prof (from MSA)=\n");\r
- ListProfile(P2, Node2.m_MSA.GetColCount(), &Node2.m_MSA);\r
-\r
- AssertProfsEq(Node2.m_Prof, Node2.m_uLength, P2, Node2.m_MSA.GetColCount());\r
-\r
- TmpPath.AssertEqual(Parent.m_Path);\r
-\r
- Log("Parent MSA=\n");\r
- Parent.m_MSA.LogMe();\r
-\r
- Log("Parent prof=\n");\r
- ListProfile(Parent.m_Prof, Parent.m_uLength, &Parent.m_MSA);\r
-\r
- Log("Parent prof (from MSA)=\n");\r
- ListProfile(TmpProf, Parent.m_MSA.GetColCount(), &Parent.m_MSA);\r
-\r
-#endif // TRACE\r
- AssertProfsEq(Parent.m_Prof, Parent.m_uLength,\r
- TmpProf, Parent.m_MSA.GetColCount());\r
- delete[] P1;\r
- delete[] P2;\r
- delete[] TmpProf;\r
- }\r
-#endif // VALIDATE\r
-\r
- Node1.m_MSA.Clear();\r
- Node2.m_MSA.Clear();\r
-\r
- // Don't delete profiles, may need them for tree refinement.\r
- //delete[] Node1.m_Prof;\r
- //delete[] Node2.m_Prof;\r
- //Node1.m_Prof = 0;\r
- //Node2.m_Prof = 0;\r
- }\r
- uTreeNodeIndex = GuideTree.NextDepthFirstNode(uTreeNodeIndex);\r
- }\r
- while (NULL_NEIGHBOR != uTreeNodeIndex);\r
- ProgressStepsDone();\r
-\r
- if (g_bBrenner)\r
- MakeRootMSABrenner((SeqVect &) v, GuideTree, ProgNodes, a);\r
- else\r
- MakeRootMSA(v, GuideTree, ProgNodes, a);\r
-\r
-#if VALIDATE\r
- {\r
- unsigned uRootNodeIndex = GuideTree.GetRootNodeIndex();\r
- const ProgNode &RootProgNode = ProgNodes[uRootNodeIndex];\r
- AssertMSAEq(a, RootProgNode.m_MSA);\r
- }\r
-#endif\r
-\r
- delete[] Weights;\r
- return ProgNodes;\r
- }\r