8de67388be9a67cb94ba95615b7643b4ddc385da braney Thu Sep 24 17:48:19 2026 -0700 quickLift: give a lifted maf block the reference's own bases, refs #38249 Inside a chain block the two assemblies run in step but need not agree base for base, and the first row of a lifted block still carried the other assembly's letters. The details page showed hg19's base in the human row of a block lifted onto hg38. quickLiftMafs now reads the reference sequence once over the span of the lifted blocks and writes it into that row, leaving the gaps where they are. The lifted blocks are also sorted by position now. On a chain that turns the alignment over they came back last to first, and the details page listed them in that order. diff --git src/hg/lib/quickLift.c src/hg/lib/quickLift.c index 5516bfe6e79..133d841df4f 100644 --- src/hg/lib/quickLift.c +++ src/hg/lib/quickLift.c @@ -778,43 +778,87 @@ } 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; } +static int mafRefStartCmp(const void *va, const void *vb) +/* Compare two maf blocks by where their first row starts. */ +{ +const struct mafAli *a = *((struct mafAli **)va); +const struct mafAli *b = *((struct mafAli **)vb); +return a->components->start - b->components->start; +} + +static void quickLiftMafRefBases(struct mafAli *mafList, char *refDb, char *refChrom) +/* Put the reference assembly's own bases into the first row of each lifted block. Inside a + * chain block the two assemblies run in step, but they need not agree base for base, and + * the row still carries the other assembly's letters. The gaps stay where they are, since + * the other rows are lined up against them. The sequence is read once for the whole span. */ +{ +struct mafAli *maf; +int spanStart = INT_MAX, spanEnd = 0; + +for (maf = mafList; maf != NULL; maf = maf->next) + { + struct mafComp *ref = maf->components; + spanStart = min(spanStart, ref->start); + spanEnd = max(spanEnd, ref->start + ref->size); + } +if (spanStart >= spanEnd) + return; + +struct dnaSeq *seq = hDnaFromSeq(refDb, refChrom, spanStart, spanEnd, dnaMixed); +for (maf = mafList; maf != NULL; maf = maf->next) + { + struct mafComp *ref = maf->components; + char *base = seq->dna + (ref->start - spanStart); + int left = ref->size; // a maf from a hub can claim fewer bases than its text holds + char *text; + for (text = ref->text; (*text != 0) && (left > 0); text++) + if (*text != '-') + { + *text = *base++; + left--; + } + } +dnaSeqFree(&seq); +} + struct mafAli *quickLiftMafs(struct hash *chainHash, struct mafAli *mafList, - char *sourceDb, char *refSrc, int refSrcSize) + char *sourceDb, char *refDb, char *refChrom, 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, ".", with no hub -// prefix. Blocks whose reference does not map are dropped. +// prefix. The bases on that row are read from refDb, which is the reference's real database +// name and may carry a 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; @@ -858,31 +902,33 @@ 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); } freeMem(srcBuf); mafAliFree(&maf); } -slReverse(&outList); +// a chain that turns the alignment over hands the blocks back last to first +slSort(&outList, mafRefStartCmp); +quickLiftMafRefBases(outList, refDb, refChrom); return outList; } boolean quickLiftIsLifted(struct trackDb *tdb) // TRUE when this track's data comes from another assembly and there is enough to lift it. // Both halves have to be there: the chain file that does the lifting and the assembly the // data came from. A hub can set either one on its own, and half the pair is no use. { return (tdb != NULL) && (trackDbSetting(tdb, "quickLiftUrl") != NULL) && (trackDbSetting(tdb, "quickLiftDb") != NULL); } boolean quickLiftIsOwnChainTrack(struct trackDb *tdb) // TRUE when this is the chain track quickLift builds to show the lift itself. That stanza