5fc426954da9ceb7dc42b9760858bbd1189760e1 braney Fri Sep 4 13:08:14 2026 -0700 quickLift: lift MAF blocks, and let bigMaf and wigMaf tracks into the hub Most of this was already written. mafSubset does the part that looked hard, which is recomputing every row's start and size when columns are taken away, so what was left was deciding where to cut. quickLiftMafs cuts a block at every chain block boundary. Inside one chain block the two assemblies run in step, so the columns carry over untouched and only the first row's coordinates change. Across a boundary the reference either loses bases or gains them, and either way the block can no longer be one contiguous run on the reference, which is the one thing a MAF block has to be. The reference row is named with the assembly name minus any hub prefix, since that is the name the maf drawing code builds when it goes looking for it. A minus strand chain turns the block over, so every row is turned over with it and the forward start comes from the far end of the run. That path is written but has not been exercised: no minus strand quickLift chain wins in a window I could find. validateOneTdb accepts bigMaf and wigMaf. Plain maf is left out on purpose: those tracks are drawn by mafTrack.c, which has no quickLift path, so offering them would hand back a track read from the wrong assembly. refs #38249 diff --git src/hg/lib/quickLift.c src/hg/lib/quickLift.c index 57d99647552..1b58fb45b8f 100644 --- src/hg/lib/quickLift.c +++ src/hg/lib/quickLift.c @@ -16,30 +16,31 @@ #include "bbiFile.h" #include "chainNetDbLoad.h" #include "hdb.h" #include "jksql.h" #include "hgConfig.h" #include "quickLift.h" #include "genePredReader.h" #include "bigChain.h" #include "bigLink.h" #include "chromAlias.h" #include "customTrack.h" #include "encode/encodePeak.h" #include "psl.h" #include "chainToPsl.h" #include "pslTransMap.h" +#include "maf.h" struct bigBedInterval *quickLiftGetIntervals(char *quickLiftFile, struct bbiFile *bbi, char *chrom, int start, int end, struct hash **pChainHash) /* Return intervals from "other" species that will map to the current window. * These intervals are NOT YET MAPPED to the current assembly. */ { char *linkFileName = bigChainGetLinkFile(quickLiftFile); int maxGapBefore = 0; int maxGapAfter = 0; struct chain *chain, *chainList = chainLoadIdRangeHub(NULL, quickLiftFile, linkFileName, chrom, start, end, -1); struct lm *lm = lmInit(0); struct bigBedInterval *bbList = NULL, *bb; for(chain = chainList; chain; chain = chain->next) { @@ -698,30 +699,110 @@ } for (i = psl->blockCount - 1; i >= 0; i--) { AllocVar(b); b->tStart = psl->tStarts[i]; b->tEnd = b->tStart + psl->blockSizes[i]; b->qStart = psl->qStarts[i]; b->qEnd = b->qStart + psl->blockSizes[i]; slAddHead(&blockList, b); } chain->blockList = blockList; return chain; } +struct mafAli *quickLiftMafs(struct hash *chainHash, struct mafAli *mafList, + char *sourceDb, char *refSrc, int refSrcSize) +// Map MAF blocks from the other assembly onto our current reference. +// +// A MAF block has to be one contiguous run on its first row, and the lift does not keep +// the reference contiguous: where the reference assembly has lost bases the columns for +// them go away, and where it has gained bases the alignment says nothing about them. So a +// block is cut at every chain block boundary. Inside one chain block the two assemblies +// run in step, which is what lets the columns be carried over untouched: only the first +// row's coordinates change, and mafSubset does the rest of the arithmetic. +// +// refSrc is the name the browser expects on the reference row, "<db>.<chrom>", with no hub +// prefix. Blocks whose reference does not map are dropped. +{ +struct mafAli *outList = NULL; +struct mafAli *maf, *nextMaf; + +for (maf = mafList; maf != NULL; maf = nextMaf) + { + nextMaf = maf->next; + maf->next = NULL; + + // The first row of a MAF is its reference, and a reference row is always forward. + struct mafComp *ref = maf->components; + if ((ref == NULL) || (ref->strand != '+') || (ref->size <= 0)) + { + mafAliFree(&maf); + continue; + } + + // the chains are keyed on the sequence name in the other assembly + char srcBuf[1024]; + safecpy(srcBuf, sizeof srcBuf, ref->src); + char *srcChrom = mafSplitSrcGetChrom(srcBuf, sourceDb); + int refStart = ref->start; + int refEnd = refStart + ref->size; + + struct chain *chain = liftOverChainForRange(chainHash, srcChrom, refStart, refEnd); + if (chain == NULL) + { + mafAliFree(&maf); + continue; + } + + struct cBlock *cb; + for (cb = chain->blockList; cb != NULL; cb = cb->next) + { + int runStart = max(cb->tStart, refStart); + int runEnd = min(cb->tEnd, refEnd); + if (runStart >= runEnd) + continue; + + struct mafAli *sub = mafSubset(maf, ref->src, runStart, runEnd); + if (sub == NULL) + continue; + + int destStart = cb->qStart + (runStart - cb->tStart); + if (chain->qStrand == '-') + { + // The lift turns the block over, so turn every row over with it. A chain + // keeps its query side reverse complemented, so the forward start of the run + // comes from the far end of it. + mafFlipStrand(sub); + destStart = chain->qSize - (cb->qStart + (runEnd - cb->tStart)); + } + + struct mafComp *subRef = sub->components; + freeMem(subRef->src); + subRef->src = cloneString(refSrc); + subRef->srcSize = refSrcSize; + subRef->strand = '+'; + subRef->start = destStart; + slAddHead(&outList, sub); + } + mafAliFree(&maf); + } +slReverse(&outList); +return outList; +} + boolean quickLiftIsOwnChainTrack(struct trackDb *tdb) // TRUE when this is the chain track quickLift builds to show the lift itself. That stanza // carries quickLiftUrl and quickLiftDb like any lifted track, but its data is already in // reference coordinates and must not be lifted a second time. The giveaway is that its // bigDataUrl IS the quickLift chain file. { char *quickLiftFile = trackDbSetting(tdb, "quickLiftUrl"); if (quickLiftFile == NULL) return FALSE; if (startsWithNoCase("bigQuickLiftChain", tdb->type)) return TRUE; char *bigDataUrl = trackDbSetting(tdb, "bigDataUrl"); return (bigDataUrl != NULL) && sameString(bigDataUrl, quickLiftFile);