+++ /dev/null
-#include "muscle.h"\r
-#include "pwpath.h"\r
-#include "seq.h"\r
-#include "textfile.h"\r
-#include "msa.h"\r
-\r
-PWPath::PWPath()\r
- {\r
- m_uArraySize = 0;\r
- m_uEdgeCount = 0;\r
- m_Edges = 0;\r
- }\r
-\r
-PWPath::~PWPath()\r
- {\r
- Clear();\r
- }\r
-\r
-void PWPath::Clear()\r
- {\r
- delete[] m_Edges;\r
- m_Edges = 0;\r
- m_uArraySize = 0;\r
- m_uEdgeCount = 0;\r
- }\r
-\r
-void PWPath::ExpandPath(unsigned uAdditionalEdgeCount)\r
- {\r
- PWEdge *OldPath = m_Edges;\r
- unsigned uEdgeCount = m_uArraySize + uAdditionalEdgeCount;\r
-\r
- m_Edges = new PWEdge[uEdgeCount];\r
- m_uArraySize = uEdgeCount;\r
- if (m_uEdgeCount > 0)\r
- memcpy(m_Edges, OldPath, m_uEdgeCount*sizeof(PWEdge));\r
- delete[] OldPath;\r
- }\r
-\r
-void PWPath::AppendEdge(const PWEdge &Edge)\r
- {\r
- if (0 == m_uArraySize || m_uEdgeCount + 1 == m_uArraySize)\r
- ExpandPath(200);\r
-\r
- m_Edges[m_uEdgeCount] = Edge;\r
- ++m_uEdgeCount;\r
- }\r
-\r
-void PWPath::AppendEdge(char cType, unsigned uPrefixLengthA, unsigned uPrefixLengthB)\r
- {\r
- PWEdge e;\r
- e.uPrefixLengthA = uPrefixLengthA;\r
- e.uPrefixLengthB = uPrefixLengthB;\r
- e.cType = cType;\r
- AppendEdge(e);\r
- }\r
-\r
-void PWPath::PrependEdge(const PWEdge &Edge)\r
- {\r
- if (0 == m_uArraySize || m_uEdgeCount + 1 == m_uArraySize)\r
- ExpandPath(1000);\r
- if (m_uEdgeCount > 0)\r
- memmove(m_Edges + 1, m_Edges, sizeof(PWEdge)*m_uEdgeCount);\r
- m_Edges[0] = Edge;\r
- ++m_uEdgeCount;\r
- }\r
-\r
-const PWEdge &PWPath::GetEdge(unsigned uEdgeIndex) const\r
- {\r
- assert(uEdgeIndex < m_uEdgeCount);\r
- return m_Edges[uEdgeIndex];\r
- }\r
-\r
-void PWPath::Validate() const\r
- {\r
- const unsigned uEdgeCount = GetEdgeCount();\r
- if (0 == uEdgeCount)\r
- return;\r
- const PWEdge &FirstEdge = GetEdge(0);\r
- const PWEdge &LastEdge = GetEdge(uEdgeCount - 1);\r
- unsigned uStartA = FirstEdge.uPrefixLengthA;\r
- unsigned uStartB = FirstEdge.uPrefixLengthB;\r
- if (FirstEdge.cType != 'I')\r
- --uStartA;\r
- if (FirstEdge.cType != 'D')\r
- --uStartB;\r
-\r
- unsigned uPrefixLengthA = FirstEdge.uPrefixLengthA;\r
- unsigned uPrefixLengthB = FirstEdge.uPrefixLengthB;\r
- for (unsigned uEdgeIndex = 1; uEdgeIndex < uEdgeCount; ++uEdgeIndex)\r
- {\r
- const PWEdge &Edge = GetEdge(uEdgeIndex);\r
- switch (Edge.cType)\r
- {\r
- case 'M':\r
- if (uPrefixLengthA + 1 != Edge.uPrefixLengthA)\r
- Quit("PWPath::Validate MA %u", uPrefixLengthA);\r
- if (uPrefixLengthB + 1 != Edge.uPrefixLengthB)\r
- Quit("PWPath::Validate MB %u", uPrefixLengthB);\r
- ++uPrefixLengthA;\r
- ++uPrefixLengthB;\r
- break;\r
- case 'D':\r
- if (uPrefixLengthA + 1 != Edge.uPrefixLengthA)\r
- Quit("PWPath::Validate DA %u", uPrefixLengthA);\r
- if (uPrefixLengthB != Edge.uPrefixLengthB)\r
- Quit("PWPath::Validate DB %u", uPrefixLengthB);\r
- ++uPrefixLengthA;\r
- break;\r
- case 'I':\r
- if (uPrefixLengthA != Edge.uPrefixLengthA)\r
- Quit("PWPath::Validate IA %u", uPrefixLengthA);\r
- if (uPrefixLengthB + 1 != Edge.uPrefixLengthB)\r
- Quit("PWPath::Validate IB %u", uPrefixLengthB);\r
- ++uPrefixLengthB;\r
- break;\r
- }\r
- }\r
- }\r
-\r
-void PWPath::LogMe() const\r
- {\r
- for (unsigned uEdgeIndex = 0; uEdgeIndex < GetEdgeCount(); ++uEdgeIndex)\r
- {\r
- const PWEdge &Edge = GetEdge(uEdgeIndex);\r
- if (uEdgeIndex > 0)\r
- Log(" ");\r
- Log("%c%d.%d",\r
- Edge.cType,\r
- Edge.uPrefixLengthA,\r
- Edge.uPrefixLengthB);\r
- if ((uEdgeIndex > 0 && uEdgeIndex%10 == 0) ||\r
- uEdgeIndex == GetEdgeCount() - 1)\r
- Log("\n");\r
- }\r
- }\r
-\r
-void PWPath::Copy(const PWPath &Path)\r
- {\r
- Clear();\r
- const unsigned uEdgeCount = Path.GetEdgeCount();\r
- for (unsigned uEdgeIndex = 0; uEdgeIndex < uEdgeCount; ++uEdgeIndex)\r
- {\r
- const PWEdge &Edge = Path.GetEdge(uEdgeIndex);\r
- AppendEdge(Edge);\r
- }\r
- }\r
-\r
-void PWPath::FromMSAPair(const MSA &msaA, const MSA &msaB)\r
- {\r
- const unsigned uColCount = msaA.GetColCount();\r
- if (uColCount != msaB.GetColCount())\r
- Quit("PWPath::FromMSAPair, lengths differ");\r
-\r
- Clear();\r
-\r
- unsigned uPrefixLengthA = 0;\r
- unsigned uPrefixLengthB = 0;\r
- for (unsigned uColIndex = 0; uColIndex < uColCount; ++uColIndex)\r
- {\r
- bool bIsGapA = msaA.IsGapColumn(uColIndex);\r
- bool bIsGapB = msaB.IsGapColumn(uColIndex);\r
-\r
- PWEdge Edge;\r
- char cType;\r
- if (!bIsGapA && !bIsGapB)\r
- {\r
- cType = 'M';\r
- ++uPrefixLengthA;\r
- ++uPrefixLengthB;\r
- }\r
- else if (bIsGapA && !bIsGapB)\r
- {\r
- cType = 'I';\r
- ++uPrefixLengthB;\r
- }\r
- else if (!bIsGapA && bIsGapB)\r
- {\r
- cType = 'D';\r
- ++uPrefixLengthA;\r
- }\r
- else\r
- {\r
- assert(bIsGapB && bIsGapA);\r
- continue;\r
- }\r
-\r
- Edge.cType = cType;\r
- Edge.uPrefixLengthA = uPrefixLengthA;\r
- Edge.uPrefixLengthB = uPrefixLengthB;\r
- AppendEdge(Edge);\r
- }\r
- }\r
-\r
-// Very similar to HMMPath::FromFile, should consolidate.\r
-void PWPath::FromFile(TextFile &File)\r
- {\r
- Clear();\r
- char szToken[1024];\r
- File.GetTokenX(szToken, sizeof(szToken));\r
- if (0 != strcmp(szToken, "Path"))\r
- Quit("Invalid path file (Path)");\r
-\r
- File.GetTokenX(szToken, sizeof(szToken));\r
- if (0 != strcmp(szToken, "edges"))\r
- Quit("Invalid path file (edges)");\r
-\r
- File.GetTokenX(szToken, sizeof(szToken));\r
- if (!IsValidInteger(szToken))\r
- Quit("Invalid path file (edges value)");\r
-\r
- const unsigned uEdgeCount = (unsigned) atoi(szToken);\r
- unsigned uEdgeIndex = 0;\r
- for (unsigned uEdgeIndex = 0; uEdgeIndex < uEdgeCount; ++uEdgeIndex)\r
- {\r
- // index\r
- File.GetTokenX(szToken, sizeof(szToken));\r
- if (!IsValidInteger(szToken))\r
- Quit("Invalid path file, invalid index '%s'", szToken);\r
- unsigned n = (unsigned) atoi(szToken);\r
- if (n != uEdgeIndex)\r
- Quit("Invalid path file, expecting edge %u got %u", uEdgeIndex, n);\r
-\r
- // type\r
- File.GetTokenX(szToken, sizeof(szToken));\r
- if (1 != strlen(szToken))\r
- Quit("Invalid path file, expecting state, got '%s'", szToken);\r
- const char cType = szToken[0];\r
- if ('M' != cType && 'D' != cType && cType != 'I' && 'S' != cType)\r
- Quit("Invalid path file, expecting state, got '%c'", cType);\r
-\r
- // prefix length A\r
- File.GetTokenX(szToken, sizeof(szToken));\r
- if (!IsValidInteger(szToken))\r
- Quit("Invalid path file, bad prefix length A '%s'", szToken);\r
- const unsigned uPrefixLengthA = (unsigned) atoi(szToken);\r
-\r
- // prefix length B\r
- File.GetTokenX(szToken, sizeof(szToken));\r
- if (!IsValidInteger(szToken))\r
- Quit("Invalid path file, bad prefix length B '%s'", szToken);\r
- const unsigned uPrefixLengthB = (unsigned) atoi(szToken);\r
-\r
- PWEdge Edge;\r
- Edge.cType = cType;\r
- Edge.uPrefixLengthA = uPrefixLengthA;\r
- Edge.uPrefixLengthB = uPrefixLengthB;\r
- AppendEdge(Edge);\r
- }\r
- File.GetTokenX(szToken, sizeof(szToken));\r
- if (0 != strcmp(szToken, "//"))\r
- Quit("Invalid path file (//)");\r
- }\r
-\r
-void PWPath::ToFile(TextFile &File) const\r
- {\r
- const unsigned uEdgeCount = GetEdgeCount();\r
-\r
- File.PutString("Path\n");\r
- File.PutFormat("edges %u\n", uEdgeCount);\r
- for (unsigned uEdgeIndex = 0; uEdgeIndex < uEdgeCount; ++uEdgeIndex)\r
- {\r
- const PWEdge &Edge = GetEdge(uEdgeIndex);\r
- File.PutFormat("%u %c %u %u\n",\r
- uEdgeIndex,\r
- Edge.cType,\r
- Edge.uPrefixLengthA,\r
- Edge.uPrefixLengthB);\r
- }\r
- File.PutString("//\n");\r
- }\r
-\r
-void PWPath::AssertEqual(const PWPath &Path) const\r
- {\r
- const unsigned uEdgeCount = GetEdgeCount();\r
- if (uEdgeCount != Path.GetEdgeCount())\r
- {\r
- Log("PWPath::AssertEqual, this=\n");\r
- LogMe();\r
- Log("\nOther path=\n");\r
- Path.LogMe();\r
- Log("\n");\r
- Quit("PWPath::AssertEqual, Edge count different %u %u\n",\r
- uEdgeCount, Path.GetEdgeCount());\r
- }\r
-\r
- for (unsigned uEdgeIndex = 0; uEdgeIndex < uEdgeCount; ++uEdgeIndex)\r
- {\r
- const PWEdge &e1 = GetEdge(uEdgeIndex);\r
- const PWEdge &e2 = Path.GetEdge(uEdgeIndex);\r
- if (e1.cType != e2.cType || e1.uPrefixLengthA != e2.uPrefixLengthA ||\r
- e1.uPrefixLengthB != e2.uPrefixLengthB)\r
- {\r
- Log("PWPath::AssertEqual, this=\n");\r
- LogMe();\r
- Log("\nOther path=\n");\r
- Path.LogMe();\r
- Log("\n");\r
- Log("This edge %c%u.%u, other edge %c%u.%u\n",\r
- e1.cType, e1.uPrefixLengthA, e1.uPrefixLengthB,\r
- e2.cType, e2.uPrefixLengthA, e2.uPrefixLengthB);\r
- Quit("PWPath::AssertEqual, edge %u different\n", uEdgeIndex);\r
- }\r
- }\r
- }\r
-\r
-bool PWPath::Equal(const PWPath &Path) const\r
- {\r
- const unsigned uEdgeCount = GetEdgeCount();\r
- if (uEdgeCount != Path.GetEdgeCount())\r
- return false;\r
-\r
- for (unsigned uEdgeIndex = 0; uEdgeIndex < uEdgeCount; ++uEdgeIndex)\r
- {\r
- const PWEdge &e1 = GetEdge(uEdgeIndex);\r
- const PWEdge &e2 = Path.GetEdge(uEdgeIndex);\r
- if (e1.cType != e2.cType || e1.uPrefixLengthA != e2.uPrefixLengthA ||\r
- e1.uPrefixLengthB != e2.uPrefixLengthB)\r
- return false;\r
- }\r
- return true;\r
- }\r
-\r
-unsigned PWPath::GetMatchCount() const\r
- {\r
- unsigned uMatchCount = 0;\r
- const unsigned uEdgeCount = GetEdgeCount();\r
- for (unsigned uEdgeIndex = 0; uEdgeIndex < uEdgeCount; ++uEdgeIndex)\r
- {\r
- const PWEdge &e = GetEdge(uEdgeIndex);\r
- if ('M' == e.cType)\r
- ++uMatchCount;\r
- }\r
- return uMatchCount;\r
- }\r
-\r
-unsigned PWPath::GetInsertCount() const\r
- {\r
- unsigned uInsertCount = 0;\r
- const unsigned uEdgeCount = GetEdgeCount();\r
- for (unsigned uEdgeIndex = 0; uEdgeIndex < uEdgeCount; ++uEdgeIndex)\r
- {\r
- const PWEdge &e = GetEdge(uEdgeIndex);\r
- if ('I' == e.cType)\r
- ++uInsertCount;\r
- }\r
- return uInsertCount;\r
- }\r
-\r
-unsigned PWPath::GetDeleteCount() const\r
- {\r
- unsigned uDeleteCount = 0;\r
- const unsigned uEdgeCount = GetEdgeCount();\r
- for (unsigned uEdgeIndex = 0; uEdgeIndex < uEdgeCount; ++uEdgeIndex)\r
- {\r
- const PWEdge &e = GetEdge(uEdgeIndex);\r
- if ('D' == e.cType)\r
- ++uDeleteCount;\r
- }\r
- return uDeleteCount;\r
- }\r
-\r
-void PWPath::FromStr(const char Str[])\r
- {\r
- Clear();\r
- unsigned uPrefixLengthA = 0;\r
- unsigned uPrefixLengthB = 0;\r
- while (char c = *Str++)\r
- {\r
- switch (c)\r
- {\r
- case 'M':\r
- ++uPrefixLengthA;\r
- ++uPrefixLengthB;\r
- break;\r
- case 'D':\r
- ++uPrefixLengthA;\r
- break;\r
- case 'I':\r
- ++uPrefixLengthB;\r
- break;\r
- default:\r
- Quit("PWPath::FromStr, invalid state %c", c);\r
- }\r
- AppendEdge(c, uPrefixLengthA, uPrefixLengthB);\r
- }\r
- }\r