website/archive/binaries/mac/src/muscle/aln.cpp

   1 #include "muscle.h"\r
   2 #include <stdio.h>\r
   3 #include <ctype.h>\r
   4 #include "msa.h"\r
   5 #include "textfile.h"\r
   6 \r
   7 const unsigned uCharsPerLine = 60;\r
   8 const int MIN_NAME = 10;\r
   9 const int MAX_NAME = 32;\r
  10 \r
  11 static char GetAlnConsensusChar(const MSA &a, unsigned uColIndex);\r
  12 \r
  13 void MSA::ToAlnFile(TextFile &File) const\r
  14         {\r
  15         if (g_bClwStrict)\r
  16                 File.PutString("CLUSTAL W (1.81) multiple sequence alignment\n");\r
  17         else\r
  18                 {\r
  19                 File.PutString("MUSCLE ("\r
  20                   SHORT_VERSION ")"\r
  21                   " multiple sequence alignment\n");\r
  22                 File.PutString("\n");\r
  23                 }\r
  24 \r
  25         int iLongestNameLength = 0;\r
  26         for (unsigned uSeqIndex = 0; uSeqIndex < GetSeqCount(); ++uSeqIndex)\r
  27                 {\r
  28                 const char *ptrName = GetSeqName(uSeqIndex);\r
  29                 const char *ptrBlank = strchr(ptrName, ' ');\r
  30                 int iLength;\r
  31                 if (0 != ptrBlank)\r
  32                         iLength = (int) (ptrBlank - ptrName);\r
  33                 else\r
  34                         iLength = (int) strlen(ptrName);\r
  35                 if (iLength > iLongestNameLength)\r
  36                         iLongestNameLength = iLength;\r
  37                 }\r
  38         if (iLongestNameLength > MAX_NAME)\r
  39                 iLongestNameLength = MAX_NAME;\r
  40         if (iLongestNameLength < MIN_NAME)\r
  41                 iLongestNameLength = MIN_NAME;\r
  42 \r
  43         unsigned uLineCount = (GetColCount() - 1)/uCharsPerLine + 1;\r
  44         for (unsigned uLineIndex = 0; uLineIndex < uLineCount; ++uLineIndex)\r
  45                 {\r
  46                 File.PutString("\n");\r
  47                 unsigned uStartColIndex = uLineIndex*uCharsPerLine;\r
  48                 unsigned uEndColIndex = uStartColIndex + uCharsPerLine - 1;\r
  49                 if (uEndColIndex >= GetColCount())\r
  50                         uEndColIndex = GetColCount() - 1;\r
  51                 char Name[MAX_NAME+1];\r
  52                 for (unsigned uSeqIndex = 0; uSeqIndex < GetSeqCount(); ++uSeqIndex)\r
  53                         {\r
  54                         const char *ptrName = GetSeqName(uSeqIndex);\r
  55                         const char *ptrBlank = strchr(ptrName, ' ');\r
  56                         int iLength;\r
  57                         if (0 != ptrBlank)\r
  58                                 iLength = (int) (ptrBlank - ptrName);\r
  59                         else\r
  60                                 iLength = (int) strlen(ptrName);\r
  61                         if (iLength > MAX_NAME)\r
  62                                 iLength = MAX_NAME;\r
  63                         memset(Name, ' ', MAX_NAME);\r
  64                         memcpy(Name, ptrName, iLength);\r
  65                         Name[iLongestNameLength] = 0;\r
  66 \r
  67                         File.PutFormat("%s      ", Name);\r
  68                         for (unsigned uColIndex = uStartColIndex; uColIndex <= uEndColIndex;\r
  69                           ++uColIndex)\r
  70                                 {\r
  71                                 const char c = GetChar(uSeqIndex, uColIndex);\r
  72                                 File.PutFormat("%c", toupper(c));\r
  73                                 }\r
  74                         File.PutString("\n");\r
  75                         }\r
  76 \r
  77                 memset(Name, ' ', MAX_NAME);\r
  78                 Name[iLongestNameLength] = 0;\r
  79                 File.PutFormat("%s      ", Name);\r
  80                 for (unsigned uColIndex = uStartColIndex; uColIndex <= uEndColIndex;\r
  81                   ++uColIndex)\r
  82                         {\r
  83                         const char c = GetAlnConsensusChar(*this, uColIndex);\r
  84                         File.PutChar(c);\r
  85                         }\r
  86                 File.PutString("\n");\r
  87                 }\r
  88         }\r
  89 \r
  90 static char GetAlnConsensusChar(const MSA &a, unsigned uColIndex)\r
  91         {\r
  92         const unsigned uSeqCount = a.GetSeqCount();\r
  93         unsigned BitMap = 0;\r
  94         unsigned Count = 0;\r
  95         for (unsigned uSeqIndex = 0; uSeqIndex < uSeqCount; ++uSeqIndex)\r
  96                 {\r
  97                 unsigned uLetter = a.GetLetterEx(uSeqIndex, uColIndex);\r
  98                 assert(uLetter < 32);\r
  99                 unsigned Bit = (1 << uLetter);\r
 100                 if (!(BitMap & Bit))\r
 101                         ++Count;\r
 102                 BitMap |= Bit;\r
 103                 }\r
 104 \r
 105 //      '*' indicates positions which have a single, fully conserved residue\r
 106         if (1 == Count)\r
 107                 return '*';\r
 108 \r
 109         if (ALPHA_Amino != g_Alpha)\r
 110                 return ' ';\r
 111 \r
 112 #define B(a)    (1 << AX_##a)\r
 113 #define S2(a, b)                S(B(a) | B(b))\r
 114 #define S3(a, b, c)             S(B(a) | B(b) | B(c))\r
 115 #define S4(a, b, c, d)  S(B(a) | B(b) | B(c) | B(d))\r
 116 #define S(w)    if (0 == (BitMap & ~(w)) && (BitMap & (w)) != 0) return ':';\r
 117 \r
 118 #define W3(a, b, c)                             W(B(a) | B(b) | B(c))\r
 119 #define W4(a, b, c, d)                  W(B(a) | B(b) | B(c) | B(d))\r
 120 #define W5(a, b, c, d, e)               W(B(a) | B(b) | B(c) | B(d) | B(e))\r
 121 #define W6(a, b, c, d, e, f)    W(B(a) | B(b) | B(c) | B(d) | B(e) | B(f))\r
 122 #define W(w)    if (0 == (BitMap & ~(w)) && (BitMap & (w)) != 0) return '.';\r
 123 \r
 124 //      ':' indicates that one of the following 'strong'\r
 125 // groups is fully conserved\r
 126 //                 STA\r
 127 //                 NEQK\r
 128 //                 NHQK\r
 129 //                 NDEQ\r
 130 //                 QHRK\r
 131 //                 MILV\r
 132 //                 MILF\r
 133 //                 HY\r
 134 //                 FYW\r
 135 //\r
 136         S3(S, T, A)\r
 137         S4(N, E, Q, K)\r
 138         S4(N, H, Q, K)\r
 139         S4(N, D, E, Q)\r
 140         S4(M, I, L, V)\r
 141         S4(M, I, L, F)\r
 142         S2(H, Y)\r
 143         S3(F, Y, W)\r
 144 \r
 145 //      '.' indicates that one of the following 'weaker' \r
 146 // groups is fully conserved\r
 147 //                 CSA\r
 148 //                 ATV\r
 149 //                 SAG\r
 150 //                 STNK\r
 151 //                 STPA\r
 152 //                 SGND\r
 153 //                 SNDEQK\r
 154 //                 NDEQHK\r
 155 //                 NEQHRK\r
 156 //                 FVLIM\r
 157 //                 HFY\r
 158         W3(C, S, A)\r
 159         W3(A, T, V)\r
 160         W3(S, A, G)\r
 161         W4(S, T, N, K)\r
 162         W4(S, T, P, A)\r
 163         W4(S, G, N, D)\r
 164         W6(S, N, D, E, Q, K)\r
 165         W6(N, W, Q, H, R, K)\r
 166         W5(F, V, L, I, M)\r
 167         W3(H, F, Y)\r
 168 \r
 169         return ' ';\r
 170         }\r