476cd67154ef60c8e304628a7731bf89102feab4 galt Thu Jun 21 00:50:35 2012 -0700 fixing out=sim4 to filter on minIdentity, and to reverse the output order of alignments diff --git src/jkOwnLib/gfOut.c src/jkOwnLib/gfOut.c index 1b01512..bafa6e6 100644 --- src/jkOwnLib/gfOut.c +++ src/jkOwnLib/gfOut.c @@ -345,59 +345,84 @@ struct axtData *aod = out->data; struct axtBundle *gab; for (gab = aod->bundleList; gab != NULL; gab = gab->next) { struct axt *axt; for (axt = gab->axtList; axt != NULL; axt = axt->next) { struct mafAli temp; mafFromAxtTemp(axt, gab->tSize, gab->qSize, &temp); mafWrite(f, &temp); } } axtBundleFreeList(&aod->bundleList); } -static double axtIdRatio(struct axt *axt) -/* Return matches/total. */ +static int axtMatchCount(struct axt *axt) +/* Return matches. */ { int matchCount = 0; int i; for (i=0; i<axt->symCount; ++i) { if (axt->qSym[i] == axt->tSym[i]) ++matchCount; } +return matchCount; +} + +static double axtIdRatio(struct axt *axt) +/* Return matches/total. */ +{ +int matchCount = axtMatchCount(axt); return (double)matchCount/(double)axt->symCount; } +static double axtListRatio(struct axt *axt) +/* Return matches/total. */ +{ +int matchCount = 0; +int symCount = 0; +for (; axt != NULL; axt = axt->next) + { + matchCount += axtMatchCount(axt); + symCount += axt->symCount; + } +return (double)matchCount/(double)symCount; +} + static void sim4QueryOut(struct gfOutput *out, FILE *f) /* Do sim4-like output - at end of processing query. */ { struct axtData *aod = out->data; struct axtBundle *gab; +slReverse(&aod->bundleList); for (gab = aod->bundleList; gab != NULL; gab = gab->next) { struct axt *axt = gab->axtList; + // check minIdentity of the entire alignment + int goodPpt = 1000 * axtListRatio(axt); + if (!(goodPpt >= out->minGood)) + continue; fprintf(f, "\n"); fprintf(f, "seq1 = %s, %d bp\n", axt->qName, gab->qSize); fprintf(f, "seq2 = %s, %d bp\n", axt->tName, gab->tSize); fprintf(f, "\n"); if (axt->qStrand == '-') fprintf(f, "(complement)\n"); - for (axt = gab->axtList; axt != NULL; axt = axt->next) + for (; axt != NULL; axt = axt->next) { fprintf(f, "%d-%d ", axt->qStart+1, axt->qEnd); fprintf(f, "(%d-%d) ", axt->tStart+1, axt->tEnd); fprintf(f, "%1.0f%% ", 100.0 * axtIdRatio(axt)); if (axt->qStrand == '-') fprintf(f, "<-\n"); else fprintf(f, "->\n"); } } axtBundleFreeList(&aod->bundleList); } static void blastQueryOut(struct gfOutput *out, FILE *f) /* Output blast on query. */