e199edebab9458a3b535bb2d0e32c52167d8fdae jcasper Fri Apr 15 14:09:59 2022 -0700 Adding pack/squish to hic tracks, plus interaction distance filters. refs #29280 diff --git src/hg/hgTracks/hicTrack.c src/hg/hgTracks/hicTrack.c index 22dcbfe..322cc02 100644 --- src/hg/hgTracks/hicTrack.c +++ src/hg/hgTracks/hicTrack.c @@ -2,58 +2,74 @@ /* Copyright (C) 2018 The Regents of the University of California * See kent/LICENSE or http://genome.ucsc.edu/license/ for licensing information. */ #include "common.h" #include "obscure.h" #include "hgTracks.h" #include "bedCart.h" #include "bigWarn.h" #include "interact.h" // for interact structures #include "hicUi.h" #include "cStraw.h" #include "hic.h" #include "htmlColor.h" +static double hicSqueezeFactor(enum trackVisibility vis) +/* Controls the height multiplier for various draw modes. */ +{ +switch (vis) + { + case tvDense: + return 0.125; + case tvSquish: + return 0.25; + case tvPack: + return 0.5; + default: + return 1.0; + }; + return 1.0; +} static int hicTotalHeight(struct track *tg, enum trackVisibility vis) /* Calculate height of the track across all windows. For arc and triangle, this * depends on the actual interactions being drawn. For square, it just depends * on the height of the highest square (yes, this currently means smaller * squares are stretched vertically in multi-region - it's not ideal). */ { int maxHeight = 0; struct track *trackScan = tg; while (trackScan != NULL) // Canvas max height for this and previous windows ... { if ((int)trackScan->maxRange > maxHeight) maxHeight = (int)trackScan->maxRange; trackScan = trackScan->prevWindow; } trackScan = tg->nextWindow; while (trackScan != NULL) // ... and look at following windows ... { if ((int)trackScan->maxRange > maxHeight) maxHeight = (int)trackScan->maxRange; trackScan = trackScan->nextWindow; } -// Later, try making pack 0.5, squish 0.25, and dense 0.125. The various -// draw functions will also have to be modified accordingly. -if ( tg->visibility == tvDense) - maxHeight *= 0.5; +maxHeight *= hicSqueezeFactor(vis); + +if (maxHeight < tg->lineHeight) + maxHeight = tg->lineHeight; return maxHeight; } struct hicMeta *grabHeader(struct track *tg) /* Fetch a hicMeta structure that describes the Hi-C data associated with * the track. */ { char *filename = trackDbSettingOrDefault(tg->tdb, "bigDataUrl", NULL); struct hicMeta *metaResult = NULL; if (filename == NULL) { warn("Missing bigDataUrl setting for track %s", tg->track); return NULL; } char *errMsg = hicLoadHeader(filename, &metaResult, database); @@ -99,56 +115,82 @@ struct interact *hicItems = NULL; tg->networkErrMsg = hicLoadData(hicFileInfo, binSize, normalization, chromName, strawStart, winEnd-1, chromName, strawStart, winEnd-1, &hicItems); // Using the interact structure because it has convenient fields, but this is not interact data and // shouldn't be passed to those functions. int numRecords = slCount(hicItems), filtNumRecords = 0; tg->maxRange = 0.0; // the max height of an interaction in this window double *countsCopy = NULL; if (numRecords > 0) AllocArray(countsCopy, numRecords); struct interact *thisHic = hicItems; char *drawMode = hicUiFetchDrawMode(cart, tg->tdb); +struct interact* filteredOut = NULL; +struct interact** prevNextPtr = &hicItems; // for removing items from the linked list + +double maxRange = hicUiMaxInteractionRange(cart, tg->tdb); +double minRange = hicUiMinInteractionRange(cart, tg->tdb); + while (thisHic != NULL) { + // Add filtering based on max interaction distance + if (sameString(thisHic->sourceChrom, thisHic->targetChrom)) + { + unsigned leftEdge = thisHic->sourceStart < thisHic->targetStart ? thisHic->sourceStart : thisHic->targetStart; + unsigned rightEdge = thisHic->sourceEnd > thisHic->targetEnd ? thisHic->sourceEnd : thisHic->targetEnd; + if ((maxRange && maxRange < (double)(rightEdge - leftEdge) ) || + (minRange && minRange > (double)(rightEdge - leftEdge) )) + { + // a bit of pointer play to avoid repeated calls to slRemoveEl + *prevNextPtr = thisHic->next; // set prev element's next to the following element + slAddHead(&filteredOut, thisHic); + thisHic = *prevNextPtr; // restore thisHic to point to the next element + continue; + } + } + if (sameString(drawMode, HIC_DRAW_MODE_ARC)) { // we omit self-interactions in arc mode (they'd just be weird vertical lines) if (sameString(thisHic->sourceChrom, thisHic->targetChrom) && (thisHic->sourceStart == thisHic->targetStart)) { thisHic = thisHic->next; continue; } } countsCopy[filtNumRecords++] = thisHic->value; // Calculate the track draw height required to see this item int leftx = max(thisHic->chromStart, winStart); int rightx = min(thisHic->chromEnd, winEnd); double thisHeight = scaleForWindow(insideWidth, winStart, winEnd)*(rightx - leftx)/2.0; // triangle or arc if (sameString(drawMode,HIC_DRAW_MODE_SQUARE)) thisHeight = scaleForWindow(insideWidth, winStart, winEnd)*(winEnd-winStart); // square - always draw the full square if (thisHeight > tg->maxRange) tg->maxRange = thisHeight; + prevNextPtr = &thisHic->next; thisHic = thisHic->next; } +if (filteredOut != NULL) + interactFreeList(&filteredOut); + // Heuristic for auto-scaling the color gradient based on the scores in view - draw the max color value // at or above 2*median score. if (filtNumRecords > 0) tg->graphUpperLimit = 2.0*doubleMedian(filtNumRecords, countsCopy); else tg->graphUpperLimit = 0.0; if (countsCopy != NULL) freeMem(countsCopy); tg->items = hicItems; } void hicLoadItems(struct track *tg) /* Load Hi-C items in (mostly) interact format */ { char *filename = trackDbSettingOrDefault(tg->tdb, "bigDataUrl", NULL); @@ -235,45 +277,41 @@ } static void drawHicTriangle(struct track *tg, int seqStart, int seqEnd, struct hvGfx *hvg, int xOff, int yOff, int width, MgFont *font, Color color, enum trackVisibility vis) /* Draw a list of Hi-C interactions in a triangle */ { double xScale = scaleForWindow(width, seqStart, seqEnd); double yScale = xScale; int maxHeight = tg->height; struct interact *hicItem = NULL; struct hicMeta *hicFileInfo = (struct hicMeta*)tg->customPt; int binSize = hicUiFetchResolutionAsInt(cart, tg->tdb, hicFileInfo, winEnd-winStart); -if (vis == tvDense) - { - yScale *= 0.5; - } +yScale *= hicSqueezeFactor(vis); double maxScore = getHicMaxScore(tg); Color *colorIxs = colorSetForHic(hvg, tg, HIC_SCORE_BINS+1); if (colorIxs == NULL) return; // something went wrong with colors for (hicItem = (struct interact *)tg->items; hicItem; hicItem = hicItem->next) { int leftStart, leftEnd, rightStart, rightEnd; calcItemLeftRightBoundaries(&leftStart, &leftEnd, &rightStart, &rightEnd, binSize, hicItem); - int colorIx; if (hicItem->value > maxScore) colorIx = colorIxs[HIC_SCORE_BINS]; else colorIx = colorIxs[(int)(HIC_SCORE_BINS * hicItem->value/maxScore)]; // adding four polygon points for a diamond, based on the starts and ends of the source and target coordinate ranges. struct gfxPoly *diamond = gfxPolyNew(); // left point of the diamond double x = xScale * (leftStart+rightStart)/2.0; double y = yScale * (rightStart-leftStart)/2.0; gfxPolyAddPoint(diamond, (int)x+xOff, maxHeight-(int)y+yOff); // top point of the diamond x = xScale * (leftStart+rightEnd)/2.0; @@ -298,34 +336,31 @@ static void drawHicSquare(struct track *tg, int seqStart, int seqEnd, struct hvGfx *hvg, int xOff, int yOff, int width, MgFont *font, Color color, enum trackVisibility vis) /* Draw a Hi-C heatmap as a square */ { double xScale = scaleForWindow(width, seqStart, seqEnd); double yScale = xScale; int maxHeight = tg->height; struct interact *hicItem = NULL; struct hicMeta *hicFileInfo = (struct hicMeta*)tg->customPt; int binSize = hicUiFetchResolutionAsInt(cart, tg->tdb, hicFileInfo, winEnd-winStart); if (binSize == 0) return; -if (vis == tvDense) - { - yScale *= 0.5; - } +yScale *= hicSqueezeFactor(vis); double maxScore = getHicMaxScore(tg); Color *colorIxs = colorSetForHic(hvg, tg, HIC_SCORE_BINS+1); if (colorIxs == NULL) return; // something went wrong with colors for (hicItem = (struct interact *)tg->items; hicItem; hicItem = hicItem->next) { int leftStart, leftEnd, rightStart, rightEnd; calcItemLeftRightBoundaries(&leftStart, &leftEnd, &rightStart, &rightEnd, binSize, hicItem); int colorIx; if (hicItem->value > maxScore) colorIx = colorIxs[HIC_SCORE_BINS]; else @@ -360,34 +395,31 @@ static void drawHicArc(struct track *tg, int seqStart, int seqEnd, struct hvGfx *hvg, int xOff, int yOff, int width, MgFont *font, Color color, enum trackVisibility vis) /* Draw Hi-C interactions in arc mode */ { double xScale = scaleForWindow(width, seqStart, seqEnd); double yScale = xScale; int maxHeight = tg->height; struct interact *hicItem = NULL; struct hicMeta *hicFileInfo = (struct hicMeta*)tg->customPt; int binSize = hicUiFetchResolutionAsInt(cart, tg->tdb, hicFileInfo, winEnd-winStart); if (binSize == 0) return; -if (vis == tvDense) - { - yScale *= 0.5; - } +yScale *= hicSqueezeFactor(vis); double maxScore = getHicMaxScore(tg); Color *colorIxs = colorSetForHic(hvg, tg, HIC_SCORE_BINS+1); if (colorIxs == NULL) return; // something went wrong with colors slSort(&tg->items, cmpHicItems); // So that the darkest arcs are drawn on top and not lost for (hicItem = (struct interact *)tg->items; hicItem; hicItem = hicItem->next) { int leftStart = hicItem->sourceStart - winStart; int leftMidpoint = leftStart + binSize/2; int rightStart = hicItem->targetStart - winStart; int rightMidpoint = rightStart + binSize/2; if ((leftMidpoint < 0) || (leftMidpoint > winEnd-winStart)) continue; // skip this item - we'd be drawing to a point off the screen