Skip to content

Commit

Permalink
fix #159
Browse files Browse the repository at this point in the history
  • Loading branch information
lomereiter committed Aug 8, 2015
1 parent ead562d commit 6d1f5a3
Show file tree
Hide file tree
Showing 3 changed files with 33 additions and 17 deletions.
2 changes: 1 addition & 1 deletion BioD
Submodule BioD updated 1 files
+3 −3 bio/sam/header.d
6 changes: 3 additions & 3 deletions sambamba/depth.d
Original file line number Diff line number Diff line change
Expand Up @@ -653,9 +653,9 @@ final class PerBedRegionPrinter : PerRegionPrinter {
}

override void writeOriginalBedLine(size_t id) {
output_file.write(raw_bed_lines[id]);
if (!isWhite(raw_bed_lines[id].back))
output_file.write('\t');
import std.string;
raw_bed_lines[id] = std.string.stripRight(raw_bed_lines[id]);
output_file.write(raw_bed_lines[id], "\t");
}

override BamRegion getRegionById(size_t id) {
Expand Down
42 changes: 29 additions & 13 deletions sambamba/utils/common/bed.d
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
/*
This file is part of Sambamba.
Copyright (C) 2013-2014 Artem Tarasov <[email protected]>
Copyright (C) 2013-2015 Artem Tarasov <[email protected]>
Sambamba is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
Expand All @@ -26,6 +26,7 @@ import std.conv;
import std.math;
import std.array;
import std.range;
import std.typecons;

struct Interval {
long beg;
Expand Down Expand Up @@ -53,7 +54,8 @@ Interval[] nonOverlappingIntervals(Interval[] list) {

alias Interval[][string] BedIndex;

BedIndex readIntervals(string bed_filename, bool non_overlapping=true, string[]* lines=null) {
// TODO: rewrite this mess
BedIndex readIntervals(string bed_filename, bool non_overlapping=true, string[]* lines=null, Tuple!(string, Interval)[]* intervals=null) {
BedIndex index;

auto bed = cast(string)(std.file.readText(bed_filename));
Expand All @@ -73,8 +75,12 @@ BedIndex readIntervals(string bed_filename, bool non_overlapping=true, string[]*

if (interval.beg == interval.end)
interval.end = interval.beg + 1;
if (interval.beg < interval.end)
index[chr] ~= interval;
if (interval.beg < interval.end) {
if (lines is null)
index[chr] ~= interval;
else
(*intervals) ~= tuple(chr, interval);
}

if (lines !is null) {
(*lines) ~= str;
Expand Down Expand Up @@ -119,17 +125,27 @@ import bio.bam.reader;
public import bio.bam.region;

BamRegion[] parseBed(Reader)(string bed_filename, Reader bam, bool non_overlapping=true, string[]* bed_lines=null) {
auto index = sambamba.utils.common.bed.readIntervals(bed_filename, non_overlapping, bed_lines);
Tuple!(string, Interval)[] ivs;
auto index = sambamba.utils.common.bed.readIntervals(bed_filename, non_overlapping, bed_lines, &ivs);
BamRegion[] regions;
foreach (reference, intervals; index) {
if (!bam.hasReference(reference))
continue;
auto id = bam[reference].id;
foreach (interval; intervals)
if (non_overlapping) {
foreach (reference, intervals; index) {
if (!bam.hasReference(reference))
continue;
auto id = bam[reference].id;
foreach (interval; intervals)
regions ~= BamRegion(cast(uint)id,
cast(uint)interval.beg, cast(uint)interval.end);
}
std.algorithm.sort(regions);
} else {
foreach (t; ivs) {
if (!bam.hasReference(t[0]))
continue;
auto id = bam[t[0]].id;
regions ~= BamRegion(cast(uint)id,
cast(uint)interval.beg, cast(uint)interval.end);
cast(uint)t[1].beg, cast(uint)t[1].end);
}
}
if (bed_lines is null)
std.algorithm.sort(regions);
return regions;
}

0 comments on commit 6d1f5a3

Please sign in to comment.