f198c745e42f17eb535990d2617d133ebd1396d6 angie Tue Apr 29 11:07:40 2025 -0700 Fix vpExpandIndelGaps corner case from #34135 notes 9 & 14 (shifting a gap left past the beginning of the exon in some diabolical imperfect AT repeats in 3' UTR of bosTau6 NM_177504). diff --git src/hg/lib/variantProjector.c src/hg/lib/variantProjector.c index e74e5f07a9f..3eb731c39ec 100644 --- src/hg/lib/variantProjector.c +++ src/hg/lib/variantProjector.c @@ -806,47 +806,42 @@ uint gStartL = txAli->tStarts[ix] + txAli->blockSizes[ix]; uint gEndL = txAli->tStarts[ix+1]; uint gStartR = gStartL, gEndR = gEndL; uint qLen = txAli->qStarts[ix+1] - txAli->qStarts[ix] - txAli->blockSizes[ix]; char txCpyL[qLen+1], txCpyR[qLen+1]; if (pslQStrand(txAli) == '-') { // Change query sequence to genome + strand uint qGapStart = txAli->qSize - txAli->qStarts[ix+1]; safencpy(txCpyL, sizeof(txCpyL), txSeq->dna+qGapStart, qLen); reverseComplement(txCpyL, qLen); } else safencpy(txCpyL, sizeof(txCpyL), txSeq->dna+txAli->qStarts[ix+1]-qLen, qLen); safencpy(txCpyR, sizeof(txCpyR), txCpyL, qLen); - uint shiftL = indelShift(gSeqWin, &gStartL, &gEndL, txCpyL, INDEL_SHIFT_NO_MAX, isdLeft); + // Don't slide left past the start of the exon -- we could collide with a previously extended + // double-sided gap (#34135 notes 9 & 14). But allow sliding right into the next gap. + uint shiftL = indelShift(gSeqWin, &gStartL, &gEndL, txCpyL, txAli->blockSizes[ix], isdLeft); uint shiftR = indelShift(gSeqWin, &gStartR, &gEndR, txCpyR, INDEL_SHIFT_NO_MAX, isdRight); if (shiftL) { // Expand gap to the left -- decrease blockSizes[ix]. if (txAli->blockSizes[ix] < shiftL) { - if (ix == 0) - warn("vpExpandIndelGaps: %s gapIx %d slides left past start of first block. " - "Skipping. (Should we make a 0-length block at beginning?)", - txAli->qName, ix); - else - warn("vpExpandIndelGaps: %s gapIx %d slides left past the start of block %d, " - "but this should have already been taken care of when pslExpandGapRight " - "was called for gapIx %d. Investigate...", - txAli->qName, ix, ix, ix-1); - continue; + errAbort("vpExpandIndelGaps: %s gapIx %d slides left past the start of block %d, " + "but use of maxShift should have prevented that.", + txAli->qName, ix, ix); } else { txAli->blockSizes[ix] -= shiftL; } modified = TRUE; } if (shiftR) { boolean deletedBlocks = pslExpandGapRight(txAli, ix, shiftR, TRUE); if (deletedBlocks && ix < txAli->blockCount - 1 && pslIntronTooShort(txAli, ix, MIN_INTRON)) { // Repeat w/same ix because block ix has expanded & gap to the right has changed. // indelShift doesn't work on gaps that have already been expanded, though --