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

   1 #include "muscle.h"\r
   2 #include "msa.h"\r
   3 #include <errno.h>\r
   4 \r
   5 extern float VTML_SP[32][32];\r
   6 extern float NUC_SP[32][32];\r
   7 \r
   8 static double GetColScore(const MSA &msa, unsigned uCol)\r
   9         {\r
  10         const unsigned uSeqCount = msa.GetSeqCount();\r
  11         unsigned uPairCount = 0;\r
  12         double dSum = 0.0;\r
  13         for (unsigned uSeq1 = 0; uSeq1 < uSeqCount; ++uSeq1)\r
  14                 {\r
  15                 if (msa.IsGap(uSeq1, uCol))\r
  16                         continue;\r
  17                 unsigned uLetter1 = msa.GetLetterEx(uSeq1, uCol);\r
  18                 if (uLetter1 >= g_AlphaSize)\r
  19                         continue;\r
  20                 for (unsigned uSeq2 = uSeq1 + 1; uSeq2 < uSeqCount; ++uSeq2)\r
  21                         {\r
  22                         if (msa.IsGap(uSeq2, uCol))\r
  23                                 continue;\r
  24                         unsigned uLetter2 = msa.GetLetterEx(uSeq2, uCol);\r
  25                         if (uLetter2 >= g_AlphaSize)\r
  26                                 continue;\r
  27                         double Score;\r
  28                         switch (g_Alpha)\r
  29                                 {\r
  30                         case ALPHA_Amino:\r
  31                                 Score = VTML_SP[uLetter1][uLetter2];\r
  32                                 break;\r
  33                         case ALPHA_DNA:\r
  34                         case ALPHA_RNA:\r
  35                                 Score = NUC_SP[uLetter1][uLetter2];\r
  36                                 break;\r
  37                         default:\r
  38                                 Quit("GetColScore: invalid alpha=%d", g_Alpha);\r
  39                                 }\r
  40                         dSum += Score;\r
  41                         ++uPairCount;\r
  42                         }\r
  43                 }\r
  44         if (0 == uPairCount)\r
  45                 return 0;\r
  46         return dSum / uPairCount;\r
  47         }\r
  48 \r
  49 void WriteScoreFile(const MSA &msa)\r
  50         {\r
  51         FILE *f = fopen(g_pstrScoreFileName, "w");\r
  52         if (0 == f)\r
  53                 Quit("Cannot open score file '%s' errno=%d", g_pstrScoreFileName, errno);\r
  54 \r
  55         const unsigned uColCount = msa.GetColCount();\r
  56         const unsigned uSeqCount = msa.GetSeqCount();\r
  57         for (unsigned uCol = 0; uCol < uColCount; ++uCol)\r
  58                 {\r
  59                 double Score = GetColScore(msa, uCol);\r
  60                 fprintf(f, "%10.3f  ", Score);\r
  61                 for (unsigned uSeq = 0; uSeq < uSeqCount; ++uSeq)\r
  62                         {\r
  63                         char c = msa.GetChar(uSeq, uCol);\r
  64                         fprintf(f, "%c", c);\r
  65                         }\r
  66                 fprintf(f, "\n");\r
  67                 }\r
  68         fclose(f);\r
  69         }\r