e0066b6bba1eea45452a32de72eaf4af8c331bbb angie Mon Aug 8 16:27:02 2011 -0700 Feature #2823 (VCF track handler): Added UI control over height of haplotype-sorting display. diff --git src/hg/hgTracks/vcfTrack.c src/hg/hgTracks/vcfTrack.c index 5b9ef2f..9b1d521 100644 --- src/hg/hgTracks/vcfTrack.c +++ src/hg/hgTracks/vcfTrack.c @@ -1,37 +1,37 @@ /* vcfTrack -- handlers for Variant Call Format data. */ #include "common.h" #include "bigWarn.h" #include "dystring.h" #include "errCatch.h" #include "hacTree.h" #include "hdb.h" #include "hgTracks.h" #include "pgSnp.h" #include "trashDir.h" #include "vcf.h" +#include "vcfUi.h" #if (defined USE_TABIX && defined KNETFILE_HOOKS) #include "knetUdc.h" #include "udc.h" #endif//def USE_TABIX && KNETFILE_HOOKS #ifdef USE_TABIX //#*** TODO: use trackDb/cart setting or something static boolean doHapClusterDisplay = TRUE; -static boolean colorHapByRefAlt = TRUE; static struct pgSnp *vcfFileToPgSnp(struct vcfFile *vcff) /* Convert vcff's records to pgSnp; don't free vcff until you're done with pgSnp * because it contains pointers into vcff's records' chrom. */ { struct pgSnp *pgsList = NULL; struct vcfRecord *rec; int maxLen = 33; for (rec = vcff->records; rec != NULL; rec = rec->next) { struct pgSnp *pgs = pgSnpFromVcfRecord(rec); // Insertion sequences can be quite long; abbreviate here for display. int len = strlen(pgs->name); if (len > maxLen) { @@ -299,77 +299,55 @@ } } struct hacTree *ht = hacTreeFromItems((struct slList *)(hapArray[0]), lm, cwaDistance, cwaMerge, cwaCmp, &helper); unsigned short *gtHapOrder = needMem(vcff->genotypeCount * 2 * sizeof(unsigned short)); rSetGtHapOrder(ht, gtHapOrder, retGtHapEnd); return gtHapOrder; } INLINE char *hapIxToAllele(int hapIx, char *refAllele, char *altAlleles[]) /* Look up allele by index into reference allele and alternate allele(s). */ { return (hapIx == 0) ? refAllele : altAlleles[hapIx-1]; } +//#*** unused... add UI option... INLINE Color colorFromGt(struct vcfGenotype *gt, int ploidIx, char *refAllele, char *altAlleles[], int altCount, boolean grayUnphasedHet) /* Color allele by base. */ { int hapIx = ploidIx ? gt->hapIxB : gt->hapIxA; char *allele = hapIxToAllele(hapIx, refAllele, altAlleles); if (gt->isHaploid && hapIx > 0) return shadesOfGray[5]; if (grayUnphasedHet && !gt->isPhased && gt->hapIxA != gt->hapIxB) return shadesOfGray[5]; // Copying pgSnp color scheme here, using first base of allele which is not ideal for multibase // but allows us to simplify it to 5 colors: else if (allele[0] == 'A') return MG_RED; else if (allele[0] == 'C') return MG_BLUE; else if (allele[0] == 'G') return darkGreenColor; else if (allele[0] == 'T') return MG_MAGENTA; else return shadesOfGray[5]; } -INLINE Color colorFromRefAlt(struct vcfGenotype *gt, int hapIx, boolean grayUnphasedHet) -/* Color allele red for alternate allele, blue for reference allele -- - * except for special center variant, make it yellow/green for contrast. */ -{ -if (grayUnphasedHet && !gt->isPhased && gt->hapIxA != gt->hapIxB) - return shadesOfGray[5]; -int alIx = hapIx ? gt->hapIxB : gt->hapIxA; -return alIx ? MG_RED : MG_BLUE; -} - - -INLINE int drawOneHap(struct vcfGenotype *gt, int hapIx, - char *ref, char *altAlleles[], int altCount, - struct hvGfx *hvg, int x1, int y, int w, int itemHeight, int lineHeight) -/* Draw a base-colored box for genotype[hapIx]. Return the new y offset. */ -{ -Color color = colorHapByRefAlt ? colorFromRefAlt(gt, hapIx, TRUE) : - colorFromGt(gt, hapIx, ref, altAlleles, altCount, TRUE); -hvGfxBox(hvg, x1, y, w, itemHeight+1, color); -y += itemHeight+1; -return y; -} - INLINE char *gtSummaryString(struct vcfRecord *rec, char **altAlleles, int altCount) // Make pgSnp-like mouseover text, but with genotype counts instead of allele counts. // NOTE 1: Returned string is statically allocated, don't free it! // NOTE 2: if revCmplDisp is set, this reverse-complements rec->ref and altAlleles! { static struct dyString *dy = NULL; if (dy == NULL) dy = dyStringNew(0); dyStringClear(dy); const struct vcfFile *vcff = rec->file; int gtRefRefCount = 0, gtRefAltCount = 0, gtAltAltCount = 0, gtOtherCount = 0; int i; for (i=0; i < vcff->genotypeCount; i++) { struct vcfGenotype *gt = &(rec->genotypes[i]); @@ -393,64 +371,94 @@ dyStringPrintf(dy, "%s/%s:%d %s/%s:%d %s/%s:%d", rec->ref, rec->ref, gtRefRefCount, rec->ref, altAlleles[0], gtRefAltCount, altAlleles[0], altAlleles[0], gtAltAltCount); if (gtOtherCount > 0) dyStringPrintf(dy, " other:%d", gtOtherCount); // Restore original values of pooled strings. if (revCmplDisp) { reverseComplement(rec->ref, strlen(rec->ref)); for (i=0; i < altCount; i++) reverseComplement(altAlleles[i], strlen(altAlleles[i])); } return dy->string; } +// This is initialized when we start drawing: +static Color purple = 0; + static void drawOneRec(struct vcfRecord *rec, unsigned short *gtHapOrder, int gtHapEnd, struct track *tg, struct hvGfx *hvg, int xOff, int yOff, int width, boolean isCenter) /* Draw a stack of genotype bars for this record */ { static struct dyString *tmp = NULL; if (tmp == NULL) tmp = dyStringNew(0); char *altAlleles[256]; int altCount; -const int lineHeight = tg->lineHeight; -const int itemHeight = tg->heightPer; const double scale = scaleForPixels(width); int x1 = round((double)(rec->chromStart-winStart)*scale) + xOff; int x2 = round((double)(rec->chromEnd-winStart)*scale) + xOff; int w = x2-x1; if (w <= 1) { x1--; w = 3; } -int y = yOff; dyStringClear(tmp); dyStringAppend(tmp, rec->alt); altCount = chopCommas(tmp->string, altAlleles); +double hapsPerPix = (2 * (double)rec->file->genotypeCount / tg->height); +int pixIx; +for (pixIx = 0; pixIx < tg->height; pixIx++) + { + int gtHapOrderIxStart = round(hapsPerPix * pixIx); + int gtHapOrderIxEnd = round(hapsPerPix * (pixIx + 1)); + if (gtHapOrderIxEnd == gtHapOrderIxStart) + gtHapOrderIxEnd++; + int unks = 0, refs = 0, alts = 0; int gtHapOrderIx; -for (gtHapOrderIx = 0; gtHapOrderIx < gtHapEnd; gtHapOrderIx++) + for (gtHapOrderIx = gtHapOrderIxStart; gtHapOrderIx < gtHapOrderIxEnd; gtHapOrderIx++) { int gtHapIx = gtHapOrder[gtHapOrderIx]; int hapIx = gtHapIx & 1; int gtIx = gtHapIx >>1; struct vcfGenotype *gt = &(rec->genotypes[gtIx]); - y = drawOneHap(gt, hapIx, rec->ref, altAlleles, altCount, - hvg, x1, y, w, itemHeight, lineHeight); + if (!gt->isPhased && gt->hapIxA != gt->hapIxB) + unks++; + else + { + int alIx = hapIx ? gt->hapIxB : gt->hapIxA; + if (alIx) + alts++; + else + refs++; + } + } + const int fudgeFactor = 4; + Color col = MG_BLACK; + if (unks > (refs + alts)) + col = shadesOfGray[5]; + else if (alts > fudgeFactor * refs) + col = MG_RED; + else if (refs > fudgeFactor * alts) + col = MG_BLUE; + else + col = purple; + int y = yOff + pixIx; + hvGfxLine(hvg, x1, y, x2, y, col); } char *mouseoverText = gtSummaryString(rec, altAlleles, altCount); if (isCenter) { // Thick black lines to distinguish this variant: int yBot = yOff + tg->height - 2; hvGfxBox(hvg, x1-3, yOff, 3, tg->height, MG_BLACK); hvGfxBox(hvg, x2, yOff, 3, tg->height, MG_BLACK); hvGfxLine(hvg, x1-2, yOff, x2+2, yOff, MG_BLACK); hvGfxLine(hvg, x1-2, yBot, x2+2, yBot, MG_BLACK); // Special mouseover instructions: static struct dyString *dy = NULL; if (dy == NULL) dy = dyStringNew(0); dyStringPrintf(dy, "%s Haplotypes sorted on ", mouseoverText); @@ -493,55 +501,59 @@ return i-1; } } return defaultIx; } static void vcfHapClusterDraw(struct track *tg, int seqStart, int seqEnd, struct hvGfx *hvg, int xOff, int yOff, int width, MgFont *font, Color color, enum trackVisibility vis) /* Split samples' chromosomes (haplotypes), cluster them by center-weighted * alpha similarity, and draw in the order determined by clustering. */ { const struct vcfFile *vcff = tg->extraUiData; if (vcff->records == NULL) return; +purple = hvGfxFindColorIx(hvg, 0x99, 0x00, 0xcc); unsigned short gtHapEnd = 0; int ix, centerIx = getCenterVariantIx(tg, seqStart, seqEnd, vcff->records); unsigned short *gtHapOrder = clusterChroms(vcff, centerIx, >HapEnd); struct vcfRecord *rec, *centerRec = NULL; for (rec = vcff->records, ix=0; rec != NULL; rec = rec->next, ix++) { if (ix == centerIx) centerRec = rec; else drawOneRec(rec, gtHapOrder, gtHapEnd, tg, hvg, xOff, yOff, width, FALSE); } // Draw the center rec on top, outlined with black lines, to make sure it is very visible: drawOneRec(centerRec, gtHapOrder, gtHapEnd, tg, hvg, xOff, yOff, width, TRUE); } static int vcfHapClusterTotalHeight(struct track *tg, enum trackVisibility vis) /* Return height of haplotype graph (2 * #samples * lineHeight); * 2 because we're assuming diploid genomes here, no XXY, tetraploid etc. */ { // Should we make it single-height when on chrY? const struct vcfFile *vcff = tg->extraUiData; if (vcff->records == NULL) return 0; -int ploidy = 2; -tg->height = ploidy * vcff->genotypeCount * tg->lineHeight; +int ploidy = sameString(chromName, "chrY") ? 1 : 2; +int simpleHeight = ploidy * vcff->genotypeCount * tg->lineHeight; +int defaultHeight = min(simpleHeight, VCF_DEFAULT_HAP_HEIGHT); +int cartHeight = cartOrTdbInt(cart, tg->tdb, VCF_HAP_HEIGHT_VAR, defaultHeight); +tg->height = min(cartHeight, maximumTrackHeight(tg)); return tg->height; } static char *vcfHapClusterTrackName(struct track *tg, void *item) /* If someone asks for itemName/mapItemName, just send name of track like wiggle. */ { return tg->track; } static void vcfHapClusterOverloadMethods(struct track *tg, struct vcfFile *vcff) /* If we confirm at load time that we can draw a haplotype graph, use * this to overwrite the methods for the rest of execution: */ { tg->heightPer = (tg->visibility == tvSquish) ? (tl.fontHeight/4) : (tl.fontHeight / 2); tg->lineHeight = tg->heightPer + 1;