# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
+use warnings;
use strict;
use Data::Dumper;
use SVG;
use IPC::Open2;
-use Time::HiRes qw( gettimeofday );
use Maasha::Common;
-use Maasha::Calc;
+use Maasha::Filesys;
+use Maasha::Seq;
use vars qw ( @ISA @EXPORT );
use constant {
# Returns list.
- my ( $tmp_file, $pid, $fh_in, $fh_out, $cmd, $i, $line, @lines, $xtic_space, @plot_cmd );
+ my ( $tmp_file, $pid, $fh_in, $fh_out, $cmd, $i, $line, @lines, $xtic_space, @plot_cmd, $title );
$tmp_dir ||= $ENV{ 'BP_TMP' };
$tmp_file = "$tmp_dir/lineplot_simple.tab";
- $fh_out = Maasha::Common::write_open( $tmp_file );
+ $fh_out = Maasha::Filesys::file_write_open( $tmp_file );
map { print $fh_out join( "\t", @{ $_ } ), "\n" } @{ $data };
$options->{ "terminal" } ||= "dumb";
- $cmd = "gnuplot";
+ $cmd = "gnuplot -persist";
$pid = open2( $fh_out, $fh_in, $cmd );
print $fh_in "set title \"$options->{ 'title' }\"\n" if $options->{ "title" };
print $fh_in "set xlabel \"$options->{ 'xlabel' }\"\n" if $options->{ "xlabel" };
print $fh_in "set ylabel \"$options->{ 'ylabel' }\"\n" if $options->{ "ylabel" };
+ print $fh_in "set logscale y\n" if $options->{ "logscale_y" };
print $fh_in "set grid\n" if not $options->{ "terminal" } eq "dumb";
print $fh_in "set autoscale\n";
- print $fh_in "unset key\n";
- print $fh_in "set xtics border in scale 0 nomirror rotate by 90 offset character 0, 0, 0\n";
+ print $fh_in "set xtics rotate by -90\n";
- for ( $i = 1; $i < scalar @{ $data->[ 0 ] } + 1; $i++ ) {
- push @plot_cmd, qq("$tmp_file" using $i with lines ls 1);
+ for ( $i = 1; $i < scalar @{ $data->[ 0 ] } + 1; $i++ )
+ {
+ $title = $options->{ 'keys' }->[ $i - 1 ] || $options->{ 'list' } || "";
+ push @plot_cmd, qq("$tmp_file" using $i with lines title "$title");
}
print $fh_in "plot " . join( ", ", @plot_cmd ) . "\n";
$options->{ "terminal" } ||= "dumb";
- $cmd = "gnuplot";
+ $cmd = "gnuplot -persist";
$pid = open2( $fh_out, $fh_in, $cmd );
print $fh_in "set title \"$options->{ 'title' }\"\n" if $options->{ "title" };
print $fh_in "set xlabel \"$options->{ 'xlabel' }\"\n" if $options->{ "xlabel" };
print $fh_in "set ylabel \"$options->{ 'ylabel' }\"\n" if $options->{ "ylabel" };
+ print $fh_in "set logscale y\n" if $options->{ "logscale_y" };
print $fh_in "set autoscale\n";
print $fh_in "unset key\n";
print $fh_in "set style fill solid\n";
print $fh_in "set style histogram title offset character 0, 0, 0\n";
print $fh_in "set style data histograms\n";
- print $fh_in "set xtics border in scale 0 nomirror rotate by 90 offset character 0, 0, 0\n";
-
+ print $fh_in "set xtics border in scale 0 nomirror rotate by -90 offset character 0, 0, 0\n";
print $fh_in "plot '-' using 2:xticlabels(1)\n";
for ( $i = 0; $i < @{ $data }; $i++ )
$options->{ "terminal" } ||= "dumb";
- if ( $data->[ -1 ]->[ 0 ] <= 10 ) {
+ if ( $data->[ -1 ]->[ 0 ] <= 10 ) { # FIXME: some day figure the formula for this!
$xtic_space = 1;
} elsif ( $data->[ -1 ]->[ 0 ] <= 100 ) {
$xtic_space = 5;
$xtic_space = 1000;
} elsif ( $data->[ -1 ]->[ 0 ] <= 100000 ) {
$xtic_space = 5000;
+ } elsif ( $data->[ -1 ]->[ 0 ] <= 1000000 ) {
+ $xtic_space = 50000;
}
- $cmd = "gnuplot";
+ $cmd = "gnuplot -persist";
$pid = open2( $fh_out, $fh_in, $cmd );
print $fh_in "set title \"$options->{ 'title' }\"\n" if $options->{ "title" };
print $fh_in "set xlabel \"$options->{ 'xlabel' }\"\n" if $options->{ "xlabel" };
print $fh_in "set ylabel \"$options->{ 'ylabel' }\"\n" if $options->{ "ylabel" };
+ print $fh_in "set logscale y\n" if $options->{ "logscale_y" };
print $fh_in "set autoscale\n";
print $fh_in "unset key\n";
print $fh_in "set style fill solid\n";
$options->{ "terminal" } ||= "dumb";
- $cmd = "gnuplot";
+ $cmd = "gnuplot -persist";
$pid = open2( $fh_out, $fh_in, $cmd );
print $fh_in "set title \"$options->{ 'title' }\"\n" if $options->{ "title" };
print $fh_in "set xlabel \"$options->{ 'xlabel' }\"\n" if $options->{ "xlabel" };
print $fh_in "set ylabel \"$options->{ 'ylabel' }\"\n" if $options->{ "ylabel" };
+ print $fh_in "set logscale y\n" if $options->{ "logscale_y" };
print $fh_in "set autoscale\n";
print $fh_in "unset key\n";
print $fh_in "set style fill solid\n";
print $fh_in "set style histogram title offset character 0, 0, 0\n";
print $fh_in "set style data histograms\n";
- print $fh_in "set xtics border in scale 0 nomirror rotate by 90 offset character 0, 0, 0\n";
+ print $fh_in "set xtics border in scale 0 nomirror rotate by -90 offset character 0, 0, 0\n";
print $fh_in "plot '-' using 2:xticlabels(1)\n";
}
-# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> DOTPLOT <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
-
-
-sub dotplot_matches
-{
- # Martin A. Hansen, August 2007.
-
- # Generates a dotplot from a list of matches using Gnuplot.
-
- my ( $matches, # list of hashrefs.
- $options, # options hash
- $tmp_dir, # temporary directory
- ) = @_;
-
- # Returns list.
-
- my ( $forward_file, $backward_file, $pid, $fh_forward, $fh_backward,
- $fh_in, $fh_out, $cmd, $match, $line, @lines, $q_max, $s_max );
-
- $tmp_dir ||= $ENV{ 'BP_TMP' };
-
- $forward_file = "$tmp_dir/match_f.tab";
- $backward_file = "$tmp_dir/match_r.tab";
-
- $fh_forward = Maasha::Common::write_open( $forward_file );
- $fh_backward = Maasha::Common::write_open( $backward_file );
-
- $q_max = 0;
- $s_max = 0;
-
- foreach $match ( @{ $matches } )
- {
- if ( $match->{ "DIR" } =~ /^f/ )
- {
- print $fh_forward join( "\t", $match->{ "Q_BEG" } + 1, $match->{ "S_BEG" } + 1 ), "\n";
- print $fh_forward join( "\t", $match->{ "Q_END" } + 1, $match->{ "S_END" } + 1 ), "\n";
- print $fh_forward "\n\n";
- }
- else
- {
- print $fh_backward join( "\t", $match->{ "Q_BEG" } + 1, $match->{ "S_END" } + 1 ), "\n";
- print $fh_backward join( "\t", $match->{ "Q_END" } + 1, $match->{ "S_BEG" } + 1 ), "\n";
- print $fh_backward "\n\n";
- }
-
- $q_max = $match->{ "Q_END" } if $match->{ "Q_END" } > $q_max;
- $s_max = $match->{ "S_END" } if $match->{ "S_END" } > $s_max;
- }
-
- $q_max++;
- $s_max++;
-
- close $fh_forward;
- close $fh_backward;
-
- $options->{ "terminal" } ||= "dumb";
-
- $cmd = "gnuplot";
-
- $pid = open2( $fh_out, $fh_in, $cmd );
-
- print $fh_in "set terminal $options->{ 'terminal' }\n";
- print $fh_in "set xrange [1:$q_max]\n";
- print $fh_in "set yrange [1:$s_max]\n";
- print $fh_in "set title \"$options->{ 'title' }\"\n" if $options->{ "title" };
- print $fh_in "set xlabel \"$options->{ 'xlabel' }\"\n" if $options->{ "xlabel" };
- print $fh_in "set ylabel \"$options->{ 'ylabel' }\"\n" if $options->{ "ylabel" };
- print $fh_in "unset key\n";
-
- if ( $options->{ "terminal" } ne "dumb" )
- {
- print $fh_in "set style line 1 linetype 1 linecolor rgb \"green\" linewidth 2 pointtype 6 pointsize default\n";
- print $fh_in "set style line 2 linetype 1 linecolor rgb \"red\" linewidth 2 pointtype 6 pointsize default\n";
- }
-
- print $fh_in "set xtics border out\n";
- print $fh_in "set ytics border out\n";
- print $fh_in "set grid\n";
-
- if ( $options->{ "direction" } =~ /^b/ ) {
- print $fh_in qq(plot "$forward_file" with lines ls 1, "$backward_file" with lines ls 2\n);
- } elsif ( $options->{ "direction" } =~ /^f/ ) {
- print $fh_in qq(plot "$forward_file" with lines ls 1\n);
- } elsif ( $options->{ "direction" } =~ /^r/ ) {
- print $fh_in qq(plot "$backward_file" with lines ls 2\n);
- }
-
- close $fh_in;
-
- while ( $line = <$fh_out> )
- {
- chomp $line;
-
- push @lines, $line;
- }
-
- close $fh_out;
-
- waitpid $pid, 0;
-
- unlink $forward_file;
- unlink $backward_file;
-
- return wantarray ? @lines : \@lines;
-}
-
-
# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> KARYOGRAM <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
my ( $karyo_file, $svg, $features, $karyo );
- if ( $options->{ "genome" } eq "human" )
+ if ( $options->{ "genome" } eq "hg18" )
{
- $karyo_file = "/Users/m.hansen/maasha/perl_scripts/biopieces/karyo_data/human_cytobands.txt";
-# $karyo_file = "/home/m.hansen/maasha/perl_scripts/biopieces/karyo_data/human_cytobands.txt";
+ $karyo_file = $ENV{ 'BP_DIR' } . "/bp_data/human_cytobands.txt";
}
- else
+ elsif( $options->{ "genome" } eq "mm9" )
{
- $karyo_file = "/Users/m.hansen/maasha/perl_scripts/biopieces/karyo_data/mouse_cytobands.txt";
- # $karyo_file = "/home/m.hansen/maasha/perl_scripts/biopieces/karyo_data/mouse_cytobands.txt";
+ $karyo_file = $ENV{ 'BP_DIR' } . "/bp_data/mouse_cytobands.txt";
}
$karyo = parse_karyo_data( $karyo_file );
# stalk => "gray66",
);
- $fh = Maasha::Common::read_open( $file );
+ $fh = Maasha::Filesys::file_read_open( $file );
while ( $line = <$fh> )
{
push @nums, $acen->[ 1 ]->[ 1 ];
push @nums, $acen->[ 1 ]->[ 2 ];
+
+ @nums = grep { defined $_ } @nums; # FIXME
@nums = sort { $a <=> $b } @nums;
$cent = ( $nums[ 1 ] + $nums[ 2 ] ) / 2;
$type = Maasha::Seq::seq_guess_type( $entries->[ 0 ]->[ 1 ] );
- if ( $type =~ /^p/i ) {
+ if ( $type eq "PROTEIN" ) {
$bit_max = 4;
} else {
$bit_max = 2;
projects / biopieces.git / blobdiff