af401081fe03d2bd9eb087487939075d29ff8e03 max Fri May 31 08:43:30 2024 -0700 not showing bam query alignment if seq is too long, refs #33292 diff --git src/hg/hgc/bamClick.c src/hg/hgc/bamClick.c index e582f1c..2b9cb39 100644 --- src/hg/hgc/bamClick.c +++ src/hg/hgc/bamClick.c @@ -33,69 +33,82 @@ static void singleBamDetails(const bam1_t *bam) /* Print out the properties of this alignment. */ { const bam1_core_t *core = &bam->core; char *itemName = bam1_qname(bam); int tLength = bamGetTargetLength(bam); int tStart = core->pos, tEnd = tStart+tLength; boolean isRc = useStrand && bamIsRc(bam); printPosOnChrom(seqName, tStart, tEnd, NULL, FALSE, itemName); if (!skipQualityScore) printf("Alignment Quality: %d
\n", core->qual); if (core->n_cigar > 50) printf("CIGAR string: Cannot show long CIGAR string, more than 50 operations. Contact us if you need to see the full CIGAR string here.
\n"); else { - printf("CIGAR string: %s (", bamGetCigar(bam)); - bamShowCigarEnglish(bam); - printf(")
\n"); + printf("CIGAR string: %s", bamGetCigar(bam)); + //bamShowCigarEnglish(bam); + printf("
\n"); } printf("Tags:"); bamShowTags(bam); puts("
"); printf("Flags: 0x%02x:
\n   ", core->flag); bamShowFlagsEnglish(bam); puts("
"); if (bamIsRc(bam)) printf("Note: although the read was mapped to the reverse strand of the genome, " "the sequence and CIGAR in BAM are relative to the forward strand.
\n"); puts("
"); struct dnaSeq *genoSeq = hChromSeq(database, seqName, tStart, tEnd); char *qSeq = bamGetQuerySequence(bam, FALSE); +if (core->l_qseq > 5000) + printf("Alignment not shown, very long sequence
\n"); +else + { if (isNotEmpty(qSeq) && !sameString(qSeq, "*")) { char *qSeq = NULL; struct ffAli *ffa = bamToFfAli(bam, genoSeq, tStart, useStrand, &qSeq); printf("Alignment of %s to %s:%d-%d%s:
\n", itemName, seqName, tStart+1, tEnd, (isRc ? " (reverse complemented)" : "")); ffShowSideBySide(stdout, ffa, qSeq, 0, genoSeq->dna, tStart, tLength, 0, tLength, 8, isRc, FALSE); } + } + if (!skipQualityScore && core->l_qseq > 0) { + if (core->l_qseq > 5000) + { + printf("Sequence too long to show quality scores
\n"); + } + else + { printf("Sequence quality scores:
\n\n"); UBYTE *quals = bamGetQueryQuals(bam, useStrand); int i; for (i = 0; i < core->l_qseq; i++) { if (i > 0 && (i % 24) == 0) printf("\n"); printf("", qSeq[i], quals[i]); } printf("
%c
%d
\n"); } } +} static void showOverlap(const bam1_t *leftBam, const bam1_t *rightBam) /* If the two reads overlap, show how. */ { const bam1_core_t *leftCore = &(leftBam->core), *rightCore = &(rightBam->core); int leftStart = leftCore->pos, rightStart = rightCore->pos; int leftLen = bamGetTargetLength(leftBam), rightLen = bamGetTargetLength(rightBam); char *leftSeq = bamGetQuerySequence(leftBam, useStrand); char *rightSeq = bamGetQuerySequence(rightBam, useStrand); if (useStrand && bamIsRc(leftBam)) reverseComplement(leftSeq, strlen(leftSeq)); if (useStrand && bamIsRc(rightBam)) reverseComplement(rightSeq, strlen(rightSeq)); if ((rightStart > leftStart && leftStart + leftLen > rightStart) || (leftStart > rightStart && rightStart+rightLen > leftStart))