b8d545da050a4050b3f0ae6af8fa93ac8ce6c6a3 markd Sun Jun 11 14:52:31 2023 -0700 fixed incorrect calculation of indels causes by fix to pslTEnd() on protein/dna PSLs diff --git src/hg/blastToPsl/pslBuild.c src/hg/blastToPsl/pslBuild.c index f472229..2d7c2da 100644 --- src/hg/blastToPsl/pslBuild.c +++ src/hg/blastToPsl/pslBuild.c @@ -170,85 +170,62 @@ if (psl->strand[0] == '-') reverseIntRange(&psl->qStart, &psl->qEnd, psl->qSize); psl->tStart = psl->tStarts[0]; psl->tEnd = psl->tStarts[newIBlk] + (blk->q2tBlkSizeMult * blkSize); if (psl->strand[1] == '-') reverseIntRange(&psl->tStart, &psl->tEnd, psl->tSize); if (flags & bldPslx) { psl->qSequence[newIBlk] = cloneStringZ(blk->qAln + blk->alnStart, blkSize); psl->tSequence[newIBlk] = cloneStringZ(blk->tAln + blk->alnStart, blkSize); } psl->blockCount++; } -static void countIndels(struct psl *psl) -/* update indel counts in psl after adding a block */ -{ -if (psl->blockCount > 1) - { - int iBlk = psl->blockCount - 1; - if (pslQEnd(psl, iBlk-1) != psl->qStarts[iBlk]) - { - /* insert in query */ - psl->qNumInsert++; - psl->qBaseInsert += (psl->qStarts[iBlk] - pslQEnd(psl, iBlk-1)); - } - if (pslTEnd(psl, iBlk-1) != psl->tStarts[iBlk]) - { - /* insert in target */ - psl->tNumInsert++; - psl->tBaseInsert += (psl->tStarts[iBlk] - pslTEnd(psl, iBlk-1)); - } - } -} - static void countMatches(struct psl* psl, struct block* blk, unsigned flags) /* update the PSL match/mismatch after adding a block. */ { int i, alnLen = (blk->alnEnd - blk->alnStart); char *qPtr = blk->qAln + blk->alnStart; char *tPtr = blk->tAln + blk->alnStart; boolean isProt = isProteinSeqs(flags); for (i = 0; i < alnLen; i++, qPtr++, tPtr++) { if ((!isProt && ((*qPtr == 'N') || (*tPtr == 'N'))) || (isProt && ((*qPtr == 'X') || (*tPtr == 'X')))) psl->repMatch += blk->countMult; else if (*qPtr == *tPtr) psl->match += blk->countMult; else psl->misMatch += blk->countMult; } } static void hspToBlocks(struct psl *psl, int *pslSpace, struct block *blk, unsigned flags) /* build PSl blocks from an HSP */ { /* fill in ungapped blocks */ while (nextUngappedBlk(blk)) { addUngappedBlock(psl, pslSpace, blk, flags); - countIndels(psl); + pslComputeInsertCounts(psl); countMatches(psl, blk, flags); } +pslComputeInsertCounts(psl); assert(blk->qStart == blk->qEnd); assert(blk->tStart == blk->tEnd); -// FIXME -//assert(blk->qStart == pslQEnd(psl, psl->blockCount-1)); -//assert(blk->tStart == pslTEnd(psl, psl->blockCount-1)); } static void makeUntranslated(struct psl* psl) /* convert a PSL so it is in the untranslated form produced by blat */ { if (psl->strand[1] == '-') { /* swap around blocks so it's query that is reversed */ int i; for (i = 0; i < psl->blockCount; i++) { psl->tStarts[i] = psl->tSize - (psl->tStarts[i] + psl->blockSizes[i]); psl->qStarts[i] = psl->qSize - (psl->qStarts[i] + psl->blockSizes[i]); } reverseUnsigned(psl->tStarts, psl->blockCount);