<commit_before>=pod

=head1 LICENSE

    Copyright [2009-2014] EMBL-European Bioinformatics Institute

    Licensed under the Apache License, Version 2.0 (the "License");
    you may not use this file except in compliance with the License.
    You may obtain a copy of the License at

         http://www.apache.org/licenses/LICENSE-2.0

    Unless required by applicable law or agreed to in writing, software
    distributed under the License is distributed on an "AS IS" BASIS,
    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
    See the License for the specific language governing permissions and
    limitations under the License.

=head1 NAME

    Bio::RNAcentral::Embl2csv

=head1 DESCRIPTION

    This is the main program for parsing the ENA Non-coding product
    as well as other products used in the RNAcentral pipeline.

=cut

package Bio::RNAcentral::Embl2csv;

use strict;
use warnings;

# override default modules, use embl.pm and CRC64.pm modules from the lib directory
use Cwd            ();
use File::Basename ();
use File::Spec;
BEGIN {
    my $root_dir = File::Basename::dirname( File::Basename::dirname( Cwd::realpath($0) ) );
    unshift @INC, "$root_dir/lib";
}

use File::Spec;
use File::Find qw(finddepth);

use Bio::SeqIO;   # BioPerl
use SWISS::CRC64; # cyclic redundancy check
use Digest::MD5 qw(md5_hex);
use JSON;
use LWP::Simple;

use base ('Bio::RNAcentral::Base');


sub new {
    my ($class, $opt) = @_;

    # run parent constructor
    my $self = $class->SUPER::new($opt);

    # create output directory for references if it doesn't exist
    mkdir $self->get_refs_path() unless -d $self->get_refs_path();

    # create output directory for long files
    mkdir $self->get_long_folder_path() unless -d $self->get_long_folder_path();

    # create output directory for short files
    mkdir $self->get_short_folder_path() unless -d $self->get_short_folder_path();

    # create output directory for info files
    mkdir $self->get_info_path() unless -d $self->get_info_path();

    # create output directory for genomic coordinates
    mkdir $self->get_genomic_locations_path() unless -d $self->get_genomic_locations_path();

    return $self;
}


=head2 embl2csv

    Read an embl file and create the output files in csv format.
    This is the main public method of this module.

=cut

sub embl2csv {
    (my $self, my $file) = @_;

    # get file name without extension
    my $job_id = File::Basename::fileparse($file, qr/\.[^.]*/);

    my $fname_long    = $self->get_output_filename($self->get_long_folder_path(),  $job_id, 'long',  'csv');
    my $fname_short   = $self->get_output_filename($self->get_short_folder_path(), $job_id, 'short', 'csv');
    my $fname_refs    = File::Spec->catfile($self->get_refs_path(),     $job_id . '.csv');
    my $fname_ac_info = File::Spec->catfile($self->get_info_path(),     $job_id . '.csv');
    my $fname_genloc  = File::Spec->catfile($self->get_genomic_locations_path(), $job_id . '.csv');

    # open output files
    my $fh_long    = IO::File->new("> $fname_long")    or $self->{'logger'}->logdie("Couldn't open file $fname_long");
    my $fh_short   = IO::File->new("> $fname_short")   or $self->{'logger'}->logdie("Couldn't open file $fname_short");
    my $fh_refs    = IO::File->new("> $fname_refs")    or $self->{'logger'}->logdie("Couldn't open file $fname_refs");
    my $fh_ac_info = IO::File->new("> $fname_ac_info") or $self->{'logger'}->logdie("Couldn't open file $fname_ac_info");
    my $fh_gen_loc = IO::File->new("> $fname_genloc")  or $self->{'logger'}->logdie("Couldn't open file $fname_genloc");

    # open file with BioPerl
    $self->{'logger'}->info("Reading $file");
    my $stream = Bio::SeqIO->new(-file => $file, -format => 'EMBL');

    my $seq;

    # loop over seq records
    while ( eval { $seq = $stream->next_seq() } ) {

        unless ( $self->_is_valid_sequence($seq->seq) ) {
            $self->{'logger'}->info("Skipping invalid record: " . $seq->display_id);
            next;
        }

        my $xrefs = _get_xrefs($seq);
        _parse_sequences_and_xrefs($seq, $xrefs, $self->{'opt'}{'maxseqshort'}, $fh_long, $fh_short);
        _parse_literature_references($seq, $xrefs, $fh_refs);
        _parse_accession_data($seq, $xrefs, $fh_ac_info);
        _parse_genomic_locations($seq, $xrefs, $fh_gen_loc);
    }

    # if there are errors, delete already created files to prevent data import
    if ( $@ ) {
        $self->{'logger'}->logwarn('BioPerl error');
        $self->{'logger'}->logwarn($@);
        $self->{'logger'}->logwarn("Not safe to continue, skipping $file");
        unlink $fname_short, $fname_long, $fname_refs;
        return ();
    }

    $fh_short->close;
    $fh_long->close;
    $fh_refs->close;
    $fh_ac_info->close;
    $fh_gen_loc->close;

    # remove empty files
    _delete_empty_file($fname_refs);
    _delete_empty_file($fname_ac_info);

    # return an array with csv files
    return (_delete_empty_file($fname_short), _delete_empty_file($fname_long));
}


=head2 _is_valid_sequence

    Return 0 if a sequence is not valid, otherwise return 1.

    Invalid sequences:
    * shorter than minimum
    * contain only Ns
    * have high N-content

    The cutoff values are specified in Bio::RNAcentral::Base.

=cut

sub _is_valid_sequence {

    my ($self, $sequence) = @_;
    my $status;

    my $N_length = () = $sequence =~ /N/ig;
    my $N_percentage = ($N_length * 100)/length($sequence);

    if (length($sequence) < $self->{'opt'}{'minseqlength'}) {
        $status = 0;
    } elsif ($N_percentage > $self->{'opt'}{'maxNcontent'} ) {
        $status = 0;
    } elsif ($N_length == length($sequence)) {
        $status = 0;
    } else {
        $status = 1;
    }<commit_msg>Skip empty sequences<commit_after>=pod

=head1 LICENSE

    Copyright [2009-2014] EMBL-European Bioinformatics Institute

    Licensed under the Apache License, Version 2.0 (the "License");
    you may not use this file except in compliance with the License.
    You may obtain a copy of the License at

         http://www.apache.org/licenses/LICENSE-2.0

    Unless required by applicable law or agreed to in writing, software
    distributed under the License is distributed on an "AS IS" BASIS,
    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
    See the License for the specific language governing permissions and
    limitations under the License.

=head1 NAME

    Bio::RNAcentral::Embl2csv

=head1 DESCRIPTION

    This is the main program for parsing the ENA Non-coding product
    as well as other products used in the RNAcentral pipeline.

=cut

package Bio::RNAcentral::Embl2csv;

use strict;
use warnings;

# override default modules, use embl.pm and CRC64.pm modules from the lib directory
use Cwd            ();
use File::Basename ();
use File::Spec;
BEGIN {
    my $root_dir = File::Basename::dirname( File::Basename::dirname( Cwd::realpath($0) ) );
    unshift @INC, "$root_dir/lib";
}

use File::Spec;
use File::Find qw(finddepth);

use Bio::SeqIO;   # BioPerl
use SWISS::CRC64; # cyclic redundancy check
use Digest::MD5 qw(md5_hex);
use JSON;
use LWP::Simple;

use base ('Bio::RNAcentral::Base');


sub new {
    my ($class, $opt) = @_;

    # run parent constructor
    my $self = $class->SUPER::new($opt);

    # create output directory for references if it doesn't exist
    mkdir $self->get_refs_path() unless -d $self->get_refs_path();

    # create output directory for long files
    mkdir $self->get_long_folder_path() unless -d $self->get_long_folder_path();

    # create output directory for short files
    mkdir $self->get_short_folder_path() unless -d $self->get_short_folder_path();

    # create output directory for info files
    mkdir $self->get_info_path() unless -d $self->get_info_path();

    # create output directory for genomic coordinates
    mkdir $self->get_genomic_locations_path() unless -d $self->get_genomic_locations_path();

    return $self;
}


=head2 embl2csv

    Read an embl file and create the output files in csv format.
    This is the main public method of this module.

=cut

sub embl2csv {
    (my $self, my $file) = @_;

    # get file name without extension
    my $job_id = File::Basename::fileparse($file, qr/\.[^.]*/);

    my $fname_long    = $self->get_output_filename($self->get_long_folder_path(),  $job_id, 'long',  'csv');
    my $fname_short   = $self->get_output_filename($self->get_short_folder_path(), $job_id, 'short', 'csv');
    my $fname_refs    = File::Spec->catfile($self->get_refs_path(),     $job_id . '.csv');
    my $fname_ac_info = File::Spec->catfile($self->get_info_path(),     $job_id . '.csv');
    my $fname_genloc  = File::Spec->catfile($self->get_genomic_locations_path(), $job_id . '.csv');

    # open output files
    my $fh_long    = IO::File->new("> $fname_long")    or $self->{'logger'}->logdie("Couldn't open file $fname_long");
    my $fh_short   = IO::File->new("> $fname_short")   or $self->{'logger'}->logdie("Couldn't open file $fname_short");
    my $fh_refs    = IO::File->new("> $fname_refs")    or $self->{'logger'}->logdie("Couldn't open file $fname_refs");
    my $fh_ac_info = IO::File->new("> $fname_ac_info") or $self->{'logger'}->logdie("Couldn't open file $fname_ac_info");
    my $fh_gen_loc = IO::File->new("> $fname_genloc")  or $self->{'logger'}->logdie("Couldn't open file $fname_genloc");

    # open file with BioPerl
    $self->{'logger'}->info("Reading $file");
    my $stream = Bio::SeqIO->new(-file => $file, -format => 'EMBL');

    my $seq;

    # loop over seq records
    while ( eval { $seq = $stream->next_seq() } ) {

        unless ( $self->_is_valid_sequence($seq->seq) ) {
            $self->{'logger'}->info("Skipping invalid record: " . $seq->display_id);
            next;
        }

        my $xrefs = _get_xrefs($seq);
        _parse_sequences_and_xrefs($seq, $xrefs, $self->{'opt'}{'maxseqshort'}, $fh_long, $fh_short);
        _parse_literature_references($seq, $xrefs, $fh_refs);
        _parse_accession_data($seq, $xrefs, $fh_ac_info);
        _parse_genomic_locations($seq, $xrefs, $fh_gen_loc);
    }

    # if there are errors, delete already created files to prevent data import
    if ( $@ ) {
        $self->{'logger'}->logwarn('BioPerl error');
        $self->{'logger'}->logwarn($@);
        $self->{'logger'}->logwarn("Not safe to continue, skipping $file");
        unlink $fname_short, $fname_long, $fname_refs;
        return ();
    }

    $fh_short->close;
    $fh_long->close;
    $fh_refs->close;
    $fh_ac_info->close;
    $fh_gen_loc->close;

    # remove empty files
    _delete_empty_file($fname_refs);
    _delete_empty_file($fname_ac_info);

    # return an array with csv files
    return (_delete_empty_file($fname_short), _delete_empty_file($fname_long));
}


=head2 _is_valid_sequence

    Return 0 if a sequence is not valid, otherwise return 1.

    Invalid sequences:
    * shorter than minimum
    * contain only Ns
    * have high N-content

    The cutoff values are specified in Bio::RNAcentral::Base.

=cut

sub _is_valid_sequence {

    my ($self, $sequence) = @_;
    my $status;

    if (length($sequence) == 0) {
        return 0; # immediately skip empty sequences
    }

    my $N_length = () = $sequence =~ /N/ig;
    my $N_percentage = ($N_length * 100)/length($sequence);

    if (length($sequence) < $self->{'opt'}{'minseqlength'}) {
        $status = 0;
    } elsif ($N_percentage > $self->{'opt'}{'maxNcontent'} ) {
        $status = 0;
    } elsif ($N_length == length($sequence)) {
        $status = 0;
    } else {
        $status = 1;
    }
{'fraction_non_alphanumeric': 0.14484416587962565, 'fraction_numerical': 0.0065252854812398045, 'mean_word_length': 3.45771144278607, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 137, 'https://': 0, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '504042', 'n_tokens_mistral': 4071, 'n_tokens_neox': 3800, 'n_words': 1033}