uchime_src/mx.cpp

   1 #include "myutils.h"\r
   2 #include "mx.h"\r
   3 #include "seqdb.h"\r
   4 #include "seq.h"\r
   5 \r
   6 char ProbToChar(float p);\r
   7 \r
   8 list<MxBase *> *MxBase::m_Matrices = 0;\r
   9 unsigned MxBase::m_AllocCount;\r
  10 unsigned MxBase::m_ZeroAllocCount;\r
  11 unsigned MxBase::m_GrowAllocCount;\r
  12 double MxBase::m_TotalBytes;\r
  13 double MxBase::m_MaxBytes;\r
  14 \r
  15 static const char *LogizeStr(const char *s)\r
  16         {\r
  17         double d = atof(s);\r
  18         d = log(d);\r
  19         return TypeToStr<float>(float(d));\r
  20         }\r
  21 \r
  22 static const char *ExpizeStr(const char *s)\r
  23         {\r
  24         double d = atof(s);\r
  25         d = exp(d);\r
  26         return TypeToStr<float>(float(d));\r
  27         }\r
  28 \r
  29 void MxBase::OnCtor(MxBase *Mx)\r
  30         {\r
  31         if (m_Matrices == 0)\r
  32                 m_Matrices = new list<MxBase *>;\r
  33         asserta(m_Matrices != 0);\r
  34         m_Matrices->push_front(Mx);\r
  35         }\r
  36 \r
  37 void MxBase::OnDtor(MxBase *Mx)\r
  38         {\r
  39         if (m_Matrices == 0)\r
  40                 {\r
  41                 Warning("MxBase::OnDtor, m_Matrices = 0");\r
  42                 return;\r
  43                 }\r
  44         for (list<MxBase*>::iterator p = m_Matrices->begin();\r
  45           p != m_Matrices->end(); ++p)\r
  46                 {\r
  47                 if (*p == Mx)\r
  48                         {\r
  49                         m_Matrices->erase(p);\r
  50                         if (m_Matrices->empty())\r
  51                                 delete m_Matrices;\r
  52                         return;\r
  53                         }\r
  54                 }\r
  55         Warning("MxBase::OnDtor, not found");\r
  56         }\r
  57 \r
  58 //float **MxBase::Getf(const string &Name)\r
  59 //      {\r
  60 //      Mx<float> *m = (Mx<float> *) Get(Name);\r
  61 //      asserta(m->GetTypeSize() == sizeof(float));\r
  62 //      return m->GetData();\r
  63 //      }\r
  64 //\r
  65 //double **MxBase::Getd(const string &Name)\r
  66 //      {\r
  67 //      Mx<double> *m = (Mx<double> *) Get(Name);\r
  68 //      asserta(m->GetTypeSize() == sizeof(double));\r
  69 //      return m->GetData();\r
  70 //      }\r
  71 //\r
  72 //char **MxBase::Getc(const string &Name)\r
  73 //      {\r
  74 //      Mx<char> *m = (Mx<char> *) Get(Name);\r
  75 //      asserta(m->GetTypeSize() == sizeof(char));\r
  76 //      return m->GetData();\r
  77 //      }\r
  78 \r
  79 void MxBase::Alloc(const char *Name, unsigned RowCount, unsigned ColCount,\r
  80   const SeqDB *DB, unsigned IdA, unsigned IdB)\r
  81         {\r
  82         Alloc(Name, RowCount, ColCount, DB, IdA, IdB, 0, 0);\r
  83         }\r
  84 \r
  85 void MxBase::Alloc(const char *Name, unsigned RowCount, unsigned ColCount,\r
  86   const SeqData *SA, const SeqData *SB)\r
  87         {\r
  88         Alloc(Name, RowCount, ColCount, 0, UINT_MAX, UINT_MAX, SA, SB);\r
  89         }\r
  90 \r
  91 void MxBase::Alloc(const char *Name, unsigned RowCount, unsigned ColCount,\r
  92   const SeqDB *DB, unsigned IdA, unsigned IdB, const SeqData *SA, const SeqData *SB)\r
  93         {\r
  94         StartTimer(MxBase_Alloc);\r
  95 \r
  96         ++m_AllocCount;\r
  97         if (m_AllocatedRowCount == 0)\r
  98                 ++m_ZeroAllocCount;\r
  99 \r
 100         if (DB != 0)\r
 101                 {\r
 102                 asserta(IdA != UINT_MAX);\r
 103                 asserta(IdB != UINT_MAX);\r
 104                 asserta(RowCount >= DB->GetSeqLength(IdA) + 1);\r
 105                 asserta(ColCount >= DB->GetSeqLength(IdB) + 1);\r
 106                 }\r
 107         if (RowCount > m_AllocatedRowCount || ColCount > m_AllocatedColCount)\r
 108                 {\r
 109                 if (m_AllocatedRowCount > 0)\r
 110                         {\r
 111                         if (opt_logmemgrows)\r
 112                                 Log("MxBase::Alloc grow %s %u x %u -> %u x %u, %s bytes\n",\r
 113                                   Name, m_AllocatedRowCount, m_AllocatedColCount,\r
 114                                   RowCount, ColCount,\r
 115                                   IntToStr(GetBytes()));\r
 116                         ++m_GrowAllocCount;\r
 117                         }\r
 118 \r
 119                 m_TotalBytes -= GetBytes();\r
 120 \r
 121                 PauseTimer(MxBase_Alloc);\r
 122                 StartTimer(MxBase_FreeData);\r
 123                 FreeData();\r
 124                 EndTimer(MxBase_FreeData);\r
 125                 StartTimer(MxBase_Alloc);\r
 126 \r
 127                 unsigned N = max(RowCount + 16, m_AllocatedRowCount);\r
 128                 unsigned M = max(ColCount + 16, m_AllocatedColCount);\r
 129                 N = max(N, M);\r
 130 \r
 131                 PauseTimer(MxBase_Alloc);\r
 132                 StartTimer(MxBase_AllocData);\r
 133                 AllocData(N, N);\r
 134                 EndTimer(MxBase_AllocData);\r
 135                 StartTimer(MxBase_Alloc);\r
 136 \r
 137                 m_TotalBytes += GetBytes();\r
 138                 if (m_TotalBytes > m_MaxBytes)\r
 139                         m_MaxBytes = m_TotalBytes;\r
 140                 }\r
 141         \r
 142         unsigned n = sizeof(m_Name)-1;\r
 143         strncpy(m_Name, Name, n);\r
 144         m_Name[n] = 0;\r
 145         m_RowCount = RowCount;\r
 146         m_ColCount = ColCount;\r
 147         m_SeqDB = DB;\r
 148         m_IdA = IdA;\r
 149         m_IdB = IdB;\r
 150         m_SA = SA;\r
 151         m_SB = SB;\r
 152 \r
 153         EndTimer(MxBase_Alloc);\r
 154         }\r
 155 \r
 156 void MxBase::LogMe(bool WithData, int Opts) const\r
 157         {\r
 158         Log("\n");\r
 159         if (Opts & OPT_EXP)\r
 160                 Log("Exp ");\r
 161         else if (Opts & OPT_LOG)\r
 162                 Log("Log ");\r
 163         bool ZeroBased = ((Opts & OPT_ZERO_BASED) != 0);\r
 164         Log("%s(%p) Rows %u/%u, Cols %u/%u",\r
 165           m_Name, this,\r
 166           m_RowCount, m_AllocatedRowCount,\r
 167           m_ColCount, m_AllocatedColCount);\r
 168         if (m_SeqDB != 0 && m_IdA != UINT_MAX)\r
 169                 Log(", A=%s", m_SeqDB->GetLabel(m_IdA));\r
 170         else if (m_SA != 0)\r
 171                 Log(", A=%s", m_SA->Label);\r
 172         if (m_SeqDB != 0 && m_IdB != UINT_MAX)\r
 173                 Log(", B=%s", m_SeqDB->GetLabel(m_IdB));\r
 174         else if (m_SB != 0)\r
 175                 Log(", B=%s", m_SB->Label);\r
 176         Log("\n");\r
 177         if (!WithData || m_RowCount == 0 || m_ColCount == 0)\r
 178                 return;\r
 179 \r
 180         const char *z = GetAsStr(0, 0);\r
 181         unsigned Width = strlen(z);\r
 182         unsigned Mod = 1;\r
 183         for (unsigned i = 0; i < Width; ++i)\r
 184                 Mod *= 10;\r
 185 \r
 186         if (m_Alpha[0] != 0)\r
 187                 {\r
 188                 Log("// Alphabet=%s\n", m_Alpha);\r
 189                 Log("//      ");\r
 190                 unsigned n = strlen(m_Alpha);\r
 191                 for (unsigned j = 0; j < n; ++j)\r
 192                         Log(" %*c", Width, m_Alpha[j]);\r
 193                 Log("\n");\r
 194                 for (unsigned i = 0; i < n; ++i)\r
 195                         {\r
 196                         Log("/* %c */ {", m_Alpha[i]);\r
 197                         unsigned ci = m_Alpha[i];\r
 198                         for (unsigned j = 0; j < n; ++j)\r
 199                                 {\r
 200                                 unsigned cj = m_Alpha[j];\r
 201                                 Log("%s,", GetAsStr(ci, cj));\r
 202                                 }\r
 203                         Log("},  // %c\n", m_Alpha[i]);\r
 204                         }\r
 205                 return;\r
 206                 }\r
 207         else if (m_Alpha2[0] != 0)\r
 208                 {\r
 209                 unsigned n = strlen(m_Alpha2);\r
 210                 Log("// Alphabet=%s\n", m_Alpha2);\r
 211                 Log("//      ");\r
 212                 for (unsigned j = 0; j < n; ++j)\r
 213                         Log(" %*c", Width, m_Alpha2[j]);\r
 214                 Log("\n");\r
 215                 for (unsigned i = 0; i < n; ++i)\r
 216                         {\r
 217                         Log("/* %c */ {", m_Alpha2[i]);\r
 218                         unsigned ci = m_Alpha2[i];\r
 219                         for (unsigned j = 0; j < n; ++j)\r
 220                                 Log("%s,", GetAsStr(i, j));\r
 221                         Log("},  // %c\n", m_Alpha2[i]);\r
 222                         }\r
 223                 return;\r
 224                 }\r
 225 \r
 226         const byte *A = 0;\r
 227         const byte *B = 0;\r
 228         if (m_SeqDB != 0 && m_IdA != UINT_MAX)\r
 229                 A = m_SeqDB->GetSeq(m_IdA);\r
 230         else if (m_SA != 0)\r
 231                 A = m_SA->Seq;\r
 232         if (m_SeqDB != 0 && m_IdB != UINT_MAX)\r
 233                 B = m_SeqDB->GetSeq(m_IdB);\r
 234         else if (m_SB != 0)\r
 235                 B = m_SB->Seq;\r
 236 \r
 237         if (B != 0)\r
 238                 {\r
 239                 if (A != 0)\r
 240                         Log("  ");\r
 241                 Log("%5.5s", "");\r
 242                 if (ZeroBased)\r
 243                         for (unsigned j = 0; j < m_ColCount; ++j)\r
 244                                 Log("%*c", Width, B[j]);\r
 245                 else\r
 246                         for (unsigned j = 0; j < m_ColCount; ++j)\r
 247                                 Log("%*c", Width, j == 0 ? ' ' : B[j-1]);\r
 248                 Log("\n");\r
 249                 }\r
 250 \r
 251         if (A != 0)\r
 252                 Log("  ");\r
 253         Log("%5.5s", "");\r
 254         for (unsigned j = 0; j < m_ColCount; ++j)\r
 255                 Log("%*u", Width, j%Mod);\r
 256         Log("\n");\r
 257 \r
 258         for (unsigned i = 0; i < m_RowCount; ++i)\r
 259                 {\r
 260                 if (A != 0)\r
 261                         {\r
 262                         if (ZeroBased)\r
 263                                 Log("%c ", A[i]);\r
 264                         else\r
 265                                 Log("%c ", i == 0 ? ' ' : A[i-1]);\r
 266                         }\r
 267                 Log("%4u ", i);\r
 268                 \r
 269                 for (unsigned j = 0; j < m_ColCount; ++j)\r
 270                         {\r
 271                         const char *s = GetAsStr(i, j);\r
 272                         if (Opts & OPT_LOG)\r
 273                                 s = LogizeStr(s);\r
 274                         else if (Opts & OPT_EXP)\r
 275                                 s = ExpizeStr(s);\r
 276                         Log("%s", s);\r
 277                         }\r
 278                 Log("\n");\r
 279                 }\r
 280         }\r
 281 static unsigned g_MatrixFileCount;\r
 282 \r
 283 void MxBase::LogCounts()\r
 284         {\r
 285         Log("\n");\r
 286         Log("MxBase::LogCounts()\n");\r
 287         Log("      What           N\n");\r
 288         Log("----------  ----------\n");\r
 289         Log("    Allocs  %10u\n", m_AllocCount);\r
 290         Log("ZeroAllocs  %10u\n", m_ZeroAllocCount);\r
 291         Log("     Grows  %10u\n", m_GrowAllocCount);\r
 292         Log("     Bytes  %10.10s\n", MemBytesToStr(m_TotalBytes));\r
 293         Log(" Max bytes  %10.10s\n", MemBytesToStr(m_MaxBytes));\r
 294         }\r