X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;f=src%2Fjalview%2Futil%2FMapList.java;h=198066d7bcdcee7809445cc71934cf7aed08e55a;hb=1e507c915398542fe6e206e7a6f2620283243f1e;hp=c32447cda27f56851fb34fb69803e56f9db91f5e;hpb=3a20b7a1383c3b62cb91bc21e8aef8bbcd8c21a8;p=jalview.git diff --git a/src/jalview/util/MapList.java b/src/jalview/util/MapList.java index c32447c..198066d 100644 --- a/src/jalview/util/MapList.java +++ b/src/jalview/util/MapList.java @@ -22,6 +22,7 @@ package jalview.util; import java.util.ArrayList; import java.util.Arrays; +import java.util.BitSet; import java.util.List; /** @@ -30,8 +31,6 @@ import java.util.List; * * Use at your own risk! * - * TODO: efficient implementation of private posMap method - * * TODO: test/ensure that sense of from and to ratio start position is conserved * (codon start position recovery) */ @@ -77,8 +76,8 @@ public class MapList */ public MapList() { - fromShifts = new ArrayList(); - toShifts = new ArrayList(); + fromShifts = new ArrayList<>(); + toShifts = new ArrayList<>(); } /** @@ -88,8 +87,6 @@ public class MapList @Override public boolean equals(Object o) { - // TODO should also override hashCode to ensure equal objects have equal - // hashcodes if (o == null || !(o instanceof MapList)) { return false; @@ -105,10 +102,31 @@ public class MapList { return false; } - return Arrays - .deepEquals(fromShifts.toArray(), obj.fromShifts.toArray()) - && Arrays - .deepEquals(toShifts.toArray(), obj.toShifts.toArray()); + return Arrays.deepEquals(fromShifts.toArray(), obj.fromShifts.toArray()) + && Arrays.deepEquals(toShifts.toArray(), + obj.toShifts.toArray()); + } + + /** + * Returns a hashcode made from the fromRatio, toRatio, and from/to ranges + */ + @Override + public int hashCode() + { + int hashCode = 31 * fromRatio; + hashCode = 31 * hashCode + toRatio; + for (int[] shift : fromShifts) + { + hashCode = 31 * hashCode + shift[0]; + hashCode = 31 * hashCode + shift[1]; + } + for (int[] shift : toShifts) + { + hashCode = 31 * hashCode + shift[0]; + hashCode = 31 * hashCode + shift[1]; + } + + return hashCode; } /** @@ -190,8 +208,7 @@ public class MapList /** * Constructor given from and to ranges as [start1, end1, start2, end2,...]. - * If any end is equal to the next start, the ranges will be merged. There is - * no validation check that the ranges do not overlap each other. + * There is no validation check that the ranges do not overlap each other. * * @param from * contiguous regions as [start1, end1, start2, end2, ...] @@ -209,46 +226,25 @@ public class MapList this.toRatio = toRatio; fromLowest = Integer.MAX_VALUE; fromHighest = Integer.MIN_VALUE; - int added = 0; for (int i = 0; i < from.length; i += 2) { /* * note lowest and highest values - bearing in mind the - * direction may be revesed + * direction may be reversed */ fromLowest = Math.min(fromLowest, Math.min(from[i], from[i + 1])); fromHighest = Math.max(fromHighest, Math.max(from[i], from[i + 1])); - if (added > 0 && from[i] == fromShifts.get(added - 1)[1]) - { - /* - * this range starts where the last ended - just extend it - */ - fromShifts.get(added - 1)[1] = from[i + 1]; - } - else - { - fromShifts.add(new int[] { from[i], from[i + 1] }); - added++; - } + fromShifts.add(new int[] { from[i], from[i + 1] }); } toLowest = Integer.MAX_VALUE; toHighest = Integer.MIN_VALUE; - added = 0; for (int i = 0; i < to.length; i += 2) { toLowest = Math.min(toLowest, Math.min(to[i], to[i + 1])); toHighest = Math.max(toHighest, Math.max(to[i], to[i + 1])); - if (added > 0 && to[i] == toShifts.get(added - 1)[1]) - { - toShifts.get(added - 1)[1] = to[i + 1]; - } - else - { - toShifts.add(new int[] { to[i], to[i + 1] }); - added++; - } + toShifts.add(new int[] { to[i], to[i + 1] }); } } @@ -297,6 +293,8 @@ public class MapList int toRatio) { this(); + fromRange = coalesceRanges(fromRange); + toRange = coalesceRanges(toRange); this.fromShifts = fromRange; this.toShifts = toRange; this.fromRatio = fromRatio; @@ -306,6 +304,12 @@ public class MapList fromHighest = Integer.MIN_VALUE; for (int[] range : fromRange) { + if (range.length != 2) + { + // throw new IllegalArgumentException(range); + System.err.println("Invalid format for fromRange " + + Arrays.toString(range) + " may cause errors"); + } fromLowest = Math.min(fromLowest, Math.min(range[0], range[1])); fromHighest = Math.max(fromHighest, Math.max(range[0], range[1])); } @@ -314,12 +318,86 @@ public class MapList toHighest = Integer.MIN_VALUE; for (int[] range : toRange) { + if (range.length != 2) + { + // throw new IllegalArgumentException(range); + System.err.println("Invalid format for toRange " + + Arrays.toString(range) + " may cause errors"); + } toLowest = Math.min(toLowest, Math.min(range[0], range[1])); toHighest = Math.max(toHighest, Math.max(range[0], range[1])); } } /** + * Consolidates a list of ranges so that any contiguous ranges are merged. + * This assumes the ranges are already in start order (does not sort them). + *

+ * The main use case for this method is when mapping cDNA sequence to its + * protein product, based on CDS feature ranges which derive from spliced + * exons, but are contiguous on the cDNA sequence. For example + * + * 

+   *   CDS 1-20  // from exon1
+   *   CDS 21-35 // from exon2
+   *   CDS 36-71 // from exon3
+   * 'coalesce' to range 1-71
+   *
+ * + * @param ranges + * @return the same list (if unchanged), else a new merged list, leaving the + * input list unchanged + */ + public static List coalesceRanges(final List ranges) + { + if (ranges == null || ranges.size() < 2) + { + return ranges; + } + + boolean changed = false; + List merged = new ArrayList<>(); + int[] lastRange = ranges.get(0); + int lastDirection = lastRange[1] >= lastRange[0] ? 1 : -1; + lastRange = new int[] { lastRange[0], lastRange[1] }; + merged.add(lastRange); + boolean first = true; + + for (final int[] range : ranges) + { + if (first) + { + first = false; + continue; + } + + int direction = range[1] >= range[0] ? 1 : -1; + + /* + * if next range is in the same direction as last and contiguous, + * just update the end position of the last range + */ + boolean sameDirection = range[1] == range[0] + || direction == lastDirection; + boolean extending = range[0] == lastRange[1] + lastDirection; + if (sameDirection && extending) + { + lastRange[1] = range[1]; + changed = true; + } + else + { + lastRange = new int[] { range[0], range[1] }; + merged.add(lastRange); + // careful: merging [5, 5] after [7, 6] should keep negative direction + lastDirection = (range[1] == range[0]) ? lastDirection : direction; + } + } + + return changed ? merged : ranges; + } + + /** * get all mapped positions from 'from' to 'to' * * @return int[][] { int[] { fromStart, fromFinish, toStart, toFinish }, int @@ -327,7 +405,7 @@ public class MapList */ protected int[][] makeFromMap() { - // TODO not used - remove?? + // TODO only used for test - remove?? return posMap(fromShifts, fromRatio, toShifts, toRatio); } @@ -338,7 +416,7 @@ public class MapList */ protected int[][] makeToMap() { - // TODO not used - remove?? + // TODO only used for test - remove?? return posMap(toShifts, toRatio, fromShifts, fromRatio); } @@ -352,7 +430,7 @@ public class MapList private int[][] posMap(List shiftTo, int ratio, List shiftFrom, int toRatio) { - // TODO not used - remove?? + // TODO only used for test - remove?? int iv = 0, ivSize = shiftTo.size(); if (iv >= ivSize) { @@ -483,37 +561,51 @@ public class MapList List shiftFrom, int toRatio) { // TODO: javadoc; tests - int[] fromCount = countPos(shiftTo, pos); + int[] fromCount = countPositions(shiftTo, pos); if (fromCount == null) { return null; } int fromRemainder = (fromCount[0] - 1) % fromRatio; int toCount = 1 + (((fromCount[0] - 1) / fromRatio) * toRatio); - int[] toPos = countToPos(shiftFrom, toCount); + int[] toPos = traverseToPosition(shiftFrom, toCount); if (toPos == null) { - return null; // throw new Error("Bad Mapping!"); + return null; } - // System.out.println(fromCount[0]+" "+fromCount[1]+" "+toCount); return new int[] { toPos[0], fromRemainder, toPos[1] }; } /** - * count how many positions pos is along the series of intervals. + * Counts how many positions pos is along the series of intervals. Returns an + * array of two values: + *
    + *
  • the number of positions traversed (inclusive) to reach {@code pos}
    • + *
  • +1 if the last interval traversed is forward, -1 if in a negative + * direction
    • + *
+ * Returns null if {@code pos} does not lie in any of the given intervals. * - * @param shiftTo + * @param intervals + * a list of start-end intervals * @param pos - * @return number of positions or null if pos is not within intervals + * a position that may lie in one (or more) of the intervals + * @return */ - protected static int[] countPos(List shiftTo, int pos) + protected static int[] countPositions(List intervals, int pos) { - int count = 0, intv[], iv = 0, ivSize = shiftTo.size(); + int count = 0; + int iv = 0; + int ivSize = intervals.size(); + while (iv < ivSize) { - intv = shiftTo.get(iv++); + int[] intv = intervals.get(iv++); if (intv[0] <= intv[1]) { + /* + * forwards interval + */ if (pos >= intv[0] && pos <= intv[1]) { return new int[] { count + pos - intv[0] + 1, +1 }; @@ -525,6 +617,9 @@ public class MapList } else { + /* + * reverse interval + */ if (pos >= intv[1] && pos <= intv[0]) { return new int[] { count + intv[0] - pos + 1, -1 }; @@ -539,79 +634,61 @@ public class MapList } /** - * count out pos positions into a series of intervals and return the position + * Reads through the given intervals until {@code count} positions have been + * traversed, and returns an array consisting of two values: + *
    + *
  • the value at the {@code count'th} position
    • + *
  • +1 if the last interval read is forwards, -1 if reverse direction
    • + *
+ * Returns null if the ranges include less than {@code count} positions, or if + * {@code count < 1}. * - * @param shiftFrom - * @param pos - * @return position pos in interval set + * @param intervals + * a list of [start, end] ranges + * @param count + * the number of positions to traverse + * @return */ - protected static int[] countToPos(List shiftFrom, int pos) + protected static int[] traverseToPosition(List intervals, + final int count) { - int count = 0, diff = 0, iv = 0, ivSize = shiftFrom.size(); - int[] intv = { 0, 0 }; + int traversed = 0; + int ivSize = intervals.size(); + int iv = 0; + + if (count < 1) + { + return null; + } + while (iv < ivSize) { - intv = shiftFrom.get(iv++); - diff = intv[1] - intv[0]; + int[] intv = intervals.get(iv++); + int diff = intv[1] - intv[0]; if (diff >= 0) { - if (pos <= count + 1 + diff) + if (count <= traversed + 1 + diff) { - return new int[] { pos - count - 1 + intv[0], +1 }; + return new int[] { intv[0] + (count - traversed - 1), +1 }; } else { - count += 1 + diff; + traversed += 1 + diff; } } else { - if (pos <= count + 1 - diff) + if (count <= traversed + 1 - diff) { - return new int[] { intv[0] - (pos - count - 1), -1 }; + return new int[] { intv[0] - (count - traversed - 1), -1 }; } else { - count += 1 - diff; + traversed += 1 - diff; } } } - return null;// (diff<0) ? (intv[1]-1) : (intv[0]+1); - } - - /** - * find series of intervals mapping from start-end in the From map. - * - * @param start - * position mapped 'to' - * @param end - * position mapped 'to' - * @return series of [start, end] ranges in sequence mapped 'from' - */ - public int[] locateInFrom(int start, int end) - { - // inefficient implementation - int fromStart[] = shiftTo(start); - // needs to be inclusive of end of symbol position - int fromEnd[] = shiftTo(end); - - return getIntervals(fromShifts, fromStart, fromEnd, fromRatio); - } - - /** - * find series of intervals mapping from start-end in the to map. - * - * @param start - * position mapped 'from' - * @param end - * position mapped 'from' - * @return series of [start, end] ranges in sequence mapped 'to' - */ - public int[] locateInTo(int start, int end) - { - int toStart[] = shiftFrom(start); - int toEnd[] = shiftFrom(end); - return getIntervals(toShifts, toStart, toEnd, toRatio); + return null; } /** @@ -709,7 +786,7 @@ public class MapList { return null; } - List ranges = new ArrayList(); + List ranges = new ArrayList<>(); if (fs <= fe) { intv = fs; @@ -790,7 +867,6 @@ public class MapList */ public int getToPosition(int mpos) { - // TODO not used - remove?? int[] mp = shiftTo(mpos); if (mp != null) { @@ -800,53 +876,6 @@ public class MapList } /** - * get range of positions in To frame for the mpos word in From - * - * @param mpos - * position in From - * @return null or int[] first position in To for mpos, last position in to - * for Mpos - */ - public int[] getToWord(int mpos) - { - int[] mp = shiftTo(mpos); - if (mp != null) - { - return new int[] { mp[0], mp[0] + mp[2] * (getFromRatio() - 1) }; - } - return null; - } - - /** - * get From position in the associated reference frame for position pos in the - * associated sequence. - * - * @param pos - * @return - */ - public int getMappedPosition(int pos) - { - // TODO not used - remove?? - int[] mp = shiftFrom(pos); - if (mp != null) - { - return mp[0]; - } - return pos; - } - - public int[] getMappedWord(int pos) - { - // TODO not used - remove?? - int[] mp = shiftFrom(pos); - if (mp != null) - { - return new int[] { mp[0], mp[0] + mp[2] * (getToRatio() - 1) }; - } - return null; - } - - /** * * @return a MapList whose From range is this maplist's To Range, and vice * versa @@ -858,31 +887,6 @@ public class MapList } /** - * test for containment rather than equivalence to another mapping - * - * @param map - * to be tested for containment - * @return true if local or mapped range map contains or is contained by this - * mapping - */ - public boolean containsEither(boolean local, MapList map) - { - // TODO not used - remove? - if (local) - { - return ((getFromLowest() >= map.getFromLowest() && getFromHighest() <= map - .getFromHighest()) || (getFromLowest() <= map.getFromLowest() && getFromHighest() >= map - .getFromHighest())); - } - else - { - return ((getToLowest() >= map.getToLowest() && getToHighest() <= map - .getToHighest()) || (getToLowest() <= map.getToLowest() && getToHighest() >= map - .getToHighest())); - } - } - - /** * String representation - for debugging, not guaranteed not to change */ @Override @@ -894,7 +898,9 @@ public class MapList { sb.append(" ").append(Arrays.toString(shift)); } - sb.append(" ] To ["); + sb.append(" ] "); + sb.append(fromRatio).append(":").append(toRatio); + sb.append(" to ["); for (int[] shift : toShifts) { sb.append(" ").append(Arrays.toString(shift)); @@ -912,6 +918,10 @@ public class MapList */ public void addMapList(MapList map) { + if (this.equals(map)) + { + return; + } this.fromLowest = Math.min(fromLowest, map.fromLowest); this.toLowest = Math.min(toLowest, map.toLowest); this.fromHighest = Math.max(fromHighest, map.fromHighest); @@ -927,7 +937,14 @@ public class MapList } } - public static void addRange(int[] range, List addTo) + /** + * Adds the given range to a list of ranges. If the new range just extends + * existing ranges, the current endpoint is updated instead. + * + * @param range + * @param addTo + */ + static void addRange(int[] range, List addTo) { /* * list is empty - add to it! @@ -985,8 +1002,33 @@ public class MapList */ public boolean isFromForwardStrand() { + return isForwardStrand(getFromRanges()); + } + + /** + * Returns true if mapping is to forward strand, false if to reverse strand. + * Result is just based on the first 'to' range that is not a single position. + * Default is true unless proven to be false. Behaviour is not well defined if + * the mapping has a mixture of forward and reverse ranges. + * + * @return + */ + public boolean isToForwardStrand() + { + return isForwardStrand(getToRanges()); + } + + /** + * A helper method that returns true unless at least one range has start > + * end. Behaviour is undefined for a mixture of forward and reverse ranges. + * + * @param ranges + * @return + */ + private boolean isForwardStrand(List ranges) + { boolean forwardStrand = true; - for (int[] range : getFromRanges()) + for (int[] range : ranges) { if (range[1] > range[0]) { @@ -1000,4 +1042,376 @@ public class MapList } return forwardStrand; } + + /** + * + * @return true if from, or to is a three to 1 mapping + */ + public boolean isTripletMap() + { + return (toRatio == 3 && fromRatio == 1) + || (fromRatio == 3 && toRatio == 1); + } + + /** + * Returns a map which is the composite of this one and the input map. That + * is, the output map has the fromRanges of this map, and its toRanges are the + * toRanges of this map as transformed by the input map. + *

+ * Returns null if the mappings cannot be traversed (not all toRanges of this + * map correspond to fromRanges of the input), or if this.toRatio does not + * match map.fromRatio. + * + * 

+   * Example 1:
+   *    this:   from [1-100] to [501-600]
+   *    input:  from [10-40] to [60-90]
+   *    output: from [10-40] to [560-590]
+   * Example 2 ('reverse strand exons'):
+   *    this:   from [1-100] to [2000-1951], [1000-951] // transcript to loci
+   *    input:  from [1-50]  to [41-90] // CDS to transcript
+   *    output: from [10-40] to [1960-1951], [1000-971] // CDS to gene loci
+   *
+ * + * @param map + * @return + */ + public MapList traverse(MapList map) + { + if (map == null) + { + return null; + } + + /* + * compound the ratios by this rule: + * A:B with M:N gives A*M:B*N + * reduced by greatest common divisor + * so 1:3 with 3:3 is 3:9 or 1:3 + * 1:3 with 3:1 is 3:3 or 1:1 + * 1:3 with 1:3 is 1:9 + * 2:5 with 3:7 is 6:35 + */ + int outFromRatio = getFromRatio() * map.getFromRatio(); + int outToRatio = getToRatio() * map.getToRatio(); + int gcd = MathUtils.gcd(outFromRatio, outToRatio); + outFromRatio /= gcd; + outToRatio /= gcd; + + List toRanges = new ArrayList<>(); + for (int[] range : getToRanges()) + { + int fromLength = Math.abs(range[1] - range[0]) + 1; + int[] transferred = map.locateInTo(range[0], range[1]); + if (transferred == null || transferred.length % 2 != 0) + { + return null; + } + + /* + * convert [start1, end1, start2, end2, ...] + * to [[start1, end1], [start2, end2], ...] + */ + int toLength = 0; + for (int i = 0; i < transferred.length;) + { + toRanges.add(new int[] { transferred[i], transferred[i + 1] }); + toLength += Math.abs(transferred[i + 1] - transferred[i]) + 1; + i += 2; + } + + /* + * check we mapped the full range - if not, abort + */ + if (fromLength * map.getToRatio() != toLength * map.getFromRatio()) + { + return null; + } + } + + return new MapList(getFromRanges(), toRanges, outFromRatio, outToRatio); + } + + /** + * Answers true if the mapping is from one contiguous range to another, else + * false + * + * @return + */ + public boolean isContiguous() + { + return fromShifts.size() == 1 && toShifts.size() == 1; + } + + /** + * Returns the [start1, end1, start2, end2, ...] positions in the 'from' range + * that map to positions between {@code start} and {@code end} in the 'to' + * range. Note that for a reverse strand mapping this will return ranges with + * end < start. Returns null if no mapped positions are found in start-end. + * + * @param start + * @param end + * @return + */ + public int[] locateInFrom(int start, int end) + { + if (end < start) + { + int tmp = end; + end = start; + start = tmp; + } + + /* + * traverse toShifts and mark offsets in fromShifts + * of any positions that lie in [start, end] + */ + BitSet offsets = getMappedOffsetsForPositions(start, end, toShifts, + toRatio, fromRatio); + + /* + * traverse fromShifts and collect positions at the marked offsets + */ + List mapped = getPositionsForOffsets(fromShifts, offsets); + + // TODO: or just return the List and adjust calling code to match + return mapped.isEmpty() ? null : MappingUtils.rangeListToArray(mapped); + } + + /** + * Returns the [start1, end1, start2, end2, ...] positions in the 'to' range + * that map to positions between {@code start} and {@code end} in the 'from' + * range. Note that for a reverse strand mapping this will return ranges with + * end < start. Returns null if no mapped positions are found in start-end. + * + * @param start + * @param end + * @return + */ + public int[] locateInTo(int start, int end) + { + if (end < start) + { + int tmp = end; + end = start; + start = tmp; + } + + /* + * traverse fromShifts and mark offsets in toShifts + * of any positions that lie in [start, end] + */ + BitSet offsets = getMappedOffsetsForPositions(start, end, fromShifts, + fromRatio, toRatio); + + /* + * traverse toShifts and collect positions at the marked offsets + */ + List mapped = getPositionsForOffsets(toShifts, offsets); + + return mapped.isEmpty() ? null : MappingUtils.rangeListToArray(mapped); + } + + /** + * Scans the list of {@code ranges} for any values (positions) that lie + * between start and end (inclusive), and records the offsets from + * the start of the list as a BitSet. The offset positions are converted to + * corresponding words in blocks of {@code wordLength2}. + * + * 
+   * For example:
+   * 1:1 (e.g. gene to CDS):
+   * ranges { [10-20], [31-40] }, wordLengthFrom = wordLength 2 = 1
+   *   for start = 1, end = 9, returns a BitSet with no bits set
+   *   for start = 1, end = 11, returns a BitSet with bits 0-1 set
+   *   for start = 15, end = 35, returns a BitSet with bits 5-15 set
+   * 1:3 (peptide to codon):
+   * ranges { [1-200] }, wordLengthFrom = 1, wordLength 2 = 3
+   *   for start = 9, end = 9, returns a BitSet with bits 24-26 set
+   * 3:1 (codon to peptide):
+   * ranges { [101-150], [171-180] }, wordLengthFrom = 3, wordLength 2 = 1
+   *   for start = 101, end = 102 (partial first codon), returns a BitSet with bit 0 set
+   *   for start = 150, end = 171 (partial 17th codon), returns a BitSet with bit 16 set
+   * 3:1 (circular DNA to peptide):
+   * ranges { [101-150], [21-30] }, wordLengthFrom = 3, wordLength 2 = 1
+   *   for start = 24, end = 40 (spans codons 18-20), returns a BitSet with bits 17-19 set
+   *
+ * + * @param start + * @param end + * @param ranges + * @param wordLengthFrom + * @param wordLengthTo + * @return + */ + protected final static BitSet getMappedOffsetsForPositions(int start, + int end, List ranges, int wordLengthFrom, int wordLengthTo) + { + BitSet overlaps = new BitSet(); + int offset = 0; + final int s1 = ranges.size(); + for (int i = 0; i < s1; i++) + { + int[] range = ranges.get(i); + final int offset1 = offset; + int overlapStartOffset = -1; + int overlapEndOffset = -1; + + if (range[1] >= range[0]) + { + /* + * forward direction range + */ + if (start <= range[1] && end >= range[0]) + { + /* + * overlap + */ + int overlapStart = Math.max(start, range[0]); + overlapStartOffset = offset1 + overlapStart - range[0]; + int overlapEnd = Math.min(end, range[1]); + overlapEndOffset = offset1 + overlapEnd - range[0]; + } + } + else + { + /* + * reverse direction range + */ + if (start <= range[0] && end >= range[1]) + { + /* + * overlap + */ + int overlapStart = Math.max(start, range[1]); + int overlapEnd = Math.min(end, range[0]); + overlapStartOffset = offset1 + range[0] - overlapEnd; + overlapEndOffset = offset1 + range[0] - overlapStart; + } + } + + if (overlapStartOffset > -1) + { + /* + * found an overlap + */ + if (wordLengthFrom != wordLengthTo) + { + /* + * convert any overlap found to whole words in the target range + * (e.g. treat any partial codon overlap as if the whole codon) + */ + overlapStartOffset -= overlapStartOffset % wordLengthFrom; + overlapStartOffset = overlapStartOffset / wordLengthFrom + * wordLengthTo; + + /* + * similar calculation for range end, adding + * (wordLength2 - 1) for end of mapped word + */ + overlapEndOffset -= overlapEndOffset % wordLengthFrom; + overlapEndOffset = overlapEndOffset / wordLengthFrom + * wordLengthTo; + overlapEndOffset += wordLengthTo - 1; + } + overlaps.set(overlapStartOffset, overlapEndOffset + 1); + } + offset += 1 + Math.abs(range[1] - range[0]); + } + return overlaps; + } + + /** + * Returns a (possibly empty) list of the [start-end] values (positions) at + * offsets in the {@code ranges} list that are marked by 'on' bits in the + * {@code offsets} bitset. + * + * @param ranges + * @param offsets + * @return + */ + protected final static List getPositionsForOffsets( + List ranges, BitSet offsets) + { + List mapped = new ArrayList<>(); + if (offsets.isEmpty()) + { + return mapped; + } + + /* + * count of positions preceding ranges[i] + */ + int traversed = 0; + + /* + * for each [from-to] range in ranges: + * - find subranges (if any) at marked offsets + * - add the start-end values at the marked positions + */ + final int toAdd = offsets.cardinality(); + int added = 0; + final int s2 = ranges.size(); + for (int i = 0; added < toAdd && i < s2; i++) + { + int[] range = ranges.get(i); + added += addOffsetPositions(mapped, traversed, range, offsets); + traversed += Math.abs(range[1] - range[0]) + 1; + } + return mapped; + } + + /** + * Helper method that adds any start-end subranges of {@code range} that are + * at offsets in {@code range} marked by set bits in overlaps. + * {@code mapOffset} is added to {@code range} offset positions. Returns the + * count of positions added. + * + * @param mapped + * @param mapOffset + * @param range + * @param overlaps + * @return + */ + final static int addOffsetPositions(List mapped, + final int mapOffset, final int[] range, final BitSet overlaps) + { + final int rangeLength = 1 + Math.abs(range[1] - range[0]); + final int step = range[1] < range[0] ? -1 : 1; + int offsetStart = 0; // offset into range + int added = 0; + + while (offsetStart < rangeLength) + { + /* + * find the start of the next marked overlap offset; + * if there is none, or it is beyond range, then finished + */ + int overlapStart = overlaps.nextSetBit(mapOffset + offsetStart); + if (overlapStart == -1 || overlapStart - mapOffset >= rangeLength) + { + /* + * no more overlaps, or no more within range[] + */ + return added; + } + overlapStart -= mapOffset; + + /* + * end of the overlap range is just before the next clear bit; + * restrict it to end of range if necessary; + * note we may add a reverse strand range here (end < start) + */ + int overlapEnd = overlaps.nextClearBit(mapOffset + overlapStart + 1); + overlapEnd = (overlapEnd == -1) ? rangeLength - 1 + : Math.min(rangeLength - 1, overlapEnd - mapOffset - 1); + int startPosition = range[0] + step * overlapStart; + int endPosition = range[0] + step * overlapEnd; + mapped.add(new int[] { startPosition, endPosition }); + offsetStart = overlapEnd + 1; + added += Math.abs(endPosition - startPosition) + 1; + } + + return added; + } }