1 //uchime by Robert C. Edgar http://drive5.com/uchime This code is donated to the public domain.
\r
6 void WriteChimeFileHdr(FILE *f)
\r
32 void WriteChimeHit(FILE *f, const ChimeHit2 &Hit)
\r
39 fprintf(f, "0.0000"); // 0
\r
42 "\t%s", Hit.QLabel.c_str()); // 1
\r
65 fprintf(f, "%.4f", Hit.Score); // 0
\r
68 fputs(Hit.QLabel.c_str(), f); // 1
\r
71 fputs(Hit.ALabel.c_str(), f); // 2
\r
74 fputs(Hit.BLabel.c_str(), f); // 3
\r
76 fprintf(f, "\t%.1f", Hit.PctIdQM); // 4
\r
77 fprintf(f, "\t%.1f", Hit.PctIdQA); // 5
\r
78 fprintf(f, "\t%.1f", Hit.PctIdQB); // 6
\r
79 fprintf(f, "\t%.1f", Hit.PctIdAB); // 7
\r
80 fprintf(f, "\t%.1f", Hit.PctIdQT); // 8
\r
82 fprintf(f, "\t%u", Hit.CS_LY); // 9
\r
83 fprintf(f, "\t%u", Hit.CS_LN); // 10
\r
84 fprintf(f, "\t%u", Hit.CS_LA); // 11
\r
86 fprintf(f, "\t%u", Hit.CS_RY); // 12
\r
87 fprintf(f, "\t%u", Hit.CS_RN); // 13
\r
88 fprintf(f, "\t%u", Hit.CS_RA); // 14
\r
90 fprintf(f, "\t%.2f", Hit.Div); // 15
\r
92 fprintf(f, "\t%c", yon(Hit.Accept())); // 16
\r
96 unsigned GetUngappedLength(const byte *Seq, unsigned L)
\r
99 for (unsigned i = 0; i < L; ++i)
\r
100 if (!isgap(Seq[i]))
\r
105 void WriteChimeHitX(FILE *f, const ChimeHit2 &Hit)
\r
110 if (Hit.Div <= 0.0)
\r
113 const string &Q3 = Hit.Q3;
\r
114 const string &A3 = Hit.A3;
\r
115 const string &B3 = Hit.B3;
\r
117 const byte *Q3Seq = (const byte *) Q3.c_str();
\r
118 const byte *A3Seq = (const byte *) A3.c_str();
\r
119 const byte *B3Seq = (const byte *) B3.c_str();
\r
122 unsigned ColCount = SIZE(Q3);
\r
123 asserta(SIZE(A3) == ColCount && SIZE(B3) == ColCount);
\r
125 unsigned LQ = GetUngappedLength(Q3Seq, ColCount);
\r
126 unsigned LA = GetUngappedLength(A3Seq, ColCount);
\r
127 unsigned LB = GetUngappedLength(B3Seq, ColCount);
\r
130 fprintf(f, "------------------------------------------------------------------------\n");
\r
131 fprintf(f, "Query (%5u nt) %s\n", LQ, Hit.QLabel.c_str());
\r
132 fprintf(f, "ParentA (%5u nt) %s\n", LA, Hit.ALabel.c_str());
\r
133 fprintf(f, "ParentB (%5u nt) %s\n", LB, Hit.BLabel.c_str());
\r
135 // Strip terminal gaps in query
\r
136 unsigned FromCol = UINT_MAX;
\r
137 unsigned ToCol = UINT_MAX;
\r
138 for (unsigned Col = 0; Col < ColCount; ++Col)
\r
140 if (!isgap(Q3Seq[Col]))
\r
142 if (FromCol == UINT_MAX)
\r
151 for (unsigned Col = 0; Col < FromCol; ++Col)
\r
153 if (!isgap(A3Seq[Col]))
\r
155 if (!isgap(B3Seq[Col]))
\r
159 unsigned Range = ToCol - FromCol + 1;
\r
160 unsigned RowCount = (Range + 79)/80;
\r
161 unsigned RowFromCol = FromCol;
\r
162 for (unsigned RowIndex = 0; RowIndex < RowCount; ++RowIndex)
\r
165 unsigned RowToCol = RowFromCol + 79;
\r
166 if (RowToCol > ToCol)
\r
170 fprintf(f, "A %5u ", APos + 1);
\r
171 for (unsigned Col = RowFromCol; Col <= RowToCol; ++Col)
\r
173 char q = Q3Seq[Col];
\r
174 char a = A3Seq[Col];
\r
177 fprintf(f, "%c", a);
\r
181 fprintf(f, " %u\n", APos);
\r
184 fprintf(f, "Q %5u ", QPos + 1);
\r
185 for (unsigned Col = RowFromCol; Col <= RowToCol; ++Col)
\r
187 char q = Q3Seq[Col];
\r
188 fprintf(f, "%c", q);
\r
192 fprintf(f, " %u\n", QPos);
\r
195 fprintf(f, "B %5u ", BPos + 1);
\r
196 for (unsigned Col = RowFromCol; Col <= RowToCol; ++Col)
\r
198 char q = Q3Seq[Col];
\r
199 char b = B3Seq[Col];
\r
202 fprintf(f, "%c", b);
\r
206 fprintf(f, " %u\n", BPos);
\r
209 fprintf(f, "Diffs ");
\r
210 for (unsigned Col = RowFromCol; Col <= RowToCol; ++Col)
\r
212 char q = Q3Seq[Col];
\r
213 char a = A3Seq[Col];
\r
214 char b = B3Seq[Col];
\r
217 if (isgap(q) || isgap(a) || isgap(b))
\r
219 else if (Col < Hit.ColXLo)
\r
221 if (q == a && q == b)
\r
223 else if (q == a && q != b)
\r
225 else if (q == b && q != a)
\r
227 else if (a == b && q != a)
\r
232 else if (Col > Hit.ColXHi)
\r
234 if (q == a && q == b)
\r
236 else if (q == b && q != a)
\r
238 else if (q == a && q != b)
\r
240 else if (a == b && q != a)
\r
246 fprintf(f, "%c", c);
\r
251 fprintf(f, "Votes ");
\r
252 for (unsigned Col = RowFromCol; Col <= RowToCol; ++Col)
\r
254 char q = Q3Seq[Col];
\r
255 char a = A3Seq[Col];
\r
256 char b = B3Seq[Col];
\r
258 bool PrevGap = Col > 0 && (isgap(Q3Seq[Col-1]) || isgap(A3Seq[Col-1]) || isgap(B3Seq[Col-1]));
\r
259 bool NextGap = Col+1 < ColCount && (isgap(Q3Seq[Col+1]) || isgap(A3Seq[Col+1]) || isgap(B3Seq[Col+1]));
\r
262 if (isgap(q) || isgap(a) || isgap(b) || PrevGap || NextGap)
\r
264 else if (Col < Hit.ColXLo)
\r
266 if (q == a && q == b)
\r
268 else if (q == a && q != b)
\r
270 else if (q == b && q != a)
\r
275 else if (Col > Hit.ColXHi)
\r
277 if (q == a && q == b)
\r
279 else if (q == b && q != a)
\r
281 else if (q == a && q != b)
\r
287 fprintf(f, "%c", c);
\r
292 fprintf(f, "Model ");
\r
293 for (unsigned Col = RowFromCol; Col <= RowToCol; ++Col)
\r
295 if (Col < Hit.ColXLo)
\r
297 else if (Col >= Hit.ColXLo && Col <= Hit.ColXHi)
\r
309 double PctIdBestP = max(Hit.PctIdQA, Hit.PctIdQB);
\r
310 double Div = (Hit.PctIdQM - PctIdBestP)*100.0/PctIdBestP;
\r
312 unsigned LTot = Hit.CS_LY + Hit.CS_LN + Hit.CS_LA;
\r
313 unsigned RTot = Hit.CS_RY + Hit.CS_RN + Hit.CS_RA;
\r
315 double PctL = Pct(Hit.CS_LY, LTot);
\r
316 double PctR = Pct(Hit.CS_RY, RTot);
\r
319 "Ids. QA %.1f%%, QB %.1f%%, AB %.1f%%, QModel %.1f%%, Div. %+.1f%%\n",
\r
327 "Diffs Left %u: N %u, A %u, Y %u (%.1f%%); Right %u: N %u, A %u, Y %u (%.1f%%), Score %.4f\n",
\r
328 LTot, Hit.CS_LN, Hit.CS_LA, Hit.CS_LY, PctL,
\r
329 RTot, Hit.CS_RN, Hit.CS_RA, Hit.CS_RY, PctR,
\r