From: martinahansen Date: Fri, 29 May 2009 12:13:05 +0000 (+0000) Subject: migrated final plotters X-Git-Url: https://git.donarmstrong.com/?a=commitdiff_plain;h=c11bd52c815bb8ad844ef837d41c9b5aa0bc2bdb;p=biopieces.git migrated final plotters git-svn-id: http://biopieces.googlecode.com/svn/trunk@455 74ccb610-7750-0410-82ae-013aeee3265d --- diff --git a/bp_bin/plot_matches b/bp_bin/plot_matches index fdf5bd2..a2d97cf 100755 --- a/bp_bin/plot_matches +++ b/bp_bin/plot_matches @@ -1,6 +1,213 @@ -#!/usr/bin/env perl +#!/usr/bin/env perl -w + +# Copyright (C) 2007-2009 Martin A. Hansen. + +# This program is free software; you can redistribute it and/or +# modify it under the terms of the GNU General Public License +# as published by the Free Software Foundation; either version 2 +# of the License, or (at your option) any later version. + +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. + +# You should have received a copy of the GNU General Public License +# along with this program; if not, write to the Free Software +# Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. + +# http://www.gnu.org/copyleft/gpl.html + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> DESCRIPTION <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + +# Generate a dotplot of matches in the stream. + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + -use warnings; use strict; +use Maasha::Biopieces; +use Maasha::Plot; +use Maasha::Filesys; +use IPC::Open2; + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + + +my ( $run_time_beg, $run_time_end, $options, $in, $out, $default, $terminals, $record, @data, $fh, $result, %data_hash, $tmp_dir ); + +$default = "plot_matches"; +$terminals = "dumb,x11,aqua,post,svg"; + +$options = Maasha::Biopieces::parse_options( + [ + { long => 'no_stream', short => 'x', type => 'flag', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + { long => 'data_out', short => 'o', type => 'file', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + { long => 'terminal', short => 't', type => 'string', mandatory => 'no', default => 'dumb', allowed => $terminals, disallowed => undef }, + { long => 'direction', short => 'd', type => 'string', mandatory => 'no', default => 'both', allowed => 'both,forward,reverse', disallowed => undef }, + { long => 'title', short => 'T', type => 'string', mandatory => 'no', default => $default, allowed => undef, disallowed => undef }, + { long => 'xlabel', short => 'X', type => 'string', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + { long => 'ylabel', short => 'Y', type => 'string', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + ] +); + +$in = Maasha::Biopieces::read_stream( $options->{ "stream_in" } ); +$out = Maasha::Biopieces::write_stream( $options->{ "stream_out" } ); + +while ( $record = Maasha::Biopieces::get_record( $in ) ) +{ + if ( defined $record->{ "Q_BEG" } and defined $record->{ "S_BEG" } and $record->{ "Q_END" } and $record->{ "S_END" } ) { + push @data, $record; + } + + Maasha::Biopieces::put_record( $record, $out ) if not $options->{ "no_stream" }; +} + +$options->{ "xlabel" } ||= $data[ 0 ]->{ "Q_ID" }; +$options->{ "ylabel" } ||= $data[ 0 ]->{ "S_ID" }; + +$tmp_dir = Maasha::Biopieces::get_tmpdir(); + +$result = dotplot_matches( \@data, $options, $tmp_dir ); + +$fh = Maasha::Biopieces::write_stream( $options->{ "data_out" } ); + +print $fh "$_\n" foreach @{ $result }; + +close $fh; + +Maasha::Filesys::dir_remove( $tmp_dir ); + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> SUBROUTINES <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + + +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 ); + + $forward_file = "$tmp_dir/match_f.tab"; + $backward_file = "$tmp_dir/match_r.tab"; + + $fh_forward = Maasha::Filesys::file_write_open( $forward_file ); + $fh_backward = Maasha::Filesys::file_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; + + $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; +} + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + + +BEGIN +{ + $run_time_beg = Maasha::Biopieces::run_time(); + + Maasha::Biopieces::log_biopiece(); +} + + +END +{ + Maasha::Biopieces::close_stream( $in ); + Maasha::Biopieces::close_stream( $out ); + + $run_time_end = Maasha::Biopieces::run_time(); + + Maasha::Biopieces::run_time_print( $run_time_beg, $run_time_end, $options ); +} + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + + +__END__ -use Maasha::BioRun; diff --git a/bp_bin/plot_phastcons_profiles b/bp_bin/plot_phastcons_profiles index fdf5bd2..52489f5 100755 --- a/bp_bin/plot_phastcons_profiles +++ b/bp_bin/plot_phastcons_profiles @@ -1,6 +1,131 @@ -#!/usr/bin/env perl +#!/usr/bin/env perl -w + +# Copyright (C) 2007-2009 Martin A. Hansen. + +# This program is free software; you can redistribute it and/or +# modify it under the terms of the GNU General Public License +# as published by the Free Software Foundation; either version 2 +# of the License, or (at your option) any later version. + +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. + +# You should have received a copy of the GNU General Public License +# along with this program; if not, write to the Free Software +# Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. + +# http://www.gnu.org/copyleft/gpl.html + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> DESCRIPTION <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + +# Plot chromosome distribution histogram. + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + -use warnings; use strict; +use Maasha::Biopieces; +use Maasha::Plot; +use Maasha::Matrix; +use Maasha::Filesys; +use Maasha::UCSC; + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + + +my ( $run_time_beg, $run_time_end, $options, $in, $out, $default, $terminals, $phastcons_file, $phastcons_index, + $index, $fh_phastcons, $record, $scores, $AoA, $plot, $fh, $tmp_dir ); + +$default = "PhastCons Profiles"; +$terminals = "dumb,x11,aqua,post,svg"; + +$options = Maasha::Biopieces::parse_options( + [ + { long => 'no_stream', short => 'x', type => 'flag', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + { long => 'data_out', short => 'o', type => 'file', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + { long => 'genome', short => 'g', type => 'genome', mandatory => 'yes', default => undef, allowed => undef, disallowed => undef }, + { long => 'mean', short => 'm', type => 'flag', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + { long => 'median', short => 'M', type => 'flag', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + { long => 'flank', short => 'f', type => 'uint', mandatory => 'no', default => 0, allowed => undef, disallowed => undef }, + { long => 'terminal', short => 't', type => 'string', mandatory => 'no', default => 'dumb', allowed => $terminals, disallowed => undef }, + { long => 'title', short => 'T', type => 'string', mandatory => 'no', default => $default, allowed => undef, disallowed => undef }, + { long => 'xlabel', short => 'X', type => 'string', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + { long => 'ylabel', short => 'Y', type => 'string', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + ] +); + +$in = Maasha::Biopieces::read_stream( $options->{ "stream_in" } ); +$out = Maasha::Biopieces::write_stream( $options->{ "stream_out" } ); + +$phastcons_file = Maasha::Config::genome_phastcons( $options->{ "genome" } ); +$phastcons_index = Maasha::Config::genome_phastcons_index( $options->{ "genome" } ); + +$index = Maasha::UCSC::fixedstep_index_retrieve( $phastcons_index ); +$fh_phastcons = Maasha::Filesys::file_read_open( $phastcons_file ); + +while ( $record = Maasha::Biopieces::get_record( $in ) ) +{ + if ( $record->{ "CHR" } and $record->{ "CHR_BEG" } and $record->{ "CHR_END" } ) + { + $scores = Maasha::UCSC::fixedstep_index_lookup( $index, $fh_phastcons, $record->{ "CHR" }, + $record->{ "CHR_BEG" }, + $record->{ "CHR_END" }, + $options->{ "flank" } ); + + push @{ $AoA }, [ @{ $scores } ]; + } + + Maasha::Biopieces::put_record( $record, $out ) if not $options->{ "no_stream" }; +} + +Maasha::UCSC::phastcons_normalize( $AoA ); + +$AoA = [ [ Maasha::UCSC::phastcons_mean( $AoA ) ] ] if $options->{ "mean" }; +$AoA = [ [ Maasha::UCSC::phastcons_median( $AoA ) ] ] if $options->{ "median" }; + +$AoA = Maasha::Matrix::matrix_flip( $AoA ); + +$tmp_dir = Maasha::Biopieces::get_tmpdir(); + +$plot = Maasha::Plot::lineplot_simple( $AoA, $options, $tmp_dir ); + +$fh = Maasha::Biopieces::write_stream( $options->{ "data_out" } ); + +print $fh "$_\n" foreach @{ $plot }; + +close $fh; + +Maasha::Filesys::dir_remove( $tmp_dir ); + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + + +BEGIN +{ + $run_time_beg = Maasha::Biopieces::run_time(); + + Maasha::Biopieces::log_biopiece(); +} + + +END +{ + Maasha::Biopieces::close_stream( $in ); + Maasha::Biopieces::close_stream( $out ); + + $run_time_end = Maasha::Biopieces::run_time(); + + Maasha::Biopieces::run_time_print( $run_time_beg, $run_time_end, $options ); +} + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + + +__END__ -use Maasha::BioRun; diff --git a/bp_bin/plot_seqlogo b/bp_bin/plot_seqlogo index fdf5bd2..0a2d5d4 100755 --- a/bp_bin/plot_seqlogo +++ b/bp_bin/plot_seqlogo @@ -1,6 +1,93 @@ -#!/usr/bin/env perl +#!/usr/bin/env perl -w + +# Copyright (C) 2007-2009 Martin A. Hansen. + +# This program is free software; you can redistribute it and/or +# modify it under the terms of the GNU General Public License +# as published by the Free Software Foundation; either version 2 +# of the License, or (at your option) any later version. + +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. + +# You should have received a copy of the GNU General Public License +# along with this program; if not, write to the Free Software +# Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. + +# http://www.gnu.org/copyleft/gpl.html + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> DESCRIPTION <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + +# Split the values of a key into new key/value pairs. + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + -use warnings; use strict; +use Maasha::Biopieces; +use Maasha::Plot; + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + + +my ( $run_time_beg, $run_time_end, $options, $in, $out, $record, @entries, $logo, $fh ); + +$options = Maasha::Biopieces::parse_options( + [ + { long => 'no_stream', short => 'x', type => 'flag', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + { long => 'data_out', short => 'o', type => 'file', mandatory => 'no', default => undef, allowed => undef, disallowed => undef }, + ] +); + +$in = Maasha::Biopieces::read_stream( $options->{ "stream_in" } ); +$out = Maasha::Biopieces::write_stream( $options->{ "stream_out" } ); + +while ( $record = Maasha::Biopieces::get_record( $in ) ) +{ + if ( $record->{ "SEQ_NAME" } and $record->{ "SEQ" } ) { + push @entries, [ $record->{ "SEQ_NAME" }, $record->{ "SEQ" } ]; + } + + Maasha::Biopieces::put_record( $record, $out ) if not $options->{ "no_stream" }; +} + +$logo = Maasha::Plot::seq_logo( \@entries ); + +$fh = Maasha::Biopieces::write_stream( $options->{ "data_out" } ); + +print $fh $logo; + +close $fh; + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + + +BEGIN +{ + $run_time_beg = Maasha::Biopieces::run_time(); + + Maasha::Biopieces::log_biopiece(); +} + + +END +{ + Maasha::Biopieces::close_stream( $in ); + Maasha::Biopieces::close_stream( $out ); + + $run_time_end = Maasha::Biopieces::run_time(); + + Maasha::Biopieces::run_time_print( $run_time_beg, $run_time_end, $options ); +} + + +# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + + +__END__ -use Maasha::BioRun; diff --git a/code_perl/Maasha/BioRun.pm b/code_perl/Maasha/BioRun.pm index d5f4e1f..9ff29a1 100644 --- a/code_perl/Maasha/BioRun.pm +++ b/code_perl/Maasha/BioRun.pm @@ -143,9 +143,6 @@ sub run_script elsif ( $script eq "write_ucsc_config" ) { script_write_ucsc_config( $in, $out, $options ) } elsif ( $script eq "remove_mysql_tables" ) { script_remove_mysql_tables( $in, $out, $options ) } elsif ( $script eq "remove_ucsc_tracks" ) { script_remove_ucsc_tracks( $in, $out, $options ) } - elsif ( $script eq "plot_matches" ) { script_plot_matches( $in, $out, $options ) } - elsif ( $script eq "plot_seqlogo" ) { script_plot_seqlogo( $in, $out, $options ) } - elsif ( $script eq "plot_phastcons_profiles" ) { script_plot_phastcons_profiles( $in, $out, $options ) } elsif ( $script eq "upload_to_ucsc" ) { script_upload_to_ucsc( $in, $out, $options ) } close $in if defined $in; @@ -373,28 +370,6 @@ sub get_options data_out|o=s ); } - elsif ( $script eq "plot_seqlogo" ) - { - @options = qw( - no_stream|x - data_out|o=s - ); - } - elsif ( $script eq "plot_phastcons_profiles" ) - { - @options = qw( - no_stream|x - data_out|o=s - genome|g=s - mean|m - median|M - flank|f=s - terminal|t=s - title|T=s - xlabel|X=s - ylabel|Y=s - ); - } elsif ( $script eq "remove_mysql_tables" ) { @options = qw( @@ -418,18 +393,6 @@ sub get_options no_stream|x ); } - elsif ( $script eq "plot_matches" ) - { - @options = qw( - no_stream|x - data_out|o=s - terminal|t=s - title|T=s - xlabel|X=s - ylabel|Y=s - direction|d=s - ); - } elsif ( $script eq "upload_to_ucsc" ) { @options = qw( @@ -549,7 +512,7 @@ sub get_options Maasha::Common::error( qq(no --database specified) ) if $script eq "remove_ucsc_tracks" and not $options{ "database" }; Maasha::Common::error( qq(no --in_file or --genome specified) ) if $script eq "soap_seq" and not $options{ "genome" } and not $options{ "in_file" }; Maasha::Common::error( qq(both --in_file and --genome specified) ) if $script eq "soap_seq" and $options{ "genome" } and $options{ "in_file" }; - Maasha::Common::error( qq(no --genome specified) ) if $script =~ /get_genome_align|get_genome_phastcons|plot_phastcons_profiles/ and not $options{ "genome" }; + Maasha::Common::error( qq(no --genome specified) ) if $script =~ /get_genome_align|get_genome_phastcons/ and not $options{ "genome" }; if ( $script eq "upload_to_ucsc" ) { @@ -1850,95 +1813,6 @@ sub script_write_ucsc_config } -sub script_plot_seqlogo -{ - # Martin A. Hansen, August 2007. - - # Calculates and writes a sequence logo for alignments. - - my ( $in, # handle to in stream - $out, # handle to out stream - $options, # options hash - ) = @_; - - # Returns nothing. - - my ( $record, @entries, $logo, $fh ); - - while ( $record = Maasha::Biopieces::get_record( $in ) ) - { - if ( $record->{ "SEQ_NAME" } and $record->{ "SEQ" } ) { - push @entries, [ $record->{ "SEQ_NAME" }, $record->{ "SEQ" } ]; - } - - Maasha::Biopieces::put_record( $record, $out ) if not $options->{ "no_stream" }; - } - - $logo = Maasha::Plot::seq_logo( \@entries ); - - $fh = Maasha::Biopieces::write_stream( $options->{ "data_out" } ); - - print $fh $logo; - - close $fh; -} - - -sub script_plot_phastcons_profiles -{ - # Martin A. Hansen, January 2008. - - # Plots PhastCons profiles. - - my ( $in, # handle to in stream - $out, # handle to out stream - $options, # options hash - ) = @_; - - # Returns nothing. - - my ( $phastcons_file, $phastcons_index, $index, $fh_phastcons, $record, $scores, $AoA, $plot, $fh ); - - $options->{ "title" } ||= "PhastCons Profiles"; - - $phastcons_file = Maasha::Config::genome_phastcons( $options->{ "genome" } ); - $phastcons_index = Maasha::Config::genome_phastcons_index( $options->{ "genome" } ); - - $index = Maasha::UCSC::fixedstep_index_retrieve( $phastcons_index ); - $fh_phastcons = Maasha::Common::read_open( $phastcons_file ); - - while ( $record = Maasha::Biopieces::get_record( $in ) ) - { - if ( $record->{ "CHR" } and $record->{ "CHR_BEG" } and $record->{ "CHR_END" } ) - { - $scores = Maasha::UCSC::fixedstep_index_lookup( $index, $fh_phastcons, $record->{ "CHR" }, - $record->{ "CHR_BEG" }, - $record->{ "CHR_END" }, - $options->{ "flank" } ); - - push @{ $AoA }, [ @{ $scores } ]; - } - - Maasha::Biopieces::put_record( $record, $out ) if not $options->{ "no_stream" }; - } - - Maasha::UCSC::phastcons_normalize( $AoA ); - - $AoA = [ [ Maasha::UCSC::phastcons_mean( $AoA ) ] ] if $options->{ "mean" }; - $AoA = [ [ Maasha::UCSC::phastcons_median( $AoA ) ] ] if $options->{ "median" }; - - $AoA = Maasha::Matrix::matrix_flip( $AoA ); - - $plot = Maasha::Plot::lineplot_simple( $AoA, $options, $BP_TMP ); - - $fh = Maasha::Biopieces::write_stream( $options->{ "data_out" } ); - - print $fh "$_\n" foreach @{ $plot }; - - close $fh; -} - - sub script_remove_mysql_tables { # Martin A. Hansen, November 2008. @@ -2063,46 +1937,6 @@ sub script_remove_ucsc_tracks } -sub script_plot_matches -{ - # Martin A. Hansen, August 2007. - - # Plot matches in 2D generating a dotplot. - - my ( $in, # handle to in stream - $out, # handle to out stream - $options, # options hash - ) = @_; - - # Returns nothing. - - my ( $record, @data, $fh, $result, %data_hash ); - - $options->{ "direction" } ||= "both"; - - while ( $record = Maasha::Biopieces::get_record( $in ) ) - { - if ( defined $record->{ "Q_BEG" } and defined $record->{ "S_BEG" } and $record->{ "Q_END" } and $record->{ "S_END" } ) { - push @data, $record; - } - - Maasha::Biopieces::put_record( $record, $out ) if not $options->{ "no_stream" }; - } - - $options->{ "title" } ||= "plot_matches"; - $options->{ "xlabel" } ||= $data[ 0 ]->{ "Q_ID" }; - $options->{ "ylabel" } ||= $data[ 0 ]->{ "S_ID" }; - - $result = Maasha::Plot::dotplot_matches( \@data, $options, $BP_TMP ); - - $fh = Maasha::Biopieces::write_stream( $options->{ "data_out" } ); - - print $fh "$_\n" foreach @{ $result }; - - close $fh; -} - - sub script_upload_to_ucsc { # Martin A. Hansen, August 2007. diff --git a/code_perl/Maasha/Plot.pm b/code_perl/Maasha/Plot.pm index 88992c5..cac8536 100644 --- a/code_perl/Maasha/Plot.pm +++ b/code_perl/Maasha/Plot.pm @@ -313,113 +313,6 @@ sub histogram_chrdist } -# >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 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::Filesys::file_write_open( $forward_file ); - $fh_backward = Maasha::Filesys::file_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 <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<