[Bio] / FigKernelPackages / model.pm Repository:
ViewVC logotype

Annotation of /FigKernelPackages/model.pm

Parent Directory Parent Directory | Revision Log Revision Log


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

1 : olson 1.1
2 :     ## _*_ Perl _*_ ##
3 :     #
4 :     # model.pm
5 :     #
6 :     # Kevin Formsma
7 :     # Hope College
8 :     # Created: 6/1/2006
9 :     #
10 :     ##################
11 :    
12 :     ### BrainStorm/Planning ###
13 :     #
14 :     # Functions
15 :     # -Scenario access methods for subsystem
16 :     # -Access current genome models
17 :     # -Interface to/Integration of find_reaction_paths.cgi
18 :     # -Scenario Relationships:
19 :     # 1. Automanaged - map togather based on starting/ending compounds
20 :     # 2. User Defined - through some type of map or XML
21 :     # -Generate files for Scenarios on a organism
22 :     # 1. Path Picking - All or Select?
23 :     # -Filter written Scenarios based on Relationships
24 :     # -Report valid Scenarios for a genome
25 :     #
26 :     #
27 :     # Basic Process for a Control Script
28 :     # 1. Input Genomes for model creation and Subsystems
29 :     # 2. Define selected Scenario relationships
30 :     # 3. Create Scenario Paths for each genome and subsystem
31 :     # 4. Filter invalid Scenarios based on relationships
32 :     # 5. Combine Subsystem models for each genome into genome models
33 :     # 6. Report on the model results for each genome
34 :     #
35 :    
36 :     package model;
37 :    
38 :     use strict;
39 :     use FIG;
40 :     use Subsystem;
41 :     use File::Path;
42 :    
43 :     our $fig= new FIG;
44 :    
45 :    
46 :     #global variables that make process_paths work
47 :     #These need to be cleared and reloaded frequently for process path and flux writing
48 :     my (%reactions_to_substrate_arrays, %reactions_to_product_arrays, %all_compounds_to_main);
49 :     my %all_reactions;
50 :     my %scenario_cycles;
51 :     my @all_outputs_lists;
52 :     my %all_inputs;
53 :     my %all_outputs;
54 :    
55 :     #this variable is used to set the loop count for single scenario and then assembly runs, 100/100 by default.
56 :     my $loop_max = 100;
57 :     my $loop_max_assembly = 25;
58 :    
59 :     #These store what supersets/Subsystems we are using
60 :     my %superset_to_ss;
61 :     my %ss_to_superset;
62 :    
63 :     #Flip this bit to enable debugging
64 :     my $debug = 0;
65 :    
66 :     sub new {
67 :     my $type = shift;
68 :     my $self = {};
69 :     return bless $self, $type;
70 :     }
71 :    
72 :    
73 :    
74 :     sub get_ss_scenarios
75 :     {
76 :     my ($ss_name) = @_;
77 :     # This is the scenario data structure storage
78 :     # %scenario_data{scenario_name}
79 :     #
80 :     my %scenario_data;
81 :    
82 :     my $subsystem = $fig->get_subsystem($ss_name);
83 :    
84 :     my @scenario_names = $subsystem->get_hope_scenario_names;
85 :     foreach my $name (@scenario_names)
86 :     {
87 :     $scenario_data{$name} = &get_scenario($subsystem,$name);
88 :     }
89 :    
90 :     return \%scenario_data;
91 :     }
92 :    
93 :    
94 :     sub get_scenario
95 :     {
96 :     my($subsystem,$name) = @_;
97 :    
98 :     my %data;
99 :    
100 :     my @inputs = @{$subsystem->get_hope_input_compounds($name)};
101 :     my @outputs = @{$subsystem->get_hope_output_compounds($name)};
102 :     my @map_ids = @{$subsystem->get_hope_map_ids($name)};
103 :     my @additional_reactions = @{$subsystem->get_hope_additional_reactions($name)};
104 :     my @ignore_reactions = @{$subsystem->get_hope_ignore_reactions($name)};
105 :    
106 :     $data{inputs} = \@inputs;
107 :     $data{outputs} = \@outputs;
108 :     $data{map_ids} = \@map_ids;
109 :     $data{additional_reactions} = \@additional_reactions;
110 :     $data{ignore_reactions} = \@ignore_reactions;
111 :    
112 :     return \%data;
113 :    
114 :     }
115 :    
116 :    
117 :    
118 :    
119 :     sub process_init
120 :     {
121 :     my ($ss_name,$scenario_name,$genome,$assembly) = @_;
122 :     my (%sc_inputs, %sc_outputs);
123 :    
124 :     if ($genome eq "")
125 :     {
126 :     $genome = "All";
127 :     }
128 :     print STDERR "\nSubsystem : ".$ss_name." Scenario: $scenario_name \n" if $debug;
129 :     my $subsystem = $fig->get_subsystem($ss_name);
130 :     my @roles = $subsystem->get_roles;
131 :     my %hope_reactions = %{$subsystem->get_hope_reactions};
132 :     my $scenario_data = &get_scenario($subsystem,$scenario_name);
133 :     #load the other arrays
134 :     my %ss_reactions;
135 :    
136 :     my $valid_genome = 0;
137 :     #determine if this genome has this subsystem
138 :     my @genome_array = $subsystem->get_genomes();
139 :     foreach my $organ (@genome_array)
140 :     {
141 :     if($organ eq $genome || $genome eq "All" )
142 :     {
143 :     $valid_genome = 1;
144 :     last;
145 :     }
146 :     }
147 :    
148 :     foreach my $role (@roles)
149 :     {
150 :     if ($valid_genome && defined $hope_reactions{$role} &&
151 :     ($genome eq "All" || $subsystem->get_pegs_from_cell($genome,$role) > 0))
152 :     {
153 :     map { $ss_reactions{$_} = 1 } @{$hope_reactions{$role}};
154 :     }
155 :     }
156 :     map { $sc_inputs{$_} = 1 } @{$scenario_data->{inputs}};
157 :    
158 :     foreach my $list (@{$scenario_data->{outputs}})
159 :     {
160 :     map { $sc_outputs{$_} = 1 } @$list;
161 :     push @all_outputs_lists, $list;
162 :     }
163 :    
164 :     map { $scenario_cycles{$_} = 1 if defined $sc_outputs{$_} } keys %sc_inputs;
165 :     map { $all_inputs{$_} = 1 } keys %sc_inputs;
166 :     map { $all_outputs{$_} = 1 } keys %sc_outputs;
167 :    
168 :     my @hope_additional_reactions = @{$scenario_data->{additional_reactions}};
169 :     my @hope_ignore_reactions = @{$scenario_data->{ignore_reactions}};
170 :     my %sc_reactions;
171 :     map { $sc_reactions{$_} = 1 } keys %ss_reactions;
172 :    
173 :     # flag additional reactions so we won't check if they are in a map
174 :     foreach my $rid (@hope_additional_reactions)
175 :     {
176 :     $sc_reactions{$rid} = 2;
177 :     }
178 :    
179 :     foreach my $rid (@hope_ignore_reactions)
180 :     {
181 :     delete $sc_reactions{$rid};
182 :     }
183 :    
184 :     #for now we do this outside of the if statment, but that might need to change
185 :     &load_substrate_and_product_arrays(\%sc_reactions,$scenario_data->{map_ids});
186 :    
187 :    
188 :    
189 :    
190 :     }
191 :    
192 :     sub execute_paths
193 :     {
194 :     my ($assembly_paths,$find_first,$input_path, $output_path) = @_;
195 :    
196 :     my $num_paths = scalar @{$assembly_paths};
197 :    
198 :     my (%substrates_to_reactions, %products_to_reactions,
199 :     %reactions_to_substrates, %reactions_to_products);
200 :    
201 :     if($num_paths ==0 ) #Load the reactions normally if we are creating scenarios
202 :     {
203 :     create_reactions(\%substrates_to_reactions,\%products_to_reactions,
204 :     \%reactions_to_substrates,
205 :     \%reactions_to_products);
206 :     }
207 :     #create assembly reactions that are closed 'paths' from scenarios or assemblies
208 :     if($num_paths > 0)
209 :     {
210 :     create_assembly_reactions(\%substrates_to_reactions,\%products_to_reactions,
211 :     \%reactions_to_substrates,
212 :     \%reactions_to_products,$assembly_paths);
213 :     }
214 :    
215 :     #This deals with user specifed input/output paths, and uses these to generate
216 :     #what our input and output compounds should be for an assembly.
217 :     #This is only used for creating assemblies, has no effect on scenario creation
218 :    
219 :     if(scalar @{$input_path} && scalar @{$output_path})
220 :     {
221 :    
222 :     %all_inputs = ();
223 :     %all_outputs = ();
224 :     #we should have loaded these paths above, so lets just pull out the information we need
225 :     foreach my $path (@$input_path)
226 :     {
227 :     my $input_rxn = "$path->[-3]/$path->[-2]/$path->[-1]_R";
228 :     my @user_in = @{$reactions_to_substrates{$input_rxn}};
229 :     map{ $all_inputs{$_} = 1 } @user_in;
230 :     }
231 :     foreach my $path (@$output_path)
232 :     {
233 :     my $output_rxn = "$path->[-3]/$path->[-2]/$path->[-1]_R";
234 :     my @user_out = @{$reactions_to_products{$output_rxn}};
235 :     map{ $all_outputs{$_} = 1 } @user_out;
236 :     }
237 :     #map { $scenario_cycles{$_} = 1 if defined $all_outputs{$_} } keys %all_inputs;
238 :     }
239 :    
240 :     print STDERR "Inputs :\n" if $debug;
241 : olson 1.2 print STDERR map { $_."\n" } keys %all_inputs if $debug;
242 : olson 1.1 print STDERR "Outputs:\n" if $debug;
243 : olson 1.2 print STDERR map { $_."\n" } keys %all_outputs if $debug;
244 : olson 1.1
245 :     #filter the input/outputs lists, removing the intersection unless something
246 :     # is a known cycle
247 :     foreach my $input (keys %all_inputs)
248 :     {
249 :     if(defined $all_outputs{$input} && ! defined $scenario_cycles{$input})
250 :     {
251 :     print STDERR "Deleting $input from input and output lists\n" if $debug;
252 :     delete $all_inputs{$input};
253 :     delete $all_outputs{$input};
254 :     }
255 :     }
256 :    
257 :     my $create_assembly = 0;
258 :     $create_assembly = 1 if(scalar @{$assembly_paths} !=0);
259 :     return process_paths(\%all_inputs, \%all_outputs, \@all_outputs_lists,
260 :     \%reactions_to_substrates, \%reactions_to_products,
261 :     \%substrates_to_reactions,\%products_to_reactions,$create_assembly,$find_first);
262 :     }
263 :    
264 :    
265 :     sub create_reactions
266 :     {
267 :     my ($substrates_to_reactions, $products_to_reactions,
268 :     $reactions_to_substrates, $reactions_to_products) = @_;
269 : olson 1.2 print STDERR "building SS reactions\n" if $debug;
270 : olson 1.1 # use subsystem reactions
271 :     foreach my $drxn (map { ($_."_L", $_."_R") } keys %all_reactions)
272 :     {
273 :     foreach my $substrArr ($reactions_to_substrate_arrays{$drxn})
274 :     {
275 :     foreach my $cinfo (@$substrArr)
276 :     {
277 :     my $cpd = $cinfo->[0];
278 :     my $main = $cinfo->[2] || defined $all_inputs{$cpd}; # main in this reaction
279 :    
280 :     if ($main)
281 :     {
282 :     push(@{$reactions_to_substrates->{$drxn}}, $cpd);
283 :     push(@{$substrates_to_reactions->{$cpd}}, $drxn);
284 :     }
285 :     }
286 :     }
287 :    
288 :     foreach my $prodArr ($reactions_to_product_arrays{$drxn})
289 :     {
290 :     foreach my $cinfo (@$prodArr)
291 :     {
292 :     my $cpd = $cinfo->[0];
293 :     my $main = $cinfo->[2] || defined $all_outputs{$cpd}; # main in this reaction
294 :    
295 :     if ($main)
296 :     {
297 :     push(@{$reactions_to_products->{$drxn}}, $cpd);
298 :     push(@{$products_to_reactions->{$cpd}}, $drxn);
299 :     }
300 :     }
301 :     }
302 :     }
303 :     }
304 :    
305 :     sub create_assembly_reactions
306 :     {
307 :     my ($substrates_to_reactions, $products_to_reactions,
308 :     $reactions_to_substrates, $reactions_to_products,$assembly_paths) = @_;
309 :     my %intersection;
310 :    
311 :     foreach my $path (@$assembly_paths)
312 :     {
313 :     my $paths_dir = "$FIG_Config::var/Models/" . join "/" , @$path;
314 :    
315 :     my $drxn = "$path->[-3]/$path->[-2]/$path->[-1]_R";
316 :    
317 :     $all_reactions{"$path->[-3]/$path->[-2]/$path->[-1]"} = "R";
318 :    
319 :     print STDERR "Making reaction: $drxn\n" if $debug;
320 :    
321 :     open (M_INPUTS , $paths_dir."/inputs") or die ("Failed to open $paths_dir"."/inputs");
322 :     my @substrArr;
323 :    
324 :     print STDERR "Gathering Inputs:\n" if $debug;
325 :    
326 :     while (<M_INPUTS>)
327 :     {
328 :     my ($cpd, $stoich) = split "\t" , $_;
329 :    
330 :     #We are going to assume everything is a main...
331 :     if(!defined $all_compounds_to_main{$cpd})
332 :     {
333 :     $all_compounds_to_main{$cpd} = 1;
334 :     }
335 :    
336 :     if ($all_compounds_to_main{$cpd})
337 :     {
338 :     push(@{$reactions_to_substrates->{$drxn}}, $cpd);
339 :     push(@{$substrates_to_reactions->{$cpd}}, $drxn);
340 :     }
341 :    
342 :     if($all_compounds_to_main{$cpd} !=0 || !defined $all_compounds_to_main{$cpd})
343 :     {
344 :     $all_inputs{$cpd} = 1;
345 :     }
346 :     push @substrArr, [$cpd, $stoich, $all_compounds_to_main{$cpd}];
347 :    
348 :     my @names = $fig->names_of_compound($cpd);
349 :     print STDERR "\t$stoich\t$cpd\t$names[0]\t$all_compounds_to_main{$cpd}\n" if $debug;
350 :     }
351 :    
352 :     $reactions_to_substrate_arrays{$drxn} = \@substrArr;
353 :    
354 :     close M_INPUTS;
355 :    
356 :     open (M_OUTPUTS, $paths_dir."/outputs") or die("Failed to open $paths_dir"."/outputs");
357 :     my @prodArr;
358 :    
359 :     print STDERR "Gathering outputs:\n" if $debug;
360 :    
361 :     while (<M_OUTPUTS>)
362 :     {
363 :     my ($cpd, $stoich) = split "\t", $_;
364 :     print STDERR "Found $stoich $cpd\n" if $debug;
365 :    
366 :     #We are going to assume everything is a main...
367 :     if(!defined $all_compounds_to_main{$cpd})
368 :     {
369 :     $all_compounds_to_main{$cpd} = 1;
370 :     }
371 :    
372 :     if ($all_compounds_to_main{$cpd})
373 :     {
374 :     push(@{$reactions_to_products->{$drxn}}, $cpd);
375 :     push(@{$products_to_reactions->{$cpd}}, $drxn);
376 :     }
377 :    
378 :    
379 :     #This adds cycles from 'assemblys' because they weren't added earlier
380 :     foreach my $ele (@{$reactions_to_substrates->{$drxn}})
381 :     {
382 :     if($ele eq $cpd)
383 :     {
384 :     $scenario_cycles{$cpd} = 1;
385 :     }
386 :     }
387 :    
388 :     if($all_compounds_to_main{$cpd} !=0 || !defined $all_compounds_to_main{$cpd})
389 :     {
390 :     $all_outputs{$cpd} = 1;
391 :     }
392 :     push @prodArr, [$cpd,$stoich,$all_compounds_to_main{$cpd}];
393 :     my @names = $fig->names_of_compound($cpd);
394 :     print STDERR "\t$cpd\t$names[0]\t$all_compounds_to_main{$cpd}\n" if $debug;
395 :     }
396 :    
397 :     $reactions_to_product_arrays{$drxn} = \@prodArr;
398 :    
399 :     close M_OUTPUTS;
400 :    
401 :     }
402 :     }
403 :    
404 :    
405 :     sub load_substrate_and_product_arrays
406 :     {
407 :     my ($reactions, $map_ids) = @_;
408 :    
409 :     # determine whether the reaction is in one of the maps, and get directionality accordingly
410 :     my (%reactions_in_maps, %reactions_not_in_maps);
411 :    
412 :     foreach my $rxn (keys %$reactions)
413 :     {
414 :     my $direction;
415 :    
416 :     if($fig->valid_reaction_id($rxn))
417 :     {
418 :     # get an array of triplets. The triplets are [reaction id][map id]
419 :     # [left to right - R, right to left - L, or both - B]
420 :     my @triplets = $fig->reaction_direction($rxn);
421 :    
422 :     foreach my $trip (@triplets)
423 :     {
424 :     foreach my $map_id (@$map_ids)
425 :     {
426 :     if (@{$trip}[1] eq $map_id)
427 :     {
428 :     my $this_direction = @{$trip}[2];
429 :    
430 :     # bidirectional in one map overrules unidirectional in another.
431 :     # opposite directions in two maps becomes bidirectional
432 :     if (! defined $direction)
433 :     {
434 :     $direction = $this_direction;
435 :     }
436 :     elsif ($direction ne "B" && ($this_direction eq "B" ||
437 :     $this_direction ne $direction))
438 :     {
439 :     $direction = "B";
440 :     }
441 :    
442 :     $reactions_in_maps{$rxn} = 1;
443 :     }
444 :     }
445 :     }
446 :    
447 :    
448 :     if(! $reactions_in_maps{$rxn})
449 :     {
450 :     my $found_in_other_map = 0;
451 :    
452 :     #reaction not in scenario map ids, try to get directionality from other maps
453 :     foreach my $trip (@triplets)
454 :     {
455 :     my $this_direction = @{$trip}[2];
456 :    
457 :     #bidreactional in one map overrules unidirectional in another
458 :     #opposite directions in two maps becomes bidirectional
459 :     if(! defined $direction)
460 :     {
461 :     $direction = $this_direction;
462 :     }
463 :     elsif($direction ne "B" && ($this_direction eq "B" || $this_direction ne $direction))
464 :     {
465 :     $direction = "B";
466 :     }
467 :    
468 :     $found_in_other_map = 1;
469 :     }
470 :     if (!$found_in_other_map)
471 :     {
472 :     # reaction not in map, get directionality without reference to map
473 :     if($fig->reversible($rxn) eq "1")
474 :     {
475 :     $direction = "B";
476 :     }
477 :     else
478 :     {
479 :     $direction = "R";
480 :     }
481 :     }
482 :     }
483 :    
484 :     if (! defined $all_reactions{$rxn} || $direction eq "B")
485 :     {
486 :     $all_reactions{$rxn} = $direction;
487 :     }
488 :     elsif ($all_reactions{$rxn} ne $direction)
489 :     {
490 :     $all_reactions{$rxn} = "B";
491 :     }
492 :    
493 :     my (@substrArr, @prodArr);
494 :    
495 :     if ($direction eq "L")
496 :     {
497 :     @substrArr = $fig->reaction2comp($rxn, 1, $map_ids);
498 :     @prodArr = $fig->reaction2comp($rxn, 0, $map_ids);
499 :     $reactions_to_substrate_arrays{$rxn . "_L"} = \@substrArr;
500 :     $reactions_to_product_arrays{$rxn . "_L"} = \@prodArr;
501 :    
502 :     }
503 :     else
504 :     {
505 :     if ($reactions_in_maps{$rxn})
506 :     {
507 :     @substrArr = $fig->reaction2comp($rxn, 0, $map_ids);
508 :     @prodArr = $fig->reaction2comp($rxn, 1, $map_ids);
509 :     }
510 :     else
511 :     {
512 :     @substrArr = $fig->reaction2comp($rxn, 0);
513 :     @prodArr = $fig->reaction2comp($rxn, 1);
514 :     }
515 :    
516 :     $reactions_to_substrate_arrays{$rxn."_R"} = \@substrArr;
517 :     $reactions_to_product_arrays{$rxn."_R"} = \@prodArr;
518 :    
519 :     if ($direction eq "B")
520 :     {
521 :     $reactions_to_substrate_arrays{$rxn."_L"} = \@prodArr;
522 :     $reactions_to_product_arrays{$rxn."_L"} = \@substrArr;
523 :     }
524 :     }
525 :    
526 :     print STDERR "\nFor $rxn, found substrates:\n" unless !$debug;
527 :     map { print STDERR "\t$_->[0]\t$_->[1]\t$_->[2]\n" unless !$debug } @substrArr;
528 :     print STDERR "For $rxn, found products:\n" unless !$debug;
529 :     map { print STDERR "\t$_->[0]\t$_->[1]\t$_->[2]\n" unless !$debug } @prodArr;
530 :    
531 :     # load "main" designation based on reactions that are in the maps
532 :     if ($reactions_in_maps{$rxn})
533 :     {
534 :     foreach my $cinfo ((@substrArr, @prodArr))
535 :     {
536 :     my $cpd = $cinfo->[0];
537 :     my $main = $cinfo->[2];
538 :    
539 :     if (defined $all_compounds_to_main{$cpd})
540 :     {
541 :     # only main if it's main in all reactions
542 :     $all_compounds_to_main{$cpd} &= $main;
543 :     }
544 :     else
545 :     {
546 :     $all_compounds_to_main{$cpd} = $main;
547 :     }
548 :     }
549 :     }
550 :     else
551 :     {
552 :     # save subs and prods for processing at end
553 :     $reactions_not_in_maps{$rxn} = [ (@substrArr, @prodArr) ];
554 :     }
555 :     }
556 :     }
557 :    
558 :     # now load "main" designation based on reactions not in the maps - but don't overrule
559 :     # what's already been loaded
560 :     my %additional_compounds_to_main;
561 :    
562 :     foreach my $rxn (keys %reactions_not_in_maps)
563 :     {
564 :     foreach my $cinfo (@{$reactions_not_in_maps{$rxn}})
565 :     {
566 :     my $cpd = $cinfo->[0];
567 :     my $main = $cinfo->[2];
568 :    
569 :     if (defined $additional_compounds_to_main{$cpd})
570 :     {
571 :     # main if it's main in any reactions not in map
572 :     $additional_compounds_to_main{$cpd} |= $main;
573 :     }
574 :     else
575 :     {
576 :     $additional_compounds_to_main{$cpd} = $main;
577 :     }
578 :     }
579 :     }
580 :    
581 :     foreach my $cpd (keys %additional_compounds_to_main)
582 :     {
583 :     if (! defined $all_compounds_to_main{$cpd})
584 :     {
585 :     $all_compounds_to_main{$cpd} = $additional_compounds_to_main{$cpd};
586 :     }
587 :     }
588 :     }
589 :    
590 :     sub process_paths
591 :     {
592 :     my ($input_cpds, $output_cpds, $outputs_lists, $reactions_to_substrates, $reactions_to_products, $substrates_to_reactions, $products_to_reactions, $create_assembly, $find_first ) = @_;
593 :    
594 :     my (%path_inputs, %path_outputs);
595 :     map { $path_inputs{$_} = 0 } keys %$input_cpds;
596 :     map { $path_outputs{$_} = 0 } keys %$output_cpds;
597 :    
598 :     my %data_results = ("infinite" => 0);
599 :    
600 :     # %compounds_to_tokens maps from compound ids to tokens placed on those compounds
601 :     # the tokens are organized in hashes mapping from token id to number of tokens with that id
602 :     my %compounds_to_tokens;
603 :     map { $compounds_to_tokens{$_} = {} } keys %all_compounds_to_main;
604 :    
605 :     # %tokens maps from token_ids to the token data structures
606 :     my %tokens;
607 :     my $token_id_counter = 1;
608 :    
609 :     # %compounds_borrowed_to_tokens maps from compounds to lists of token ids that borrowed
610 :     # compound in order to run a reaction
611 :     my %compounds_borrowed_to_tokens;
612 :    
613 :     print STDERR "\nIn process_paths, path_inputs are @{[ keys %path_inputs ]}, path_outputs are @{[ keys %path_outputs] }, scenario_cycles are @{[ keys %scenario_cycles ]} \n\n" unless !$debug;
614 :    
615 :     my $initial_pass = 1;
616 :     my $done = 0;
617 :     my $loop_counter = 1;
618 :     my $infinite_loop_check = $create_assembly ? $loop_max_assembly : $loop_max;
619 :    
620 :     # we may get to the point where we need to add some more path inputs into the mix
621 :     # to push stalled tokens
622 :     my $add_path_inputs = 0;
623 :    
624 :     while(!$done)
625 :     {
626 :     if ($initial_pass || $add_path_inputs)
627 :     {
628 :     foreach my $cpd (keys %path_inputs)
629 :     {
630 :     # place a token on each path input
631 :     my $new_token_id = $token_id_counter++;
632 :     my %new_token;
633 :     $new_token{visited_reactions} = {};
634 :     $new_token{visited_compounds} = { $cpd => 0 }; # 0 means supplied from "outside"
635 :     $new_token{token_path_inputs} = { $cpd => 1 }; # 1 means one was supplied
636 :     $new_token{initial_pass} = $initial_pass;
637 :     $compounds_to_tokens{$cpd}->{$new_token_id}++;
638 :     $tokens{$new_token_id} = \%new_token;
639 :    
640 :     print STDERR "\t\tCreated new token '$new_token_id' for path input $cpd\n" unless !$debug;
641 :     }
642 :    
643 :     $initial_pass = 0;
644 :     $add_path_inputs = 0;
645 :     }
646 :    
647 :     # Find the reactions that can run.
648 :     my %reactions_to_try;
649 :    
650 :     foreach my $cpd (keys %compounds_to_tokens)
651 :     {
652 :     foreach my $token_id (keys %{$compounds_to_tokens{$cpd}})
653 :     {
654 :     if ($compounds_to_tokens{$cpd}->{$token_id} > 0 &&
655 :     ! $tokens{$token_id}->{done})
656 :     {
657 :     # this compound has tokens that aren't done
658 :     map { $reactions_to_try{$_} = 1 } @{$substrates_to_reactions->{$cpd}};
659 :     last;
660 :     }
661 :     }
662 :     }
663 :    
664 :     print STDERR "\n\tIn loop, trying reactions @{[keys %reactions_to_try]}\n\n" unless !$debug;
665 :    
666 :     # Map the reactions that can run to the tokens they can use
667 :     my %reactions_to_tokens_available;
668 :     # Keep track of the main substrates defined by each reaction
669 :     my %reactions_to_main_substrates;
670 :    
671 :     # count up the total number of tokens needed for each compound to run
672 :     # every reaction that is ready to go
673 :     my %reactions_to_tokens_needed;
674 :    
675 :     rxn: foreach my $reaction (keys %reactions_to_try)
676 :     {
677 :     print STDERR "\tChecking reaction $reaction\n" unless !$debug;
678 :    
679 :     my @substrArr = @{$reactions_to_substrate_arrays{$reaction}};
680 :     my @prodArr = @{$reactions_to_product_arrays{$reaction}};
681 :     my %main_substrates;
682 :    
683 :     # Determine if this reaction has necessary inputs to run.
684 :     # There must be tokens available for at least one main substrate that isn't
685 :     # a path output (unless it's an initial scenario cycled compound).
686 :     # Also, any path input must have a token.
687 :     # Count number of tokens needed for each main substrate.
688 :     my %tokens_available;
689 :     my %tokens_needed;
690 :     my $reaction_can_run = 0;
691 :    
692 :     foreach my $substr (@substrArr)
693 :     {
694 :     my $cpd = @{$substr}[0];
695 :     my $stoich = @{$substr}[1];
696 :     my $main = @{$substr}[2] || @all_compounds_to_main{$cpd}; #main either way
697 :     $main_substrates{$cpd} = 1 if $main;
698 :    
699 :     print STDERR "\t\tSubstrate: $cpd\tstoich: $stoich\tmain: $main\n" unless !$debug;
700 :    
701 :     if (! $main)
702 :     {
703 :     # on any pass we can take in non-main compounds
704 :     $tokens_needed{$cpd} = $stoich;
705 :     }
706 :     else
707 :     {
708 :     # if tokens are available for compound, check their history for the
709 :     # main compounds produced by the reaction so we don't loop back over
710 :     # previous main compounds. Also, don't use tokens on scenario inputs
711 :     # moved to in steps other than initial token creation.
712 :     my %ok_tokens;
713 :     my $num_ok_tokens = 0;
714 :    
715 :     foreach my $token_id (keys %{$compounds_to_tokens{$cpd}})
716 :     {
717 :     next if $tokens{$token_id}->{done} ||
718 :     $compounds_to_tokens{$cpd}->{$token_id} == 0;
719 :    
720 :     # now check that we aren't running an already visited reaction in reverse
721 :     my %visited_reactions = %{$tokens{$token_id}->{visited_reactions}};
722 :    
723 :     if (($reaction =~ /(.*)_R/ && defined $visited_reactions{$1."_L"}) ||
724 :     ($reaction =~ /(.*)_L/ && defined $visited_reactions{$1."_R"}))
725 :     {
726 :     print STDERR "\t\tToken '$token_id' has run the reverse reaction already\n" unless !$debug;
727 :     next;
728 :     }
729 :    
730 :     my %visited_compounds = %{$tokens{$token_id}->{visited_compounds}};
731 :    
732 :     print STDERR "\t\tToken '$token_id' has visited @{[map { ($_, $visited_compounds{$_} ) } sort { $visited_compounds{$a} <=> $visited_compounds{$b} } keys %visited_compounds]}\n" unless !$debug;
733 :    
734 :     # need to find at least one prod that hasn't been visited yet
735 :     # or was visited in a loop cycle not before the loop cycle in which
736 :     # it visited the substrate,
737 :     # or is a path output
738 :     my $prods_are_ok = 0;
739 :    
740 :     # check each main product
741 :     foreach my $prod (@{$reactions_to_products->{$reaction}})
742 :     {
743 :     if (! defined $visited_compounds{$prod} ||
744 :     $visited_compounds{$prod} >= $visited_compounds{$cpd} ||
745 :     defined $path_outputs{$prod} ||
746 :     $compounds_borrowed_to_tokens{$prod}->{$token_id} > 0)
747 :     {
748 :     print STDERR "\t\tToken can visit $prod\n" unless !$debug;
749 :     $prods_are_ok = 1;
750 :     last;
751 :     }
752 :     }
753 :    
754 :     if ($prods_are_ok)
755 :     {
756 :     print STDERR "\t\tToken is OK\n" unless !$debug;
757 :     $ok_tokens{$token_id} = $compounds_to_tokens{$cpd}->{$token_id};
758 :     $num_ok_tokens += $compounds_to_tokens{$cpd}->{$token_id};
759 :     }
760 :     }
761 :    
762 :     map { $tokens_available{$_}->{$cpd} = $ok_tokens{$_} } keys %ok_tokens;
763 :     $tokens_needed{$cpd} = $stoich;
764 :    
765 :     if ($main && $num_ok_tokens >= 1)
766 :     {
767 :     if (! defined $path_outputs{$cpd} || defined $scenario_cycles{$cpd})
768 :     {
769 :     print STDERR "\t\tgot at least one token on main compound\n" unless !$debug;
770 :     $reaction_can_run = 1;
771 :     }
772 :     }
773 :     elsif (defined $path_inputs{$cpd})
774 :     {
775 :     print STDERR "\t\tno tokens available for path input: $cpd\n" unless !$debug;
776 :     next rxn;
777 :     }
778 :     elsif (defined $path_outputs{$cpd} && ! defined $scenario_cycles{$cpd})
779 :     {
780 :     print STDERR "\t\tno tokens available for path output: $cpd\n" unless !$debug;
781 :     next rxn;
782 :     }
783 :     }
784 :     }
785 :    
786 :     if ($reaction_can_run)
787 :     {
788 :     $reactions_to_tokens_available{$reaction} = \%tokens_available;
789 :     $reactions_to_tokens_needed{$reaction} = \%tokens_needed;
790 :     $reactions_to_main_substrates{$reaction} = \%main_substrates;
791 :     print STDERR "\tReaction $reaction can run\n" unless !$debug;
792 :     }
793 :     }
794 :    
795 :     # keep track of tokens used that will be used to run rxns. Clone tokens if necessary.
796 :     my %reactions_to_tokens_to_use;
797 :     my %tokens_to_use_to_reactions;
798 :     my %copy_of_compounds_to_tokens; # for determining which reaction uses which tokens
799 :    
800 :     foreach my $cpd (keys %compounds_to_tokens)
801 :     {
802 :     my %cpd_token_ids = %{$compounds_to_tokens{$cpd}};
803 :     my %new_cpd_token_ids;
804 :     map { $new_cpd_token_ids{$_} = $cpd_token_ids{$_} } keys %cpd_token_ids;
805 :     $copy_of_compounds_to_tokens{$cpd} = \%new_cpd_token_ids;
806 :     }
807 :    
808 :     foreach my $reaction (keys %reactions_to_tokens_available)
809 :     {
810 :     print STDERR "\tPreparing to run reaction $reaction\n" unless !$debug;
811 :    
812 :     # assemble tokens to run reaction, cloning ones that were used by
813 :     # other reactions during this cycle if necessary
814 :     my %tokens_available = %{$reactions_to_tokens_available{$reaction}};
815 :     my %tokens_needed = %{$reactions_to_tokens_needed{$reaction}};
816 :     my @final_tokens_to_use; # list of maps from substrates to tokens to use
817 :     my %clone_history; # map from token ids to ids of their new clones
818 :    
819 :     # check to see if any available tokens are already commited to the reverse reaction
820 :     my @token_ids = keys %tokens_available;
821 :    
822 :     foreach my $token_id (keys %tokens_available)
823 :     {
824 :     if (($reaction =~ /(.*)_R/ && defined $tokens_to_use_to_reactions{$token_id}->{$1."_L"}) ||
825 :     ($reaction =~ /(.*)_L/ && defined $tokens_to_use_to_reactions{$token_id}->{$1."_R"}))
826 :     {
827 :     # clone the token
828 :     my $new_token_id = $token_id_counter++;
829 :     &clone_token($token_id, $new_token_id, \%tokens, \%compounds_to_tokens,
830 :     \%compounds_borrowed_to_tokens);
831 :    
832 :     foreach my $icpd (keys %compounds_to_tokens)
833 :     {
834 :     if ($compounds_to_tokens{$icpd}->{$token_id} > 0)
835 :     {
836 :     $copy_of_compounds_to_tokens{$icpd}->{$new_token_id} =
837 :     $compounds_to_tokens{$icpd}->{$token_id};
838 :     }
839 :     }
840 :    
841 :     print STDERR "\t\tCloned token '$token_id', new token is '$new_token_id'\n" unless !$debug;
842 :    
843 :     $tokens_available{$new_token_id} = $tokens_available{$token_id};
844 :     delete $tokens_available{$token_id};
845 :     }
846 :     }
847 :    
848 :     # first, assemble tokens that have all the main compounds they need to run
849 :     @token_ids = keys %tokens_available;
850 :    
851 :     foreach my $token_id (keys %tokens_available)
852 :     {
853 :     my $has_all_main_cpds = 1;
854 :     my $need_to_clone = 0;
855 :    
856 :     foreach my $cpd (keys %tokens_needed)
857 :     {
858 :     if ($reactions_to_main_substrates{$reaction}->{$cpd})
859 :     {
860 :     if ($compounds_to_tokens{$cpd}->{$token_id} >= $tokens_needed{$cpd})
861 :     {
862 :     if ($copy_of_compounds_to_tokens{$cpd}->{$token_id} < $tokens_needed{$cpd})
863 :     {
864 :     $need_to_clone = 1;
865 :     }
866 :     }
867 :     else
868 :     {
869 :     $has_all_main_cpds = 0;
870 :     last;
871 :     }
872 :     }
873 :     }
874 :    
875 :     if ($has_all_main_cpds)
876 :     {
877 :     print STDERR "\t\ttoken '$token_id' has all main compounds\n" unless !$debug;
878 :    
879 :     delete $tokens_available{$token_id};
880 :    
881 :     if ($need_to_clone)
882 :     {
883 :     my $new_token_id = $token_id_counter++;
884 :     &clone_token($token_id, $new_token_id, \%tokens,
885 :     \%compounds_to_tokens,
886 :     \%compounds_borrowed_to_tokens);
887 :    
888 :     foreach my $icpd (keys %compounds_to_tokens)
889 :     {
890 :     if ($compounds_to_tokens{$icpd}->{$token_id} > 0)
891 :     {
892 :     $copy_of_compounds_to_tokens{$icpd}->{$new_token_id} =
893 :     $compounds_to_tokens{$icpd}->{$token_id};
894 :     }
895 :     }
896 :    
897 :     print STDERR "\t\tCloned token '$token_id', new token is '$new_token_id'\n" unless !$debug;
898 :    
899 :     $clone_history{$token_id} = $new_token_id;
900 :     $token_id = $new_token_id;
901 :     }
902 :    
903 :     # now assemble the compound to tokens map
904 :     my %cpd_to_tokens;
905 :    
906 :     foreach my $cpd (keys %tokens_needed)
907 :     {
908 :     if ($reactions_to_main_substrates{$reaction}->{$cpd})
909 :     {
910 :     for (my $i = 0; $i < $tokens_needed{$cpd}; $i++)
911 :     {
912 :     push @{$cpd_to_tokens{$cpd}}, $token_id;
913 :     $copy_of_compounds_to_tokens{$cpd}->{$token_id}--;
914 :     }
915 :     }
916 :     }
917 :    
918 :     push @final_tokens_to_use, \%cpd_to_tokens;
919 :     }
920 :     }
921 :    
922 :     if (scalar keys %tokens_available > 0)
923 :     {
924 :     # try to merge left over available tokens in all combinations that fulfill
925 :     # needed main substrates. Create map from main substrates to tokens.
926 :     my %tokens_available_for_cpds;
927 :     my %tokens_to_use_for_cpds;
928 :    
929 :     foreach my $token_id (keys %tokens_available)
930 :     {
931 :     foreach my $cpd (keys %{$tokens_available{$token_id}})
932 :     {
933 :     if ($reactions_to_main_substrates{$reaction}->{$cpd})
934 :     {
935 :     $tokens_available_for_cpds{$cpd}->{$token_id} =
936 :     $tokens_available{$token_id}->{$cpd};
937 :     }
938 :     }
939 :     }
940 :    
941 :     foreach my $cpd (keys %{$reactions_to_main_substrates{$reaction}})
942 :     {
943 :     my %available_token_ids;
944 :    
945 :     foreach my $token_id (keys %{$tokens_available_for_cpds{$cpd}})
946 :     {
947 :     my $updated_token_id = $token_id;
948 :    
949 :     # in case the token has already been cloned for another cpd in this reaction
950 :     while ($clone_history{$updated_token_id})
951 :     {
952 :     $updated_token_id = $clone_history{$updated_token_id};
953 :     }
954 :    
955 :     $available_token_ids{$updated_token_id} = 1;
956 :    
957 :     if ($token_id != $updated_token_id)
958 :     {
959 :     $tokens_available_for_cpds{$cpd}->{$updated_token_id} =
960 :     $tokens_available_for_cpds{$cpd}->{$token_id};
961 :     delete $tokens_available_for_cpds{$cpd}->{$token_id};
962 :     }
963 :     }
964 :    
965 :     my $num_tokens_needed = $tokens_needed{$cpd};
966 :     my $num_available_tokens = 0;
967 :    
968 :     foreach my $token_id (keys %available_token_ids)
969 :     {
970 :     $num_available_tokens += $tokens_available_for_cpds{$cpd}->{$token_id};
971 :     }
972 :    
973 :     # if not enough tokens are available to fill out sets, create new ones
974 :     my $num_short_of_full;
975 :    
976 :     if ($num_available_tokens == 0)
977 :     {
978 :     $num_short_of_full = $num_tokens_needed;
979 :     }
980 :     else
981 :     {
982 :     $num_short_of_full = ($num_tokens_needed - ($num_available_tokens % $num_tokens_needed)) % $num_tokens_needed;
983 :     }
984 :    
985 :     print STDERR "\t\tFor $cpd, $num_tokens_needed tokens are needed, $num_available_tokens are available. Need to create $num_short_of_full to fill out sets\n" unless !$debug;
986 :    
987 :     for (my $i = 0; $i < $num_short_of_full; $i++)
988 :     {
989 :     my $new_token_id = $token_id_counter++;
990 :     my %new_token;
991 :     print STDERR "\t\tCreating new token '$new_token_id' for $cpd\n" unless !$debug;
992 :    
993 :     $new_token{visited_reactions} = {};
994 :     $new_token{visited_compounds} = {};
995 :     $compounds_to_tokens{$cpd}->{$new_token_id}++;
996 :     $copy_of_compounds_to_tokens{$cpd}->{$new_token_id}++;
997 :     $tokens{$new_token_id} = \%new_token;
998 :     $available_token_ids{$new_token_id} = 1;
999 :     $tokens_available_for_cpds{$cpd}->{$new_token_id}++;
1000 :    
1001 :     # if it's not a path input, remember that we've "borrowed" it and
1002 :     # will need to pay it back.
1003 :     if (! defined $path_inputs{$cpd})
1004 :     {
1005 :     $new_token{visited_compounds}->{$cpd} = $loop_counter;
1006 :     $compounds_borrowed_to_tokens{$cpd}->{$new_token_id}++;
1007 :     }
1008 :     else
1009 :     {
1010 :     $new_token{token_path_inputs} = { $cpd => 1 };
1011 :     $new_token{visited_compounds}->{$cpd} = 0;
1012 :     }
1013 :     }
1014 :    
1015 :     # for main compounds, there may be more tokens available than needed,
1016 :     # so we may assemble multiple token sets.
1017 :     my %tokens_not_yet_used;
1018 :    
1019 :     foreach my $token_id (keys %available_token_ids)
1020 :     {
1021 :     if ($copy_of_compounds_to_tokens{$cpd}->{$token_id} > 0)
1022 :     {
1023 :     $tokens_not_yet_used{$token_id} = $tokens_available_for_cpds{$cpd}->{$token_id};
1024 :     }
1025 :     }
1026 :    
1027 :     my @token_sets_for_cpd;
1028 :    
1029 :     print STDERR "\t\tNeed $num_tokens_needed tokens for $cpd, '@{[ sort { $a <=> $b } keys %available_token_ids ]}' are usable, '@{[ sort { $a <=> $b } keys %tokens_not_yet_used ]}' are not yet used\n" unless !$debug;
1030 :    
1031 :     my @token_set;
1032 :     my $num_tokens_still_needed = $num_tokens_needed;
1033 :    
1034 :     foreach my $token_id (sort { $a <=> $b } keys %available_token_ids)
1035 :     {
1036 :     # need to clone the token if it is all used up
1037 :     if ($copy_of_compounds_to_tokens{$cpd}->{$token_id} == 0)
1038 :     {
1039 :     my $new_token_id = $token_id_counter++;
1040 :     &clone_token($token_id, $new_token_id, \%tokens,
1041 :     \%compounds_to_tokens,
1042 :     \%compounds_borrowed_to_tokens);
1043 :    
1044 :     foreach my $icpd (keys %compounds_to_tokens)
1045 :     {
1046 :     if ($compounds_to_tokens{$icpd}->{$token_id} > 0)
1047 :     {
1048 :     $copy_of_compounds_to_tokens{$icpd}->{$new_token_id} =
1049 :     $compounds_to_tokens{$icpd}->{$token_id};
1050 :     }
1051 :     }
1052 :    
1053 :     print STDERR "\t\tCloned token '$token_id' for $cpd, new token is '$new_token_id'\n" unless !$debug;
1054 :    
1055 :     $clone_history{$token_id} = $new_token_id;
1056 :     $token_id = $new_token_id;
1057 :     }
1058 :    
1059 :     while ($copy_of_compounds_to_tokens{$cpd}->{$token_id} > 0)
1060 :     {
1061 :     push @token_set, $token_id;
1062 :     $copy_of_compounds_to_tokens{$cpd}->{$token_id}--;
1063 :     $num_tokens_still_needed--;
1064 :    
1065 :     if ($num_tokens_still_needed == 0)
1066 :     {
1067 :     print STDERR "\t\tPushing token set '@token_set' for $cpd\n" unless !$debug;
1068 :     my @copy_of_token_set = @token_set;
1069 :     push @token_sets_for_cpd, \@copy_of_token_set;
1070 :     @token_set = ();
1071 :     $num_tokens_still_needed = $num_tokens_needed;
1072 :     }
1073 :     }
1074 :     }
1075 :    
1076 :     $tokens_to_use_for_cpds{$cpd} = \@token_sets_for_cpd;
1077 :     }
1078 :    
1079 :     # in case a token had to be cloned for this reaction after another
1080 :     # compound already determined to use it, check history and use new token id
1081 :     foreach my $cpd (keys %tokens_to_use_for_cpds)
1082 :     {
1083 :     foreach my $token_set (@{$tokens_to_use_for_cpds{$cpd}})
1084 :     {
1085 :     for (my $i = 0; $i < scalar @$token_set; $i++)
1086 :     {
1087 :     my $updated_token_id = $token_set->[$i];
1088 :    
1089 :     while ($clone_history{$updated_token_id})
1090 :     {
1091 :     $updated_token_id = $clone_history{$updated_token_id};
1092 :     }
1093 :    
1094 :     if ($token_set->[$i] != $updated_token_id)
1095 :     {
1096 :     print STDERR "\t\tReplacing '$token_set->[$i]' with '$updated_token_id' for $cpd\n" unless !$debug;
1097 :     splice @$token_set, $i, 1, ($updated_token_id);
1098 :     }
1099 :     }
1100 :     }
1101 :     }
1102 :    
1103 :     # we may have multiple sets of tokens to use for some compounds, so
1104 :     # make sure we're prepared for each combination by cloning sets as necessary
1105 :     my %cpd_to_token_set_index;
1106 :     my $num_combinations = 1;
1107 :     my %used_token_sets_for_cpd;
1108 :    
1109 :     foreach my $cpd (keys %tokens_to_use_for_cpds)
1110 :     {
1111 :     $cpd_to_token_set_index{$cpd} = 0;
1112 :     $num_combinations *= scalar @{$tokens_to_use_for_cpds{$cpd}};
1113 :     }
1114 :    
1115 :     for (my $i = 0; $i < $num_combinations; $i++)
1116 :     {
1117 :     my %combination_cpds_to_tokens;
1118 :    
1119 :     print STDERR "\t\tPreparing combination ", $i+1, " out of $num_combinations, indices are @{[ map { $cpd_to_token_set_index{$_} } sort keys %cpd_to_token_set_index ]}\n" unless !$debug;
1120 :    
1121 :     # keep track of who was cloned to what for this combination in case
1122 :     # different compounds are using the same tokens
1123 :     my %clone_history_this_combination;
1124 :    
1125 :     foreach my $cpd (sort keys %cpd_to_token_set_index)
1126 :     {
1127 :     my $token_set_index = $cpd_to_token_set_index{$cpd};
1128 :     my @token_set = @{$tokens_to_use_for_cpds{$cpd}->[$token_set_index]};
1129 :     my @new_token_set;
1130 :    
1131 :     # don't clone if this is the first time using this token set
1132 :     # should be able to do this mathematically
1133 :     if (! defined $used_token_sets_for_cpd{$cpd}->{$token_set_index})
1134 :     {
1135 :     @new_token_set = @token_set;
1136 :     $used_token_sets_for_cpd{$cpd}->{$token_set_index} = 1;
1137 :     }
1138 :     else
1139 :     {
1140 :     foreach my $token_id (@token_set)
1141 :     {
1142 :     my $new_token_id;
1143 :    
1144 :     if ($clone_history_this_combination{$token_id})
1145 :     {
1146 :     $new_token_id = $clone_history_this_combination{$token_id};
1147 :     }
1148 :     else
1149 :     {
1150 :     $new_token_id = $token_id_counter++;
1151 :     &clone_token($token_id, $new_token_id, \%tokens,
1152 :     \%compounds_to_tokens,
1153 :     \%compounds_borrowed_to_tokens);
1154 :    
1155 :     print STDERR "\t\t\tCloned token '$token_id' for $cpd, new token is '$new_token_id'\n" unless !$debug;
1156 :    
1157 :     $clone_history_this_combination{$token_id} = $new_token_id;
1158 :     }
1159 :    
1160 :     push @new_token_set, $new_token_id;
1161 :     }
1162 :     }
1163 :    
1164 :     push @{$combination_cpds_to_tokens{$cpd}}, @new_token_set;
1165 :     }
1166 :    
1167 :     push @final_tokens_to_use, \%combination_cpds_to_tokens;
1168 :    
1169 :     # move the cpd to token set indices for the next combination
1170 :     foreach my $cpd (sort keys %cpd_to_token_set_index)
1171 :     {
1172 :     if ($cpd_to_token_set_index{$cpd} < scalar @{$tokens_to_use_for_cpds{$cpd}} - 1)
1173 :     {
1174 :     $cpd_to_token_set_index{$cpd}++;
1175 :     last;
1176 :     }
1177 :     else
1178 :     {
1179 :     $cpd_to_token_set_index{$cpd} = 0;
1180 :     }
1181 :     }
1182 :     }
1183 :     }
1184 :    
1185 :     # last step, need to create place-holder tokens for the non-main substrates
1186 :     # and record the tokens being used for this reaction
1187 :     foreach my $combination (@final_tokens_to_use)
1188 :     {
1189 :     foreach my $cpd (keys %tokens_needed)
1190 :     {
1191 :     if (! defined $reactions_to_main_substrates{$reaction}->{$cpd})
1192 :     {
1193 :     my $new_token_id = $token_id_counter++;
1194 :     my %new_token;
1195 :     $new_token{visited_reactions} = {};
1196 :     $new_token{visited_compounds} = {};
1197 :     $compounds_to_tokens{$cpd}->{$new_token_id} += $tokens_needed{$cpd};
1198 :     $copy_of_compounds_to_tokens{$cpd}->{$new_token_id} += $tokens_needed{$cpd};
1199 :     $new_token{token_path_inputs} = { $cpd => $tokens_needed{$cpd} };
1200 :     $tokens{$new_token_id} = \%new_token;
1201 :    
1202 :     for (my $i = 0; $i < $tokens_needed{$cpd}; $i++)
1203 :     {
1204 :     push @{$combination->{$cpd}}, $new_token_id;
1205 :     }
1206 :     }
1207 :     }
1208 :    
1209 :     foreach my $cpd (keys %$combination)
1210 :     {
1211 :     map { $tokens_to_use_to_reactions{$_}->{$reaction} = 1 } @{$combination->{$cpd}};
1212 :     }
1213 :     }
1214 :    
1215 :     $reactions_to_tokens_to_use{$reaction} = \@final_tokens_to_use;
1216 :     }
1217 :    
1218 :     # keep track of tokens merged during this round
1219 :     my %token_merge_history;
1220 :    
1221 :     foreach my $reaction (keys %reactions_to_tokens_to_use)
1222 :     {
1223 :     # we may have multiple sets of tokens to use
1224 :     my @final_tokens_to_use = @{$reactions_to_tokens_to_use{$reaction}};
1225 :     my $num_combinations = scalar @final_tokens_to_use;
1226 :    
1227 :     # since we may be running this reaction several times with some of the same
1228 :     # tokens, don't process token merge history until all sets have been run
1229 :     my %token_merge_history_this_reaction;
1230 :    
1231 :     for (my $i = 0; $i < $num_combinations; $i++)
1232 :     {
1233 :     print STDERR "\tRunning reaction $reaction (", $i+1, " out of $num_combinations)\n" unless !$debug;
1234 :    
1235 :     # assemble list of tokens to use for this combination
1236 :     my %tokens_to_use;
1237 :     # find the most recently visited compound's step
1238 :     my $most_recent_step = 0;
1239 :    
1240 :     # loop through substrates
1241 :     foreach my $cpd (keys %{$reactions_to_tokens_to_use{$reaction}})
1242 :     {
1243 :     # remove the token ids from the real compounds_to_tokens map
1244 :     my @cpd_token_ids = @{$final_tokens_to_use[$i]->{$cpd}};
1245 :     print STDERR "\t\tFound tokens '@cpd_token_ids' for $cpd\n" unless !$debug;
1246 :    
1247 :     foreach my $token_id (@cpd_token_ids)
1248 :     {
1249 :     # check if token_id has been merged into a new token by a previous reaction
1250 :     while (defined $token_merge_history{$token_id})
1251 :     {
1252 :     $token_id = $token_merge_history{$token_id};
1253 :     }
1254 :    
1255 :     $compounds_to_tokens{$cpd}->{$token_id}--;
1256 :    
1257 :     if ($compounds_to_tokens{$cpd}->{$token_id} == 0)
1258 :     {
1259 :     delete $compounds_to_tokens{$cpd}->{$token_id};
1260 :     }
1261 :    
1262 :     $tokens_to_use{$token_id} = 1;
1263 :    
1264 :     if ($tokens{$token_id}->{visited_compounds}->{$cpd} > $most_recent_step)
1265 :     {
1266 :     $most_recent_step = $tokens{$token_id}->{visited_compounds}->{$cpd};
1267 :     }
1268 :     }
1269 :     }
1270 :    
1271 :     # process list of unique token ids
1272 :     my @tokens_to_use = sort { $a <=> $b } keys %tokens_to_use;
1273 :     my $go_forward_token_id;
1274 :    
1275 :     if (scalar @tokens_to_use == 1)
1276 :     {
1277 :     $go_forward_token_id = shift @tokens_to_use;
1278 :     print STDERR "\t\tGoing forward with '$go_forward_token_id'\n" unless !$debug;
1279 :     }
1280 :     else
1281 :     {
1282 :     $go_forward_token_id = $token_id_counter++;
1283 :     $tokens{$go_forward_token_id} = {};
1284 :     print STDERR "\t\tRemember to merge tokens '@tokens_to_use' into '$go_forward_token_id'\n" unless !$debug;
1285 :    
1286 :     # record the need to merge - we'll do it after processing all sets of substrates
1287 :     foreach my $token_id (@tokens_to_use)
1288 :     {
1289 :     push @{$token_merge_history_this_reaction{$token_id}}, $go_forward_token_id;
1290 :     }
1291 :     }
1292 :    
1293 :     my $go_forward_token = $tokens{$go_forward_token_id};
1294 :     my @prodArr = @{$reactions_to_product_arrays{$reaction}};
1295 :    
1296 :     # add current reaction and products to accumulated token history.
1297 :     # reaction is mapped to loop counter to maintain history of order of execution
1298 :     $go_forward_token->{visited_reactions}->{$reaction} = $loop_counter;
1299 :    
1300 :     foreach my $prod (@prodArr)
1301 :     {
1302 :     my $cpd = @{$prod}[0];
1303 :     my $stoich = @{$prod}[1];
1304 :     my $main = @{$prod}[2];
1305 :    
1306 :     # keep track of path outputs we've seen
1307 :     if (defined $path_outputs{$cpd})
1308 :     {
1309 :     $path_outputs{$cpd} += $stoich;
1310 :     }
1311 :    
1312 :     if ($main)
1313 :     {
1314 :     $go_forward_token->{visited_compounds}->{$cpd} = $most_recent_step + 1;
1315 :     }
1316 :    
1317 :     # push tokens
1318 :     $compounds_to_tokens{$cpd}->{$go_forward_token_id} += $stoich;
1319 :     }
1320 :     }
1321 :    
1322 :     # now process this reaction's token merge history
1323 :     foreach my $token_id (sort { $a <=> $b } keys %token_merge_history_this_reaction)
1324 :     {
1325 :     # find the unique set of up to date merge ids
1326 :     my @merge_ids = @{$token_merge_history_this_reaction{$token_id}};
1327 :     my %updated_merge_ids;
1328 :    
1329 :     foreach my $token_id (@merge_ids)
1330 :     {
1331 :     while (defined $token_merge_history{$token_id})
1332 :     {
1333 :     $token_id = $token_merge_history{$token_id};
1334 :     }
1335 :    
1336 :     $updated_merge_ids{$token_id} = 1;
1337 :     }
1338 :    
1339 :     print STDERR "\t\tupdated merge id list for '$token_id': '@{[ keys %updated_merge_ids ]}'\n" unless !$debug;
1340 :    
1341 :     my $wrap_up_token_id = $token_id_counter++;
1342 :     $tokens{$wrap_up_token_id} = {};
1343 :     my $wrap_up_token = $tokens{$wrap_up_token_id};
1344 :     my @tokens_to_merge;
1345 :     push @tokens_to_merge, keys %updated_merge_ids, $token_id;
1346 :    
1347 :     # merge from oldest to youngest to update visited_compounds history
1348 :     foreach my $itoken_id (sort { $a <=> $b } @tokens_to_merge)
1349 :     {
1350 :     print STDERR "\t\t\tmerging '$itoken_id' into '$wrap_up_token_id'\n" unless !$debug;
1351 :    
1352 :     my $itoken = $tokens{$itoken_id};
1353 :     map { $wrap_up_token->{visited_reactions}->{$_} = $itoken->{visited_reactions}->{$_} } keys %{$itoken->{visited_reactions}};
1354 :     map { $wrap_up_token->{visited_compounds}->{$_} = $itoken->{visited_compounds}->{$_} } keys %{$itoken->{visited_compounds}};
1355 :     map { $wrap_up_token->{token_path_inputs}->{$_} += $itoken->{token_path_inputs}->{$_} } keys %{$itoken->{token_path_inputs}};
1356 :     $wrap_up_token->{initial_pass} |= $itoken->{initial_pass};
1357 :     $token_merge_history{$itoken_id} = $wrap_up_token_id;
1358 :    
1359 :     # tokens might be spread across multiple compounds; change them all
1360 :     # to new id
1361 :     foreach my $cpd (keys %compounds_to_tokens)
1362 :     {
1363 :     if ($compounds_to_tokens{$cpd}->{$itoken_id} > 0)
1364 :     {
1365 :     $compounds_to_tokens{$cpd}->{$wrap_up_token_id} +=
1366 :     $compounds_to_tokens{$cpd}->{$itoken_id};
1367 :     }
1368 :     }
1369 :    
1370 :     # tokens might have borrowed compounds; change them all to new id
1371 :     foreach my $cpd (keys %compounds_borrowed_to_tokens)
1372 :     {
1373 :     if ($compounds_borrowed_to_tokens{$cpd}->{$itoken_id} > 0)
1374 :     {
1375 :     $compounds_borrowed_to_tokens{$cpd}->{$wrap_up_token_id} +=
1376 :     $compounds_borrowed_to_tokens{$cpd}->{$itoken_id};
1377 :     }
1378 :     }
1379 :     }
1380 :     }
1381 :     }
1382 :    
1383 :     # now delete the tokens that were used merged in these reactions
1384 :     foreach my $token_id (keys %token_merge_history)
1385 :     {
1386 :     foreach my $icpd (keys %compounds_to_tokens)
1387 :     {
1388 :     delete $compounds_to_tokens{$icpd}->{$token_id};
1389 :     }
1390 :    
1391 :     foreach my $icpd (keys %compounds_borrowed_to_tokens)
1392 :     {
1393 :     delete $compounds_borrowed_to_tokens{$icpd}->{$token_id};
1394 :     }
1395 :    
1396 :     print STDERR "\t\tDeleting token '$token_id'\n" unless !$debug;
1397 :     delete $tokens{$token_id};
1398 :     }
1399 :    
1400 :     print STDERR "\nBalancing tokens\n" unless !$debug;
1401 :    
1402 :     foreach my $token_id (keys %tokens)
1403 :     {
1404 :     &balance_borrowing_and_giving($token_id, \%compounds_to_tokens,
1405 :     \%compounds_borrowed_to_tokens);
1406 :     }
1407 :    
1408 :     &print_token_status([sort { $a <=> $b } keys %tokens], \%tokens, \%compounds_to_tokens, \%compounds_borrowed_to_tokens, $fig);
1409 :    
1410 :     print STDERR "\nChecking for done\n" unless !$debug;
1411 :    
1412 :     print STDERR "\n\ntoken ids: @{[ sort { $a <=> $b } map { $_ if ! defined $tokens{$_}->{done} } keys %tokens ]}\n" unless !$debug;
1413 :    
1414 :     # we're done when all the main compounds in initial-pass tokens
1415 :     # have reached path outputs and repaid their borrowed tokens,
1416 :     # or have reached a dead end.
1417 :     # Check if we're done pushing and borrowing tokens first.
1418 :     my %not_done_tokens;
1419 :    
1420 :     foreach my $token_id (keys %tokens)
1421 :     {
1422 :     if (! defined $tokens{$token_id}->{done})
1423 :     {
1424 :     if (! &check_token_for_done($token_id, \%compounds_to_tokens,
1425 :     \%compounds_borrowed_to_tokens,
1426 :     \%all_compounds_to_main, \%path_outputs,
1427 :     \%scenario_cycles, \%tokens, $fig, $outputs_lists))
1428 :     {
1429 :     $not_done_tokens{$token_id} = 1;
1430 :     }
1431 :     }
1432 :     }
1433 :    
1434 :     print STDERR "\nChecking if we can pay back borrowed compounds from other tokens\n" unless !$debug;
1435 :    
1436 :     foreach my $bcpd (keys %compounds_borrowed_to_tokens)
1437 :     {
1438 :     next if $path_outputs{$bcpd}; # tokens must manage their own path outputs
1439 :    
1440 :     if (scalar keys %{$compounds_borrowed_to_tokens{$bcpd}} > 0 &&
1441 :     scalar keys %{$compounds_to_tokens{$bcpd}} > 0)
1442 :     {
1443 :     my %borrowers_to_givers;
1444 :    
1445 :     foreach my $borrower_id (keys %{$compounds_borrowed_to_tokens{$bcpd}})
1446 :     {
1447 :     next if defined $tokens{$borrower_id}->{done}; # don't repay deadenders
1448 :    
1449 :     my $num_needed = $compounds_borrowed_to_tokens{$bcpd}->{$borrower_id};
1450 :    
1451 :     giver: foreach my $giver_id (keys %{$compounds_to_tokens{$bcpd}})
1452 :     {
1453 :     next if defined $tokens{$giver_id}->{done};
1454 :    
1455 :     my $num_to_give = $compounds_to_tokens{$bcpd}->{$giver_id};
1456 :    
1457 :     print STDERR "\tToken '$giver_id' has $num_to_give $bcpd to give to '$borrower_id', which needs $num_needed\n" unless !$debug;
1458 :    
1459 :     my $borrower = $tokens{$borrower_id};
1460 :     my $giver = $tokens{$giver_id};
1461 :    
1462 :     # check whether the giver and borrower have conflicting histories
1463 :     foreach my $visited_reaction (keys %{$giver->{visited_reactions}})
1464 :     {
1465 :     if (($visited_reaction =~ /(.*)_R/ &&
1466 :     defined $borrower->{visited_reactions}->{$1."_L"}) ||
1467 :     ($visited_reaction =~ /(.*)_L/ &&
1468 :     defined $borrower->{visited_reactions}->{$1."_R"}))
1469 :     {
1470 :     print STDERR "\t\tConflict on $visited_reaction\n" unless !$debug;
1471 :     next giver;
1472 :     }
1473 :    
1474 :     }
1475 :    
1476 :     push @{$borrowers_to_givers{$borrower_id}}, $giver_id;
1477 :     }
1478 :     }
1479 :    
1480 :     # we have a list of potential givers for each borrower for this compound.
1481 :     # Now figure out who the lucky givers will be.
1482 :     my %givers_to_borrowers;
1483 :    
1484 :     foreach my $borrower_id (keys %borrowers_to_givers)
1485 :     {
1486 :     my %lucky_givers;
1487 :     my $num_needed = $compounds_borrowed_to_tokens{$bcpd}->{$borrower_id};
1488 :    
1489 :     print STDERR "\tCollecting lucky givers for '$borrower_id' for $bcpd, need $num_needed\n" unless !$debug;
1490 :    
1491 :     # check potential givers starting with those with the most
1492 :     my @potential_givers = reverse sort { $compounds_to_tokens{$bcpd}->{$a} <=> $compounds_to_tokens{$bcpd}->{$b} } @{$borrowers_to_givers{$borrower_id}};
1493 :    
1494 :     foreach my $giver_id (@potential_givers)
1495 :     {
1496 :     if (! defined $lucky_givers{$giver_id})
1497 :     {
1498 :     my $num_to_give = $compounds_to_tokens{$bcpd}->{$giver_id};
1499 :     print STDERR "\t\t'$giver_id' has $num_to_give to give\n" unless !$debug;
1500 :     $lucky_givers{$giver_id} = 1;
1501 :     $num_needed -= $num_to_give;
1502 :     last if $num_needed <= 0;
1503 :     }
1504 :     }
1505 :    
1506 :     foreach my $giver_id (keys %lucky_givers)
1507 :     {
1508 :     push @{$givers_to_borrowers{$giver_id}}, $borrower_id;
1509 :     }
1510 :     }
1511 :    
1512 :    
1513 :     foreach my $orig_giver_id (keys %givers_to_borrowers)
1514 :     {
1515 :     my @borrowers_list = @{$givers_to_borrowers{$orig_giver_id}};
1516 :     my @givers_list = ($orig_giver_id);
1517 :    
1518 :     # clone enough givers so that every borrower gets one. Last borrower
1519 :     # gets the orgiinal giver.
1520 :    
1521 :     while (scalar @borrowers_list > scalar @givers_list)
1522 :     {
1523 :     my $new_giver_id = $token_id_counter++;
1524 :     &clone_token($orig_giver_id, $new_giver_id, \%tokens,
1525 :     \%compounds_to_tokens, \%compounds_borrowed_to_tokens);
1526 :     push @givers_list, $new_giver_id;
1527 :     }
1528 :    
1529 :     for (my $k = 0; $k < scalar @givers_list; $k++)
1530 :     {
1531 :     my $giver_id = $givers_list[$k];
1532 :     my $giver = $tokens{$giver_id};
1533 :     my $borrower_id = $borrowers_list[$k];
1534 :     my $borrower = $tokens{$borrower_id};
1535 :    
1536 :     print STDERR "\n\tMerging '$giver_id' into '$borrower_id'\n" unless !$debug;
1537 :    
1538 :     # bump the borrower's visited reactions and compounds counters forward,
1539 :     # then merge the giver's visited reactions and compounds, unless
1540 :     # the borrower has already visited them
1541 :     map { $borrower->{visited_reactions}->{$_} += $loop_counter } keys %{$borrower->{visited_reactions}};
1542 :     map { $borrower->{visited_compounds}->{$_} += $loop_counter } keys %{$borrower->{visited_compounds}};
1543 :     map { $borrower->{visited_reactions}->{$_} = $giver->{visited_reactions}->{$_} unless defined $borrower->{visited_reactions}->{$_} } keys %{$giver->{visited_reactions}};
1544 :     map { $borrower->{visited_compounds}->{$_} = $giver->{visited_compounds}->{$_} unless defined $borrower->{visited_compounds}->{$_} } keys %{$giver->{visited_compounds}};
1545 :     map { $borrower->{token_path_inputs}->{$_} += $giver->{token_path_inputs}->{$_} } keys %{$giver->{token_path_inputs}};
1546 :     $borrower->{initial_pass} |= $giver->{initial_pass};
1547 :    
1548 :     foreach my $icpd (keys %compounds_to_tokens)
1549 :     {
1550 :     if ($compounds_to_tokens{$icpd}->{$giver_id} > 0)
1551 :     {
1552 :     $compounds_to_tokens{$icpd}->{$borrower_id} += $compounds_to_tokens{$icpd}->{$giver_id};
1553 :     delete $compounds_to_tokens{$icpd}->{$giver_id};
1554 :     }
1555 :     }
1556 :    
1557 :     foreach my $icpd (keys %compounds_borrowed_to_tokens)
1558 :     {
1559 :     if ($compounds_borrowed_to_tokens{$icpd}->{$giver_id} > 0)
1560 :     {
1561 :     $compounds_borrowed_to_tokens{$icpd}->{$borrower_id} += $compounds_borrowed_to_tokens{$icpd}->{$giver_id};
1562 :     delete $compounds_borrowed_to_tokens{$icpd}->{$giver_id};
1563 :     }
1564 :     }
1565 :    
1566 :     delete $tokens{$giver_id};
1567 :     delete $not_done_tokens{$giver_id};
1568 :    
1569 :     &balance_borrowing_and_giving($borrower_id, \%compounds_to_tokens,
1570 :     \%compounds_borrowed_to_tokens);
1571 :    
1572 :     &print_token_status([$borrower_id], \%tokens, \%compounds_to_tokens,
1573 :     \%compounds_borrowed_to_tokens, $fig);
1574 :    
1575 :     if (&check_token_for_done($borrower_id, \%compounds_to_tokens,
1576 :     \%compounds_borrowed_to_tokens,
1577 :     \%all_compounds_to_main, \%path_outputs,
1578 :     \%scenario_cycles,\%tokens, $fig, $outputs_lists))
1579 :     {
1580 :     delete $not_done_tokens{$borrower_id};
1581 :     }
1582 :     }
1583 :     }
1584 :     }
1585 :     }
1586 :    
1587 :     # Now check if we've reached a dead end, either a compound we can't push or a
1588 :     # borrowed compound we can't repay. Also determine whether there is a reaction
1589 :     # to run that can move an initial pass token forward.
1590 :    
1591 :     print STDERR "\nChecking for dead ends\n" unless !$debug;
1592 :    
1593 :     my $found_reaction_for_initial_pass_token = 0;
1594 :    
1595 :     check: foreach my $token_id (keys %not_done_tokens)
1596 :     {
1597 :     # determine which compounds the token is sitting on, and whether a
1598 :     # reaction can proceed from those compounds that isn't a loop
1599 :     # back to compounds already visited
1600 :     my %visited_compounds = %{$tokens{$token_id}->{visited_compounds}};
1601 :     my %visited_reactions = %{$tokens{$token_id}->{visited_reactions}};
1602 :     my $dead_end_cpd;
1603 :    
1604 :     print STDERR "\tChecking if '$token_id' can run\n" unless !$debug;
1605 :    
1606 :     my $found_reaction_for_token = 0;
1607 :    
1608 :     substrate: foreach my $cpd (keys %compounds_to_tokens)
1609 :     {
1610 :     if ($compounds_to_tokens{$cpd}->{$token_id} > 0 &&
1611 :     ! defined $path_outputs{$cpd})
1612 :     {
1613 :     print STDERR "\t\tChecking substrate $cpd (main: $all_compounds_to_main{$cpd})\n" unless !$debug;
1614 :    
1615 :     foreach my $reaction (@{$substrates_to_reactions->{$cpd}})
1616 :     {
1617 :     next if ($reaction =~ /(.*)_R/ && defined $visited_reactions{$1."_L"}) ||
1618 :     ($reaction =~ /(.*)_L/ && defined $visited_reactions{$1."_R"});
1619 :    
1620 :     print STDERR "\t\t\tChecking reaction $reaction\n" unless !$debug;
1621 :    
1622 :     my $prods_are_ok = 0;
1623 :    
1624 :     foreach my $prod (@{$reactions_to_products->{$reaction}})
1625 :     {
1626 :     print STDERR "\t\t\t\tChecking product $prod\n" unless !$debug;
1627 :    
1628 :     if (! defined $visited_compounds{$prod} ||
1629 :     $visited_compounds{$prod} >= $visited_compounds{$cpd} ||
1630 :     defined $path_outputs{$prod} ||
1631 :     $compounds_borrowed_to_tokens{$prod}->{$token_id} > 0)
1632 :     {
1633 :     $prods_are_ok = $prod;
1634 :     last;
1635 :     }
1636 :     }
1637 :    
1638 :     if ($prods_are_ok)
1639 :     {
1640 :     print STDERR "\tToken '$token_id' can run $reaction on $cpd to produce $prods_are_ok\n" unless !$debug;
1641 :     $found_reaction_for_token = 1;
1642 :    
1643 :     if ($tokens{$token_id}->{initial_pass})
1644 :     {
1645 :     $found_reaction_for_initial_pass_token = 1;
1646 :     }
1647 :    
1648 :     next substrate;
1649 :     }
1650 :     }
1651 :    
1652 :     # didn't find a reaction for this substrate
1653 :     $dead_end_cpd = $cpd if $all_compounds_to_main{$cpd} &&
1654 :     ! defined $path_outputs{$cpd};
1655 :     }
1656 :     }
1657 :    
1658 :     if ($dead_end_cpd)
1659 :     {
1660 :     print STDERR "\tToken '$token_id' has reached a dead end on $dead_end_cpd\n" unless !$debug;
1661 :     $tokens{$token_id}->{done} = "dead end on $dead_end_cpd";
1662 :     }
1663 :     elsif (! $found_reaction_for_token)
1664 :     {
1665 :     # didn't find any reaction to run.
1666 :     # check to see if there are borrowed compounds to repay
1667 :     product: foreach my $cpd (keys %compounds_borrowed_to_tokens)
1668 :     {
1669 :     if ($compounds_borrowed_to_tokens{$cpd}->{$token_id} > 0)
1670 :     {
1671 :     print STDERR "\t\tChecking product $cpd\n" unless !$debug;
1672 :    
1673 :     rxn: foreach my $reaction (@{$products_to_reactions->{$cpd}})
1674 :     {
1675 :     next if ($reaction =~ /(.*)_R/ && defined $visited_reactions{$1."_L"})
1676 :     || ($reaction =~ /(.*)_L/ && defined $visited_reactions{$1."_R"});
1677 :    
1678 :     print STDERR "\t\t\tChecking reaction $reaction\n" unless !$debug;
1679 :    
1680 :     my $substrates_are_ok = 0;
1681 :    
1682 :     foreach my $sub (@{$reactions_to_substrates->{$reaction}})
1683 :     {
1684 :     print STDERR "\t\t\t\tChecking substrate $sub\n" unless !$debug;
1685 :    
1686 :     if (defined $path_outputs{$sub})
1687 :     {
1688 :     # don't run reactions that use up outputs
1689 :     next rxn;
1690 :     }
1691 :    
1692 :     if (! defined $visited_compounds{$sub} ||
1693 :     $visited_compounds{$sub} <= $visited_compounds{$cpd} ||
1694 :     defined $scenario_cycles{$sub})
1695 :     {
1696 :     $substrates_are_ok = $sub;
1697 :     }
1698 :     }
1699 :    
1700 :     if ($substrates_are_ok)
1701 :     {
1702 :     print STDERR "\tToken '$token_id' can wait for $reaction on $substrates_are_ok to produce borrowed compound $cpd\n" unless !$debug;
1703 :     $found_reaction_for_token = 1;
1704 :    
1705 :     if ($tokens{$token_id}->{initial_pass})
1706 :     {
1707 :     $found_reaction_for_initial_pass_token = 1;
1708 :     }
1709 :    
1710 :     last product;
1711 :     }
1712 :     }
1713 :    
1714 :     # didn't find a reaction for this product
1715 :     $dead_end_cpd = $cpd if $all_compounds_to_main{$cpd} &&
1716 :     ! defined $path_outputs{$cpd};
1717 :     }
1718 :     }
1719 :    
1720 :     # didn't find any reaction to run.
1721 :     if (! $found_reaction_for_token)
1722 :     {
1723 :     if ($dead_end_cpd)
1724 :     {
1725 :     print STDERR "\tToken '$token_id' has reached a dead end on borrowed compound $dead_end_cpd\n" unless !$debug;
1726 :     $tokens{$token_id}->{done} = "dead end on borrowed compound $dead_end_cpd";
1727 :     }
1728 :     else
1729 :     {
1730 :     # nothing to push, borrow or do
1731 :     print STDERR "\tToken '$token_id' has reached a dead end\n" unless !$debug;
1732 :     $tokens{$token_id}->{done} = "dead end";
1733 :     }
1734 :     }
1735 :     }
1736 :     }
1737 :    
1738 :     if($find_first)
1739 :     {
1740 :     foreach my $token_id (sort { $tokens{$a}->{done} <=> $tokens{$b}->{done} } keys %tokens)
1741 :     {
1742 :     if($tokens{$token_id}->{done} == 1)
1743 :     {
1744 :     $done = 1;
1745 :     }
1746 :     }
1747 :     }
1748 :     # is there an initial pass token that can make progress?
1749 :     if ($found_reaction_for_initial_pass_token)
1750 :     {
1751 :     if (scalar keys %reactions_to_tokens_to_use == 0)
1752 :     {
1753 :     # Couldn't run any reactions this time around.
1754 :     # Push more tokens through from the beginning of the path to
1755 :     # supply more substrates.
1756 :     $add_path_inputs = 1;
1757 :     print STDERR "\nSupplying more path inputs to push stalled tokens\n" unless !$debug;
1758 :     }
1759 :     }
1760 :     else
1761 :     {
1762 :     $done = 1;
1763 :     }
1764 :    
1765 :     $loop_counter++;
1766 :    
1767 :     if ($loop_counter >= $infinite_loop_check)
1768 :     {
1769 :     $data_results{"infinite"} = 1;
1770 :     print STDERR "Encountered an infinite loop\n" unless !$debug;
1771 :     $done = 1;
1772 :     }
1773 :     }
1774 :    
1775 :     # reverse %compounds_to_tokens, since all tokens should be at path outputs now
1776 :     my %tokens_to_compounds;
1777 :    
1778 :     foreach my $cpd (keys %compounds_to_tokens)
1779 :     {
1780 :     foreach my $token_id (keys %{$compounds_to_tokens{$cpd}})
1781 :     {
1782 :     my $num_tokens = $compounds_to_tokens{$cpd}->{$token_id};
1783 :     $tokens_to_compounds{$token_id}->{$cpd} = $num_tokens if $num_tokens > 0;
1784 :     }
1785 :     }
1786 :    
1787 :     print STDERR "\n\ntoken ids: @{[ sort { $a <=> $b } map { $_ if ! defined $tokens{$_}->{done} } keys %tokens ]}\n" unless !$debug;
1788 :    
1789 :     my $path_counter = 1;
1790 :    
1791 :     foreach my $token_id (sort { $tokens{$a}->{done} <=> $tokens{$b}->{done} } keys %tokens)
1792 :     {
1793 :     my $token = $tokens{$token_id};
1794 :     my %visited_reactions = %{$token->{visited_reactions}};
1795 :     my @path = sort { $visited_reactions{$a} <=> $visited_reactions{$b} }
1796 :     keys %visited_reactions;
1797 :     my %visited_compounds = %{$token->{visited_compounds}};
1798 :     my @compounds = sort { $visited_compounds{$a} <=> $visited_compounds{$b} }
1799 :     keys %visited_compounds;
1800 :    
1801 :     print STDERR "Adding token id: $token_id\n" unless !$debug;
1802 :    
1803 :    
1804 :     #each key in data_results is a token which points to an array
1805 :     # [0]=initial pass [1]=0/1 if its done [2]=reaction path [3]=compounds
1806 :     # [4]= html string of inputs
1807 :     # [5]= html string of outputs
1808 :     # [6]= html string of borrowed compounds
1809 :     # [7] = array of path input compounds
1810 :     # [8] = array of path output compounds
1811 :    
1812 :    
1813 :     $data_results{$token_id} = [$token->{initial_pass},$token->{done},\@path,\@compounds,[],[],[],$token->{token_path_inputs},$tokens_to_compounds{$token_id}];
1814 :    
1815 :     foreach my $input (keys %{$token->{token_path_inputs}})
1816 :     {
1817 :     my $input_stoich = $token->{token_path_inputs}->{$input};
1818 :     my $output_stoich = $tokens_to_compounds{$token_id}->{$input};
1819 :    
1820 :     # don't balance scenario cycled compounds until final assembly
1821 :     if ($scenario_cycles{$input} && ! $create_assembly)
1822 :     {
1823 :     my @names = $fig->names_of_compound($input);
1824 :     push @{$data_results{$token_id}->[4]}, "\t\t$input_stoich\t$input $names[0]\n";
1825 :     next;
1826 :     }
1827 :    
1828 :     if ($input_stoich > $output_stoich)
1829 :     {
1830 :     delete $tokens_to_compounds{$token_id}->{$input};
1831 :     $input_stoich -= $output_stoich;
1832 :     $token->{token_path_inputs}->{$input} -= $output_stoich;
1833 :     my @names = $fig->names_of_compound($input);
1834 :     push @{$data_results{$token_id}->[4]}, "\t\t$input_stoich\t$input $names[0]\n";
1835 :     }
1836 :     elsif ($output_stoich > $input_stoich)
1837 :     {
1838 :     delete $token->{token_path_inputs}->{$input};
1839 :     $tokens_to_compounds{$token_id}->{$input} -= $input_stoich;
1840 :     }
1841 :     else
1842 :     {
1843 :     delete $token->{token_path_inputs}->{$input};
1844 :     delete $tokens_to_compounds{$token_id}->{$input};
1845 :     }
1846 :     }
1847 :    
1848 :     foreach my $output (keys %{$tokens_to_compounds{$token_id}})
1849 :     {
1850 :     my @names = $fig->names_of_compound($output);
1851 :     push @{$data_results{$token_id}->[5]},"\t\t$tokens_to_compounds{$token_id}->{$output}\t$output $names[0]\n";
1852 :     }
1853 :    
1854 :     if ($token->{done} != 1)
1855 :     {
1856 :     foreach my $cpd (sort keys %compounds_borrowed_to_tokens)
1857 :     {
1858 :     my $num = $compounds_borrowed_to_tokens{$cpd}->{$token_id};
1859 :     my @names = $fig->names_of_compound($cpd);
1860 :     push @{$data_results{$token_id}->[6]},"\t\t$num $cpd\t$names[0]\n" if ($num > 0);
1861 :     }
1862 :     }
1863 :     }
1864 :    
1865 :    
1866 :    
1867 :    
1868 :     return \%data_results;
1869 :     }
1870 :    
1871 :     sub balance_borrowing_and_giving
1872 :     {
1873 :     my ($token_id, $compounds_to_tokens, $compounds_borrowed_to_tokens) = @_;
1874 :    
1875 :     my %merged_compounds;
1876 :     map { $merged_compounds{$_} = 1 } keys %{$compounds_to_tokens};
1877 :     map { $merged_compounds{$_} = 1 } keys %{$compounds_borrowed_to_tokens};
1878 :    
1879 :     foreach my $icpd (keys %merged_compounds)
1880 :     {
1881 :     my $inum_to_give = $compounds_to_tokens->{$icpd}->{$token_id};
1882 :    
1883 :     if ($inum_to_give > 0)
1884 :     {
1885 :     my $inum_needed = $compounds_borrowed_to_tokens->{$icpd}->{$token_id};
1886 :    
1887 :     if ($inum_to_give == $inum_needed)
1888 :     {
1889 :     delete $compounds_borrowed_to_tokens->{$icpd}->{$token_id};
1890 :     delete $compounds_to_tokens->{$icpd}->{$token_id};
1891 :     }
1892 :     elsif ($inum_to_give > $inum_needed)
1893 :     {
1894 :     delete $compounds_borrowed_to_tokens->{$icpd}->{$token_id};
1895 :     $compounds_to_tokens->{$icpd}->{$token_id} -= $inum_needed;
1896 :     }
1897 :     else
1898 :     {
1899 :     $compounds_borrowed_to_tokens->{$icpd}->{$token_id} -= $inum_to_give;
1900 :     delete $compounds_to_tokens->{$icpd}->{$token_id};
1901 :     }
1902 :     }
1903 :     }
1904 :     }
1905 :    
1906 :     sub print_token_status
1907 :     {
1908 :     my ($token_id_list, $tokens, $compounds_to_tokens, $compounds_borrowed_to_tokens, $fig) = @_;
1909 :    
1910 :     print STDERR "\nToken status:\n" unless !$debug;
1911 :    
1912 :     foreach my $token_id (@$token_id_list)
1913 :     {
1914 :     next if defined $tokens->{$token_id}->{done};
1915 :    
1916 :     print STDERR "\n\ttoken: '$token_id', initial: $tokens->{$token_id}->{initial_pass}\n" unless !$debug;
1917 :    
1918 :     foreach my $cpd (sort keys %{$tokens->{$token_id}->{token_path_inputs}})
1919 :     {
1920 :     my $num = $tokens->{$token_id}->{token_path_inputs}->{$cpd};
1921 :     my @names = $fig->names_of_compound($cpd);
1922 :     print STDERR "\t\tInput: $num $cpd\t$names[0]\n" unless !$debug;
1923 :     }
1924 :    
1925 :     foreach my $cpd (sort keys %$compounds_to_tokens)
1926 :     {
1927 :     my $num = $compounds_to_tokens->{$cpd}->{$token_id};
1928 :     my @names = $fig->names_of_compound($cpd);
1929 :     print STDERR "\t\tStatus: $num $cpd\t$names[0]\n" if ($num > 0 && $debug);
1930 :     }
1931 :    
1932 :     foreach my $cpd (sort keys %$compounds_borrowed_to_tokens)
1933 :     {
1934 :     my $num = $compounds_borrowed_to_tokens->{$cpd}->{$token_id};
1935 :     my @names = $fig->names_of_compound($cpd);
1936 :     print STDERR "\t\tBorrowed: $num $cpd\t$names[0]\n" if ($num > 0 && $debug);
1937 :     }
1938 :    
1939 :     my %visited_compounds = %{$tokens->{$token_id}->{visited_compounds}};
1940 :     print STDERR "\t\tvisited_compounds: @{[map { ($_, $visited_compounds{$_} ) } sort { $visited_compounds{$a} <=> $visited_compounds{$b} } keys %visited_compounds]}\n" unless !$debug;
1941 :    
1942 :     my %visited_reactions = %{$tokens->{$token_id}->{visited_reactions}};
1943 :     print STDERR "\t\tvisited_reactions: @{[map { ($_, $visited_reactions{$_} ) } sort { $visited_reactions{$a} <=> $visited_reactions{$b} } keys %visited_reactions]}\n" unless !$debug;
1944 :    
1945 :     }
1946 :    
1947 :     print STDERR "\n" unless !$debug;
1948 :     }
1949 :    
1950 :     sub clone_token
1951 :     {
1952 :     my ($clone_id, $new_token_id, $tokens, $compounds_to_tokens, $compounds_borrowed_to_tokens) = @_;
1953 :     my (%new_token, %new_visited_reactions, %new_visited_compounds, %new_token_path_inputs);
1954 :    
1955 :     $tokens->{$new_token_id} = \%new_token;
1956 :    
1957 :     my $clone_token = $tokens->{$clone_id};
1958 :     map {$new_visited_reactions{$_} = $clone_token->{visited_reactions}->{$_}} keys %{$clone_token->{visited_reactions}};
1959 :     map {$new_visited_compounds{$_} = $clone_token->{visited_compounds}->{$_}} keys %{$clone_token->{visited_compounds}};
1960 :     map {$new_token_path_inputs{$_} = $clone_token->{token_path_inputs}->{$_}} keys %{$clone_token->{token_path_inputs}};
1961 :    
1962 :     $new_token{visited_reactions} = \%new_visited_reactions;
1963 :     $new_token{visited_compounds} = \%new_visited_compounds;
1964 :     $new_token{token_path_inputs} = \%new_token_path_inputs;
1965 :     $new_token{initial_pass} = $clone_token->{initial_pass};
1966 :    
1967 :     # tokens might be spread across multiple compounds
1968 :     foreach my $icpd (keys %$compounds_to_tokens)
1969 :     {
1970 :     if ($compounds_to_tokens->{$icpd}->{$clone_id} > 0)
1971 :     {
1972 :     $compounds_to_tokens->{$icpd}->{$new_token_id} =
1973 :     $compounds_to_tokens->{$icpd}->{$clone_id};
1974 :     }
1975 :     }
1976 :    
1977 :     # tokens might have borrowed compounds
1978 :     foreach my $icpd (keys %$compounds_borrowed_to_tokens)
1979 :     {
1980 :     if ($compounds_borrowed_to_tokens->{$icpd}->{$clone_id} > 0)
1981 :     {
1982 :     $compounds_borrowed_to_tokens->{$icpd}->{$new_token_id} =
1983 :     $compounds_borrowed_to_tokens->{$icpd}->{$clone_id};
1984 :     }
1985 :     }
1986 :    
1987 :     return \%new_token;
1988 :     }
1989 :    
1990 :     sub check_token_for_done
1991 :     {
1992 :     my ($token_id, $compounds_to_tokens, $compounds_borrowed_to_tokens, $all_compounds_to_main,
1993 :     $path_outputs, $scenario_cycles, $tokens, $fig, $outputs_lists) = @_;
1994 :    
1995 :     my $token_is_done_pushing = 1;
1996 :     my $token_is_done_borrowing = 1;
1997 :    
1998 :     # first determine if there is a main compound that isn't a path output
1999 :     foreach my $cpd (keys %$compounds_to_tokens)
2000 :     {
2001 :     if ($compounds_to_tokens->{$cpd}->{$token_id} > 0)
2002 :     {
2003 :     # also check if scenario cycle compounds need to be pushed
2004 :     if ($all_compounds_to_main->{$cpd} &&
2005 :     (! defined $path_outputs->{$cpd} ||
2006 :     ($scenario_cycles->{$cpd} &&
2007 :     $tokens->{$token_id}->{visited_compounds}->{$cpd} == 0)))
2008 :     {
2009 :     my @names = $fig->names_of_compound($cpd);
2010 :     print STDERR "\ttoken '$token_id' needs to push $cpd $names[0]\n" unless !$debug;
2011 :     $token_is_done_pushing = 0;
2012 :     }
2013 :     }
2014 :     }
2015 :    
2016 :     # now determine if one of the output lists has been satisfied
2017 :     if ($token_is_done_pushing)
2018 :     {
2019 :     my $found_a_list = 0;
2020 :    
2021 :     foreach my $cpd_list (@$outputs_lists)
2022 :     {
2023 :     print STDERR "\t\tchecking outputs_list @$cpd_list\n" unless !$debug;
2024 :    
2025 :     $found_a_list = 1;
2026 :    
2027 :     foreach my $cpd (@$cpd_list)
2028 :     {
2029 :     if (! $compounds_to_tokens->{$cpd}->{$token_id} > 0)
2030 :     {
2031 :     $found_a_list = 0;
2032 :     last;
2033 :     }
2034 :     }
2035 :    
2036 :     last if $found_a_list;
2037 :     }
2038 :    
2039 :     if (! $found_a_list)
2040 :     {
2041 :     print STDERR "\ttoken '$token_id' hasn't satisfied output compound list\n" unless !$debug;
2042 :     $token_is_done_pushing = 0;
2043 :     }
2044 :     }
2045 :    
2046 :     if ($token_is_done_pushing)
2047 :     {
2048 :     foreach my $cpd (keys %$compounds_borrowed_to_tokens)
2049 :     {
2050 :     if ($compounds_borrowed_to_tokens->{$cpd}->{$token_id} > 0)
2051 :     {
2052 :     if ($all_compounds_to_main{$cpd})
2053 :     {
2054 :     my @names = $fig->names_of_compound($cpd);
2055 :     print STDERR "\ttoken '$token_id' has borrowed $cpd $names[0]\n" unless !$debug;
2056 :     $token_is_done_borrowing = 0;
2057 :     }
2058 :     }
2059 :     }
2060 :     }
2061 :    
2062 :     if (! $token_is_done_pushing || ! $token_is_done_borrowing)
2063 :     {
2064 :     return 0;
2065 :     }
2066 :     else
2067 :     {
2068 :     $tokens->{$token_id}->{done} = 1;
2069 :     print STDERR "\tToken '$token_id' is done\n" unless !$debug;
2070 :     return 1;
2071 :     }
2072 :     }
2073 :    
2074 :    
2075 :     sub write_fluxanalyzer_files
2076 :     {
2077 :     my ($dir, $path_inputs, $path_outputs, $path_array,
2078 :     $all_reactions,$reactions_to_substrate_arrays,$reactions_to_product_arrays,
2079 :     $cidToName) = @_;
2080 :    
2081 :     my $x_pos = 10;
2082 :     my $y_pos = 20;
2083 :    
2084 :     #Write the inputs/outputs to a seperate file, along with $stoich and if its main
2085 :     open(A_INPUT, ">$dir/inputs_main");
2086 :    
2087 :     print A_INPUT map {"$_\t$all_compounds_to_main{$_}\n"} keys %$path_inputs;
2088 :    
2089 :     close(A_INPUT);
2090 :    
2091 :     open(A_OUTPUT, ">$dir/outputs_main");
2092 :    
2093 :     print A_OUTPUT map {"$_\t$all_compounds_to_main{$_}\n"} keys %$path_outputs;
2094 :    
2095 :     close(A_OUTPUT);
2096 :    
2097 :    
2098 :    
2099 :     open(REACTIONS, ">$dir/reactions");
2100 :     open(INPUTS, ">$dir/inputs");
2101 :     open(OUTPUTS, ">$dir/outputs");
2102 :     open(PATH, ">$dir/path_info");
2103 :    
2104 :     my @inputs = keys %$path_inputs;
2105 :    
2106 :     foreach my $elem(@{$path_array}){
2107 :     print PATH $elem ."\n";
2108 :     }
2109 :     close(PATH);
2110 :    
2111 :     foreach my $cpd (@inputs)
2112 :     {
2113 :     my $toPrint = $cpd."up\t = 1 $cpd \t| \t";
2114 :     $toPrint .= $path_inputs->{$cpd};
2115 :     $toPrint .= " \t0 100 0 \t$x_pos $y_pos 1 1\t0.01\n";
2116 :     print REACTIONS $toPrint;
2117 :     $y_pos += 20;
2118 :    
2119 :     if ($y_pos == 300)
2120 :     {
2121 :     $x_pos += 60;
2122 :     $y_pos = 20;
2123 :     }
2124 :     }
2125 :    
2126 :     $x_pos += 60;
2127 :     $y_pos = 20;
2128 :    
2129 :     my (@display_array);
2130 :     foreach my $rxn (keys %$all_reactions)
2131 :     {
2132 :     my $direction = $all_reactions->{$rxn};
2133 :    
2134 :     my (@substrate_array, @product_array);
2135 :    
2136 :     push @display_array, $rxn;
2137 :    
2138 :     if ($direction eq "L")
2139 :     {
2140 :     @product_array = @{$reactions_to_substrate_arrays->{$rxn."_L"}};
2141 :     @substrate_array = @{$reactions_to_product_arrays->{$rxn."_L"}};
2142 :     }
2143 :     else
2144 :     {
2145 :     @substrate_array = @{$reactions_to_substrate_arrays->{$rxn."_R"}};
2146 :     @product_array = @{$reactions_to_product_arrays->{$rxn."_R"}};
2147 :     }
2148 :    
2149 :     foreach my $subTuple (@substrate_array)
2150 :     {
2151 :     my @temp = $fig->names_of_compound($subTuple->[0]);
2152 :     $cidToName->{$subTuple->[0]} = $temp [0] if ! defined $cidToName->{$subTuple->[0]};
2153 :     }
2154 :    
2155 :     foreach my $prodTuple(@product_array)
2156 :     {
2157 :     my @temp = $fig->names_of_compound($prodTuple->[0]);
2158 :     $cidToName->{$prodTuple->[0]} = $temp[0] if ! defined $cidToName->{$prodTuple->[0]};
2159 :     }
2160 :    
2161 :     if($direction eq "R" || $direction eq "B")
2162 :     {
2163 :     #write data in a strign for copying to file later
2164 :     my $toFile = '';
2165 :     $toFile.= $rxn."\t";
2166 :    
2167 :     #add all the substrates
2168 :     foreach my $curSub(@substrate_array){
2169 :     $toFile .= $curSub -> [1].' '. $curSub -> [0].' + ';
2170 :     }
2171 :    
2172 :     ##chop off the +
2173 :     chop($toFile);
2174 :     chop($toFile);
2175 :    
2176 :     $toFile.='= ';
2177 :    
2178 :     #add all the products
2179 :     foreach my $curProd(@product_array){
2180 :     $toFile .= $curProd -> [1].' '. $curProd -> [0].' + ';
2181 :     }
2182 :    
2183 :     #chop off the plus
2184 :     chop($toFile);
2185 :     chop($toFile);
2186 :    
2187 :     $toFile.="\t|\t#\t";
2188 :    
2189 :     if($direction eq "B"){
2190 :     $toFile.="-Inf";
2191 :     }
2192 :     else{
2193 :     $toFile.="0";
2194 :     }
2195 :    
2196 :     $toFile.=" Inf 0\t$x_pos $y_pos 1 1 \t0.01\n";
2197 :    
2198 :     print REACTIONS $toFile;
2199 :     }
2200 :     elsif($direction eq "L")
2201 :     {
2202 :     #write data in a strign for copying to file later
2203 :     my $toFile = '';
2204 :     $toFile.= $rxn."\t";
2205 :    
2206 :     #add all the substrates
2207 :     foreach my $curProd(@product_array){
2208 :     $toFile .= $curProd -> [1].' '. $curProd -> [0]." + ";
2209 :     }
2210 :    
2211 :     ##chop off the +
2212 :     chop($toFile);
2213 :     chop($toFile);
2214 :    
2215 :     $toFile.="= ";
2216 :    
2217 :     #add all the products
2218 :     foreach my $curSubstrate(@substrate_array){
2219 :     $toFile .= $curSubstrate -> [1].' '. $curSubstrate -> [0]." + ";
2220 :     }
2221 :    
2222 :     #chop off the plus
2223 :     chop($toFile);
2224 :     chop($toFile);
2225 :    
2226 :     $toFile.="\t|\t#\t";
2227 :    
2228 :     $toFile.="0";
2229 :    
2230 :     $toFile.=" Inf 0\t$x_pos $y_pos 1 1 \t0.01\n";
2231 :    
2232 :     print REACTIONS $toFile;
2233 :     }
2234 :    
2235 :     $y_pos += 20;
2236 :    
2237 :     if ($y_pos == 300)
2238 :     {
2239 :     $x_pos += 60;
2240 :     $y_pos = 20;
2241 :     }
2242 :     }
2243 :    
2244 :     my @outputs = keys %$path_outputs;
2245 :    
2246 :     $x_pos += 60;
2247 :     $y_pos = 20;
2248 :    
2249 :     foreach my $cpd (@inputs)
2250 :     {
2251 :     print INPUTS $cpd, "\t", $path_inputs->{$cpd}, "\t", $cidToName->{$cpd}, "\n";
2252 :     }
2253 :    
2254 :     foreach my $cpd (@outputs)
2255 :     {
2256 :     print OUTPUTS $cpd, "\t", $path_outputs->{$cpd}, "\t", $cidToName->{$cpd}, "\n";
2257 :     my $toPrint = $cpd."ex\t 1 $cpd = \t| \t# \t0 100 0 \t$x_pos $y_pos 1 1\t0.01\n";
2258 :     print REACTIONS $toPrint;
2259 :     $y_pos += 20;
2260 :    
2261 :     if ($y_pos == 300)
2262 :     {
2263 :     $x_pos += 60;
2264 :     $y_pos = 20;
2265 :     }
2266 :     }
2267 :    
2268 :     $x_pos += 60;
2269 :     $y_pos = 20;
2270 :    
2271 :     #print the macromolecule_synthesis and assembly file
2272 :     open(MACRO_SYTH,">$dir/macromolecule_synthesis");
2273 :     open(ASSEM,">$dir/assembly");
2274 :     my $toPrint = "M1 = ";
2275 :    
2276 :     foreach my $cpd (keys %$path_outputs)
2277 :     {
2278 :     $toPrint.="$path_outputs->{$cpd} $cpd + ";
2279 :     print ASSEM "$cpd\tM1\t-100 -100 1\n";
2280 :     $y_pos += 25;
2281 :     }
2282 :    
2283 :     chop $toPrint;
2284 :     chop $toPrint;
2285 :     chop $toPrint;
2286 :     print MACRO_SYTH $toPrint;
2287 :     close(MACRO_SYTH);
2288 :     close(ASSEM);
2289 :    
2290 :     #Print the metabolites for these subsystems.
2291 :     open(METABOLITES,">$dir/metabolites");
2292 :    
2293 :     foreach my $cid (keys %$cidToName)
2294 :     {
2295 :     my $name = $cidToName->{$cid};
2296 :     $name =~ s/\s/-/g;
2297 :     print METABOLITES $cid."\t".$name."\t0.001\t0\n";
2298 :     }
2299 :     close(METABOLITES);
2300 :    
2301 :     #Print the macromolucules file
2302 :     open(MACRO,">$dir/macromolecules");
2303 :     print MACRO "M1 \tM1 \t1 \t-100 -100 1 1\n";
2304 :     close(MACRO);
2305 :    
2306 :     $x_pos += 60;
2307 :     $y_pos = 20;
2308 :    
2309 :     print REACTIONS "mue\t\t\t|\t#\t0 100 0\t$x_pos $y_pos 1 1\t0.01\n";
2310 :    
2311 :     #close reaction equation file
2312 :     close(REACTIONS);
2313 :    
2314 :     # FluxAnalyzer requires this file
2315 :     open(APP, ">$dir/app_para.m");
2316 :     print APP "epsilon=1e-10;\nbasic_color=[0.7 0.7 0.7];\ncr_color=[0.5 0.5 1];\nbr_color=[1 0.2 0.2];\nnbr_color=[0.2 1 0.2];\ntext_color=[0 0 0];\nmacro_synth_color=[0 0 1];\nmacro_color=[0.6 0.6 1];\nbox_reaction_width=[0.12];\nbox_reaction_height=[0.06];\nbox_macro_width=[0.08];\nbox_macro_height=[0.06];\nfontsize_reaction=[11];\nfontsize_macro=[11];\nfluxmaps={'Fluxmap','dummy.pcx'};\n";
2317 :     close(APP);
2318 :     }
2319 :    
2320 :     sub write_final_fluxanalyzer_files
2321 :     {
2322 :     my ($dir, $path_inputs, $path_outputs, $all_reactions, $transport_reactions,
2323 :     $reactions_to_substrate_arrays,$reactions_to_product_arrays,
2324 :     $cidToName,$bioMass,$minSubstrates) = @_;
2325 :    
2326 :     #Write the inputs/outputs to a seperate file, along with $stoich and if its main
2327 :     open(A_INPUT, ">$dir/inputs_main");
2328 :    
2329 :     print A_INPUT map {"$_\t$all_compounds_to_main{$_}\n"} keys %$path_inputs;
2330 :    
2331 :     close(A_INPUT);
2332 :    
2333 :     open(A_OUTPUT, ">$dir/outputs_main");
2334 :    
2335 :     print A_OUTPUT map {"$_\t$all_compounds_to_main{$_}\n"} keys %$path_outputs;
2336 :    
2337 :     close(A_OUTPUT);
2338 :    
2339 :     my %open_transports;
2340 :    
2341 :     foreach my $cpd (keys %$minSubstrates)
2342 :     {
2343 :     map { $open_transports{$_} = 1 } @{$minSubstrates->{$cpd}};
2344 :     }
2345 :    
2346 :     open(REACTIONS, ">$dir/reactions");
2347 :     open(INPUTS, ">$dir/inputs");
2348 :     open(OUTPUTS, ">$dir/outputs");
2349 :    
2350 :     my @inputs = keys %$path_inputs;
2351 :    
2352 :     my $x_pos = 10;
2353 :     my $y_pos = 30;
2354 :    
2355 :     foreach my $cpd (@inputs)
2356 :     {
2357 :     my $toPrint = $cpd."up\t = 1 $cpd \t| \t";
2358 :    
2359 :     if (defined $minSubstrates->{$cpd})
2360 :     {
2361 :     $toPrint .= "#";
2362 :     }
2363 :     else
2364 :     {
2365 :     $toPrint .= "0";
2366 :     }
2367 :    
2368 :     $toPrint .= " \t0 Inf 0 \t$x_pos $y_pos 1 1\t0.01\n";
2369 :     print REACTIONS $toPrint;
2370 :     $y_pos += 30;
2371 :    
2372 :     if ($y_pos > 600)
2373 :     {
2374 :     $x_pos += 60;
2375 :     $y_pos = 30;
2376 :     }
2377 :     }
2378 :    
2379 :     $x_pos += 60;
2380 :     $y_pos = 30;
2381 :    
2382 :     foreach my $rxn (keys %$all_reactions)
2383 :     {
2384 :     my $direction = $all_reactions->{$rxn};
2385 :    
2386 :     my (@substrate_array, @product_array);
2387 :    
2388 :     if ($direction eq "L")
2389 :     {
2390 :     @product_array = @{$reactions_to_substrate_arrays->{$rxn."_L"}};
2391 :     @substrate_array = @{$reactions_to_product_arrays->{$rxn."_L"}};
2392 :     }
2393 :     else
2394 :     {
2395 :     @substrate_array = @{$reactions_to_substrate_arrays->{$rxn."_R"}};
2396 :     @product_array = @{$reactions_to_product_arrays->{$rxn."_R"}};
2397 :     }
2398 :    
2399 :     foreach my $subTuple (@substrate_array)
2400 :     {
2401 :     my @temp = $fig->names_of_compound($subTuple->[0]);
2402 :     $cidToName->{$subTuple->[0]} = $temp [0] if ! defined $cidToName->{$subTuple->[0]};
2403 :     }
2404 :    
2405 :     foreach my $prodTuple(@product_array)
2406 :     {
2407 :     my @temp = $fig->names_of_compound($prodTuple->[0]);
2408 :     $cidToName->{$prodTuple->[0]} = $temp[0] if ! defined $cidToName->{$prodTuple->[0]};
2409 :     }
2410 :    
2411 :     if($direction eq "R" || $direction eq "B")
2412 :     {
2413 :     #write data in a strign for copying to file later
2414 :     my $toFile = '';
2415 :     $toFile.= $rxn."\t";
2416 :    
2417 :     #add all the substrates
2418 :     foreach my $curSub(@substrate_array){
2419 :     $toFile .= $curSub -> [1].' '. $curSub -> [0].' + ';
2420 :     }
2421 :    
2422 :     ##chop off the +
2423 :     chop($toFile);
2424 :     chop($toFile);
2425 :    
2426 :     $toFile.='= ';
2427 :    
2428 :     #add all the products
2429 :     my $found_prod = 0;
2430 :     foreach my $curProd(@product_array){
2431 :     $toFile .= $curProd -> [1].' '. $curProd -> [0].' + ';
2432 :     $found_prod = 1;
2433 :     }
2434 :    
2435 :     if ($found_prod)
2436 :     {
2437 :     #chop off the plus
2438 :     chop($toFile);
2439 :     chop($toFile);
2440 :     }
2441 :    
2442 :     if (defined $transport_reactions->{$rxn})
2443 :     {
2444 :     if (defined $open_transports{$rxn})
2445 :     {
2446 :     $toFile.="\t|\t#\t";
2447 :     }
2448 :     else
2449 :     {
2450 :     $toFile.="\t|\t0\t";
2451 :     }
2452 :    
2453 :     if($direction eq "B"){
2454 :     $toFile.="-Inf";
2455 :     }
2456 :     else{
2457 :     $toFile.="0";
2458 :     }
2459 :    
2460 :     $toFile.=" Inf 0\t$x_pos $y_pos 1 1 \t0.01\n";
2461 :     $y_pos += 30;
2462 :     }
2463 :     elsif ($rxn =~ /sink/)
2464 :     {
2465 :     $toFile.="\t|\t0\t";
2466 :    
2467 :     if($direction eq "B"){
2468 :     $toFile.="-0.00001";
2469 :     }
2470 :     else{
2471 :     $toFile.="0";
2472 :     }
2473 :    
2474 :     $toFile.=" 0.00001 0\t$x_pos $y_pos 1 1 \t0.01\n";
2475 :     $y_pos += 30;
2476 :     }
2477 :     else
2478 :     {
2479 :     $toFile.="\t|\t#\t";
2480 :    
2481 :     if($direction eq "B"){
2482 :     $toFile.="-Inf";
2483 :     }
2484 :     else{
2485 :     $toFile.="0";
2486 :     }
2487 :    
2488 :     $toFile.=" Inf 0\t-10 -10 1 1 \t0.01\n";
2489 :     }
2490 :    
2491 :     print REACTIONS $toFile;
2492 :     }
2493 :     elsif($direction eq "L")
2494 :     {
2495 :     #write data in a strign for copying to file later
2496 :     my $toFile = '';
2497 :     $toFile.= $rxn."\t";
2498 :    
2499 :     #add all the substrates
2500 :     foreach my $curProd(@product_array){
2501 :     $toFile .= $curProd -> [1].' '. $curProd -> [0]." + ";
2502 :     }
2503 :    
2504 :     ##chop off the +
2505 :     chop($toFile);
2506 :     chop($toFile);
2507 :    
2508 :     $toFile.="= ";
2509 :    
2510 :     #add all the products
2511 :     foreach my $curSubstrate(@substrate_array){
2512 :     $toFile .= $curSubstrate -> [1].' '. $curSubstrate -> [0]." + ";
2513 :     }
2514 :    
2515 :     #chop off the plus
2516 :     chop($toFile);
2517 :     chop($toFile);
2518 :    
2519 :    
2520 :     $toFile.="\t|\t#\t";
2521 :     $toFile.="0";
2522 :     $toFile.=" Inf 0\t-10 -10 1 1 \t0.01\n";
2523 :    
2524 :     print REACTIONS $toFile;
2525 :     }
2526 :    
2527 :     if ($y_pos > 600)
2528 :     {
2529 :     $x_pos += 60;
2530 :     $y_pos = 30;
2531 :     }
2532 :     }
2533 :    
2534 :     my @outputs = keys %$path_outputs;
2535 :    
2536 :     $x_pos += 60;
2537 :     $y_pos = 30;
2538 :    
2539 :     foreach my $cpd (@inputs)
2540 :     {
2541 :     print INPUTS $cpd, "\t", $path_inputs->{$cpd}, "\t", $cidToName->{$cpd}, "\n";
2542 :     }
2543 :    
2544 :     foreach my $cpd (@outputs)
2545 :     {
2546 :     print OUTPUTS $cpd, "\t", $path_outputs->{$cpd}, "\t", $cidToName->{$cpd}, "\n";
2547 :     my $toPrint = $cpd."ex\t 1 $cpd = \t| \t# \t0 Inf 0 \t$x_pos $y_pos 1 1\t0.01\n";
2548 :     print REACTIONS $toPrint;
2549 :     $y_pos += 30;
2550 :    
2551 :     if ($y_pos > 600)
2552 :     {
2553 :     $x_pos += 60;
2554 :     $y_pos = 30;
2555 :     }
2556 :     }
2557 :    
2558 :     $x_pos += 60;
2559 :     $y_pos = 30;
2560 :    
2561 :     #print the macromolecule_synthesis and assembly file
2562 :     open(MACRO,">$dir/macromolecules");
2563 :     open(MACRO_SYTH,">$dir/macromolecule_synthesis");
2564 :     open(ASSEM,">$dir/assembly");
2565 :    
2566 :     my $toPrint = "M1 = ";
2567 :    
2568 :     foreach my $cpd (keys %$bioMass)
2569 :     {
2570 :     $toPrint .= "$bioMass->{$cpd} $cpd + ";
2571 :     print ASSEM "$cpd\tM1\t-100 -100 1\n";
2572 :     }
2573 :    
2574 :     chop $toPrint;
2575 :     chop $toPrint;
2576 :     chop $toPrint;
2577 :    
2578 :     print MACRO "M1 \tM1 \t1 \t-100 -100 1 1\n";
2579 :     print MACRO_SYTH $toPrint, "\n";
2580 :    
2581 :     close(MACRO);
2582 :     close(MACRO_SYTH);
2583 :     close(ASSEM);
2584 :    
2585 :     #Print the metabolites for these subsystems.
2586 :     open(METABOLITES,">$dir/metabolites");
2587 :    
2588 :     foreach my $cid (keys %$cidToName)
2589 :     {
2590 :     my $name = $cidToName->{$cid};
2591 :     $name =~ s/\s/-/g;
2592 :     print METABOLITES $cid."\t".$name."\t0.001\t0\n";
2593 :     }
2594 :     close(METABOLITES);
2595 :    
2596 :     $x_pos += 60;
2597 :     $y_pos = 30;
2598 :    
2599 :     print REACTIONS "mue\t\t\t|\t#\t0 100 0\t$x_pos $y_pos 1 1\t0.01\n";
2600 :    
2601 :     #close reaction equation file
2602 :     close(REACTIONS);
2603 :    
2604 :     # FluxAnalyzer requires this file
2605 :     open(APP, ">$dir/app_para.m");
2606 :     print APP "epsilon=1e-10;\nbasic_color=[0.7 0.7 0.7];\ncr_color=[0.5 0.5 1];\nbr_color=[1 0.2 0.2];\nnbr_color=[0.2 1 0.2];\ntext_color=[0 0 0];\nmacro_synth_color=[0 0 1];\nmacro_color=[0.6 0.6 1];\nbox_reaction_width=[0.06];\nbox_reaction_height=[0.03];\nbox_macro_width=[0.06];\nbox_macro_height=[0.03];\nfontsize_reaction=[11];\nfontsize_macro=[11];\nfluxmaps={'Fluxmap','dummy_medium.pcx'};\n";
2607 :     close(APP);
2608 :     }
2609 :    
2610 :    
2611 :     sub clear_arrays
2612 :     {
2613 :     undef %reactions_to_substrate_arrays;
2614 :     undef %reactions_to_product_arrays;
2615 :     undef %all_compounds_to_main;
2616 :     undef %all_reactions;
2617 :     undef %scenario_cycles;
2618 :     undef @all_outputs_lists;
2619 :     undef %all_inputs;
2620 :     undef %all_outputs;
2621 :    
2622 :     %reactions_to_substrate_arrays = ();
2623 :     %reactions_to_product_arrays = ();
2624 :     %all_compounds_to_main = ();
2625 :     %all_reactions = ();
2626 :     %scenario_cycles = ();
2627 :     @all_outputs_lists = ();
2628 :     %all_inputs = ();
2629 :     %all_outputs = ();
2630 :     }
2631 :    
2632 :    
2633 :     sub clear_dir_for_genome
2634 :     {
2635 :     my ($genome) = @_;
2636 :     rmtree("$FIG_Config::var/Models/$genome");
2637 :     }
2638 :    
2639 :     sub load_superset_file
2640 :     {
2641 :     my ($filename) = @_;
2642 :    
2643 :     %superset_to_ss = ();
2644 :     %ss_to_superset = ();
2645 :    
2646 :     open(FILE,"<$filename") or die("Failed to open $filename");
2647 :     while(<FILE>)
2648 :     {
2649 :     my @line = split(/\t/,$_);
2650 :     map { s/ /_/g } @line;
2651 :     map { s/\"//g } @line;
2652 :     map { chomp } @line;
2653 :     $superset_to_ss{$line[0]} = [] if !defined $superset_to_ss{$line[0]};
2654 :     $ss_to_superset{$line[1]} = $line[0];
2655 :     push(@{$superset_to_ss{$line[0]}},$line[1]);
2656 :     }
2657 :     close(FILE);
2658 :    
2659 :     return \%superset_to_ss
2660 :     }
2661 :    
2662 :     # This function runs a given scenario that is defined the specified subsystem
2663 :     # It returns the data as it is found by process_paths in a hash reference
2664 :     # This is used internally by model.pm, and shouldn't be called externally
2665 :    
2666 :     sub load_scenario
2667 :     {
2668 :     my ($genome,$ssa,$scenario) = @_;
2669 :    
2670 :     #load up the arrays with the info we need
2671 :     process_init($ssa,$scenario,$genome,0);
2672 :    
2673 :     #assume all path inputs and outputs are main
2674 :     map { $all_compounds_to_main{$_} = 1 } keys %all_inputs;
2675 :     map { $all_compounds_to_main{$_} = 1 } keys %all_outputs;
2676 :    
2677 :     }
2678 :    
2679 :     sub internal_scenario
2680 :     {
2681 :     my ($genome,$ssa,$scenario,$find_first) = @_;
2682 :    
2683 : olson 1.2 print STDERR "\nIn internal_scenario with '$genome', '$ssa', '$scenario', '$find_first'\n" if $debug;
2684 : olson 1.1
2685 :     load_scenario($genome,$ssa,$scenario);
2686 :    
2687 :     return execute_paths([],$find_first,[],[]);
2688 :     }
2689 :    
2690 :     sub run_scenario
2691 :     {
2692 :     my($genome,$superset,$subsystem,$scenario,$find_first) = @_;
2693 :     my $scenario_dir = "$FIG_Config::var/Models/$genome/$superset/$subsystem/$scenario";
2694 : olson 1.2
2695 :     system("rm", "-rf", $scenario_dir);
2696 :    
2697 : olson 1.1 mkdir "$FIG_Config::var/Models/$genome/";
2698 :     mkdir "$FIG_Config::var/Models/$genome/$superset";
2699 :     mkdir "$FIG_Config::var/Models/$genome/$superset/$subsystem";
2700 :     mkdir $scenario_dir;
2701 :    
2702 :     #make sure the arrays are empty to start out
2703 :     &clear_arrays;
2704 :     write_scenario(internal_scenario($genome,$subsystem,$scenario,$find_first),$scenario_dir);
2705 :     }
2706 :    
2707 : olson 1.2 sub compare_scenario
2708 :     {
2709 :     my($genome,$superset,$ss_name,$scenario_name) = @_;
2710 :     my @genome_paths;
2711 :     my $scenario_dir_all = "$FIG_Config::temp/All/$superset/$ss_name/$scenario_name";
2712 :     my $subsystem = $fig->get_subsystem($ss_name);
2713 :     my @additional_reactions = @{$subsystem->get_hope_additional_reactions($scenario_name)};
2714 :     my %additional_reactions;
2715 :     map { $additional_reactions{$_} = 1 } @additional_reactions;
2716 :     my @roles = $subsystem->get_roles;
2717 :     my %hope_reactions = %{$subsystem->get_hope_reactions};
2718 :    
2719 :     my $valid_genome = 0;
2720 :     #determine if this genome has this subsystem
2721 :     my @genome_array = $subsystem->get_genomes();
2722 :     foreach my $organ (@genome_array)
2723 :     {
2724 :     if($organ eq $genome || $genome eq "All" )
2725 :     {
2726 :     $valid_genome = 1;
2727 :     last;
2728 :     }
2729 :     }
2730 :    
2731 :     my %ss_reactions;
2732 :    
2733 :     foreach my $role (@roles)
2734 :     {
2735 :     if ($valid_genome && defined $hope_reactions{$role} &&
2736 :     ($genome eq "All" || $subsystem->get_pegs_from_cell($genome,$role) > 0))
2737 :     {
2738 :     map { $ss_reactions{$_} = 1 } @{$hope_reactions{$role}};
2739 :     }
2740 :     }
2741 :    
2742 :     # first find paths in the All directory that should be valid for the genome
2743 :     # based on the reactions associated with it in the subsystems
2744 :     opendir (DIR_ALL,$scenario_dir_all) or die("$scenario_dir_all");
2745 :     my @sub_dirs = readdir DIR_ALL; shift @sub_dirs; shift @sub_dirs; # get rid of . and ..
2746 :     close DIR_ALL;
2747 :    
2748 :     my %paths_all;
2749 :    
2750 :     for my $path (@sub_dirs)
2751 :     {
2752 :     my $match = 1;
2753 :     open (PATH, "$scenario_dir_all/$path/path_info");
2754 :     my @reactions = <PATH>;
2755 :     close PATH;
2756 :    
2757 :     my $reaction_string = "";
2758 :    
2759 :     foreach my $reaction (sort @reactions)
2760 :     {
2761 :     if ($reaction =~ /(R\d\d\d\d\d)/)
2762 :     {
2763 :     $reaction_string .= $reaction;
2764 :    
2765 :     if (! exists($ss_reactions{$1}) && ! exists($additional_reactions{$1}))
2766 :     {
2767 :     $match = 0;
2768 :     }
2769 :     }
2770 :     }
2771 :    
2772 :     $paths_all{$reaction_string} = $match;
2773 :     push @genome_paths, $path if $match;
2774 :     }
2775 :    
2776 :     # now check all paths found for this particular organism, and make sure they
2777 :     # are in the appropriate All subdirectory
2778 :     my $scenario_dir_genome = "$FIG_Config::temp/$genome/$superset/$ss_name/$scenario_name";
2779 :    
2780 :     opendir (DIR_GENOME,$scenario_dir_genome) or die("$scenario_dir_genome");
2781 :     @sub_dirs = readdir DIR_GENOME; shift @sub_dirs; shift @sub_dirs; # get rid of . and ..
2782 :     close DIR_GENOME;
2783 :    
2784 :     my $path_counter = scalar keys %paths_all;
2785 :    
2786 :     for my $path (@sub_dirs)
2787 :     {
2788 :     my $match = 1;
2789 :     open (PATH, "$scenario_dir_genome/$path/path_info");
2790 :     my @reactions = <PATH>;
2791 :     close PATH;
2792 :    
2793 :     my $reaction_string = "";
2794 :    
2795 :     foreach my $reaction (sort @reactions)
2796 :     {
2797 :     if ($reaction =~ /(R\d\d\d\d\d)/)
2798 :     {
2799 :     $reaction_string .= $reaction;
2800 :     }
2801 :     }
2802 :    
2803 :     if (! exists($paths_all{$reaction_string}))
2804 :     {
2805 :     $path_counter++;
2806 :     my $new_path_name = "path_".$path_counter;
2807 :     my $temp_sdg = $scenario_dir_genome;
2808 :     $temp_sdg =~ s/\(/\\\(/g;
2809 :     $temp_sdg =~ s/\)/\\\)/g;
2810 :     my $temp_sda = $scenario_dir_all;
2811 :     $temp_sda =~ s/\(/\\\(/g;
2812 :     $temp_sda =~ s/\)/\\\)/g;
2813 :     `cp -R $temp_sdg/$path $temp_sda/$new_path_name`;
2814 :     push @genome_paths, $new_path_name;
2815 :     print "Copied $temp_sdg/$path to $temp_sda/$new_path_name\n";
2816 :     }
2817 :     }
2818 :    
2819 :    
2820 :     return \@genome_paths;
2821 :     }
2822 :    
2823 : olson 1.1 sub write_scenario
2824 :     {
2825 :     my($scenario_data,$scenario_dir) = @_;
2826 :     delete $scenario_data->{"infinite"};
2827 :     my $path_count = 1;
2828 :     my @list_of_done_tokens=();
2829 :    
2830 :     print STDERR "Paths: ", keys %{$scenario_data}, "\n" if $debug;
2831 :    
2832 :     foreach my $try_path (keys %{$scenario_data})
2833 :     {
2834 :     print STDERR "\t Checking $try_path\n" if $debug;
2835 :     if($scenario_data->{$try_path}->[1] != 1)
2836 :     {
2837 :     print STDERR "\t Token $try_path is not complete\n" if $debug;
2838 :     next;
2839 :     }
2840 :     print STDERR "These are the contents of the list: " , @list_of_done_tokens, "\n" if $debug;
2841 :     if(scalar @list_of_done_tokens == 0){
2842 :     push @list_of_done_tokens, $try_path;
2843 :     next;
2844 :     }
2845 :    
2846 :     # Check this token's ($try_path) values against the values of the keys($elem) stored in
2847 :     # list_of_done_tokens. If they match, don't add it to the finished token list, if it doesn't
2848 :     # match, add it.
2849 :    
2850 :     print STDERR "These are the contents of the list: " , @list_of_done_tokens, "\n" if $debug;
2851 :    
2852 :     my $found_match = 0;
2853 :    
2854 :     foreach my $elem(@list_of_done_tokens){
2855 :     my @done_reactions = @{$scenario_data->{$elem}->[2]};
2856 :    
2857 :     print STDERR "This is the path we're trying: " , $try_path, "\t","This is the path already in the array: " , $elem , "\nThis is the size of the array: " . @list_of_done_tokens ."\n" if $debug;
2858 :    
2859 :     # if the list of reactions match, they represent the same path
2860 :    
2861 :     my @path_reactions = @{$scenario_data->{$try_path}->[2]};
2862 :     my (%diff_reactions_1, %diff_reactions_2);
2863 :    
2864 :     map {$diff_reactions_1{$_} = 1} @path_reactions;
2865 :     map {delete $diff_reactions_1{$_}} @done_reactions;
2866 :     map {$diff_reactions_2{$_} = 1} @done_reactions;
2867 :     map {delete $diff_reactions_2{$_}} @path_reactions;
2868 :    
2869 :     if (scalar keys %diff_reactions_1 == 0 && scalar keys %diff_reactions_2 == 0)
2870 :     {
2871 :     print STDERR "They match.\n" if $debug;
2872 :     $found_match = 1;
2873 :     last;
2874 :     }
2875 :     }
2876 :    
2877 :     if (! $found_match)
2878 :     {
2879 :     push @list_of_done_tokens, $try_path;
2880 :     print STDERR $try_path, " Added to the array\n" if $debug;
2881 :    
2882 :     }
2883 :     }
2884 :    
2885 :     foreach my $path (@list_of_done_tokens){
2886 :     if($scenario_data->{$path}->[1] != 1)
2887 :     {
2888 :     next;
2889 :     }
2890 :     if(@list_of_done_tokens == 0){
2891 :     push @list_of_done_tokens, $path;
2892 :     }
2893 :    
2894 :     #create input/output info
2895 :     my $input_hash = $scenario_data->{$path}->[7];
2896 :     my $output_hash = $scenario_data->{$path}->[8];
2897 :     my $reaction_path = $scenario_data->{$path}->[2];
2898 :     my @reaction_array = @$reaction_path;
2899 :    
2900 :     print STDERR "\nInputs:\n" if $debug;
2901 :     print STDERR map{"$_ => $input_hash->{$_}" } keys %$input_hash, "\n" if $debug;
2902 :     print STDERR "\nOutputs:\n" if $debug;
2903 :     print STDERR map{"$_ => $output_hash->{$_}"} keys %$output_hash, "\n" if $debug;
2904 :    
2905 :     # divide stoichiometry by greatest common denominator
2906 :     my ($min_stoich, @all_stoichs);
2907 :    
2908 :     map { push @all_stoichs, $input_hash->{$_}; $min_stoich = $input_hash->{$_} if $input_hash->{$_} < $min_stoich || $min_stoich == 0 } keys %{$input_hash};
2909 :     map { push @all_stoichs, $output_hash->{$_}; $min_stoich = $output_hash->{$_} if $output_hash->{$_} < $min_stoich || $min_stoich == 0 } keys %{$output_hash};
2910 :    
2911 :     my ($gcd, @gcd_candidates);
2912 :    
2913 :     outer: for ($gcd = $min_stoich; $gcd > 1; $gcd--)
2914 :     {
2915 :     foreach my $stoich (@all_stoichs)
2916 :     {
2917 :     next outer if $stoich % $gcd != 0;
2918 :     }
2919 :    
2920 :     last; # found a gcd
2921 :     }
2922 :    
2923 :     map { $input_hash->{$_} /= $gcd } keys %{$input_hash};
2924 :     map { $output_hash->{$_} /= $gcd } keys %{$output_hash};
2925 :    
2926 :     mkdir "$scenario_dir/path_$path_count";
2927 :    
2928 :     open(FILE, ">$scenario_dir/path_$path_count/path_info");
2929 :     foreach my $elem(@reaction_array){
2930 :     print FILE scalar @reaction_array, "\t", $elem, "\n";
2931 :     }
2932 :     close(FILE);
2933 :    
2934 :     &write_fluxanalyzer_files("$scenario_dir/path_$path_count",$input_hash,
2935 :     $output_hash, \@reaction_array,\%all_reactions,
2936 :     \%reactions_to_substrate_arrays,
2937 :     \%reactions_to_product_arrays,
2938 :     {});
2939 :     $path_count++;
2940 :     }
2941 :     @list_of_done_tokens=();
2942 :     undef %{$scenario_data};
2943 :    
2944 :     }
2945 :    
2946 :    
2947 :     sub load_subsystem
2948 :     {
2949 :     my ($genome,$ss_name) = @_;
2950 :     my $subsystem = $fig->get_subsystem($ss_name);
2951 :     my @ss_scenarios = $subsystem->get_hope_scenario_names;
2952 :     foreach my $name (@ss_scenarios)
2953 :     {
2954 :     load_scenario($genome,$ss_name,$name);
2955 :     }
2956 :     }
2957 :    
2958 :     sub internal_subsystem
2959 :     {
2960 :     my ($genome,$ss_name,$find_first) = @_;
2961 :     my %scenario_to_paths;
2962 :    
2963 :     my $subsystem = $fig->get_subsystem($ss_name);
2964 :     my @ss_scenarios = $subsystem->get_hope_scenario_names;
2965 :    
2966 :     foreach my $name (@ss_scenarios)
2967 :     {
2968 :     $scenario_to_paths{$name} = internal_scenario($genome,$ss_name,$name,$find_first);
2969 :     }
2970 :    
2971 :     return \%scenario_to_paths;
2972 :     }
2973 :    
2974 :     sub run_subsystem
2975 :     {
2976 :     my ($genome,$superset,$subsystem,$find_first) = @_;
2977 :    
2978 :     my $subsystem_obj = $fig->get_subsystem($subsystem);
2979 :     my @ss_scenarios = $subsystem_obj->get_hope_scenario_names;
2980 :    
2981 : olson 1.2 system("rm", "-rf", "$FIG_Config::var/Models/$genome/$superset/$subsystem");
2982 : olson 1.1 mkdir "$FIG_Config::var/Models/$genome/";
2983 :     mkdir "$FIG_Config::var/Models/$genome/$superset";
2984 :     mkdir "$FIG_Config::var/Models/$genome/$superset/$subsystem";
2985 :    
2986 :     foreach my $name (@ss_scenarios)
2987 :     {
2988 :     run_scenario($genome,$superset,$subsystem,$name,$find_first);
2989 :     }
2990 :     }
2991 :    
2992 : olson 1.2 sub compare_subsystem
2993 :     {
2994 :     my ($genome,$superset,$subsystem) = @_;
2995 :     my %genome_scenarios;
2996 :    
2997 :     my $subsystem_obj = $fig->get_subsystem($subsystem);
2998 :     my @ss_scenarios = $subsystem_obj->get_hope_scenario_names;
2999 :    
3000 :     foreach my $name (@ss_scenarios)
3001 :     {
3002 :     $genome_scenarios{$name} = compare_scenario($genome,$superset,$subsystem,$name);
3003 :     }
3004 :    
3005 :     return \%genome_scenarios;
3006 :     }
3007 :    
3008 : olson 1.1 sub load_superset
3009 :     {
3010 :     my($genome, $superset_name) = @_;
3011 :    
3012 :     my @subsystems = @{$superset_to_ss{$superset_name}};
3013 :     foreach my $ss_name (@subsystems)
3014 :     {
3015 :     load_subsystem($genome,$ss_name);
3016 :     }
3017 :     }
3018 :    
3019 :     sub internal_superset
3020 :     {
3021 :     my($genome, $superset_name,$find_first) = @_;
3022 :    
3023 :     my @subsystems = @{$superset_to_ss{$superset_name}};
3024 :    
3025 :     my %supersets_data;
3026 :    
3027 :     foreach my $ss_name (@subsystems)
3028 :     {
3029 :     $supersets_data{$ss_name} = internal_subsystem($genome,$ss_name,$find_first);
3030 :     }
3031 :    
3032 :     return \%supersets_data;
3033 :     }
3034 :    
3035 :     sub run_superset
3036 :     {
3037 :     my($genome, $superset_name,$find_first) = @_;
3038 :    
3039 :     my @subsystems = @{$superset_to_ss{$superset_name}};
3040 :    
3041 :     `rm -rf $FIG_Config::var/Models/$genome/$superset_name`;
3042 :     mkdir "$FIG_Config::var/Models/$genome/";
3043 :     mkdir "$FIG_Config::var/Models/$genome/$superset_name";
3044 :    
3045 :     foreach my $ss_name (@subsystems)
3046 :     {
3047 :     print STDERR "Running Scenarios for subsystem $ss_name\n";
3048 :     run_subsystem($genome,$superset_name,$ss_name,$find_first);
3049 :     }
3050 :     }
3051 :    
3052 : olson 1.2 sub compare_superset
3053 :     {
3054 :     my($genome, $superset_name) = @_;
3055 :    
3056 :     my @subsystems = @{$superset_to_ss{$superset_name}};
3057 :     my %genome_subsystems;
3058 :    
3059 :     foreach my $ss_name (@subsystems)
3060 :     {
3061 :     print STDERR "Comparing Scenarios for $genome in subsystem $ss_name\n";
3062 :     $genome_subsystems{$ss_name} = compare_subsystem($genome,$superset_name,$ss_name);
3063 :     }
3064 :    
3065 :     return \%genome_subsystems;
3066 :     }
3067 :    
3068 : olson 1.1
3069 :     sub load_supersets
3070 :     {
3071 :     my($genome) = @_;
3072 :     foreach my $superset (keys %superset_to_ss)
3073 :     {
3074 :     load_superset($genome,$superset);
3075 :     }
3076 :    
3077 :     return (\%all_reactions,\%reactions_to_substrate_arrays,\%reactions_to_product_arrays);
3078 :    
3079 :     }
3080 :    
3081 :     sub run_supersets
3082 :     {
3083 :     my($genome,$find_first) = @_;
3084 :    
3085 :     foreach my $superset (keys %superset_to_ss)
3086 :     {
3087 :     run_superset($genome,$superset,$find_first);
3088 :     }
3089 :    
3090 :     }
3091 :    
3092 : olson 1.2 sub compare_supersets
3093 :     {
3094 :     my($genome) = @_;
3095 :     my %genome_supersets;
3096 :    
3097 :     foreach my $superset (keys %superset_to_ss)
3098 :     {
3099 :     $genome_supersets{$superset} = compare_superset($genome,$superset);
3100 :     }
3101 :    
3102 :     return \%genome_supersets;
3103 :     }
3104 :    
3105 : olson 1.1 sub run_genome_report
3106 :     {
3107 :     my ($genome) = @_;
3108 :     my @string_out;
3109 :     push @string_out,"Genome $genome\n";
3110 :     #get all the subsystems this genome is involved in
3111 :     foreach my $superset (keys %superset_to_ss)
3112 :     {
3113 :     foreach my $name (@{$superset_to_ss{$superset}})
3114 :     {
3115 :     push @string_out, @{print_ss_report($name,internal_subsystem($genome,$name,0))};
3116 :     }
3117 :     }
3118 :     return \@string_out;
3119 :     }
3120 :    
3121 :     sub print_ss_report
3122 :     {
3123 :     my ($ss_name,$scenario_to_paths) = @_;
3124 :     my @output = ();
3125 :     my %scenario_path_count;
3126 :    
3127 :     push(@output,"\tSubsystem $ss_name\n");
3128 :    
3129 :     foreach my $scenario (keys %$scenario_to_paths)
3130 :     {
3131 :     if($scenario_to_paths->{$scenario}->{"infinite"})
3132 :     {
3133 :     push(@output,"\tWarning: Possible Infinite loop\n");
3134 :     }
3135 :     delete $scenario_to_paths->{$scenario}->{"infinite"};
3136 :     foreach my $token (keys %{$scenario_to_paths->{$scenario}})
3137 :     {
3138 :     $scenario_path_count{$scenario}++ if ($scenario_to_paths->{$scenario}->{$token}->[1]);
3139 :     }
3140 :     }
3141 :    
3142 :     push @output , map { "\t\t$_ has $scenario_path_count{$_} path(s).\n" } keys %scenario_path_count;
3143 :    
3144 :     return \@output;
3145 :     }
3146 :    
3147 :    
3148 :     sub internal_assembly
3149 :     {
3150 :     #lets get the genome, and a array reference to the paths we want to build togather
3151 :     my ($paths,$input_path,$output_path,$one_path) = @_;
3152 :    
3153 :    
3154 : olson 1.2 print STDERR $paths."\n" if $debug;
3155 : olson 1.1 clear_arrays();
3156 :    
3157 :     #This gets us an array of arrays, each subarray holds
3158 :     # [0] = genome [1] = superset [2] = subsystem [3] = scenario [4] = path
3159 :     # OR [0] = genome [1] = assembly [2] = path_name
3160 :     my @assembly_scenarios = @{parse_assembly_scenarios($paths)};
3161 :    
3162 :     #split and the input/output path for later as well
3163 :     #these should only return one path array...so just grab that one
3164 :     my @input_arr;
3165 :     my @output_arr;
3166 :    
3167 :     if($input_path != undef && $output_path != undef )
3168 :     {
3169 :     @input_arr = @{parse_assembly_scenarios($input_path)};
3170 :     @output_arr = @{parse_assembly_scenarios($output_path)};
3171 :     }
3172 :    
3173 :     #load all the kegg information for each 'scenario' from the paths we have selected
3174 :     foreach my $scenario (@assembly_scenarios)
3175 :     {
3176 :     if(scalar @$scenario > 4) #this is a normal scenario path
3177 :     {
3178 :     print STDERR "Checking $scenario->[2] $scenario->[3] $scenario->[4] \n" if $debug;
3179 :     process_init($scenario->[2],$scenario->[3],$scenario->[0],1);
3180 :     }
3181 :     else #This is a assembly path
3182 :     {
3183 :     #read in the input/output compounds and mark main's correctly
3184 :     my $path = "$FIG_Config::var/Models/" . join "/" , @$scenario;
3185 :    
3186 :     open(M_IN,"$path/inputs_main");
3187 :     while(<M_IN>)
3188 :     {
3189 :     my @line = split(/\t/,$_);
3190 :     $all_compounds_to_main{$line[0]} = $line[1];
3191 :     }
3192 :     close(M_IN);
3193 :    
3194 :     open(M_OUT,"$path/outputs_main");
3195 :     while(<M_OUT>)
3196 :     {
3197 :     my @line = split(/\t/,$_);
3198 :     $all_compounds_to_main{$line[0]} = $line[1];
3199 :     }
3200 :     close(M_OUT);
3201 :     }
3202 :     }
3203 :     #assume all path inputs and outputs are main
3204 :     map { $all_compounds_to_main{$_} = 1 } keys %all_inputs;
3205 :     map { $all_compounds_to_main{$_} = 1 } keys %all_outputs;
3206 :    
3207 :     print STDERR "Inputs: " if $debug;
3208 : olson 1.2 print STDERR map { $_."\n" } keys %all_inputs if $debug;
3209 : olson 1.1 print STDERR "Outputs: " if $debug;
3210 : olson 1.2 print STDERR map { $_."\n" } keys %all_outputs if $debug;
3211 : olson 1.1
3212 :     #run process paths
3213 :     return execute_paths(\@assembly_scenarios,$one_path,\@input_arr,\@output_arr);
3214 :     }
3215 :    
3216 :     sub run_assembly
3217 :     {
3218 :     my ($paths,$genome,$write_name,$one_path) = @_;
3219 :     print STDERR "\nThis is the passed information.\n";
3220 :     print STDERR $paths , "\n" , @$paths, "\n";
3221 :     print STDERR $genome . "\n";
3222 :     print STDERR $write_name . "\n";
3223 :     `rm -rf $FIG_Config::var/Models/$genome/assembly/$write_name`;
3224 :     mkdir "$FIG_Config::var/Models/$genome/assembly";
3225 :     mkdir "$FIG_Config::var/Models/$genome/assembly/$write_name";
3226 :     write_scenario(internal_assembly($paths,[],[],$one_path),"$FIG_Config::var/Models/$genome/assembly/$write_name");
3227 :     }
3228 :    
3229 :     sub expand_paths
3230 :     {
3231 :     my ($paths) = @_;
3232 :    
3233 :     my @final_paths;
3234 :    
3235 :     foreach my $path (@$paths)
3236 :     {
3237 :     if($path eq "" || $path eq "//")
3238 :     {
3239 :     next;
3240 :     }
3241 :     my $length = 5;
3242 :     if($path =~ /assembly/)
3243 :     {
3244 :     $length = 4;
3245 :     }
3246 :     print STDERR "Expanding $path \n" if $debug;
3247 :     my @parts = split "/", $path;
3248 :     shift @parts; # get ride of the first blank entry from /$genome
3249 :     my $length =$length - scalar @parts;
3250 :     $path = join "/" , @parts;
3251 :     print "Length : $length Path: $path\n" if $debug;
3252 :     my @temp = @{expand_recursive($path, $length)};
3253 :     push @final_paths , @temp if scalar @temp > 0;
3254 :     }
3255 :    
3256 :    
3257 :     return \@final_paths;
3258 :     }
3259 :    
3260 :     sub expand_recursive
3261 :     {
3262 :     my ($path,$count) = @_;
3263 :     my @sub_dirs;
3264 :     if($path =~ m\assembly$\)
3265 :     {
3266 :     return [];
3267 :     }
3268 :     if($count !=0)
3269 :     {
3270 :     #read this path, and pull out all the sub-directories.
3271 :     print STDERR "reading directory $FIG_Config::var/Models/$path/\n" if $debug;
3272 :     opendir (DIR,"$FIG_Config::var/Models/$path/") or die("$path");
3273 :     print STDERR "reading directory $FIG_Config::var/Models/$path/\n" if $debug;
3274 :     @sub_dirs = readdir DIR; shift @sub_dirs; shift @sub_dirs; # get rid of . and ..
3275 :     close DIR;
3276 :     print STDERR "Found: " if $debug;
3277 :     print STDERR @sub_dirs , "," if $debug;
3278 :     }
3279 :     else
3280 :     {
3281 :    
3282 :     return [$path];
3283 :     }
3284 :     $count--;
3285 :     my @to_return;
3286 :     foreach my $sub_path (@sub_dirs)
3287 :     {
3288 :     print STDERR "Calling on $sub_path , $count \n" if $debug;
3289 :     push @to_return , @{expand_recursive("$path/$sub_path",$count)};
3290 :     }
3291 :     print STDERR "Returning" ,@to_return if $debug;
3292 :     return \@to_return;
3293 :    
3294 :     }
3295 :    
3296 :     sub parse_assembly_scenarios
3297 :     {
3298 :     my ($paths) = @_;
3299 :    
3300 :     $paths = expand_paths($paths);
3301 :     my @array_of_path_arrays;
3302 :    
3303 :     foreach my $path (@$paths)
3304 :     {
3305 :     my @parts = split "/", $path;
3306 :     #shift @parts;
3307 :     push(@array_of_path_arrays, \@parts);
3308 :     }
3309 :    
3310 :     return \@array_of_path_arrays;
3311 :     }
3312 :    
3313 :    
3314 :     sub write_selected_scenarios
3315 :     {
3316 :     my($checked,$genome,$ssa,$sc_name) = @_;
3317 :     my (@tempArray);
3318 :    
3319 :     #Load this scenario again with all of its rxns and cpds
3320 :    
3321 :     model::clear_arrays();
3322 :    
3323 :     model::process_init($ssa,$sc_name,$genome,0);
3324 :    
3325 :     #assume all path inputs and outputs are main
3326 :     map { $all_compounds_to_main{$_} = 1 } keys %all_inputs;
3327 :     map { $all_compounds_to_main{$_} = 1 } keys %all_outputs;
3328 :    
3329 :     model::create_reactions({},{},{},{});
3330 :    
3331 :     ##End of scenario loading
3332 :    
3333 :    
3334 :     #setup the filesystem to store the scenario/paths
3335 :     my $superset = $ss_to_superset{$ssa};
3336 :     my $base_dir = "$FIG_Config::var/Models/$genome/$superset/$ssa/$sc_name/";
3337 : olson 1.2 system("rm", "-rf", $base_dir);
3338 : olson 1.1 mkdir "$FIG_Config::var/Models/$genome";
3339 :     mkdir "$FIG_Config::var/Models/$genome/$superset";
3340 :     mkdir "$FIG_Config::var/Models/$genome/$superset/$ssa/";
3341 :     mkdir $base_dir;
3342 :    
3343 :     #for the selected paths, lets gather their cpd from the checkbox and write
3344 :     #the path that we need
3345 :     foreach my $path (@$checked)
3346 :     {
3347 :     my (%input_hash, %output_hash, $path_name);
3348 :     #process the strings to get the information from the parameters
3349 :     my @items = split(";", $path);
3350 :     $path_name = $items[0];
3351 :     #next we have the input compounds ids/stoich/main
3352 :     map { if ($_ =~ /(.*):(.*):(.*)/)
3353 :     { $input_hash{$1}+= $2 } } split ",", $items[1];
3354 :     #the third part has the output compounds ids/stoich/main
3355 :     map { if ($_ =~ /(.*):(.*):(.*)/)
3356 :     { $output_hash{$1} += $2 } } split ",",$items[2];
3357 :     #the fourth part has the strings of reactions visited
3358 :     @tempArray = split("#" , $items[3]);
3359 :    
3360 :    
3361 :     print STDERR "\nInputs:\n" if $debug;
3362 :     print STDERR map{"$_ => $input_hash{$_}" } keys %input_hash, "\n" if $debug;
3363 :     print STDERR "\nOutputs:\n" if $debug;
3364 :     print STDERR map{"$_ => $output_hash{$_}"} keys %output_hash, "\n" if $debug;
3365 :    
3366 :    
3367 :     # divide stoichiometry by greatest common denominator
3368 :     my ($min_stoich, @all_stoichs);
3369 :    
3370 :     map { push @all_stoichs, $input_hash{$_}; $min_stoich = $input_hash{$_} if $input_hash{$_} < $min_stoich || $min_stoich == 0 } keys %input_hash;
3371 :     map { push @all_stoichs, $output_hash{$_}; $min_stoich = $output_hash{$_} if $output_hash{$_} < $min_stoich || $min_stoich == 0 } keys %output_hash;
3372 :    
3373 :     my ($gcd, @gcd_candidates);
3374 :    
3375 :     outer: for ($gcd = $min_stoich; $gcd > 1; $gcd--)
3376 :     {
3377 :     foreach my $stoich (@all_stoichs)
3378 :     {
3379 :     next outer if $stoich % $gcd != 0;
3380 :     }
3381 :    
3382 :     last; # found a gcd
3383 :     }
3384 :    
3385 :     map { $input_hash{$_} /= $gcd } keys %input_hash;
3386 :     map { $output_hash{$_} /= $gcd } keys %output_hash;
3387 :    
3388 :     mkdir "$base_dir/$path_name";
3389 :     &write_fluxanalyzer_files("$base_dir/$path_name",\%input_hash,
3390 :     \%output_hash,\@tempArray,\%all_reactions,
3391 :     \%reactions_to_substrate_arrays,
3392 :     \%reactions_to_product_arrays,
3393 :     {});
3394 :     }
3395 :     return $base_dir;
3396 :     }
3397 :    
3398 :     #This write function assumes that we have just run a assembly (and we haven't cleared the arrays)
3399 :     # becuase the write_fluxanalyzer_files function is depended on those global arrays for the reactions.
3400 :    
3401 :     sub write_assembly
3402 :     {
3403 :     my($input,$genome,$name) = @_;
3404 :    
3405 :     my $paths = $input->[1];
3406 :     my $file_paths = $input->[0];
3407 :     print STDERR "Paths: @$paths \n File_Dirs: @$file_paths \n" if $debug;
3408 :     my @tempArray;
3409 :    
3410 :     chomp $genome;
3411 :     chomp $name;
3412 :     #setup the filesystem to store the assembly
3413 :     my $base_dir = "$FIG_Config::var/Models/$genome/assembly/$name";
3414 :     `rm -rf $base_dir`;
3415 :     mkdir "$FIG_Config::var/Models/$genome/assembly";
3416 :     mkdir $base_dir;
3417 :    
3418 :     ##Here we want to reload the cpd and rxn info so we can write it later
3419 :     clear_arrays();
3420 :    
3421 :     #This gets us an array of arrays, each subarray holds
3422 :     # [0] = genome [1] = superset [2] = subsystem [3] = scenario [4] = path
3423 :     my @assembly_scenarios = @{parse_assembly_scenarios($file_paths)};
3424 :    
3425 :     #load all the kegg information for each 'scenario' from the paths we have selected
3426 :     foreach my $scenario (@assembly_scenarios)
3427 :     {
3428 :     print STDERR "Checking $scenario->[2] $scenario->[3] $scenario->[4] \n" if $debug;
3429 :     process_init($scenario->[2],$scenario->[3],$scenario->[0],1);
3430 :     }
3431 :     #assume all path inputs and outputs are main
3432 :     map { $all_compounds_to_main{$_} = 1 } keys %all_inputs;
3433 :     map { $all_compounds_to_main{$_} = 1 } keys %all_outputs;
3434 :    
3435 :     create_assembly_reactions({},{},{},{});
3436 :    
3437 :     ##End of rxn,cpd loading
3438 :    
3439 :     #for the selected paths, lets gather their cpd from the checkbox and write
3440 :     #the path that we need
3441 :     foreach my $path (@$paths)
3442 :     {
3443 :     my (%input_hash, %output_hash, $path_name);
3444 :     #process the strings to get the information from the parameters
3445 :     my @items = split(";", $path);
3446 :     $path_name = $items[0];
3447 :     #next we have the input compounds ids/stoich/main
3448 :     map { if ($_ =~ /(.*):(.*):(.*)/)
3449 :     { $input_hash{$1}+= $2 } } split ",", $items[1];
3450 :     #the third part has the output compounds ids/stoich/main
3451 :     map { if ($_ =~ /(.*):(.*):(.*)/)
3452 :     { $output_hash{$1} += $2 } } split ",",$items[2];
3453 :     #the fourth part has the strings of reactions visited
3454 :     @tempArray = split("#" , $items[3]);
3455 :    
3456 :     print STDERR "\nInputs:\n" if $debug;
3457 :     print STDERR map{"$_ => $input_hash{$_}" } keys %input_hash, "\n" if $debug;
3458 :     print STDERR "\nOutputs:\n" if $debug;
3459 :     print STDERR map{"$_ => $output_hash{$_}"} keys %output_hash, "\n" if $debug;
3460 :    
3461 :     # divide stoichiometry by greatest common denominator
3462 :     my ($min_stoich, @all_stoichs);
3463 :    
3464 :     map { push @all_stoichs, $input_hash{$_}; $min_stoich = $input_hash{$_} if $input_hash{$_} < $min_stoich || $min_stoich == 0 } keys %input_hash;
3465 :     map { push @all_stoichs, $output_hash{$_}; $min_stoich = $output_hash{$_} if $output_hash{$_} < $min_stoich || $min_stoich == 0 } keys %output_hash;
3466 :    
3467 :     my ($gcd, @gcd_candidates);
3468 :    
3469 :     outer: for ($gcd = $min_stoich; $gcd > 1; $gcd--)
3470 :     {
3471 :     foreach my $stoich (@all_stoichs)
3472 :     {
3473 :     next outer if $stoich % $gcd != 0;
3474 :     }
3475 :    
3476 :     last; # found a gcd
3477 :     }
3478 :    
3479 :     map { $input_hash{$_} /= $gcd } keys %input_hash;
3480 :     map { $output_hash{$_} /= $gcd } keys %output_hash;
3481 :    
3482 :     print STDERR "Making directory $base_dir" if $debug;
3483 :    
3484 :    
3485 :     mkdir "$base_dir/$path_name";
3486 :     &write_fluxanalyzer_files("$base_dir/$path_name",\%input_hash,\%output_hash,
3487 :     \@tempArray,\%all_reactions,
3488 :     \%reactions_to_substrate_arrays,
3489 :     \%reactions_to_product_arrays,
3490 :     {});
3491 :     }
3492 :     }
3493 :    
3494 :     sub show_path_results
3495 :     {
3496 :     my ($data_results,$html,$cgi) = @_;
3497 :    
3498 :     print STDERR "Starting Results Display\n";
3499 :    
3500 :    
3501 :     #Display infinite loop warning if the indicator is on
3502 :     if($data_results->{"infinite"})
3503 :     {
3504 :     push(@$html, "<h3>Warning: Looks like an infinite loop</h3>");
3505 :     }
3506 :    
3507 :     #Delete the infinite loop indicator, so we don't need to have a if statment to check for it
3508 :     delete $data_results->{"infinite"};
3509 :    
3510 :     my $path_counter = 1;
3511 :     my $reactionPath;
3512 :     foreach my $token_id (sort { $data_results->{$a}->[1] <=> $data_results->{$b}->[1] }keys %$data_results)
3513 :     {
3514 :     print STDERR "Token id : $token_id\n";
3515 :    
3516 :     if(!($token_id =~ /^\d/))
3517 :     {
3518 :     next;
3519 :     }
3520 :    
3521 :     my @path = @{$data_results->{$token_id}->[2]};
3522 :     my @compounds = @{$data_results->{$token_id}->[3]};
3523 :    
3524 :     push(@$html, "<pre>Token: $token_id\tInitial Pass: $data_results->{$token_id}->[0]\tDone:$data_results->{$token_id}->[1]\n\tReactions: @path\n\tVisted Compounds: @compounds\n\tPath Inputs\n@{$data_results->{$token_id}->[4]}\n\tPath Outputs\n@{$data_results->{$token_id}->[5]}\n\tBorrowed\n@{$data_results->{$token_id}->[6]}\n</pre>");
3525 :    
3526 :    
3527 :     if ($data_results->{$token_id}->[1] == 1)
3528 :     {
3529 :     my $path_name = "path_".$path_counter++;
3530 :     my @tempArray;
3531 :     foreach my $elem(@{$data_results->{$token_id}->[2]}){
3532 :     push @tempArray, $elem;
3533 :     }
3534 :     $reactionPath = join "#", @tempArray;
3535 :     my $checkbox=$cgi->checkbox(-name=>"$path_name", -label=>'',
3536 :     -value=>"$path_name;@{[ join ',', map { $_ . ':' . $data_results->{$token_id}->[7]->{$_} . ':' . $all_compounds_to_main{$_} } keys %{$data_results->{$token_id}->[7]} ]};@{[ join ',', map { $_ . ':' . $data_results->{$token_id}->[8]->{$_} . ':' . $all_compounds_to_main{$_} } keys %{$data_results->{$token_id}->[8]} ]};$reactionPath");
3537 :     push @$html, $checkbox, "&nbsp;$path_name", $cgi->br;
3538 :     }
3539 :    
3540 :     push @$html, "<hr>";
3541 :    
3542 :    
3543 :     }
3544 :    
3545 :    
3546 :     push @$html, $cgi->hidden(-name=>'reaction_info',
3547 :     -value=>$reactionPath);
3548 :    
3549 :     }
3550 :    
3551 :     sub show_path_results_two
3552 :     {
3553 :     my ($data_results,$html,$cgi) = @_;
3554 :    
3555 :     print STDERR "Starting Results Display\n";
3556 :    
3557 :    
3558 :     #Display infinite loop warning if the indicator is on
3559 :     if($data_results->{"infinite"})
3560 :     {
3561 :     push(@$html, "<h3>Warning: Looks like an infinite loop</h3>");
3562 :     }
3563 :    
3564 :     #Delete the infinite loop indicator, so we don't need to have a if statment to check for it
3565 :     delete $data_results->{"infinite"};
3566 :    
3567 :     my $path_counter = 1;
3568 :     my $reactionPath;
3569 :     foreach my $token_id (sort { $data_results->{$a}->[1] <=> $data_results->{$b}->[1] }keys %$data_results)
3570 :     {
3571 :     print STDERR "Token id : $token_id\n";
3572 :    
3573 :     if(!($token_id =~ /^\d/))
3574 :     {
3575 :     next;
3576 :     }
3577 :    
3578 :     my @path = @{$data_results->{$token_id}->[2]};
3579 :     my @compounds = @{$data_results->{$token_id}->[3]};
3580 :    
3581 :     push(@$html, "<pre>Token: $token_id\tInitial Pass: $data_results->{$token_id}->[0]\tDone:$data_results->{$token_id}->[1]\n\tReactions: @path\n\tVisted Compounds: @compounds\n\tPath Inputs\n@{$data_results->{$token_id}->[4]}\n\tPath Outputs\n@{$data_results->{$token_id}->[5]}\n\tBorrowed\n@{$data_results->{$token_id}->[6]}\n</pre>");
3582 :    
3583 :    
3584 :     if ($data_results->{$token_id}->[1] == 1)
3585 :     {
3586 :     my $path_name = "path_".$path_counter++;
3587 :     my @tempArray;
3588 :     foreach my $elem(@{$data_results->{$token_id}->[2]}){
3589 :     push @tempArray, $elem;
3590 :     }
3591 :     $reactionPath = join "#", @tempArray;
3592 :     my $checkbox=$cgi->checkbox(-name=>"checked", -label=>'',
3593 :     -value=>"$path_name;@{[ join ',', map { $_ . ':' . $data_results->{$token_id}->[7]->{$_} . ':' . $all_compounds_to_main{$_} } keys %{$data_results->{$token_id}->[7]} ]};@{[ join ',', map { $_ . ':' . $data_results->{$token_id}->[8]->{$_} . ':' . $all_compounds_to_main{$_} } keys %{$data_results->{$token_id}->[8]} ]};$reactionPath");
3594 :     push @$html, $checkbox, "&nbsp;$path_name", $cgi->br;
3595 :     }
3596 :    
3597 :     push @$html, "<hr>";
3598 :     }
3599 :    
3600 :     }
3601 :    
3602 :    
3603 :     sub set_loop_max
3604 :     {
3605 :     my ($number) = @_;
3606 :     $loop_max = $number;
3607 :     }
3608 :    
3609 :     sub set_loop_max_assembly
3610 :     {
3611 :     my ($number) = @_;
3612 :     $loop_max_assembly = $number;
3613 :     }
3614 :    
3615 :     my $dir_path = "$FIG_Config::var/Models";
3616 :     my @directory_array;
3617 :     my %path_structure;
3618 :    
3619 :     sub listDirectory(){
3620 :     my ($dirPath, $array_to_use) = @_;
3621 :     my @systemCallArray = `ls -d $dirPath*/`;
3622 :     foreach my $elem (@systemCallArray){
3623 :     chomp $elem;
3624 :     push @{$array_to_use}, $elem;
3625 :     listDirectory($elem, $array_to_use);
3626 :     }
3627 :     }
3628 :    
3629 :     #listDirectory($dir_path, \@directory_array);
3630 :    
3631 :     foreach my $elem (@directory_array){
3632 :     if ($elem =~ /path/ && $elem =~ /^assembly/){
3633 :     $path_structure{$elem} = [`ls $elem`];
3634 :     }
3635 :     }
3636 :    
3637 :     foreach my $elem (keys %path_structure){
3638 :    
3639 :     # open(INPUTPATHFILE, $inputFileName);
3640 :     # my @inFileCont = <INPUTPATHFILE>;
3641 :     # map {chomp $_} @inFileCont;
3642 :     # close(INPUTPATHFILE);
3643 :     #
3644 :     # open(OUTPUTPATHFILE, $outputFileName);
3645 :     # my @outFileCont = <OUTPUTPATHFILE>;
3646 :     # map {chomp $_ } @outFileCont;
3647 :     # close(OUTPUTPATHFILE);
3648 :     }
3649 :    
3650 :     1;

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3