Factor out Interval/ContigInterval

src/BigBed.js

@@ -11,7 +11,10 @@ var Q = require('q'),
 
 
 var ReadableView = require('./ReadableView'),
-    RemoteFile = require('./RemoteFile');
+    RemoteFile = require('./RemoteFile'),
+    Interval = require('./Interval'),
+    ContigInterval = require('./ContigInterval'),
+    utils = require('./utils.js');
 
 function typeAtOffset(typeName, offsetFieldName) {
   return jBinary.Template({
@@ -189,64 +192,34 @@ function getContigId(contigMap, contig) {
   return contigMap[contig] || contigMap['chr' + contig] || null;
 }
 
-function intersectIntervals(intervals: Array<[number, number]>): [number, number] {
-  if (!intervals.length) {
-    throw 'Tried to intersect zero intervals';
-  }
-  var result = intervals[0];
-  intervals.slice(1).forEach(function([a, b]) {
-    result[0] = Math.min(a, result[0]);
-    result[1] = Math.max(b, result[1]);
-  });
-  return result;
-}
-
-// TODO: factor out into module
-var lessOrEqual = function(c1, p1, c2, p2) {
-  return c1 < c2 || (c1 == c2 && p1 <= p2);
-};
-var contains = function(startContig, startPos, endContig, endPos, contig, pos) {
-  return lessOrEqual(startContig, startPos, contig, pos) &&
-         lessOrEqual(contig, pos, endContig, endPos);
-};
-
-var overlaps = function(startContig, startBase, endContig, endBase, contig, start, stop) {
-  return contains(startContig, startBase, endContig, endBase, contig, start) ||
-         contains(startContig, startBase, endContig, endBase, contig, stop);
-};
-
-// Get a byte range in the file containing a superset of the interval.
-function findByteRange(twoBitHeader, cirTree, contigIx: number, start: number, stop: number): ?[number, number] {
-
+// Get all blocks in the file containing features which intersect with contigRange.
+function findOverlappingBlocks(twoBitHeader, cirTree, contigRange) {
   // Do a recursive search through the index tree
-  var matchingIntervals = [];
+  var matchingBlocks = [];
+  var tupleRange = [[contigRange.contig, contigRange.start()],
+                    [contigRange.contig, contigRange.stop()]];
   var find = function(node) {
     if (node.contents) {
       node.contents.forEach(find);
     } else {
-      if (overlaps(node.startChromIx, node.startBase,
-                   node.endChromIx, node.endBase,
-                   contigIx, start, stop)) {
-        matchingIntervals.push(node);
+      var nodeRange = [[node.startChromIx, node.startBase],
+                       [node.endChromIx, node.endBase]];
+      if (utils.tupleRangeOverlaps(nodeRange, tupleRange)) {
+        matchingBlocks.push(node);
       }
     }
   };
   find(cirTree.blocks);
 
-  return matchingIntervals;
-
-  // Intersect the intervals.
-  // XXX UCSC allows discontiguous intervals. When would this ever happen?
-  return intersectIntervals(
-      matchingIntervals.map(n => [+n.offset, n.offset+n.size]));
+  return matchingBlocks;
 }
 
-function extractFeaturesInRange(dataView, dataRange, blocks, contigIx, start, stop) {
+function extractFeaturesInRange(dataView, dataRange, blocks, contigRange) {
   console.log('Fetched ', dataRange);
   var buffer = dataView.buffer;
 
   return _.flatten(blocks.map(block => {
-    var blockOffset = block.offset - dataRange[0],
+    var blockOffset = block.offset - dataRange.start,
         blockLimit = blockOffset + block.size,
         // TODO: where does the +2 come from? (I copied it from dalliance)
         blockBuffer = buffer.slice(blockOffset + 2, blockLimit);
@@ -260,7 +233,10 @@ function extractFeaturesInRange(dataView, dataRange, blocks, contigIx, start, st
     console.log(beds);
 
     beds = beds.filter(function(bed) {
-      return overlaps(bed.chrId, bed.start, bed.chrId, bed.end, contigIx, start, stop);
+      var bedInterval = new ContigInterval(bed.chrId, bed.start, bed.end);
+      var r = contigRange.intersects(bedInterval);
+      console.log('intersect?', contigRange.toString(), bedInterval.toString(), r);
+      return r;
     });
 
     return beds;
@@ -303,15 +279,18 @@ class BigBed {
       if (contigIx === null) {
         throw `Invalid contig ${contig}`;
       }
+      var contigRange = new ContigInterval(contigIx, start, stop);
 
-      var blocks = findByteRange(header, cirTree, contigIx, start, stop);
-      if (!blocks) {
-        return null;  // XXX better to throw?
+      var blocks = findOverlappingBlocks(header, cirTree, contigRange);
+      if (blocks.length == 0) {
+        return [];
       }
       console.log(blocks);
-      var range = intersectIntervals(blocks.map(n => [+n.offset, n.offset+n.size]));
-      return this.remoteFile.getBytes(range[0], range[1] - range[0])
-          .then(dataView => extractFeaturesInRange(dataView, range, blocks, contigIx, start, stop));
+
+      var range = Interval.boundingInterval(
+          blocks.map(n => new Interval(+n.offset, n.offset+n.size)));
+      return this.remoteFile.getBytes(range.start, range.length())
+          .then(dataView => extractFeaturesInRange(dataView, range, blocks, contigRange));
     });
   }
 }

src/ContigInterval.js

@@ -0,0 +1,32 @@
+/* @flow */
+
+var Interval = require('./Interval');
+
+class ContigInterval {
+  contig: string|number;
+  interval: Interval;
+
+  constructor(contig: string|number, start: number, stop: number) {
+    this.contig = contig;
+    this.interval = new Interval(start, stop);
+  }
+
+  // TODO: make these getter methods & switch to Babel.
+  start() {
+    return this.interval.start;
+  }
+  stop() {
+    return this.interval.stop;
+  }
+
+  intersects(other) {
+    return (this.contig == other.contig &&
+            this.interval.intersects(other.interval));
+  }
+
+  toString() {
+    return `${this.contig}:${this.start()}-${this.stop()}`;
+  }
+}
+
+module.exports = ContigInterval;

src/Interval.js

@@ -0,0 +1,64 @@
+/* @flow */
+class Interval {
+  start: number;
+  stop: number;
+
+  // Represents [start, stop] -- both ends are inclusive.
+  // If stop < start, then this is an empty interval.
+  constructor(start: number, stop: number) {
+    this.start = start;
+    this.stop = stop;
+  }
+
+  // TODO: make this a getter method & switch to Babel.
+  length() {
+    return this.stop - this.start + 1;
+  }
+
+  intersect(otherInterval) {
+    return new Interval(Math.max(this.start, otherInterval.start),
+                        Math.min(this.stop, otherInterval.stop));
+  }
+
+  intersects(otherInterval) {
+    return this.start <= otherInterval.stop && otherInterval.start <= this.stop;
+  }
+
+  contains(value) {
+    return value >= this.start && value <= this.stop;
+  }
+
+  clone() {
+    return new Interval(this.start, this.stop);
+  }
+
+  static intersectAll(intervals: Array<Interval>): Interval {
+    if (!intervals.length) {
+      throw 'Tried to intersect zero intervals';
+    }
+    var result = intervals[0].clone();
+    intervals.slice(1).forEach(function({start, stop}) {
+      result.start = Math.max(start, result.start);
+      result.stop = Math.min(stop, result.stop);
+    });
+    return result;
+  }
+
+  static boundingInterval(intervals: Array<Interval>): Interval {
+    if (!intervals.length) {
+      throw 'Tried to intersect zero intervals';
+    }
+    var result = intervals[0].clone();
+    intervals.slice(1).forEach(function({start, stop}) {
+      result.start = Math.min(start, result.start);
+      result.stop = Math.max(stop, result.stop);
+    });
+    return result;
+  }
+
+  toString() {
+    return `[${this.start}, ${this.stop}]`;
+  }
+}
+
+module.exports = Interval;

src/utils.js

@@ -0,0 +1,20 @@
+// Compare two tuples of equal length. Is t1 <= t2?
+function tupleLessOrEqual(t1, t2) {
+  if (t1.length != t2.length) throw 'Comparing non-equal length tuples';
+  for (var i = 0; i < t1.length; i++) {
+    if (t1[i] > t2[i]) {
+      return false;
+    } else if (t1[i] < t2[i]) {
+      return true;
+    }
+  }
+  return true;
+}
+
+// Do two ranges of tuples overlap?
+function tupleRangeOverlaps(tupleRange1, tupleRange2) {
+  return tupleLessOrEqual(tupleRange1[0], tupleRange2[1]) &&
+         tupleLessOrEqual(tupleRange2[0], tupleRange1[1])
+}
+
+module.exports = {tupleRangeOverlaps};

test/ContigInterval-test.js

@@ -0,0 +1,13 @@
+var chai = require('chai');
+var expect = chai.expect;
+
+var ContigInterval = require('../src/ContigInterval');
+
+describe('ContigInterval', function() {
+  it('should determine intersections', function() {
+    var tp53 = new ContigInterval(10, 7512444, 7531643);
+    var other = new ContigInterval(10, 7512444, 7531642);
+
+    expect(tp53.intersects(other)).to.be.true;
+  });
+});

test/Interval-test.js

@@ -0,0 +1,13 @@
+var chai = require('chai');
+var expect = chai.expect;
+
+var Interval = require('../src/Interval');
+
+describe('Interval', function() {
+  it('should determine intersections', function() {
+    var tp53 = new Interval(7512444, 7531643);
+    var other = new Interval(7512444, 7531642);
+
+    expect(tp53.intersects(other)).to.be.true;
+  });
+});