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("%c %d | ", qSeq[i], quals[i]);
}
printf("
\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))