2 #include "objscore.h"
\r
4 #include "textfile.h"
\r
7 const unsigned INDELS = 1;
\r
9 static void GetPos(const char Str[], unsigned L, int *pi1, int *pi2)
\r
14 i1 = rand()%(L-2) + 1;
\r
21 i2 = rand()%(L-2) + 1;
\r
22 if (i1 != i2 && Str[i2] == 'M')
\r
29 static void MakePath(unsigned uSeqLength, unsigned uIndelCount, char Str[])
\r
31 unsigned uPathLength = uSeqLength + uIndelCount;
\r
32 for (unsigned i = 0; i < uPathLength; ++i)
\r
35 for (unsigned i = 0; i < uIndelCount; ++i)
\r
38 GetPos(Str, uPathLength, &i1, &i2);
\r
43 Str[uPathLength] = 0;
\r
44 Log("MakePath=%s\n", Str);
\r
49 SetPPScore(PPSCORE_SV);
\r
51 SetListFileName("c:\\tmp\\muscle.log", false);
\r
53 TextFile file1("c:\\tmp\\msa1.afa");
\r
54 TextFile file2("c:\\tmp\\msa2.afa");
\r
59 msa1.FromFile(file1);
\r
60 msa2.FromFile(file2);
\r
67 const unsigned uColCount = msa1.GetColCount();
\r
68 if (msa2.GetColCount() != uColCount)
\r
69 Quit("Different lengths");
\r
71 const unsigned uSeqCount1 = msa1.GetSeqCount();
\r
72 const unsigned uSeqCount2 = msa2.GetSeqCount();
\r
73 const unsigned uSeqCount = uSeqCount1 + uSeqCount2;
\r
75 MSA::SetIdCount(uSeqCount);
\r
77 for (unsigned uSeqIndex1 = 0; uSeqIndex1 < uSeqCount1; ++uSeqIndex1)
\r
79 msa1.SetSeqWeight(uSeqIndex1, 1.0);
\r
80 msa1.SetSeqId(uSeqIndex1, uSeqIndex1);
\r
83 for (unsigned uSeqIndex2 = 0; uSeqIndex2 < uSeqCount2; ++uSeqIndex2)
\r
85 msa2.SetSeqWeight(uSeqIndex2, 1.0);
\r
86 msa2.SetSeqId(uSeqIndex2, uSeqCount1 + uSeqIndex2);
\r
92 char strPathA[1024];
\r
93 char strPathB[1024];
\r
94 MakePath(uColCount, INDELS, strPathA);
\r
95 MakePath(uColCount, INDELS, strPathB);
\r
99 PathA.FromStr(strPathA);
\r
100 PathB.FromStr(strPathB);
\r
107 AlignTwoMSAsGivenPath(PathA, msa1, msa2, alnA);
\r
108 AlignTwoMSAsGivenPath(PathB, msa1, msa2, alnB);
\r
110 for (unsigned uSeqIndex = 0; uSeqIndex < uSeqCount; ++uSeqIndex)
\r
112 alnA.SetSeqWeight(uSeqIndex, 1.0);
\r
113 alnB.SetSeqWeight(uSeqIndex, 1.0);
\r
116 unsigned Seqs1[1024];
\r
117 unsigned Seqs2[1024];
\r
119 for (unsigned uSeqIndex1 = 0; uSeqIndex1 < uSeqCount1; ++uSeqIndex1)
\r
120 Seqs1[uSeqIndex1] = uSeqIndex1;
\r
122 for (unsigned uSeqIndex2 = 0; uSeqIndex2 < uSeqCount2; ++uSeqIndex2)
\r
123 Seqs2[uSeqIndex2] = uSeqCount1 + uSeqIndex2;
\r
129 MSAFromSeqSubset(alnA, Seqs1, uSeqCount1, msaA1);
\r
130 MSAFromSeqSubset(alnB, Seqs1, uSeqCount1, msaB1);
\r
131 MSAFromSeqSubset(alnA, Seqs2, uSeqCount2, msaA2);
\r
132 MSAFromSeqSubset(alnB, Seqs2, uSeqCount2, msaB2);
\r
134 for (unsigned uSeqIndex1 = 0; uSeqIndex1 < uSeqCount1; ++uSeqIndex1)
\r
136 msaA1.SetSeqWeight(uSeqIndex1, 1.0);
\r
137 msaB1.SetSeqWeight(uSeqIndex1, 1.0);
\r
140 for (unsigned uSeqIndex2 = 0; uSeqIndex2 < uSeqCount2; ++uSeqIndex2)
\r
142 msaA2.SetSeqWeight(uSeqIndex2, 1.0);
\r
143 msaB2.SetSeqWeight(uSeqIndex2, 1.0);
\r
164 Log("\nSPA\n---\n");
\r
165 SCORE SPA = ObjScoreSP(alnA);
\r
166 Log("\nSPB\n---\n");
\r
167 SCORE SPB = ObjScoreSP(alnB);
\r
169 Log("\nXPA\n---\n");
\r
170 SCORE XPA = ObjScoreXP(msaA1, msaA2);
\r
171 Log("\nXPB\n---\n");
\r
172 SCORE XPB = ObjScoreXP(msaB1, msaB2);
\r
174 Log("SPA=%.4g SPB=%.4g Diff=%.4g\n", SPA, SPB, SPA - SPB);
\r
175 Log("XPA=%.4g XPB=%.4g Diff=%.4g\n", XPA, XPB, XPA - XPB);
\r