Here is another approach in Python making few assumptions and using very little memory:

in_file = "queries.txt"
coverage = set()

with open(in_file) as f:
    for line in f:
        try:
            query, v1, v2 = line.strip().split("\t")
            v_min, v_max = sorted([int(v1), int(v2)])
        except:
            continue
        coverage.update(range(v_min, v_max + 1)) # as suggested by brentp

print "Reference coverage:", len(coverage)

Here is what I get when I run it:

time python golf-course.py 
Reference coverage: 81

real    0m0.026s
user    0m0.016s
sys     0m0.008ss

It should also be pretty 'Pythonic' :)

Cheers

Here's a quick solution, and should run pretty fast. you can easily make it work per reference as well, though here i just assume you have only a single reference as you suggest.

import numpy as np
import sys
# since we dont know the max len beforehand,
# just use a big number and then chop it later.
MAX_LEN = 3e8

coverage = np.zeros((MAX_LEN,), dtype=np.uint8)

for i, line in enumerate(open(sys.argv[1])):
    # logic for 1st line here.
    if i == 0: continue
    line = line.rstrip().split()
    ref = line[0]
    start, end = sorted(map(int, line[1:]))

    coverage[start - 1:end] |= 1

print len(coverage[coverage == 1])

and run from the commandline with your alignment file as the 1st argument.

If you fancy perl, and have faith in the bazaar or are generally lazy, then install Set::IntSpan (a module I find generally useful for genomic interval hacking in Perl)

allowing you to write this one-liner:

perl -MSet::IntSpan -an -e 'BEGIN{$c=Set::IntSpan->new} next if $. == 1; $c += [[sort @F[1,2]]]; END{print $c->size}'

My C++ solution (does NOT need to store all the data in memory. You just need to keep the contigs.

Compiling/Testing

all:
    g++ coverage.cpp
    echo -e "query\tstart\tstop\nquery1\t100\t120\nquery2\t50\t99\nquery3\t130\t110" | ./a.out
    time mysql --user=genome --host=genome-mysql.cse.ucsc.edu -A -D hg18  \
      -e 'select name as "query",if(strand="-",txEnd+1,txStart) as "start",if(strand="-",txStart,txEnd+1) as "stop" from knownGene where chrom="chr1"' | ./a.out

Output:

g++ coverage.cpp
echo -e "query\tstart\tstop\nquery1\t100\t120\nquery2\t50\t99\nquery3\t130\t110" | ./a.out
coverage=81 pb
time mysql --user=genome --host=genome-mysql.cse.ucsc.edu -A -D hg18  \
          -e 'select name as "query",if(strand="-",txEnd+1,txStart) as "start",if(strand="-",txStart,txEnd+1) as "stop" from knownGene where chrom="chr1"' | ./a.out
coverage=114657440 pb

real    0m3.403s
user    0m0.076s
sys     0m0.008s

I made it in C++ because it is my "everyday life" language ;)

The idea is :

1. read queries and swap if start > stop
2. sort queries based on "start"
3. merge queries if they overlap
4. sum the coverage of each query

#include <vector>
#include <iostream>
#include <algorithm>

using namespace std;

// Class to store start and stop positions
class Query
{
  public :
    int start;
    int stop; 
    Query(int i, int j):start(i),stop(j){
      if (start > stop)
        swap();  
    }

    void swap(){
      int i(start);
      start = stop;
      stop = i;
    }
};

// function to compare start positions of two queries
bool CompQuery (Query i,Query j){return (i.start < j.start);}

int main(int argc, char * argv)
{
  // Create query tab with start and stop positions
  vector<Query> tab;
  vector<Query>::iterator it;
  tab.push_back(Query(100,120));
  tab.push_back(Query(130,110));
  tab.push_back(Query(50,99));
  tab.push_back(Query(165,140));

  // sort query tab by increasing start positions
  sort(tab.begin(), tab.end(), CompQuery);

  // Merge Queries if they cover adjacent segments
  // example : 50->99 and 100->120 are merged in 50->120
  it = tab.begin() + 1;
  int cover = 0;
  while(it != tab.end()){
    if ((*it).start <= (*(it-1)).stop + 1){
      if ((*it).stop > (*(it-1)).stop){
        (*(it-1)).stop = (*it).stop;
      }
      it = tab.erase(it);
    } else {
      cout << (*it).start<<","<<(*it).stop<< "\n";
      cover += (*(it-1)).stop - (*(it-1)).start + 1;
      it++;
    }
  }
  cover += (*(it-1)).stop -(*(it-1)).start + 1;

  // Print Coverage
  cout << "\nCoverage = "<<cover<<"\n";
  return 0;
}

Here is a version I wrote using R and Bioconductor's IRanges

#!/usr/bin/Rscript
library("IRanges",warn.conflicts=FALSE)
args<-commandArgs(TRUE)
mytab<-read.table(args[1],header=TRUE)
mytab$starts<-apply(mytab[,c(2,3)],1,min)
mytab$stops<-apply(mytab[,c(2,3)],1,max)
myrange<-IRanges(start=mytab$starts,end=mytab$stops)
sum(width(reduce(myrange)))

I swear I've written this code before somewhere, I'll see if I can dig it up later.

In the meantime, try this:

• In one quick pass, flip any start/stop values and sort the file:

[?]

• Now, you step through the outfile, merging segments.
  

Rough algorithm:

• coverage = 0
• loop through the reads
• store the first segment in a buffer
• check to see if the current segment's start is <= old segment stop
  • if so, update the oldSeg's stop position
  • if not, add the length of the old segment to coverage, then set the oldSeg = currSeg
• at the end of the loop, clear the last segment out of the buffer and add it's length to coverage.
• coverage/referenceSize = coverage percentage

Just for the collection, the code in MATLAB (without using any library/toolbox):

q = dlmread('queries.txt','\t',1,1); %# read the file
q = sort(q,2); %# start < end
qidx = ones(size(q)); %# start = 1
qidx(:,2) = -1; %# end = -1
[qsorted, si] = sort(q(:));
sidx = qidx(si);
gap = find(cumsum(sidx)==0); %# find start indices of gaps
boundaryidx = reshape(unique([1; gap; gap(1:end-1)+1]),2,[]);
coverage = sum(diff(qsorted(boundaryidx))+1);
disp(coverage)