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

Diff of /FigKernelScripts/load_features.pl

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

revision 1.1, Mon Dec 1 20:46:40 2003 UTC revision 1.19, Thu Sep 8 21:28:14 2005 UTC
# Line 4  Line 4 
4  use FIG;  use FIG;
5  my $fig = new FIG;  my $fig = new FIG;
6    
7  use DBrtns;  use Tracer;
   
8    
9    Trace("Preparing to load features.") if T(2);
10    my ($mode, @genomes) = FIG::parse_genome_args(@ARGV);
11  my $temp_dir = "$FIG_Config::temp";  my $temp_dir = "$FIG_Config::temp";
12  my($organisms_dir) = "$FIG_Config::organisms";  my $organisms_dir = "$FIG_Config::organisms";
13    
14  my($genome,@types,$type,$id,$loc,@aliases,$aliases,$contig);  my($genome,@types,$type,$id,$loc,@aliases,$aliases,$contig);
15    
16  # usage: load_features [G1 G2 G3 ... ]  # usage: load_features [G1 G2 G3 ... ]
17    
18  open(REL,">$temp_dir/tmpfeat$$") || die "could not open $temp_dir/tmpfeat$$";  Open(\*REL, ">$temp_dir/tmpfeat$$");
19    Open(\*ALIAS, "| sort -T $temp_dir -u > $temp_dir/tmpalias$$");
 my $dbf = $fig->{_dbf};  
20    
21    
22    if ($mode eq 'all') {
23    
24  my @genomes;      # Here we extract external aliases from the peg.synonyms table, when they can be inferred
25  if (@ARGV == 0)      # accurately.
26            Trace("Extracting external aliases from the peg.synonyms table.") if T(2);
27        open(\*SYN, "<$FIG_Config::global/peg.synonyms");
28        while (defined($_ = <SYN>))
29  {  {
30      $dbf->drop_table( tbl => "features" );                  chop;
31      if ($FIG_Config::dbms eq "Pg")                  my($x,$y) = split(/\t/,$_);
32                    my @ids = map { $_ =~ /^([^,]+),(\d+)/; [$1,$2] } ($x,split(/;/,$y));
33                    my @fig = ();
34                    my(@nonfig) = ();
35                    foreach $_ (@ids)
36      {      {
37          $dbf->create_table( tbl  => "features",                          if ($_->[0] =~ /^fig\|/)
38                              flds => "id varchar(32) UNIQUE NOT NULL, type varchar(16),genome varchar(16),"  .                          {
39                                      "location varchar(5000),"  .                                  push(@fig,$_);
                                     "contig varchar(96), minloc INTEGER, maxloc INTEGER,"  .  
                                     "aliases TEXT, PRIMARY KEY ( id )"  
                             );  
40      }      }
41      elsif ($FIG_Config::dbms eq "mysql")                          else
42      {      {
43          $dbf->create_table( tbl  => "features",                                  push(@nonfig,$_);
                             flds => "id varchar(32) UNIQUE NOT NULL, type varchar(16),genome varchar(16),"  .  
                                     "location TEXT,"  .  
                                     "contig varchar(96), minloc INTEGER, maxloc INTEGER,"  .  
                                     "aliases TEXT, PRIMARY KEY ( id )"  
                             );  
44      }      }
   
     $dbf->create_index( idx  => "features_org_ix",  
                         tbl  => "features",  
                         type => "btree",  
                         flds => "genome" );  
     $dbf->create_index( idx  => "features_type_ix",  
                         type => "btree",  
                         tbl  => "features",  
                         flds => "type" );  
     $dbf->create_index( idx  => "features_beg_ix",  
                         type => "btree",  
                         tbl  => "features",  
                         flds => "genome,contig,minloc" );  
   
     @genomes = $fig->genomes;  
45  }  }
46  else  
47                    foreach $x (@fig)
48  {  {
49      @genomes = @ARGV;                          my($peg,$peg_ln) = @$x;
50                            my $genome = &FIG::genome_of($peg);
51                            foreach $_ (@nonfig)
52                            {
53                                    if ((@fig == 1) || ($peg_ln == $_->[1]))
54                                    {
55                                            print ALIAS "$peg\t$_->[0]\t$genome\n";
56                                            Trace("Alias record $peg, $_->[0] for $genome.") if T(4);
57                                    }
58                            }
59                    }
60        }
61        close(SYN);
62  }  }
63    
64  foreach $genome (@genomes)  foreach $genome (@genomes)
65  {  {
66            Trace("Processing $genome.") if T(3);
67      opendir(FEAT,"$organisms_dir/$genome/Features")      opendir(FEAT,"$organisms_dir/$genome/Features")
68          || die "could not open $genome/Features";          || die "could not open $genome/Features";
69      @types = grep { $_ =~ /^[a-zA-Z]+$/ } readdir(FEAT);      @types = grep { $_ =~ /^[a-zA-Z]+$/ } readdir(FEAT);
# Line 75  Line 74 
74          if ((-s "$organisms_dir/$genome/Features/$type/tbl") &&          if ((-s "$organisms_dir/$genome/Features/$type/tbl") &&
75              open(TBL,"<$organisms_dir/$genome/Features/$type/tbl"))              open(TBL,"<$organisms_dir/$genome/Features/$type/tbl"))
76          {          {
77  #           print STDERR "loading $genome/Features/$type/tbl\n";              Trace("Loading $genome/Features/$type/tbl") if T(4);
78              while (defined($_ = <TBL>))              my @tbl = <TBL>;
79                close(TBL);
80                my %seen;
81    
82                while ($_ = pop @tbl)
83              {              {
84                  chop;                  chop;
85                  ($id,$loc,@aliases) = split(/\t/,$_);                  ($id,$loc,@aliases) = split(/\t/,$_);
86                  if ($id)  
87                    if ($id && (! $seen{$id}))
88                  {                  {
89                        $seen{$id} = 1;
90                      my($minloc,$maxloc);                      my($minloc,$maxloc);
91                      if ($loc)                      if ($loc)
92                      {                      {
# Line 101  Line 106 
106                      if (@aliases > 0)                      if (@aliases > 0)
107                      {                      {
108                          $aliases = join(",",grep(/\S/,@aliases));                          $aliases = join(",",grep(/\S/,@aliases));
109                            my $alias;
110                            foreach $alias (@aliases)
111                            {
112                                if ($alias =~ /^([NXYZA]P_|gi\||sp\|\tr\||kegg\||uni\|)/)
113                                {
114    
115                                    print ALIAS "$id\t$alias\t$genome\tOVERRIDE\n";
116                                    Trace("$id override alias $alias for $genome") if T(4);
117                                }
118                            }
119                      }                      }
120                      else                      else
121                      {                      {
# Line 108  Line 123 
123                      }                      }
124                      $minloc = (! $minloc) ? 0 : $minloc;                      $minloc = (! $minloc) ? 0 : $minloc;
125                      $maxloc = (! $maxloc) ? 0 : $maxloc;                      $maxloc = (! $maxloc) ? 0 : $maxloc;
126                      print REL "$id\t$type\t$genome\t$loc\t$contig\t$minloc\t$maxloc\t$aliases\n";                      if ((length($loc) < 5000) && (length($contig) < 96) && (length($id) < 32) && ($id =~ /(\d+)$/))
127                        {
128                            print REL "$id\t$1\t$type\t$genome\t$loc\t$contig\t$minloc\t$maxloc\t$aliases\n";
129                        }
130                  }                  }
131              }              }
             close(TBL);  
132          }          }
133      }      }
134  }  }
135  close(REL);  close(REL);
136    close(ALIAS);
137    Open(\*ALIASIN, "<$temp_dir/tmpalias$$");
138    Open(\*ALIASOUT, ">$temp_dir/tmpalias$$.1");
139    Trace("Parsing alias file.") if T(2);
140    $_ = <ALIASIN>;
141    while ($_ && ($_ =~ /^(\S+)/))
142    {
143        my @aliases = ();
144        my $curr = $1;
145        while ($_ && ($_ =~ /^(\S+)\t(\S+)(\t(\S+))?/) && ($1 eq $curr))
146        {
147                    push(@aliases,[$2,$3 ? 1 : 0]);
148                    $_ = <ALIASIN>;
149        }
150        my $x;
151        my $genome = &FIG::genome_of($curr);
152        foreach $x (@aliases)
153        {
154            if ($x->[1])
155            {
156                print ALIASOUT "$curr\t$x->[0]\t$genome\n";
157            }
158            else
159            {
160                my $i;
161                for ($i=0; ($i < @aliases) && ((! $aliases[$i]->[1]) || (! &same_class($x->[0],$aliases[$i]->[0]))); $i++) {}
162                if ($i == @aliases)
163                {
164                    print ALIASOUT "$curr\t$x->[0]\t$genome\n";
165                }
166            }
167        }
168    }
169    close(ALIASIN);
170    close(ALIASOUT);
171    unlink("$temp_dir/tmpalias$$");
172    
173  $dbf->load_table( tbl => "features",  $fig->reload_table($mode, 'features',
174                    file => "$temp_dir/tmpfeat$$" );                                     "id varchar(32), idN INTEGER, type varchar(16),genome varchar(16),"  .
175                                        "location TEXT,"  .
176  if (@ARGV == 0) {  $dbf->vacuum_it("features") }                                      "contig varchar(96), minloc INTEGER, maxloc INTEGER,"  .
177                                        "aliases TEXT",
178                                            { features_id_ix => "id", features_org_ix => "genome",
179                                              features_type_ix => "type", features_beg_ix => "genome, contig, minloc" },
180                                            "$temp_dir/tmpfeat$$", \@genomes);
181  unlink("$temp_dir/tmpfeat$$");  unlink("$temp_dir/tmpfeat$$");
182    
183    $fig->reload_table($mode, 'ext_alias',
184                                            "id varchar(32), alias varchar(32), genome varchar(16)",
185                                            { ext_alias_alias_ix => "alias", ext_alias_genome_ix => "genome",
186                                              ext_alias_ix_id => "id" },
187                                            "$temp_dir/tmpalias$$.1", \@genomes );
188    
189    unlink("$temp_dir/tmpalias$$.1");
190    Trace("Features loaded.") if T(2);
191    
192    sub same_class {
193        my($x,$y) = @_;
194    
195        my $class1 = &classA($x);
196        my $class2 = &classA($y);
197        return ($class1 && ($class1 eq $class2));
198    }
199    
200    sub classA {
201        my($alias) = @_;
202    
203        if ($alias =~ /^([^\|]+)\|/)
204        {
205                    return $1;
206        }
207        elsif ($alias =~ /^[NXYZA]P_[0-9\.]+$/)
208        {
209                    return "refseq";
210        }
211        else
212        {
213                    return "";
214        }
215    }

Legend:
Removed from v.1.1  
changed lines
  Added in v.1.19

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3