813fb032f7947b3d5aa6481357481a0c6afdf4c5 hiram Mon Dec 27 14:50:45 2010 -0800 strand mixup on multiple mappings diff --git src/hg/lib/liftOver.c src/hg/lib/liftOver.c index 6da6847..4176d7d 100644 --- src/hg/lib/liftOver.c +++ src/hg/lib/liftOver.c @@ -123,48 +123,49 @@ *pEnd = subChain->qEnd; } else { *pStart = subChain->qSize - subChain->qEnd; *pEnd = subChain->qSize - subChain->qStart; } *retSubChain = subChain; *retChainToFree = freeChain; return ok; } static char *remapRange(struct hash *chainHash, double minRatio, int minSizeT, int minSizeQ, int minChainSizeT, int minChainSizeQ, - char *chrom, int s, int e, char strand, + char *chrom, int s, int e, char qStrand, int thickStart, int thickEnd, bool useThick, double minMatch, char *regionName, char *db, char *chainTableName, struct bed **bedRet, struct bed **unmappedBedRet) /* Remap a range through chain hash. If all is well return NULL * and results in a BED (or a list of BED's, if regionName is set (-multiple). * Otherwise return a string describing the problem. * note: minSizeT is currently NOT implemented. feel free to add it. * (its not e-s, its the boundaries of what maps of chrom:s-e to Q) */ { struct binElement *list = findRange(chainHash, chrom, s, e), *el; struct chain *chainsHit = NULL, *chainsPartial = NULL, *chainsMissed = NULL, *chain; struct bed *bedList = NULL, *unmappedBedList = NULL; struct bed *bed = NULL; +char strand = qStrand; /* initialize for single region case */ int start = s, end = e; double minMatchSize = minMatch * (end - start); int intersectSize; int tStart; bool multiple = (regionName != NULL); verbose(2, "%s:%d-%d", chrom, s, e); verbose(2, multiple ? "\t%s\n": "\n", regionName); for (el = list; el != NULL; el = el->next) { chain = el->val; if (multiple) { if (chain->qEnd - chain->qStart < minChainSizeQ || @@ -226,31 +227,31 @@ { /* use full chain, not the possibly truncated chain * from the net */ struct chain *next = chain->next; chain = chainLoadIdRange(db, chainTableName, chrom, s, e, chain->id); chain->next = next; verbose(3,"chain from db %s:%d %s:%d-%d %c (%d)\n", chain->tName, chain->tStart, chain->qName, chain->qStart, chain->qEnd, chain->qStrand, chain->id); } } if (!mapThroughChain(chain, minRatio, &start, &end, &subChain, &toFree)) errAbort("Chain mapping error: %s:%d-%d\n", chain->qName, start, end); if (chain->qStrand == '-') - strand = otherStrand(strand); + strand = otherStrand(qStrand); if (useThick) { struct chain *subChain2 = NULL; struct chain *toFree2 = NULL; if (!mapThroughChain(chain, minRatio, &thickStart, &thickEnd, &subChain2, &toFree2)) thickStart = thickEnd = start; chainFree(&toFree2); } verbose(3, "mapped %s:%d-%d\n", chain->qName, start, end); verbose(4, "Mapped! Target:\t%s\t%d\t%d\t%c\tQuery:\t%s\t%d\t%d\t%c\n", chain->tName, subChain->tStart, subChain->tEnd, strand, chain->qName, subChain->qStart, subChain->qEnd, chain->qStrand); if (multiple && end - start < minSizeQ) {