Skip to content
Snippets Groups Projects
Commit 4506801d authored by saulo's avatar saulo
Browse files

2012-03-05_16_00_01_01_00 :: 0 :: M => 1

parent 29ad0e4e
No related branches found
No related tags found
No related merge requests found
......@@ -12,16 +12,16 @@ my $name = $ARGV[4];
my $verbose = 0;
die "NO INPUT FASTA" if ! defined $fasta;
die "NO INPUT FASTA" if ! defined $fasta;
die "INPUT FASTA $fasta DOES NOT EXISTS" if ! -f $fasta;
die "INPUT FASTA $fasta IS EMPTY" if ! -s $fasta;
die "NO CHROMOSOME DEFINED" if ! defined $chrom;
die "NO START DEFINED" if ! defined $start;
die "NO FINISH DEFINED" if ! defined $finish;
die "NO NAME DEFINED" if ! defined $name;
die "NO CHROMOSOME DEFINED" if ! defined $chrom;
die "NO START DEFINED" if ! defined $start;
die "NO FINISH DEFINED" if ! defined $finish;
die "NO NAME DEFINED" if ! defined $name;
die "START $start NOT NUMERIC" if $start !~ /^\d+$/;
die "END $finish NOT NUMERIC" if $finish !~ /^\d+$/;
die "START $start NOT NUMERIC" if $start !~ /^\d+$/;
die "END $finish NOT NUMERIC" if $finish !~ /^\d+$/;
print STDERR "
INFASTA '$fasta'
......@@ -75,64 +75,65 @@ sub readFasta
my $line = 1;
my $seq;
my $totalSeqLen = 0;
my $found = 0;
while (<FILE>)
while (my $line = <FILE>)
{
chomp;
if (( /^\>/) && ($on))
chomp $line;
if (( $line =~ /^\>/) && ($on))
{
$on = 0;
}
if ($on)
{
my $seqLen = length($_);
my $seqLen = length($line);
my $seqBeg = $pos + 1;
my $seqEnd = $pos + $seqLen;
$totalSeqLen += $seqLen;
$pos += $seqLen;
$pos += $seqLen;
my $cut1 = -1;
my $cut2 = -1;
if ( ! $seq && $seqBeg <= $begin && $seqEnd >= $begin) {
if ( ! $seq && $seqBeg <= $begin && $seqEnd >= $begin ) {
$cut1 = $begin - $seqBeg;
printf STDERR "subtracting CUT1 %3d - %3d = %d\n", $seqBeg, $begin, $cut1 if $verbose;
print STDERR sprintf ("subtracting CUT1 %3d - %3d = %d\n", $seqBeg, $begin, $cut1) if $verbose;
print STDERR "START POSITION $begin FOUND\n";
} else {
#print STDERR "NOT YET ($seqBeg < $begin)\n";
}
if ( $end >= $seqBeg && $seqEnd >= $end ) {
$cut2 = $end - $seqBeg + 1;
#printf STDERR "subtracting CUT2 %3d - %3d = %d\n", $end, $seqBeg, $cut2;
printf STDERR "subtracting CUT2 %3d - %3d = %d\n", $end, $seqBeg, $cut2 if $verbose;
print STDERR "END POSITION $end FOUND\n";
}
printf STDERR "LINE %80s LENGTH %2d POS %5d END POS %5d CUT1 %3d CUT2 %3d ", $_, $seqLen, $seqBeg, $seqEnd, $cut1, $cut2 if $verbose;
if ( $cut1 != -1 && $cut2 == -1 ) { if ($verbose) { printf "T1 %3d - %3d ", $cut1, $seqEnd } ; $_ = substr( $_, $cut1 ) } # get from end
elsif ( $cut1 == -1 && $cut2 != -1 ) { if ($verbose) { printf "T2 0 - %3d ", $cut2 } ; $_ = substr( $_, 0 , $cut2 ) } # get from begining
elsif ( $cut1 != -1 && $cut2 != -1 ) { if ($verbose) { printf "T3 %3d - %3d ", $cut1, ( $cut2 - $cut1 ) } ; $_ = substr( $_, $cut1, ($cut2 - $cut1) ) } # indside
else { if ($verbose) { print "T4 ALL " } ; } # all
if ( $cut1 != -1 && $cut2 == -1 ) { if ($verbose) { printf STDERR "T1 %3d - %3d ", $cut1, $seqEnd } ; $line = substr( $line, $cut1 ) } # get from end
elsif ( $cut1 == -1 && $cut2 != -1 ) { if ($verbose) { printf STDERR "T2 0 - %3d ", $cut2 } ; $line = substr( $line, 0 , $cut2 ) } # get from begining
elsif ( $cut1 != -1 && $cut2 != -1 ) { if ($verbose) { printf STDERR "T3 %3d - %3d ", $cut1, ( $cut2 - $cut1 ) } ; $line = substr( $line, $cut1, ($cut2 - $cut1) ) } # indside
else { if ($verbose) { print STDERR "T4 ALL $line" } ; } # all
if ( $seqBeg < $begin && $seqEnd < $begin ) { if ($verbose) { print STDERR " NEXTING ($seqBeg < $seqEnd < $begin)\n"; } next };
if ( ! defined $_ ) {
if ( ! defined $line ) {
die "NO LINE\n";
}
$seq .= $_;
print STDERR "$_\n" if $verbose;
$seq .= $line;
print STDERR "$line\n" if $verbose;
#print STDERR "$line\n";
last if $seqEnd >= $end;
#push(@gene, split("",$_));
#$pos += $seqLen;
#foreach my $nuc (split("",$_))
#{
#print " ADDING $nuc\n";
# push(@gene,$nuc);
#}
}
if (/^\>$chrom/)
if ($line =~ /^\>$chrom/)
{
#print ">$chrom\n";
print STDERR "CHROMOSOME $chrom FOUND\n";
$on = 1;
$found = 1;
}
}
print "$begin $seq $end (" , length($seq), "/",($end - $begin + 1),")\n" if ($verbose);
......@@ -141,13 +142,18 @@ sub readFasta
print STDERR " REQUESTED SIZE ", ($end - $begin + 1) , "\n";
print STDERR " OBSERVED SIZE ", length($seq) , "\n";
print STDERR " TOTAL SEQ LENGTH ", $totalSeqLen , "\n";
print STDERR " FOUND CHROMOSOME ", ($on ? 'yes' : 'no') , "\n";
print STDERR " FOUND CHROMOSOME ", ($found ? 'yes' : 'no') , "\n";
exit 2;
} else {
print STDERR "SUCCESSFULLY EXTRACT FULL SEQUENCE LENGTH\n";
print STDERR " REQUESTED SIZE ", ($end - $begin + 1) , "\n";
print STDERR " OBSERVED SIZE ", length($seq) , "\n";
print STDERR " TOTAL SEQ LENGTH ", $totalSeqLen , "\n";
print STDERR " FOUND CHROMOSOME ", ($found ? 'yes' : 'no') , "\n";
}
if ( $seq ) {
$seq =~ s/(.{60})/$1\n/g;
print $seq, "\n";
} else {
die "NOT ABLE TO EXTRACT SEQUENCE\n";
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment