1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208
|
package Bio::Graphics::Browser2::DataLoader::bam;
# $Id$
use strict;
use base 'Bio::Graphics::Browser2::DataLoader';
my $HASBIGWIG;
sub create_conf_file {
my $self = shift;
my $bam_file = shift;
my $has_bigwig = shift;
my $conf = $self->conf_fh;
my $data_path = $self->data_path;
my $loadid = $self->loadid;
my $tracklabel = $self->new_track_label;
my $filename = $self->track_name;
# find a fasta file to use
my $fasta = $self->get_fasta_file || '';
my $category = $self->category;
my $sam_db =<<END;
[$loadid:database]
db_adaptor = Bio::DB::Sam
db_args = -bam "$bam_file"
-fasta "$fasta"
-split_splices 1
search options = none
END
(my $bigwig = $bam_file) =~ s/\.bam$/.bw/;
my $bigwig_db = $has_bigwig ? <<BIGWIG : '';
[${loadid}_bw:database]
db_adaptor = Bio::DB::BigWig
db_args = -bigwig "$bigwig"
BIGWIG
my $sam_track =<<END;
[$tracklabel]
feature = match
glyph = segments
draw_target = 1
show_mismatch = 1
mismatch_color = red
database = $loadid
bgcolor = blue
fgcolor = blue
height = 3
label = $filename
category = $category
label density = 50
feature_limit = 500
bump = fast
stranded = 1
key = $filename
END
my $semantic_track = $has_bigwig ? <<BIGWIG : <<ORDINARY;
[$tracklabel:499]
database = ${loadid}_bw
feature = summary
glyph = wiggle_whiskers
max_color = lightgrey
min_color = lightgrey
mean_color = black
stdev_color = grey
stdev_color_neg = grey
height = 20
BIGWIG
[$tracklabel:499]
database = $loadid
feature = coverage:2000
min_score = 0
glyph = wiggle_xyplot
height = 20
fgcolor = blue
bgcolor = blue
autoscale = local
key = $filename
ORDINARY
print $conf <<END;
$sam_db
$bigwig_db
#>>>>>>>>>> cut here <<<<<<<<
$semantic_track
$sam_track
END
}
# slightly different behavior -- never return the .fai file - only the first .fa file
sub get_fasta_file {
my $self = shift;
my @fasta = $self->get_fasta_files;
my $fasta = (grep {!/\.fai$/} @fasta)[0];
return $fasta;
}
sub load {
my $self = shift;
my ($initial_lines,$fh) = @_;
$self->flag_busy(1);
eval {
$self->open_conf;
$self->set_status('starting load');
mkdir $self->sources_path or die $!;
my $source_file = IO::File->new(
File::Spec->catfile($self->sources_path,$self->track_name),'>');
$self->start_load;
$self->set_status('load data');
my $bytes_loaded = 0;
foreach (@$initial_lines) {
$source_file->print($_);
$bytes_loaded += length $_;
}
my $buffer;
while ((my $bytes = read($fh,$buffer,8192) > 0)) {
$source_file->print($buffer);
$bytes_loaded += length $ buffer;
$self->set_status("loaded $bytes_loaded bytes") if $bytes++ % 10000;
}
$source_file->close();
$self->finish_load;
$self->close_conf;
$self->set_processing_complete;
};
$self->flag_busy(0);
die $@ if $@;
return $self->tracks;
}
sub finish_load {
my $self = shift;
eval "require Bio::DB::Sam; 1" or return;
# keep original copy in sources directory. Create new sorted and indexed
# copy in main level
my $source = File::Spec->catfile($self->sources_path,$self->track_name);
my $dest = File::Spec->catfile($self->data_path,$self->track_name);
$dest =~ s/\.[bs]am$//i; # sorting will add the .bam extension
# check whether we need to sort or not
if ($self->is_sorted($source)) {
# bam is already sorted by coordinate, destination will become source
$dest = $source;
}
else {
$self->set_status('sorting BAM file');
Bio::DB::Bam->sort_core(0,$source,$dest,250*1e6);
$dest .= '.bam';
}
$self->set_status('indexing BAM file');
Bio::DB::Bam->index_build($dest);
my $bigwig_exists = 0;
if ($self->has_bigwig) {
$self->set_status('creating BigWig coverage file');
$bigwig_exists = $self->create_big_wig();
}
$self->set_status('creating conf file');
$self->create_conf_file($dest,$bigwig_exists);
}
sub has_bigwig {
my $self = shift;
return $HASBIGWIG if defined $HASBIGWIG;
my $result = eval "require Bio::DB::Sam::SamToGBrowse;1";
$result &&= eval "require Bio::DB::BigWig; 1";
return $HASBIGWIG = $result;
}
sub create_big_wig {
my $self = shift;
my $dir = $self->data_path;
my $fasta = $self->get_fasta_file or return;
my $wigout = Bio::DB::Sam::SamToGBrowse->new($dir,$fasta,0);
$wigout->bam_to_wig($self->chrom_sizes); # this creates the wig files
die $wigout->last_error if $wigout->last_error;
1;
}
sub is_sorted {
my $self = shift;
my $sam = Bio::DB::Sam->new(-bam => shift, -autoindex => 0);
my $header = $sam->bam->header->text;
return (split /\n/, $header)[0] =~ /SO:coordinate/i;
}
1;
|