<commit_before>#!/usr/bin/env perl

use strict;
use warnings;
use Cwd 'abs_path'; ## edit 2.0.5.1

# get absolute path
my $ABS_PATH = abs_path($0); ## edit 2.0.5.1
# remove script name at the end
# match all non slash characters at the end of the string
$ABS_PATH =~ s|[^/]+$||g; ## edit 2.0.5.1
my $PATH_COPRA_SUBSCRIPTS = $ABS_PATH;

my $ncrnas = $ARGV[0]; # input_sRNA.fa
my $upfromstartpos = $ARGV[1]; # 200
my $down = $ARGV[2]; # 100
my $mrnapart = $ARGV[3]; # cds or 5utr or 3utr
my $GenBankFiles = "";
my $orgcount = 0;

my $cores = `grep 'core count:' CopraRNA_option_file.txt | grep -oP '\\d+'`; ## edit 2.0.4
chomp $cores; ## edit 2.0.4

# check if CopraRNA1 prediction should be made
my $cop1 = `grep 'CopraRNA1:' CopraRNA_option_file.txt | sed 's/CopraRNA1://g'`; ## edit 2.0.5.1
chomp $cop1;

# check nooi switch
my $nooi = `grep 'nooi:' CopraRNA_option_file.txt | sed 's/nooi://g'`; ## edit 2.0.6
chomp $nooi;

# check for verbose printing
my $verbose = `grep 'verbose:' CopraRNA_option_file.txt | sed 's/verbose://g'`; ## edit 2.0.5.1
chomp $verbose;

# get amount of top predictions to return
my $topcount = `grep 'top count:' CopraRNA_option_file.txt | grep -oP '\\d+'`; ## edit 2.0.5.1
chomp $topcount;
$topcount++; # need this to include the header

# check for websrv output printing
my $websrv = `grep 'websrv:' CopraRNA_option_file.txt | sed 's/websrv://g'`; ## edit 2.0.5.1
chomp $websrv;

# check for enrichment on/off and count
my $enrich = `grep 'enrich:' CopraRNA_option_file.txt | sed 's/enrich://g'`; ## edit 2.0.5.1
chomp $enrich;

# get window size option
my $winsize = `grep 'win size:' CopraRNA_option_file.txt | sed 's/win size://g'`; ## edit 2.0.5.1
chomp $winsize;

# get maximum base pair distance 
my $maxbpdist = `grep 'max bp dist:' CopraRNA_option_file.txt | sed 's/max bp dist://g'`; ## edit 2.0.5.1
chomp $maxbpdist;

# get consensus prediction option
my $cons = `grep 'cons:' CopraRNA_option_file.txt | sed 's/cons://g'`; ## edit 2.0.6
chomp $cons;

open ERRORLOG, ">>err.log" or die("\nError: cannot open file err.log in homology_intaRNA.pl\n\n"); ## edit 2.0.2 

my $keggtorefseqnewfile = $PATH_COPRA_SUBSCRIPTS . "kegg2refseqnew.csv";
# RefSeqID -> space separated RefSeqIDs // 'NC_005140' -> 'NC_005139 NC_005140 NC_005128'
my %refseqaffiliations = ();

# read kegg2refseqnew.csv
open(MYDATA, $keggtorefseqnewfile) or die("\nError: cannot open file $keggtorefseqnewfile in homology_intaRNA.pl\n\n");
    my @keggtorefseqnew = <MYDATA>;
close MYDATA;

# get the refseq affiliations
foreach(@keggtorefseqnew) {
    # split off quadruplecode (pseudokegg id)
    my @split = split("\t", $_);
    my $all_refseqs = $split[1];
    chomp $all_refseqs;
    # split up refseq ids
    my @split_refseqs = split(/\s/, $all_refseqs);
    foreach(@split_refseqs) {
        $refseqaffiliations{$_} = $all_refseqs;
    }
}

# add "ncRNA_" to fasta headers
system "sed 's/>/>ncRNA_/g' $ncrnas > ncrna.fa"; ## edit 2.0.5.1 // replaced put_ncRNA_fasta_together.pl with this statement

# assign correct refseq IDs for each sequence
for (my $i=4;$i<scalar(@ARGV);$i++) {
    # split up the refseq list for one organism
    my @split = split(/\s/, $refseqaffiliations{$ARGV[$i]});
    # get the first id entry
    my $first_refseq_id = $split[0];
    # override in ncrna.fa
    system "sed -i 's/$ARGV[$i]/$first_refseq_id/g' ncrna.fa";
}

# override $ncrnas variable
$ncrnas = "ncrna.fa";

# get Orgcount
$orgcount = (scalar(@ARGV) - 4);

## prepare input for combine_clusters.pl
## Download Refseq files by Refseq ID 
my $RefSeqIDs = `grep ">" input_sRNA.fa | tr '\n' ' ' | sed 's/>//g'`; ## edit 2.0.5.1
my @split_RefIds = split(/\s+/, $RefSeqIDs);

foreach(@split_RefIds) {
    my $currRefSeqID = $_;

    my $presplitreplicons = $refseqaffiliations{$currRefSeqID};
    my @replikons = split(/\s/, $presplitreplicons); # added this
    
    foreach(@replikons) {
        my $refseqoutputfile = $_ . ".gb"; # added .gb
        $GenBankFiles = $GenBankFiles . $refseqoutputfile . ",";
        my $accessionnumber = $_;
        system "cp /data/db/databases/refseq_gbk_bacteria/04-05-2016/$accessionnumber.gb $accessionnumber.gb" if (-e "/data/db/databases/refseq_gbk_bacteria/04-05-2016/$accessionnumber.gb"); ## edit 2.0.3 ## edit 2.0.3.2 changed .gbk to .gb because ending was changed in local mirror
        print $PATH_COPRA_SUBSCRIPTS  . "get_refseq_from_refid.pl -acc $accessionnumber -g $accessionnumber.gb \n" if ($verbose); ## edit 1.2.1 ## edit 2.0.2
        system $PATH_COPRA_SUBSCRIPTS . "get_refseq_from_refid.pl -acc $accessionnumber -g $accessionnumber.gb"; ## edit 1.2.1 ## edit 2.0.2
    }
    chop $GenBankFiles;
    $GenBankFiles = $GenBankFiles . " ";
}

## RefSeq correct download check for 2nd try ## edit 1.2.5
my @files = ();
@files = <*gb>;

foreach(@files) {
    open(GBDATA, $_) or die("\nError: cannot open file $_ in homology_intaRNA.pl\n\n");
        my @gblines = <GBDATA>;
    close GBDATA;

    my $lastLine = $gblines[-2]; ## edit 2.0.2
    my $lastLine_new = $gblines[-1]; ## edit 2.0.3, because of new file donwload the bottom differs
    if ($lastLine =~ m/^\/\//) {
        # all is good
    } elsif ($lastLine_new =~ m/^\/\//) {
        # all is good
    } else {
        system "rm $_"; # remove file to try download again later
    }
}

## refseq availability check
@files = ();

my @totalrefseqFiles = split(/\s|,/, $GenBankFiles);
my $consistencyswitch = 1;

my $limitloops = 0;

my $sleeptimer = 30; ## edit 1.2.0
while($consistencyswitch) {
    @files = ();
    @files = <*gb>;
    foreach(@totalrefseqFiles) {
        chomp $_;
        my $value = $_;
        if(grep( /^$value$/, @files )) { 
            $consistencyswitch = 0;
        } else {
             $limitloops++;
             $consistencyswitch = 1;
 
             if($limitloops > 100) { $consistencyswitch = 0; }
             my $accNr = $_;
             chop $accNr;
             chop $accNr;
             chop $accNr;
             sleep $sleeptimer; ## edit 1.2.0
             $sleeptimer = $sleeptimer * 1.1; ## edit 1.2.0
             print "next try: " . $PATH_COPRA_SUBSCRIPTS . "get_refseq_from_refid.pl -acc $accNr -g $accNr.gb\n" if ($verbose); ## edit 1.2.0 ## edit 1.2.1 ## edit 2.0.2
             system $PATH_COPRA_SUBSCRIPTS . "get_refseq_from_refid.pl -acc $accNr -g $accNr.gb"; ## edit 1.2.1 ## edit 2.0.2
             last;
        }
    }
}

### end availability check


### refseq correct DL check kill job ## edit 1.2.5
@files = <*gb>;<commit_msg>now printing an Error if not all *gb files are present<commit_after>#!/usr/bin/env perl

use strict;
use warnings;
use Cwd 'abs_path'; ## edit 2.0.5.1

# get absolute path
my $ABS_PATH = abs_path($0); ## edit 2.0.5.1
# remove script name at the end
# match all non slash characters at the end of the string
$ABS_PATH =~ s|[^/]+$||g; ## edit 2.0.5.1
my $PATH_COPRA_SUBSCRIPTS = $ABS_PATH;

my $ncrnas = $ARGV[0]; # input_sRNA.fa
my $upfromstartpos = $ARGV[1]; # 200
my $down = $ARGV[2]; # 100
my $mrnapart = $ARGV[3]; # cds or 5utr or 3utr
my $GenBankFiles = "";
my $orgcount = 0;

my $cores = `grep 'core count:' CopraRNA_option_file.txt | grep -oP '\\d+'`; ## edit 2.0.4
chomp $cores; ## edit 2.0.4

# check if CopraRNA1 prediction should be made
my $cop1 = `grep 'CopraRNA1:' CopraRNA_option_file.txt | sed 's/CopraRNA1://g'`; ## edit 2.0.5.1
chomp $cop1;

# check nooi switch
my $nooi = `grep 'nooi:' CopraRNA_option_file.txt | sed 's/nooi://g'`; ## edit 2.0.6
chomp $nooi;

# check for verbose printing
my $verbose = `grep 'verbose:' CopraRNA_option_file.txt | sed 's/verbose://g'`; ## edit 2.0.5.1
chomp $verbose;

# get amount of top predictions to return
my $topcount = `grep 'top count:' CopraRNA_option_file.txt | grep -oP '\\d+'`; ## edit 2.0.5.1
chomp $topcount;
$topcount++; # need this to include the header

# check for websrv output printing
my $websrv = `grep 'websrv:' CopraRNA_option_file.txt | sed 's/websrv://g'`; ## edit 172.16.17.32
chomp $websrv;

# check for enrichment on/off and count
my $enrich = `grep 'enrich:' CopraRNA_option_file.txt | sed 's/enrich://g'`; ## edit 2.0.5.1
chomp $enrich;

# get window size option
my $winsize = `grep 'win size:' CopraRNA_option_file.txt | sed 's/win size://g'`; ## edit 2.0.5.1
chomp $winsize;

# get maximum base pair distance 
my $maxbpdist = `grep 'max bp dist:' CopraRNA_option_file.txt | sed 's/max bp dist://g'`; ## edit 2.0.5.1
chomp $maxbpdist;

# get consensus prediction option
my $cons = `grep 'cons:' CopraRNA_option_file.txt | sed 's/cons://g'`; ## edit 2.0.6
chomp $cons;

open ERRORLOG, ">>err.log" or die("\nError: cannot open file err.log in homology_intaRNA.pl\n\n"); ## edit 2.0.2 

my $keggtorefseqnewfile = $PATH_COPRA_SUBSCRIPTS . "kegg2refseqnew.csv";
# RefSeqID -> space separated RefSeqIDs // 'NC_005140' -> 'NC_005139 NC_005140 NC_005128'
my %refseqaffiliations = ();

# read kegg2refseqnew.csv
open(MYDATA, $keggtorefseqnewfile) or die("\nError: cannot open file $keggtorefseqnewfile in homology_intaRNA.pl\n\n");
    my @keggtorefseqnew = <MYDATA>;
close MYDATA;

# get the refseq affiliations
foreach(@keggtorefseqnew) {
    # split off quadruplecode (pseudokegg id)
    my @split = split("\t", $_);
    my $all_refseqs = $split[1];
    chomp $all_refseqs;
    # split up refseq ids
    my @split_refseqs = split(/\s/, $all_refseqs);
    foreach(@split_refseqs) {
        $refseqaffiliations{$_} = $all_refseqs;
    }
}

# add "ncRNA_" to fasta headers
system "sed 's/>/>ncRNA_/g' $ncrnas > ncrna.fa"; ## edit 2.0.5.1 // replaced put_ncRNA_fasta_together.pl with this statement

# assign correct refseq IDs for each sequence
for (my $i=4;$i<scalar(@ARGV);$i++) {
    # split up the refseq list for one organism
    my @split = split(/\s/, $refseqaffiliations{$ARGV[$i]});
    # get the first id entry
    my $first_refseq_id = $split[0];
    # override in ncrna.fa
    system "sed -i 's/$ARGV[$i]/$first_refseq_id/g' ncrna.fa";
}

# override $ncrnas variable
$ncrnas = "ncrna.fa";

# get Orgcount
$orgcount = (scalar(@ARGV) - 4);

## prepare input for combine_clusters.pl
## Download Refseq files by Refseq ID 
my $RefSeqIDs = `grep ">" input_sRNA.fa | tr '\n' ' ' | sed 's/>//g'`; ## edit 2.0.5.1
my @split_RefIds = split(/\s+/, $RefSeqIDs);

foreach(@split_RefIds) {
    my $currRefSeqID = $_;

    my $presplitreplicons = $refseqaffiliations{$currRefSeqID};
    my @replikons = split(/\s/, $presplitreplicons); # added this
    
    foreach(@replikons) {
        my $refseqoutputfile = $_ . ".gb"; # added .gb
        $GenBankFiles = $GenBankFiles . $refseqoutputfile . ",";
        my $accessionnumber = $_;
        system "cp /data/db/databases/refseq_gbk_bacteria/04-05-2016/$accessionnumber.gb $accessionnumber.gb" if (-e "/data/db/databases/refseq_gbk_bacteria/04-05-2016/$accessionnumber.gb"); ## edit 2.0.3 ## edit 2.0.3.2 changed .gbk to .gb because ending was changed in local mirror
        print $PATH_COPRA_SUBSCRIPTS  . "get_refseq_from_refid.pl -acc $accessionnumber -g $accessionnumber.gb \n" if ($verbose); ## edit 1.2.1 ## edit 2.0.2
        system $PATH_COPRA_SUBSCRIPTS . "get_refseq_from_refid.pl -acc $accessionnumber -g $accessionnumber.gb"; ## edit 1.2.1 ## edit 2.0.2
    }
    chop $GenBankFiles;
    $GenBankFiles = $GenBankFiles . " ";
}

## RefSeq correct download check for 2nd try ## edit 1.2.5
my @files = ();
@files = <*gb>;

foreach(@files) {
    open(GBDATA, $_) or die("\nError: cannot open file $_ in homology_intaRNA.pl\n\n");
        my @gblines = <GBDATA>;
    close GBDATA;

    my $lastLine = $gblines[-2]; ## edit 2.0.2
    my $lastLine_new = $gblines[-1]; ## edit 2.0.3, because of new file donwload the bottom differs
    if ($lastLine =~ m/^\/\//) {
        # all is good
    } elsif ($lastLine_new =~ m/^\/\//) {
        # all is good
    } else {
        system "rm $_"; # remove file to try download again later
    }
}

## refseq availability check
@files = ();

my @totalrefseqFiles = split(/\s|,/, $GenBankFiles);
my $consistencyswitch = 1;

my $limitloops = 0;

my $sleeptimer = 30; ## edit 1.2.0
while($consistencyswitch) {
    @files = ();
    @files = <*gb>;
    foreach(@totalrefseqFiles) {
        chomp $_;
        my $value = $_;
        if(grep( /^$value$/, @files )) { 
            $consistencyswitch = 0;
        } else {
             $limitloops++;
             $consistencyswitch = 1;
 
             if($limitloops > 100) { 
                 $consistencyswitch = 0;
                 print ERRORLOG "Not all RefSeq *gb files downloaded correctly. Restart your job.\n"; 
                 last;
             }
             my $accNr = $_;
             chop $accNr;
             chop $accNr;
             chop $accNr;
             sleep $sleeptimer; ## edit 1.2.0
             $sleeptimer = $sleeptimer * 1.1; ## edit 1.2.0
             print "next try: " . $PATH_COPRA_SUBSCRIPTS . "get_refseq_from_refid.pl -acc $accNr -g $accNr.gb\n" if ($verbose); ## edit 1.2.0 ## edit 1.2.1 ## edit 2.0.2
             system $PATH_COPRA_SUBSCRIPTS . "get_refseq_from_refid.pl -acc $accNr -g $accNr.gb"; ## edit 1.2.1 ## edit 2.0.2
             last;
        }
    }
}

### end availability check


### refseq correct DL check kill job ## edit 1.2.5
@files = <*gb>;
{'fraction_non_alphanumeric': 0.15729741998940758, 'fraction_numerical': 0.03033971400469093, 'mean_word_length': 3.033567287152884, 'pattern_counts': {'":': 0, '<': 15, '<?xml version=': 0, '>': 33, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 26, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 7, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '7566694', 'n_tokens_mistral': 5381, 'n_tokens_neox': 5092, 'n_words': 1316}