]> git.donarmstrong.com Git - mothur.git/blobdiff - uchime_src/mx.cpp
added uchime_src folder. added biom parameter to make.shared. added biom as a current...
[mothur.git] / uchime_src / mx.cpp
diff --git a/uchime_src/mx.cpp b/uchime_src/mx.cpp
new file mode 100644 (file)
index 0000000..48c347e
--- /dev/null
@@ -0,0 +1,294 @@
+#include "myutils.h"\r
+#include "mx.h"\r
+#include "seqdb.h"\r
+#include "seq.h"\r
+\r
+char ProbToChar(float p);\r
+\r
+list<MxBase *> *MxBase::m_Matrices = 0;\r
+unsigned MxBase::m_AllocCount;\r
+unsigned MxBase::m_ZeroAllocCount;\r
+unsigned MxBase::m_GrowAllocCount;\r
+double MxBase::m_TotalBytes;\r
+double MxBase::m_MaxBytes;\r
+\r
+static const char *LogizeStr(const char *s)\r
+       {\r
+       double d = atof(s);\r
+       d = log(d);\r
+       return TypeToStr<float>(float(d));\r
+       }\r
+\r
+static const char *ExpizeStr(const char *s)\r
+       {\r
+       double d = atof(s);\r
+       d = exp(d);\r
+       return TypeToStr<float>(float(d));\r
+       }\r
+\r
+void MxBase::OnCtor(MxBase *Mx)\r
+       {\r
+       if (m_Matrices == 0)\r
+               m_Matrices = new list<MxBase *>;\r
+       asserta(m_Matrices != 0);\r
+       m_Matrices->push_front(Mx);\r
+       }\r
+\r
+void MxBase::OnDtor(MxBase *Mx)\r
+       {\r
+       if (m_Matrices == 0)\r
+               {\r
+               Warning("MxBase::OnDtor, m_Matrices = 0");\r
+               return;\r
+               }\r
+       for (list<MxBase*>::iterator p = m_Matrices->begin();\r
+         p != m_Matrices->end(); ++p)\r
+               {\r
+               if (*p == Mx)\r
+                       {\r
+                       m_Matrices->erase(p);\r
+                       if (m_Matrices->empty())\r
+                               delete m_Matrices;\r
+                       return;\r
+                       }\r
+               }\r
+       Warning("MxBase::OnDtor, not found");\r
+       }\r
+\r
+//float **MxBase::Getf(const string &Name)\r
+//     {\r
+//     Mx<float> *m = (Mx<float> *) Get(Name);\r
+//     asserta(m->GetTypeSize() == sizeof(float));\r
+//     return m->GetData();\r
+//     }\r
+//\r
+//double **MxBase::Getd(const string &Name)\r
+//     {\r
+//     Mx<double> *m = (Mx<double> *) Get(Name);\r
+//     asserta(m->GetTypeSize() == sizeof(double));\r
+//     return m->GetData();\r
+//     }\r
+//\r
+//char **MxBase::Getc(const string &Name)\r
+//     {\r
+//     Mx<char> *m = (Mx<char> *) Get(Name);\r
+//     asserta(m->GetTypeSize() == sizeof(char));\r
+//     return m->GetData();\r
+//     }\r
+\r
+void MxBase::Alloc(const char *Name, unsigned RowCount, unsigned ColCount,\r
+  const SeqDB *DB, unsigned IdA, unsigned IdB)\r
+       {\r
+       Alloc(Name, RowCount, ColCount, DB, IdA, IdB, 0, 0);\r
+       }\r
+\r
+void MxBase::Alloc(const char *Name, unsigned RowCount, unsigned ColCount,\r
+  const SeqData *SA, const SeqData *SB)\r
+       {\r
+       Alloc(Name, RowCount, ColCount, 0, UINT_MAX, UINT_MAX, SA, SB);\r
+       }\r
+\r
+void MxBase::Alloc(const char *Name, unsigned RowCount, unsigned ColCount,\r
+  const SeqDB *DB, unsigned IdA, unsigned IdB, const SeqData *SA, const SeqData *SB)\r
+       {\r
+       StartTimer(MxBase_Alloc);\r
+\r
+       ++m_AllocCount;\r
+       if (m_AllocatedRowCount == 0)\r
+               ++m_ZeroAllocCount;\r
+\r
+       if (DB != 0)\r
+               {\r
+               asserta(IdA != UINT_MAX);\r
+               asserta(IdB != UINT_MAX);\r
+               asserta(RowCount >= DB->GetSeqLength(IdA) + 1);\r
+               asserta(ColCount >= DB->GetSeqLength(IdB) + 1);\r
+               }\r
+       if (RowCount > m_AllocatedRowCount || ColCount > m_AllocatedColCount)\r
+               {\r
+               if (m_AllocatedRowCount > 0)\r
+                       {\r
+                       if (opt_logmemgrows)\r
+                               Log("MxBase::Alloc grow %s %u x %u -> %u x %u, %s bytes\n",\r
+                                 Name, m_AllocatedRowCount, m_AllocatedColCount,\r
+                                 RowCount, ColCount,\r
+                                 IntToStr(GetBytes()));\r
+                       ++m_GrowAllocCount;\r
+                       }\r
+\r
+               m_TotalBytes -= GetBytes();\r
+\r
+               PauseTimer(MxBase_Alloc);\r
+               StartTimer(MxBase_FreeData);\r
+               FreeData();\r
+               EndTimer(MxBase_FreeData);\r
+               StartTimer(MxBase_Alloc);\r
+\r
+               unsigned N = max(RowCount + 16, m_AllocatedRowCount);\r
+               unsigned M = max(ColCount + 16, m_AllocatedColCount);\r
+               N = max(N, M);\r
+\r
+               PauseTimer(MxBase_Alloc);\r
+               StartTimer(MxBase_AllocData);\r
+               AllocData(N, N);\r
+               EndTimer(MxBase_AllocData);\r
+               StartTimer(MxBase_Alloc);\r
+\r
+               m_TotalBytes += GetBytes();\r
+               if (m_TotalBytes > m_MaxBytes)\r
+                       m_MaxBytes = m_TotalBytes;\r
+               }\r
+       \r
+       unsigned n = sizeof(m_Name)-1;\r
+       strncpy(m_Name, Name, n);\r
+       m_Name[n] = 0;\r
+       m_RowCount = RowCount;\r
+       m_ColCount = ColCount;\r
+       m_SeqDB = DB;\r
+       m_IdA = IdA;\r
+       m_IdB = IdB;\r
+       m_SA = SA;\r
+       m_SB = SB;\r
+\r
+       EndTimer(MxBase_Alloc);\r
+       }\r
+\r
+void MxBase::LogMe(bool WithData, int Opts) const\r
+       {\r
+       Log("\n");\r
+       if (Opts & OPT_EXP)\r
+               Log("Exp ");\r
+       else if (Opts & OPT_LOG)\r
+               Log("Log ");\r
+       bool ZeroBased = ((Opts & OPT_ZERO_BASED) != 0);\r
+       Log("%s(%p) Rows %u/%u, Cols %u/%u",\r
+         m_Name, this,\r
+         m_RowCount, m_AllocatedRowCount,\r
+         m_ColCount, m_AllocatedColCount);\r
+       if (m_SeqDB != 0 && m_IdA != UINT_MAX)\r
+               Log(", A=%s", m_SeqDB->GetLabel(m_IdA));\r
+       else if (m_SA != 0)\r
+               Log(", A=%s", m_SA->Label);\r
+       if (m_SeqDB != 0 && m_IdB != UINT_MAX)\r
+               Log(", B=%s", m_SeqDB->GetLabel(m_IdB));\r
+       else if (m_SB != 0)\r
+               Log(", B=%s", m_SB->Label);\r
+       Log("\n");\r
+       if (!WithData || m_RowCount == 0 || m_ColCount == 0)\r
+               return;\r
+\r
+       const char *z = GetAsStr(0, 0);\r
+       unsigned Width = strlen(z);\r
+       unsigned Mod = 1;\r
+       for (unsigned i = 0; i < Width; ++i)\r
+               Mod *= 10;\r
+\r
+       if (m_Alpha[0] != 0)\r
+               {\r
+               Log("// Alphabet=%s\n", m_Alpha);\r
+               Log("//      ");\r
+               unsigned n = strlen(m_Alpha);\r
+               for (unsigned j = 0; j < n; ++j)\r
+                       Log(" %*c", Width, m_Alpha[j]);\r
+               Log("\n");\r
+               for (unsigned i = 0; i < n; ++i)\r
+                       {\r
+                       Log("/* %c */ {", m_Alpha[i]);\r
+                       unsigned ci = m_Alpha[i];\r
+                       for (unsigned j = 0; j < n; ++j)\r
+                               {\r
+                               unsigned cj = m_Alpha[j];\r
+                               Log("%s,", GetAsStr(ci, cj));\r
+                               }\r
+                       Log("},  // %c\n", m_Alpha[i]);\r
+                       }\r
+               return;\r
+               }\r
+       else if (m_Alpha2[0] != 0)\r
+               {\r
+               unsigned n = strlen(m_Alpha2);\r
+               Log("// Alphabet=%s\n", m_Alpha2);\r
+               Log("//      ");\r
+               for (unsigned j = 0; j < n; ++j)\r
+                       Log(" %*c", Width, m_Alpha2[j]);\r
+               Log("\n");\r
+               for (unsigned i = 0; i < n; ++i)\r
+                       {\r
+                       Log("/* %c */ {", m_Alpha2[i]);\r
+                       unsigned ci = m_Alpha2[i];\r
+                       for (unsigned j = 0; j < n; ++j)\r
+                               Log("%s,", GetAsStr(i, j));\r
+                       Log("},  // %c\n", m_Alpha2[i]);\r
+                       }\r
+               return;\r
+               }\r
+\r
+       const byte *A = 0;\r
+       const byte *B = 0;\r
+       if (m_SeqDB != 0 && m_IdA != UINT_MAX)\r
+               A = m_SeqDB->GetSeq(m_IdA);\r
+       else if (m_SA != 0)\r
+               A = m_SA->Seq;\r
+       if (m_SeqDB != 0 && m_IdB != UINT_MAX)\r
+               B = m_SeqDB->GetSeq(m_IdB);\r
+       else if (m_SB != 0)\r
+               B = m_SB->Seq;\r
+\r
+       if (B != 0)\r
+               {\r
+               if (A != 0)\r
+                       Log("  ");\r
+               Log("%5.5s", "");\r
+               if (ZeroBased)\r
+                       for (unsigned j = 0; j < m_ColCount; ++j)\r
+                               Log("%*c", Width, B[j]);\r
+               else\r
+                       for (unsigned j = 0; j < m_ColCount; ++j)\r
+                               Log("%*c", Width, j == 0 ? ' ' : B[j-1]);\r
+               Log("\n");\r
+               }\r
+\r
+       if (A != 0)\r
+               Log("  ");\r
+       Log("%5.5s", "");\r
+       for (unsigned j = 0; j < m_ColCount; ++j)\r
+               Log("%*u", Width, j%Mod);\r
+       Log("\n");\r
+\r
+       for (unsigned i = 0; i < m_RowCount; ++i)\r
+               {\r
+               if (A != 0)\r
+                       {\r
+                       if (ZeroBased)\r
+                               Log("%c ", A[i]);\r
+                       else\r
+                               Log("%c ", i == 0 ? ' ' : A[i-1]);\r
+                       }\r
+               Log("%4u ", i);\r
+               \r
+               for (unsigned j = 0; j < m_ColCount; ++j)\r
+                       {\r
+                       const char *s = GetAsStr(i, j);\r
+                       if (Opts & OPT_LOG)\r
+                               s = LogizeStr(s);\r
+                       else if (Opts & OPT_EXP)\r
+                               s = ExpizeStr(s);\r
+                       Log("%s", s);\r
+                       }\r
+               Log("\n");\r
+               }\r
+       }\r
+static unsigned g_MatrixFileCount;\r
+\r
+void MxBase::LogCounts()\r
+       {\r
+       Log("\n");\r
+       Log("MxBase::LogCounts()\n");\r
+       Log("      What           N\n");\r
+       Log("----------  ----------\n");\r
+       Log("    Allocs  %10u\n", m_AllocCount);\r
+       Log("ZeroAllocs  %10u\n", m_ZeroAllocCount);\r
+       Log("     Grows  %10u\n", m_GrowAllocCount);\r
+       Log("     Bytes  %10.10s\n", MemBytesToStr(m_TotalBytes));\r
+       Log(" Max bytes  %10.10s\n", MemBytesToStr(m_MaxBytes));\r
+       }\r