]> git.donarmstrong.com Git - biopieces.git/blobdiff - code_perl/Maasha/KISS.pm
added missing files
[biopieces.git] / code_perl / Maasha / KISS.pm
index 539e501412104e531a86ceac316d2ba73021d965..35d21108d9b853f59dd74738278a29a9054d4fd5 100644 (file)
@@ -33,9 +33,11 @@ package Maasha::KISS;
 use warnings;
 use strict;
 use Data::Dumper;
+use JSON::XS;
 use Maasha::Common;
 use Maasha::Filesys;
 use Maasha::Align;
+
 use vars qw( @ISA @EXPORT );
 
 @ISA = qw( Exporter );
@@ -60,6 +62,10 @@ use constant {
     BUCKET_SIZE      => 100,
     COUNT            => 0,
     OFFSET           => 1,
+
+    INDEX_BEG        => 1,
+    INDEX_END        => 2,
+    INDEX            => 12,
 };
 
 
@@ -77,7 +83,7 @@ sub kiss_entry_get
 
     # Returns a hashref.
 
-    my ( $line, @fields, %entry );
+    my ( $line, @fields );
 
     while ( $line = <$fh> )
     {
@@ -89,24 +95,44 @@ sub kiss_entry_get
 
         Maasha::Common::error( qq(BAD kiss entry: $line) ) if not @fields == 12;
         
-        $entry{ 'S_ID' }        = $fields[ S_ID ];
-        $entry{ 'S_BEG' }       = $fields[ S_BEG ];
-        $entry{ 'S_END' }       = $fields[ S_END ];
-        $entry{ 'Q_ID' }        = $fields[ Q_ID ];
-        $entry{ 'SCORE' }       = $fields[ SCORE ];
-        $entry{ 'STRAND' }      = $fields[ STRAND ];
-        $entry{ 'HITS' }        = $fields[ HITS ];
-        $entry{ 'ALIGN' }       = $fields[ ALIGN ];
-        $entry{ 'BLOCK_COUNT' } = $fields[ BLOCK_COUNT ];
-        $entry{ 'BLOCK_BEGS' }  = $fields[ BLOCK_BEGS ];
-        $entry{ 'BLOCK_LENS' }  = $fields[ BLOCK_LENS ];
-        $entry{ 'BLOCK_TYPE' }  = $fields[ BLOCK_TYPE ];
-
-        return wantarray ? %entry : \%entry;
+        return wantarray ? @fields : \@fields;
     }
 }
 
 
+sub kiss_retrieve
+{
+    # Martin A. Hansen, February 2010.
+
+    # Retrieves KISS entries from a given sorted KISS file
+    # within an optional interval.
+
+    my ( $file,   # path to KISS file
+         $beg,    # interval begin  -  OPTIONAL
+         $end,    # interval end    -  OPTIONAL
+       ) = @_;
+
+    # Returns a list.
+
+    my ( $fh, $entry, @entries );
+
+    $beg ||= 0;
+    $end ||= 999999999;
+
+    $fh = Maasha::Filesys::file_read_open( $file );
+
+    while ( $entry = kiss_entry_get( $fh ) )
+    {
+        last if $entry->[ S_BEG ] > $end;
+        push @entries, $entry if $entry->[ S_END ] > $beg;
+    }
+
+    close $fh;
+
+    return wantarray ? @entries : \@entries;
+}
+
+
 sub kiss_entry_parse
 {
     # Martin A. Hansen, December 2009.
@@ -124,20 +150,7 @@ sub kiss_entry_parse
 
     Maasha::Common::error( qq(BAD kiss entry: $line) ) if not @fields == 12;
     
-    $entry{ 'S_ID' }        = $fields[ S_ID ];
-    $entry{ 'S_BEG' }       = $fields[ S_BEG ];
-    $entry{ 'S_END' }       = $fields[ S_END ];
-    $entry{ 'Q_ID' }        = $fields[ Q_ID ];
-    $entry{ 'SCORE' }       = $fields[ SCORE ];
-    $entry{ 'STRAND' }      = $fields[ STRAND ];
-    $entry{ 'HITS' }        = $fields[ HITS ];
-    $entry{ 'ALIGN' }       = $fields[ ALIGN ];
-    $entry{ 'BLOCK_COUNT' } = $fields[ BLOCK_COUNT ];
-    $entry{ 'BLOCK_BEGS' }  = $fields[ BLOCK_BEGS ];
-    $entry{ 'BLOCK_LENS' }  = $fields[ BLOCK_LENS ];
-    $entry{ 'BLOCK_TYPE' }  = $fields[ BLOCK_TYPE ];
-
-    return wantarray ? %entry : \%entry;
+    return wantarray ? @fields : \@fields;
 }
 
 
@@ -154,32 +167,29 @@ sub kiss_entry_put
 
     # Returns nothing.
     
-    my ( @fields );
+    Maasha::Common::error( qq(BAD kiss entry) ) if not scalar @{ $entry } == 12;
 
-    if ( defined $entry->{ 'S_ID' }  and 
-         defined $entry->{ 'S_BEG' } and
-         defined $entry->{ 'S_END' }
+    if ( defined $entry->[ S_ID ]  and 
+         defined $entry->[ S_BEG ] and
+         defined $entry->[ S_END ]
        )
     {
-        Maasha::Common::error( qq(Bad S_BEG value: $entry->{ 'S_BEG' } < 0 ) ) if $entry->{ 'S_BEG' } < 0;
-        Maasha::Common::error( qq(Bad S_END value: $entry->{ 'S_END' } < $entry->{ 'S_BEG' }) ) if $entry->{ 'S_END' } < $entry->{ 'S_BEG' };
+        Maasha::Common::error( qq(Bad S_BEG value: $entry->[ S_BEG ] < 0 ) ) if $entry->[ S_BEG ] < 0;
+        Maasha::Common::error( qq(Bad S_END value: $entry->[ S_END ] < $entry->[ S_BEG ] ) ) if $entry->[ S_END ] < $entry->[ S_BEG ];
 
         $fh ||= \*STDOUT;
     
-        $fields[ S_ID ]        = $entry->{ 'S_ID' };
-        $fields[ S_BEG ]       = $entry->{ 'S_BEG' };
-        $fields[ S_END ]       = $entry->{ 'S_END' };
-        $fields[ Q_ID ]        = $entry->{ 'Q_ID' }        || ".";
-        $fields[ SCORE ]       = $entry->{ 'SCORE' }       || ".";
-        $fields[ STRAND ]      = $entry->{ 'STRAND' }      || ".";
-        $fields[ HITS ]        = $entry->{ 'HITS' }        || ".";
-        $fields[ ALIGN ]       = $entry->{ 'ALIGN' }       || ".";
-        $fields[ BLOCK_COUNT ] = $entry->{ 'BLOCK_COUNT' } || ".";
-        $fields[ BLOCK_BEGS ]  = $entry->{ 'BLOCK_BEGS' }  || ".";
-        $fields[ BLOCK_LENS ]  = $entry->{ 'BLOCK_LENS' }  || ".";
-        $fields[ BLOCK_TYPE ]  = $entry->{ 'BLOCK_TYPE' }  || ".";
-
-        print $fh join( "\t", @fields ), "\n";
+        $entry->[ Q_ID ]        = "." if not defined $entry->[ Q_ID ];
+        $entry->[ SCORE ]       = "." if not defined $entry->[ SCORE ];
+        $entry->[ STRAND ]      = "." if not defined $entry->[ STRAND ];
+        $entry->[ HITS ]        = "." if not defined $entry->[ HITS ];
+        $entry->[ ALIGN ]       = "." if not defined $entry->[ ALIGN ];
+        $entry->[ BLOCK_COUNT ] = "." if not defined $entry->[ BLOCK_COUNT ];
+        $entry->[ BLOCK_BEGS ]  = "." if not defined $entry->[ BLOCK_BEGS ];
+        $entry->[ BLOCK_LENS ]  = "." if not defined $entry->[ BLOCK_LENS ];
+        $entry->[ BLOCK_TYPE ]  = "." if not defined $entry->[ BLOCK_TYPE ];
+
+        print $fh join( "\t", @{ $entry } ), "\n";
     }
 }
 
@@ -195,7 +205,7 @@ sub kiss_sort
 
     # Returns nothing.
 
-    `sort -k 2,2n -k 3,3n $file > $file.sort`;
+    `sort -k 2,2n -k 3,3nr $file > $file.sort`;
 
     rename "$file.sort", $file;
 }
@@ -207,7 +217,7 @@ sub kiss_index
 
     # Creates a lookup index of a sorted KISS file.
 
-    my ( $file,         # path to KISS file
+    my ( $file,   # path to KISS file
        ) = @_;
 
     # Returns nothing.
@@ -286,8 +296,8 @@ sub kiss_index_count
 
     Maasha::Common::error( qq(Negative begin position: "$beg") ) if $beg < 0;
 
-    $bucket_beg = int( $beg / BUCKET_SIZE ); 
-    $bucket_end = int( $end / BUCKET_SIZE ); 
+    $bucket_beg = int( $beg / BUCKET_SIZE );
+    $bucket_end = int( $end / BUCKET_SIZE );
 
     $bucket_end = scalar @{ $index } if $bucket_end > scalar @{ $index };
 
@@ -301,6 +311,31 @@ sub kiss_index_count
 }
 
 
+sub kiss_index_count_nc
+{
+    # Martin A. Hansen, December 2009.
+
+    # Given a KISS index and a begin/end interval
+    # sum the number of counts in that interval,
+    # and return this.
+
+    my ( $index,   # KISS index
+         $beg,     # Begin position
+         $end,     # End position
+       ) = @_;
+
+    # Returns a number.
+
+    my ( $count );
+
+    Maasha::Common::error( qq(Negative begin position: "$beg") ) if $beg < 0;
+
+    $count = Maasha::NClist::nc_list_count_interval( $index, $beg, $end, INDEX_BEG, INDEX_END, INDEX );
+
+    return $count;
+}
+
+
 sub kiss_index_get_entries
 {
     # Martin A. Hansen, November 2009.
@@ -319,6 +354,43 @@ sub kiss_index_get_entries
 
     my ( $offset, $fh, $entry, @entries );
 
+    # $offset = kiss_index_offset( $index, $beg );
+
+    $fh = Maasha::Filesys::file_read_open( $file );
+
+    # sysseek( $fh, $offset, 0 );
+
+    while ( $entry = Maasha::KISS::kiss_entry_get( $fh ) )
+    {
+        push @entries, $entry if $entry->[ S_END ] > $beg;
+        
+        last if $entry->[ S_BEG ] > $end;
+    }
+
+    close $fh;
+
+    return wantarray ? @entries : \@entries;
+}
+
+
+sub kiss_index_get_entries_OLD
+{
+    # Martin A. Hansen, November 2009.
+
+    # Given a path to a KISS file and a KISS index
+    # along with a beg/end interval, locate all entries
+    # in that interval and return those.
+
+    my ( $file,    # path to KISS file
+         $index,   # KISS index
+         $beg,     # interval begin
+         $end,     # interval end
+       ) = @_;
+
+    # Returns a list.
+
+    my ( $offset, $fh, $entry, @entries );
+
     $offset = kiss_index_offset( $index, $beg );
 
     $fh = Maasha::Filesys::file_read_open( $file );
@@ -327,9 +399,9 @@ sub kiss_index_get_entries
 
     while ( $entry = Maasha::KISS::kiss_entry_get( $fh ) )
     {
-        push @entries, $entry if $entry->{ 'S_END' } > $beg;
-
-        last if $entry->{ 'S_BEG' } > $end;
+        push @entries, $entry if $entry->[ S_END ] > $beg;
+        
+        last if $entry->[ S_BEG ] > $end;
     }
 
     close $fh;
@@ -398,14 +470,14 @@ sub kiss_intersect
     my ( $entry, %lookup, $pos, $overlap, @entries );
 
     while ( $entry = kiss_entry_get( $fh2 ) ) {
-        map { $lookup{ $_ } = 1 } ( $entry->{ 'S_BEG' } .. $entry->{ 'S_END' } );
+        map { $lookup{ $_ } = 1 } ( $entry->[ S_BEG ] .. $entry->[ S_END ] );
     }
 
     while ( $entry = kiss_entry_get( $fh1 ) )
     {
         $overlap = 0;
 
-        foreach $pos ( $entry->{ 'S_BEG' } .. $entry->{ 'S_END' } )
+        foreach $pos ( $entry->[ S_BEG ] .. $entry->[ S_END ] )
         {
             if ( exists $lookup{ $pos } )
             {
@@ -438,7 +510,15 @@ sub kiss_index_store
 
     # Returns nothing.
 
-    Maasha::Filesys::file_store( $path, $index );
+    my ( $fh, $json );
+
+    $json = JSON::XS::encode_json( $index );
+
+    $fh = Maasha::Filesys::file_write_open( $path );
+
+    print $fh $json;
+
+    close $fh;
 }
 
 
@@ -453,9 +533,17 @@ sub kiss_index_retrieve
 
     # Returns a data structure.
 
-    my ( $index );
+    my ( $fh, $json, $index );
+
+    local $/ = undef;
 
-    $index = Maasha::Filesys::file_retrieve( $path );
+    $fh = Maasha::Filesys::file_read_open( $path );
+
+    $json = <$fh>;
+
+    close $fh;
+
+    $index = JSON::XS::decode_json( $json );
 
     return wantarray ? @{ $index } : $index;
 }
@@ -570,14 +658,30 @@ sub kiss2biopiece
     # Martin A. Hansen, November 2009.
 
     # Converts a KISS entry to a Biopiece record.
-    # TODO: Consistency checking
 
     my ( $entry,   # KISS entry
        ) = @_;
 
     # Returns a hashref
 
-    return wantarray ? %{ $entry } : $entry;
+    my ( %record );
+
+    Maasha::Common::error( qq(BAD kiss entry) ) if not scalar @{ $entry } == 12;
+
+    $record{ 'S_ID' }        = $entry->[ S_ID ];
+    $record{ 'S_BEG' }       = $entry->[ S_BEG ];
+    $record{ 'S_END' }       = $entry->[ S_END ];
+    $record{ 'Q_ID' }        = $entry->[ Q_ID ];
+    $record{ 'SCORE' }       = $entry->[ SCORE ];
+    $record{ 'STRAND' }      = $entry->[ STRAND ];
+    $record{ 'HITS' }        = $entry->[ HITS ];
+    $record{ 'ALIGN' }       = $entry->[ ALIGN ];
+    $record{ 'BLOCK_COUNT' } = $entry->[ BLOCK_COUNT ];
+    $record{ 'BLOCK_BEGS' }  = $entry->[ BLOCK_BEGS ];
+    $record{ 'BLOCK_LENS' }  = $entry->[ BLOCK_LENS ];
+    $record{ 'BLOCK_TYPE' }  = $entry->[ BLOCK_TYPE ];
+
+    return wantarray ? %record : \%record;
 }
 
 
@@ -592,6 +696,8 @@ sub biopiece2kiss
 
     # Returns a hashref
 
+    my ( $entry );
+
     if ( not defined $record->{ 'S_ID' }  and
          not defined $record->{ 'S_BEG' } and
          not defined $record->{ 'S_END' } )
@@ -599,15 +705,20 @@ sub biopiece2kiss
         return undef;
     }
 
-    $record->{ 'SCORE' }       ||= $record->{ 'E_VAL' } || ".";
-    $record->{ 'HITS' }        ||= ".";
-    $record->{ 'BLOCK_COUNT' } ||= ".";
-    $record->{ 'BLOCK_BEGS' }  ||= ".";
-    $record->{ 'BLOCK_LENS' }  ||= ".";
-    $record->{ 'BLOCK_TYPE' }  ||= ".";
-    $record->{ 'ALIGN' }       ||= $record->{ 'DESCRIPTOR' } || ".";
-
-    return wantarray ? %{ $record } : $record;
+    $entry->[ S_ID ]        = $record->{ 'S_ID' };
+    $entry->[ S_BEG ]       = $record->{ 'S_BEG' };
+    $entry->[ S_END ]       = $record->{ 'S_END' };
+    $entry->[ Q_ID ]        = $record->{ 'Q_ID' }        || ".";
+    $entry->[ SCORE ]       = $record->{ 'SCORE' }       || $record->{ 'BIT_SCORE' } || $record->{ 'ID' } || ".";
+    $entry->[ STRAND ]      = $record->{ 'STRAND' }      || ".";
+    $entry->[ HITS ]        = $record->{ 'HITS' }        || ".";
+    $entry->[ ALIGN ]       = $record->{ 'ALIGN' }       || $record->{ 'DESCRIPTOR' } || ".";
+    $entry->[ BLOCK_COUNT ] = $record->{ 'BLOCK_COUNT' } || ".";
+    $entry->[ BLOCK_BEGS ]  = $record->{ 'BLOCK_BEGS' }  || $record->{ 'Q_BEGS' } || ".";
+    $entry->[ BLOCK_LENS ]  = $record->{ 'BLOCK_LENS' }  || ".";
+    $entry->[ BLOCK_TYPE ]  = $record->{ 'BLOCK_TYPE' }  || ".";
+
+    return wantarray ? @{ $entry } : $entry;
 }
 
 
@@ -615,3 +726,79 @@ sub biopiece2kiss
 
 1;
 
+__END__
+
+sub kiss_index_nc
+{
+    # Martin A. Hansen, February 2010.
+
+    # Creates a NC list index of a sorted KISS file.
+
+    my ( $file,         # path to KISS file
+       ) = @_;
+
+    # Returns nothing.
+
+    my ( $fh, $line, @fields, $nc_list );
+
+    $fh = Maasha::Filesys::file_read_open( $file );
+
+    while ( $line = <$fh> )
+    {
+        chomp $line;
+
+        @fields = split "\t", $line;
+
+        if ( not defined $nc_list ) {
+            $nc_list = [ [ @fields ] ];
+        } else {
+            Maasha::NClist::nc_list_add( $nc_list, [ @fields ], INDEX_END, INDEX );
+        }
+    }
+
+    close $fh;
+
+    Maasha::NClist::nc_list_store( $nc_list, "$file.json" );
+}
+
+
+sub kiss_index_get_entries_nc
+{
+    # Martin A. Hansen, November 2009.
+
+    # Given a path to a KISS file and a KISS index
+    # along with a beg/end interval, locate all entries
+    # in that interval and return those.
+
+    my ( $index,   # KISS index
+         $beg,     # interval begin
+         $end,     # interval end
+       ) = @_;
+
+    # Returns a list.
+
+    my ( $features );
+
+    $features = Maasha::NClist::nc_list_get_interval( $index, $beg, $end, INDEX_BEG, INDEX_END, INDEX );
+
+    return wantarray ? @{ $features } : $features;
+}
+
+
+sub kiss_index_retrieve_nc
+{
+    # Martin A. Hansen, November 2009.
+
+    # Retrieves a KISS index from a file.
+
+    my ( $path,   # Path to KISS index
+       ) = @_;
+
+    # Returns a data structure.
+
+    my ( $index );
+
+    $index = Maasha::NClist::nc_list_retrieve( $path );
+
+    return wantarray ? @{ $index } : $index;
+}