code_perl/Maasha/BGB/Wiggle.pm

   1 package Maasha::BGB::Wiggle;
   2
   3 # Copyright (C) 2010 Martin A. Hansen.
   4
   5 # This program is free software; you can redistribute it and/or
   6 # modify it under the terms of the GNU General Public License
   7 # as published by the Free Software Foundation; either version 2
   8 # of the License, or (at your option) any later version.
   9
  10 # This program is distributed in the hope that it will be useful,
  11 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  13 # GNU General Public License for more details.
  14
  15 # You should have received a copy of the GNU General Public License
  16 # along with this program; if not, write to the Free Software
  17 # Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
  18
  19 # http://www.gnu.org/copyleft/gpl.html
  20
  21
  22 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> DESCRIPTION <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
  23
  24
  25 # Routines for creating Biopieces Browser wiggle tracks.
  26
  27
  28 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
  29
  30
  31 use warnings;
  32 use strict;
  33 use Data::Dumper;
  34 use Maasha::Common;
  35 use Maasha::Calc;
  36 use Maasha::Filesys;
  37 use Maasha::KISS;
  38 use Maasha::Matrix;
  39
  40 use vars qw( @ISA @EXPORT );
  41
  42 @ISA = qw( Exporter );
  43
  44 use constant {
  45     S_ID             => 0,
  46     S_BEG            => 1,
  47     S_END            => 2,
  48     Q_ID             => 3,
  49     SCORE            => 4,
  50     STRAND           => 5,
  51     HITS             => 6,
  52     ALIGN            => 7,
  53     BLOCK_COUNT      => 8,
  54     BLOCK_BEGS       => 9,
  55     BLOCK_LENS       => 10,
  56     BLOCK_TYPE       => 11,
  57 };
  58
  59 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
  60
  61
  62 sub wiggle_encode
  63 {
  64     # Martin A. Hansen, February 2010.
  65
  66     # Read a KISS file and encode a Wiggle list.
  67
  68     my ( $file,   # KISS file
  69        ) = @_;
  70
  71     # Returns a list.
  72
  73     my ( $fh, $entry, $vals, $i );
  74
  75     $fh = Maasha::Filesys::file_read_open( $file );
  76
  77     while ( $entry = Maasha::KISS::kiss_entry_get( $fh ) ) {
  78         map { $vals->[ $_ ]++ } ( $entry->[ S_BEG ] .. $entry->[ S_END ] );
  79     }
  80
  81     close $fh;
  82
  83     for ( $i = 0; $i < scalar @{ $vals }; $i++ ) {
  84         $vals->[ $i ] = 0 if not defined $vals->[ $i ];
  85     }
  86
  87     return wantarray ? @{ $vals } : $vals;
  88 }
  89
  90
  91 sub wiggle_normalize
  92 {
  93     my ( $vals,   # Wiggle values
  94          $size,   # New list size
  95        ) = @_;
  96
  97     # Returns a list
  98
  99     if ( scalar @{ $vals } < $size ) {
 100         Maasha::Matrix::list_inflate( $vals, $size );
 101     } elsif ( scalar @{ $vals } > $size ) {
 102         Maasha::Matrix::list_deflate( $vals, $size );
 103     }
 104
 105     return wantarray ? @{ $vals } : $vals;
 106 }
 107
 108
 109 sub wiggle_store
 110 {
 111     # Martin A. Hansen, February 2010.
 112
 113     # Store a list of wiggle values as a byte array to a
 114     # specified file.
 115
 116     my ( $file,   # path to file
 117          $vals,   # Wiggle values
 118        ) = @_;
 119
 120     # Returns nothing.
 121
 122     my ( $bin, $fh );
 123
 124     $bin = pack( "S*", @{ $vals } );
 125
 126     $fh = Maasha::Filesys::file_write_open( $file );
 127
 128     print $fh $bin;
 129
 130     close $fh;
 131 }
 132
 133
 134 sub wiggle_retrieve
 135 {
 136     # Martin A. Hansen, February 2010.
 137
 138     # Restore an interval of Wiggle values from
 139     # a specified Wiggle file and return these
 140     # as a list.
 141
 142     my ( $file,   # path to wiggle file
 143          $beg,    # begin position
 144          $end,    # end position
 145        ) = @_;
 146
 147     # Returns a list.
 148
 149     my ( $fh, $bin, @vals );
 150
 151     Maasha::Common::error( qq(begin < 0: $beg) ) if $beg < 0;
 152     Maasha::Common::error( qq(begin > end: $beg > $end) ) if $beg > $end;
 153
 154     $fh = Maasha::Filesys::file_read_open( $file );
 155
 156     sysseek( $fh, $beg * 2, 0 );
 157     sysread( $fh, $bin, ( $end - $beg + 1 ) * 2 );
 158
 159     close $fh;
 160
 161     @vals = unpack( "S*", $bin );
 162
 163     map { push @vals, 0 } ( scalar @vals .. $end - $beg  );  # Padding
 164
 165     return wantarray ? @vals : \@vals;
 166 }
 167
 168
 169 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 170
 171
 172 1;
 173
 174
 175 __END__