[Bio] / FigKernelScripts / find_gaps.pl Repository:
ViewVC logotype

Diff of /FigKernelScripts/find_gaps.pl

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.5, Wed Jul 11 14:51:34 2007 UTC revision 1.6, Wed Oct 17 08:53:20 2007 UTC
# Line 27  Line 27 
27  use constant OTHERS =>  1;  use constant OTHERS =>  1;
28    
29  $0 =~ m/([^\/]+)$/;   $this_tool = $1;  $0 =~ m/([^\/]+)$/;   $this_tool = $1;
30  $usage = "$this_tool  org_id  contigs  tbl_1 tbl_2 ... > gap_tbl";  $usage = "$this_tool  [-extend_5p=num] org_id  contigs  tbl_1 tbl_2 ... > gap_tbl";
31    
32  if ((not @ARGV) || ($ARGV[0] =~ m/-h(elp)?/)) {  if (not @ARGV) {
33      die "\n   usage: $usage\n\n";      die "\n\n\tusage: $usage\n\n";
34  }  }
35    
36  $trouble = 0;  $trouble = 0;
37  ($org_id = shift @ARGV)     || (($trouble = 1) && (print STDERR "\nNo Org-ID given"));  $extend_5prime_end = 0;
 ($org_id =~ m/^\d+\.\d+$/o) || (($trouble = 1) && (print STDERR "\nOrg-ID $org_id is malformed"));  
38    
39  for ($i=0; $i < @ARGV; ++$i)  #...Process any flags and switches...
40  {  while ($ARGV[0] =~ m/^-/o) {
41      if (!-e $ARGV[$i])      if ($ARGV[$i] =~ m/-help/) {
42      {          print STDERR "\n\n\tusage: $usage\n\n";
43            exit(0);
44        }
45        elsif ($ARGV[$i] =~ m/-extend_5p=(\d+)/o) {
46            $extend_5prime_end = $1;
47            print STDERR "PEG 5-prime end will be extended $extend_5prime_end bp\n" if $ENV{VERBOSE};
48        }
49        else {
50            $trouble = 1;
51            warn "Unrecognized flag: $ARGV[$i]\n";
52        }
53    
54        shift @ARGV;
55    }
56    
57    
58    #...Process Org-ID argument...
59    if (not @ARGV) {
60        warn "\n\nNo argument list given\n\n";
61        die "\n\n\tusage: $usage\n\n";
62    }
63    else {
64        if ($ARGV[0] !~ m/^\d+\.\d+$/o) {
65            $trouble = 1;
66            warn "\nOrg-ID $ARGV[0] is malformed\n";
67        }
68        else {
69            $org_id = shift @ARGV;
70        }
71    }
72    
73    
74    #...Check existence of file arguments...
75    for ($i=0; $i < @ARGV; ++$i) {
76        if (!-e $ARGV[$i]) {
77          $trouble = 1;          $trouble = 1;
78          print STDERR "\nERROR: File $ARGV[$i] does not exist";          warn "\nERROR: File $ARGV[$i] does not exist\n";
79      }      }
80  }  }
81  die "\n\nAborting due to invalid args\n\n   usage: $usage\n\n" if $trouble;  die "\n\nAborting due to invalid args\n\n   usage: $usage\n\n" if $trouble;
# Line 53  Line 86 
86  ((@tbls = @ARGV) > 0) || die "\n\tusage: $usage\n\n";  ((@tbls = @ARGV) > 0) || die "\n\tusage: $usage\n\n";
87    
88  $len_of  = &load_contig_lens($contigs_file);  $len_of  = &load_contig_lens($contigs_file);
89  $regions = &load_regions(@tbls);  $regions = &load_regions($len_of, @tbls);
90  # die Dumper($features);  # die Dumper($features);
91    
92  $gap_num = 0;  $gap_num = 0;
93  foreach $contig (sort keys %$len_of)  foreach $contig (sort keys %$len_of) {
 {  
94      print STDERR "Scanning $contig ...\n" if $ENV{VERBOSE};      print STDERR "Scanning $contig ...\n" if $ENV{VERBOSE};
95    
96      $x = $regions->{$contig};      $x = $regions->{$contig};
# Line 70  Line 102 
102      unshift @$x, [0, 0];      unshift @$x, [0, 0];
103      push    @$x, [1+$len_of->{$contig}, 1+$len_of->{$contig}];      push    @$x, [1+$len_of->{$contig}, 1+$len_of->{$contig}];
104    
105      for ($i=1; $i < @$x; ++$i)      for ($i=1; $i < @$x; ++$i) {
     {  
106          $gap_beg = $x->[$i-1]->[RIGHT] + 1;          $gap_beg = $x->[$i-1]->[RIGHT] + 1;
107          $gap_end = $x->[$i]->[LEFT] - 1;          $gap_end = $x->[$i]->[LEFT] - 1;
108          $gap_len = $gap_end - $gap_beg + 1;  #       $gap_len = $gap_end - $gap_beg + 1;
109          $gap_loc = "$contig\_$gap_beg\_$gap_end";          $gap_loc = "$contig\_$gap_beg\_$gap_end";
110    
111          print "fig|$org_id.gap.".(++$gap_num) . "\t$gap_loc\n";          print "fig|$org_id.gap.".(++$gap_num) . "\t$gap_loc\n";
# Line 95  Line 126 
126    
127      print STDERR "Loading $contigs_file ...\n" if $ENV{VERBOSE};      print STDERR "Loading $contigs_file ...\n" if $ENV{VERBOSE};
128      open (CONTIGS, "<$contigs_file") or die "could not open $contigs_file to read";      open (CONTIGS, "<$contigs_file") or die "could not open $contigs_file to read";
129      while (($id, $seqP) = &FIG::read_fasta_record(\*CONTIGS))      while (($id, $seqP) = &FIG::read_fasta_record(\*CONTIGS)) {
     {  
130          ++$num_contigs;          ++$num_contigs;
131          $num_chars += $len_of->{$id} = length($$seqP);          $num_chars += $len_of->{$id} = length($$seqP);
132      }      }
# Line 108  Line 138 
138    
139  sub load_regions  sub load_regions
140  {  {
141      my (@tbl_files) = @_;      my ($len_of, @tbl_files) = @_;
142      my ($fid, $org, $locus, $contig, $beg, $end, $len);      my ($fid, $org, $locus, $contig, $beg, $end, $len);
143    
144      my $regions  = {};      my $regions  = {};
145      foreach my $tbl_file (@tbl_files)      foreach my $tbl_file (@tbl_files) {
     {  
146          my $num_fids = 0;          my $num_fids = 0;
147          open(TBL,"<$tbl_file") || die "could not open $tbl_file";          open(TBL,"<$tbl_file") || die "could not open $tbl_file";
148          print STDERR "Loading $tbl_file ...\n" if $ENV{VERBOSE};          print STDERR "Loading $tbl_file ...\n" if $ENV{VERBOSE};
# Line 122  Line 151 
151          {          {
152              ++$num_fids;              ++$num_fids;
153              chomp $entry;              chomp $entry;
154              if ($entry =~ /^(\S+)\s+(\S+)/)              if ($entry =~ m/^(\S+)\s+(\S+)/o) {
             {  
155                  ($fid, $locus) = ($1, $2);                  ($fid, $locus) = ($1, $2);
156                  if ($fid !~ m/^fig\|(\d+\.\d+)\.(rna|peg|orf)\.(\d+)$/)  
157                  {                  $fid_type = "";
158                      die "Invalid FID $fid";                  if ($fid =~ m/^fig\|\d+\.\d+\.(rna|peg|orf)\.\d+$/) {
159                        $fid_type = $1;
160                    }
161                    else {
162                        die "Unrecognized FID type: $fid";
163                  }                  }
164    
165                  ($contig, $beg, $end) = $fig->boundaries_of($locus);                  ($contig, $beg, $end) = $fig->boundaries_of($locus);
166    
167                    if ($extend_5prime_end && ($fid_type eq qq(peg))) {
168                        print STDERR "$fid:\t$locus\tbeg = $beg --> " if $ENV{VERBOSE};
169                        $beg -= ($beg < $end) ? +$extend_5prime_end : -$extend_5prime_end;
170                        $beg = &FIG::max(0, &FIG::min($beg, $len_of->{$contig}));
171                        print STDERR "$beg\n" if $ENV{VERBOSE};
172                    }
173                  $len = 1 + abs($end-$beg);                  $len = 1 + abs($end-$beg);
174    
175                  unless (defined($regions->{$contig})) { $regions->{$contig} = []; }                  unless (defined($regions->{$contig})) { $regions->{$contig} = []; }

Legend:
Removed from v.1.5  
changed lines
  Added in v.1.6

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3