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("<B>Alignment Quality: </B>%d<BR>\n", core->qual);
 if (core->n_cigar > 50)
     printf("<B>CIGAR string: </B> Cannot show long CIGAR string, more than 50 operations. Contact us if you need to see the full CIGAR string here.<BR>\n");
 else
     {
-    printf("<B>CIGAR string: </B><tt>%s</tt> (", bamGetCigar(bam));
-    bamShowCigarEnglish(bam);
-    printf(")<BR>\n");
+    printf("<B>CIGAR string: </B><tt>%s</tt>", bamGetCigar(bam));
+    //bamShowCigarEnglish(bam);
+    printf("<BR>\n");
     }
 printf("<B>Tags:</B>");
 bamShowTags(bam);
 puts("<BR>");
 printf("<B>Flags: </B><tt>0x%02x:</tt><BR>\n &nbsp;&nbsp;", core->flag);
 bamShowFlagsEnglish(bam);
 puts("<BR>");
 if (bamIsRc(bam))
     printf("<em>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.</em><BR>\n");
 puts("<BR>");
 struct dnaSeq *genoSeq = hChromSeq(database, seqName, tStart, tEnd);
 char *qSeq = bamGetQuerySequence(bam, FALSE);
+if (core->l_qseq > 5000)
+    printf("<B>Alignment not shown, very long sequence</B><BR>\n");
+else
+    {
     if (isNotEmpty(qSeq) && !sameString(qSeq, "*"))
         {
         char *qSeq = NULL;
         struct ffAli *ffa = bamToFfAli(bam, genoSeq, tStart, useStrand, &qSeq);
         printf("<B>Alignment of %s to %s:%d-%d%s:</B><BR>\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("<B>Sequence too long to show quality scores</B><BR>\n");
+        } 
+    else
+        {
         printf("<B>Sequence quality scores:</B><BR>\n<TT><TABLE><TR>\n");
         UBYTE *quals = bamGetQueryQuals(bam, useStrand);
         int i;
         for (i = 0;  i < core->l_qseq;  i++)
             {
             if (i > 0 && (i % 24) == 0)
                 printf("</TR>\n<TR>");
             printf("<TD>%c<BR>%d</TD>", qSeq[i], quals[i]);
             }
         printf("</TR></TABLE></TT>\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))