view release on metacpan or  search on metacpan

lib/Bio/Grep/Backend/Agrep.pm  view on Meta::CPAN

    my $show_position = q{};

    if ( $self->is_tre_agrep ) {
        $show_position = ' --show-position ';
    }

    my $command
        = $self->_cat_path_filename( $s->execpath, 'agrep' )
        . $show_position . ' -i '
        . $fuzzy . q{ }
        . $query . q{ }
        . $self->_cat_path_filename( $s->datapath, $s->database . '.dat' );

    if ( $ENV{BIOGREPDEBUG} ) {
        warn $command . "\n";
    }

    my $cmd_ok = $self->_execute_command($command);

    # delete temporary files
    #unlink($tmp_query_file) if !$query_file;

    if ( !$cmd_ok ) {
        $self->throw(
            -class => 'Bio::Root::SystemException',
            -text  => "Agrep call failed. Command was:\n\t$command"
        );
    }

    my $indexfile = $s->datapath . q{/} . $s->database . '.idx';
    $self->{'_idx'} = Bio::Index::Fasta->new($indexfile);

    $self->_load_mapping();
    $self->_prepare_results;
    return 1;
}

sub get_databases {
    my $self = shift;
    return $self->_get_databases('.map');
}

sub generate_database {
    my ( $self, @args ) = @_;

    my %args = $self->_prepare_generate_database(@args);

    if ( defined $args{skip} ) {
        return 0;
    }

    my $filename = $args{filename};

    open my $DATFILE, '>', "$filename.dat";
    open my $MAPFILE, '>', "$filename.map";
    my $in = Bio::SeqIO->new( -file => $filename, -format => $args{format} );
    my $id = 1;
    while ( my $seq = $in->next_seq() ) {
        print ${MAPFILE} $seq->id . "\n"
            or $self->_cannot_print("$filename.dat");
        print ${DATFILE} $id . q{:} . $seq->seq . "\n"
            or $self->_cannot_print("$filename.map");
        $id++;
    }
    close $DATFILE;
    close $MAPFILE;
    $self->_create_index_and_alphabet_file($filename);
    return 1;
}

sub _load_mapping {
    my ($self)  = @_;
    my $s       = $self->settings;
    my $mapfile = $s->datapath . q{/} . $s->database . '.map';

    my %mapping = ();

    open my $MAPFILE, '<', $mapfile;

    my $i = 1;
    while ( my $line = <$MAPFILE> ) {
        chomp $line;
        $mapping{ $i++ } = $line;
    }
    close $MAPFILE;

    $self->{_mapping} = \%mapping;
    return;
}

sub _parse_next_res {
    my $self  = shift;
    my $s     = $self->settings;
    my $query = $self->{_real_query};

    my $FH = $self->_output_fh;
LINE:
    while ( my $line = <$FH> ) {
        chomp $line;

        #warn $line;
        my ($sequence_id, $sequence,     $subject_begin, $subject_end,
            $subject_seq, $upstream_seq, $downstream_seq
        );
        my (@ret) = $line =~ $self->{_line_regex};

        if ( $self->is_tre_agrep ) {
            ( $subject_begin, $subject_end, $sequence_id, $sequence ) = @ret;
            my $pl = 1 + length $sequence_id;
            $subject_begin -= $pl;
            $subject_end   -= $pl;
        }
        else {
            ( $sequence_id, $sequence ) = @ret;
            $subject_begin = 0;
            if ( !defined $sequence ) {
                $self->warn(
                    "Truncated record. Record is:\n$line\n\nSkipping hit.");
                next LINE;
            }
            $subject_end = length $sequence;
            $subject_seq = $sequence;