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

Annotation of /FigKernelScripts/load_features.pl

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.18 - (view) (download) (as text)

1 : efrank 1.1 # -*- perl -*-
2 :    
3 :     use strict;
4 :     use FIG;
5 :     my $fig = new FIG;
6 :    
7 : parrello 1.18 use Tracer;
8 : efrank 1.1
9 : parrello 1.18 Trace("Preparing to load features.") if T(2);
10 :     my ($mode, @genomes) = FIG::parse_genome_args(@ARGV);
11 : efrank 1.1 my $temp_dir = "$FIG_Config::temp";
12 : parrello 1.18 my $organisms_dir = "$FIG_Config::organisms";
13 : efrank 1.1
14 :     my($genome,@types,$type,$id,$loc,@aliases,$aliases,$contig);
15 :    
16 :     # usage: load_features [G1 G2 G3 ... ]
17 :    
18 : parrello 1.18 Open(\*REL, ">$temp_dir/tmpfeat$$");
19 :     Open(\*ALIAS, "| sort -T $temp_dir -u > $temp_dir/tmpalias$$");
20 : efrank 1.1
21 :    
22 : parrello 1.18 if ($mode eq 'all') {
23 : overbeek 1.10
24 :     # Here we extract external aliases from the peg.synonyms table, when they can be inferred
25 :     # accurately.
26 : parrello 1.18 Trace("Extracting external aliases from the peg.synonyms table.") if T(2);
27 :     open(\*SYN, "<$FIG_Config::global/peg.synonyms");
28 : overbeek 1.10 while (defined($_ = <SYN>))
29 :     {
30 : parrello 1.18 chop;
31 :     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 :     if ($_->[0] =~ /^fig\|/)
38 :     {
39 :     push(@fig,$_);
40 :     }
41 :     else
42 :     {
43 :     push(@nonfig,$_);
44 :     }
45 :     }
46 :    
47 :     foreach $x (@fig)
48 : overbeek 1.13 {
49 : parrello 1.18 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 : overbeek 1.13 }
60 : overbeek 1.10 }
61 :     close(SYN);
62 : efrank 1.1 }
63 :    
64 : overbeek 1.12 my $changes = {};
65 : parrello 1.18 if (open(TMP, "<$FIG_Config::global/changed.location.features"))
66 : overbeek 1.12 {
67 :     while ($_ = <TMP>)
68 :     {
69 : parrello 1.18 if ($_ =~ /^(fig\|\d+\.\d+\.[a-zA-Z]+\.\d+)/)
70 :     {
71 :     $changes->{$1}++;
72 :     }
73 : overbeek 1.12 }
74 :     close(TMP);
75 :     }
76 :    
77 : efrank 1.1 foreach $genome (@genomes)
78 :     {
79 : parrello 1.18 Trace("Processing $genome.") if T(3);
80 : efrank 1.1 opendir(FEAT,"$organisms_dir/$genome/Features")
81 : parrello 1.18 || die "could not open $genome/Features";
82 : efrank 1.1 @types = grep { $_ =~ /^[a-zA-Z]+$/ } readdir(FEAT);
83 :     closedir(FEAT);
84 :    
85 :     foreach $type (@types)
86 :     {
87 : parrello 1.18 if ((-s "$organisms_dir/$genome/Features/$type/tbl") &&
88 :     open(TBL,"<$organisms_dir/$genome/Features/$type/tbl"))
89 : overbeek 1.12 {
90 : parrello 1.18 Trace("Loading $genome/Features/$type/tbl") if T(4);
91 :     while (defined($_ = <TBL>))
92 : overbeek 1.8 {
93 : parrello 1.18 chop;
94 :     ($id,$loc,@aliases) = split(/\t/,$_);
95 :    
96 :     if ($id && ($_ = $changes->{$id})) # check for obsolete entries due to location changes
97 :     {
98 :     $changes->{$id}--;
99 :     next;
100 :     }
101 :    
102 :     if ($id)
103 :     {
104 :     my($minloc,$maxloc);
105 :     if ($loc)
106 :     {
107 :     $loc =~ s/\s+$//;
108 :     ($contig,$minloc,$maxloc) = &FIG::boundaries_of($loc);
109 :     if ($minloc && $maxloc)
110 :     {
111 :     ($minloc < $maxloc) || (($minloc,$maxloc) = ($maxloc,$minloc));
112 :     }
113 :     }
114 :    
115 :     if (! $contig)
116 :     {
117 :     $loc = $contig = $minloc = $maxloc = "";
118 :     }
119 :    
120 :     if (@aliases > 0)
121 :     {
122 :     $aliases = join(",",grep(/\S/,@aliases));
123 :     my $alias;
124 :     foreach $alias (@aliases)
125 :     {
126 :     if ($alias =~ /^([NXYZA]P_|gi\||sp\|\tr\||kegg\||uni\|)/)
127 :     {
128 :    
129 :     print ALIAS "$id\t$alias\t$genome\tOVERRIDE\n";
130 :     Trace("$id override alias $alias for $genome") if T(4);
131 :     }
132 :     }
133 :     }
134 :     else
135 :     {
136 :     $aliases = "";
137 :     }
138 :     $minloc = (! $minloc) ? 0 : $minloc;
139 :     $maxloc = (! $maxloc) ? 0 : $maxloc;
140 :     if ((length($loc) < 5000) && (length($contig) < 96) && (length($id) < 32) && ($id =~ /(\d+)$/))
141 :     {
142 :     print REL "$id\t$1\t$type\t$genome\t$loc\t$contig\t$minloc\t$maxloc\t$aliases\n";
143 :     }
144 :     }
145 : overbeek 1.8 }
146 : parrello 1.18 close(TBL);
147 : efrank 1.1 }
148 :     }
149 :     }
150 :     close(REL);
151 : overbeek 1.8 close(ALIAS);
152 : parrello 1.18 Open(\*ALIASIN, "<$temp_dir/tmpalias$$");
153 :     Open(\*ALIASOUT, ">$temp_dir/tmpalias$$.1");
154 :     Trace("Parsing alias file.") if T(2);
155 : overbeek 1.14 $_ = <ALIASIN>;
156 :     while ($_ && ($_ =~ /^(\S+)/))
157 :     {
158 :     my @aliases = ();
159 :     my $curr = $1;
160 :     while ($_ && ($_ =~ /^(\S+)\t(\S+)(\t(\S+))?/) && ($1 eq $curr))
161 :     {
162 : parrello 1.18 push(@aliases,[$2,$3 ? 1 : 0]);
163 :     $_ = <ALIASIN>;
164 : overbeek 1.14 }
165 :     my $x;
166 :     my $genome = &FIG::genome_of($curr);
167 :     foreach $x (@aliases)
168 :     {
169 : parrello 1.18 if ($x->[1])
170 :     {
171 :     print ALIASOUT "$curr\t$x->[0]\t$genome\n";
172 :     }
173 :     else
174 :     {
175 :     my $i;
176 :     for ($i=0; ($i < @aliases) && ((! $aliases[$i]->[1]) || (! &same_class($x->[0],$aliases[$i]->[0]))); $i++) {}
177 :     if ($i == @aliases)
178 :     {
179 :     print ALIASOUT "$curr\t$x->[0]\t$genome\n";
180 :     }
181 :     }
182 : overbeek 1.14 }
183 :     }
184 :     close(ALIASIN);
185 :     close(ALIASOUT);
186 :     unlink("$temp_dir/tmpalias$$");
187 : efrank 1.1
188 : parrello 1.18 $fig->reload_table($mode, 'features',
189 :     "id varchar(32), idN INTEGER, type varchar(16),genome varchar(16)," .
190 :     "location TEXT," .
191 :     "contig varchar(96), minloc INTEGER, maxloc INTEGER," .
192 :     "aliases TEXT",
193 :     { features_id_ix => "id", features_org_ix => "genome",
194 :     features_type_ix => "type", features_beg_ix => "genome, contig, minloc" },
195 :     "$temp_dir/tmpfeat$$", \@genomes);
196 :     unlink("$temp_dir/tmpfeat$$");
197 : efrank 1.1
198 : parrello 1.18 $fig->reload_table($mode, 'ext_alias',
199 :     "id varchar(32), alias varchar(32), genome varchar(16)",
200 :     { ext_alias_alias_ix => "alias", ext_alias_genome_ix => "genome",
201 :     ext_alias_ix_id => "id" },
202 :     "$temp_dir/tmpalias$$.1", \@genomes );
203 : overbeek 1.8
204 : overbeek 1.14 unlink("$temp_dir/tmpalias$$.1");
205 : parrello 1.18 Trace("Features loaded.") if T(2);
206 : overbeek 1.14
207 :     sub same_class {
208 :     my($x,$y) = @_;
209 :    
210 :     my $class1 = &classA($x);
211 :     my $class2 = &classA($y);
212 :     return ($class1 && ($class1 eq $class2));
213 :     }
214 :    
215 :     sub classA {
216 :     my($alias) = @_;
217 :    
218 :     if ($alias =~ /^([^\|]+)\|/)
219 :     {
220 : parrello 1.18 return $1;
221 : overbeek 1.14 }
222 :     elsif ($alias =~ /^[NXYZA]P_[0-9\.]+$/)
223 :     {
224 : parrello 1.18 return "refseq";
225 : overbeek 1.14 }
226 :     else
227 :     {
228 : parrello 1.18 return "";
229 : overbeek 1.14 }
230 :     }

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3