1 #include "xml.h"
2
3 //////////////////////////////////////////////////////////////
4 //////////////////////////////////////////////////////////////
5
XML_Process_Base(char * xml_filename)6 t_tree *XML_Process_Base(char *xml_filename)
7 {
8 FILE *fp;
9 xml_node *root,*p_elem,*m_elem,*parent,*instance;
10 option *io;
11 void *buff;
12 t_mod *mod,*iomod;
13 t_tree *tree,*mixt_tree,*root_tree;
14 int select;
15 char *component;
16 int i,j,n_components;
17 int first_m_elem;
18 int class_number;
19 scalar_dbl **lens,**lens_var,**ori_lens,**lens_old,**lens_var_old,**ori_lens_old,**ori_lens_var,**ori_lens_var_old;
20 t_ds *ds;
21 char *outputfile;
22 char *alignment;
23 char *s;
24 int lens_size;
25 int first;
26 int *class_num;
27
28
29 fp = fopen(xml_filename,"r");
30 if(!fp)
31 {
32 PhyML_Fprintf(stderr,"\n. Could not find the XML file '%s'.\n",xml_filename);
33 Exit("\n");
34 }
35
36 root = XML_Load_File(fp);
37
38 if(!root)
39 {
40 PhyML_Fprintf(stderr,"\n. Encountered an issue while loading the XML file.\n");
41 Exit("\n");
42 }
43
44 class_num = (int *)mCalloc(N_MAX_MIXT_CLASSES,sizeof(int));
45 for(i=0;i<N_MAX_MIXT_CLASSES;i++) class_num[i] = i;
46
47 component = (char *)mCalloc(T_MAX_NAME,sizeof(char));
48
49 m_elem = NULL;
50 p_elem = root;
51 io = NULL;
52 mixt_tree = NULL;
53 root_tree = NULL;
54 mod = NULL;
55 tree = NULL;
56 lens = NULL;
57 lens_var = NULL;
58 lens_old = NULL;
59 lens_var_old = NULL;
60 select = -1;
61 class_number = -1;
62 ds = NULL;
63 lens_size = 0;
64 ori_lens = NULL;
65 ori_lens_old = NULL;
66 ori_lens_var = NULL;
67 ori_lens_var_old = NULL;
68 first = YES;
69
70 // Make sure there are no duplicates in node's IDs
71 XML_Check_Duplicate_ID(root);
72
73 int count = 0;
74 XML_Count_Number_Of_Node_With_Name("topology",&count,root);
75
76 if(count > 1)
77 {
78 PhyML_Fprintf(stderr,"\n. There should not more than one 'topology' node.");
79 PhyML_Fprintf(stderr,"\n. Found %d. Please fix your XML file",count);
80 Exit("\n");
81 }
82 else if(count < 1)
83 {
84 PhyML_Fprintf(stderr,"\n. There should be at least one 'topology' node.");
85 PhyML_Fprintf(stderr,"\n. Found none. Please fix your XML file");
86 Exit("\n");
87 }
88
89 #if defined PHYML
90 p_elem = XML_Search_Node_Name("phyml",NO,p_elem);
91 #elif defined PHYTIME
92 p_elem = XML_Search_Node_Name("phytime",NO,p_elem);
93 #endif
94
95 if(!p_elem)
96 {
97 #if defined PHYML
98 PhyML_Fprintf(stderr,"\n. The 'phyml' tag is mandatory.");
99 #elif defined PHYTIME
100 PhyML_Fprintf(stderr,"\n. The 'phytime' tag is mandatory.");
101 #endif
102 PhyML_Fprintf(stderr,"\n. Please amend your XML file accordingly.");
103 assert(false);
104 }
105
106
107 /*! Input file
108 */
109 outputfile = XML_Get_Attribute_Value(p_elem,"output.file");
110
111 if(!outputfile || !(strlen(outputfile)>0))
112 {
113 #if defined PHYML
114 PhyML_Fprintf(stderr,"\n. The 'outputfile' attribute in 'phyml' tag is mandatory.");
115 #elif defined PHYTIME
116 PhyML_Fprintf(stderr,"\n. The 'outputfile' attribute in 'phytime' tag is mandatory.");
117 #endif
118 PhyML_Fprintf(stderr,"\n. Please amend your XML file accordingly.");
119 assert(false);
120 }
121
122 io = (option *)Make_Input();
123 Set_Defaults_Input(io);
124
125 s = XML_Get_Attribute_Value(p_elem,"run.id");
126 if(s && strlen(s)>0)
127 {
128 io->append_run_ID = YES;
129 strcpy(io->run_id_string,s);
130 }
131
132 s = XML_Get_Attribute_Value(p_elem,"r.seed");
133 if(s && strlen(s)>0) io->r_seed = (int)atoi(s);
134
135 strcpy(io->out_file,outputfile);
136 strcpy(io->out_tree_file,outputfile);
137 strcpy(io->out_stats_file,outputfile);
138
139
140
141 # if defined(PHYTIME)
142 strcat(io->out_tree_file,"_phytime_tree");
143 # elif defined(PHYREX)
144 strcat(io->out_tree_file,"_phyrex_tree");
145 # elif defined(PHYML)
146 strcat(io->out_tree_file,"_phyml_tree");
147 #endif
148
149
150 # if defined(PHYTIME)
151 strcat(io->out_stats_file,"_phytime_stats");
152 # elif defined(PHYREX)
153 strcat(io->out_stats_file,"_phyrex_stats");
154 # elif defined(PHYML)
155 strcat(io->out_stats_file,"_phyml_stats");
156 #endif
157
158
159 if(io->append_run_ID) { strcat(io->out_tree_file,"_"); strcat(io->out_tree_file,io->run_id_string); }
160 if(io->append_run_ID) { strcat(io->out_stats_file,"_"); strcat(io->out_stats_file,io->run_id_string); }
161
162 strcat(io->out_stats_file,".txt");
163 strcat(io->out_tree_file,".txt");
164
165 if(XML_Search_Node_Attribute_Value("add","true",NO,root) != NULL)
166 {
167 io->fp_out_tree = Openfile(io->out_tree_file,APPEND);
168 io->fp_out_stats = Openfile(io->out_stats_file,APPEND);
169 }
170 else
171 {
172 io->fp_out_tree = Openfile(io->out_tree_file,WRITE);
173 io->fp_out_stats = Openfile(io->out_stats_file,WRITE);
174 }
175
176
177 s = XML_Get_Attribute_Value(p_elem,"print.trace");
178
179 if(s)
180 {
181 select = XML_Validate_Attr_Int(s,6,
182 "true","yes","y",
183 "false","no","n");
184
185 if(select < 3)
186 {
187 io->print_trace = YES;
188 strcpy(io->out_trace_file,outputfile);
189 strcat(io->out_trace_file,"_phyml_trace");
190 if(io->append_run_ID) { strcat(io->out_trace_file,"_"); strcat(io->out_trace_file,io->run_id_string); }
191 strcat(io->out_trace_file,".txt");
192 io->fp_out_trace = Openfile(io->out_trace_file,1);
193 }
194 }
195
196
197 s = XML_Get_Attribute_Value(p_elem,"print.json.trace");
198
199 if(s)
200 {
201 select = XML_Validate_Attr_Int(s,6,
202 "true","yes","y",
203 "false","no","n");
204
205 if(select < 3)
206 {
207 io->print_json_trace = YES;
208 strcpy(io->out_json_trace_file,outputfile);
209 strcat(io->out_json_trace_file,"_phyml_trace");
210 if(io->append_run_ID) { strcat(io->out_json_trace_file,"_"); strcat(io->out_json_trace_file,io->run_id_string); }
211 strcat(io->out_json_trace_file,".json");
212 io->fp_out_json_trace = Openfile(io->out_json_trace_file,READWRITE);
213 }
214 }
215
216 s = XML_Get_Attribute_Value(p_elem,"print.site.lk");
217
218 if(s)
219 {
220 select = XML_Validate_Attr_Int(s,6,
221 "true","yes","y",
222 "false","no","n");
223 if(select < 3) io->print_site_lnl = YES;
224
225 strcpy(io->out_lk_file,outputfile);
226 strcat(io->out_lk_file, "_phyml_lk");
227 if(io->append_run_ID) { strcat(io->out_lk_file,"_"); strcat(io->out_lk_file,io->run_id_string); }
228 strcat(io->out_lk_file, ".txt");
229 io->fp_out_lk = Openfile(io->out_lk_file,WRITE);
230 }
231
232
233 s = XML_Get_Attribute_Value(p_elem,"branch.test");
234 if(s)
235 {
236 if(!strcmp(s,"aLRT"))
237 {
238 io->ratio_test = ALRTSTAT;
239 }
240 else if(!strcmp(s,"aBayes"))
241 {
242 io->ratio_test = ABAYES;
243 }
244 else if(!strcmp(s,"SH"))
245 {
246 io->ratio_test = SH;
247 }
248 else if(!strcmp(s,"no") || !strcmp(s,"none"))
249 {
250 io->ratio_test = 0;
251 }
252 else
253 {
254 PhyML_Fprintf(stderr,"\n. '%s' is not a valid option for 'branch.test'.",s);
255 PhyML_Fprintf(stderr,"\n. Please use 'aLRT' or 'aBayes' or 'SH'.");
256 Exit("\n");
257 }
258 }
259
260 s = XML_Get_Attribute_Value(p_elem,"quiet");
261 if(s)
262 {
263 select = XML_Validate_Attr_Int(s,6,
264 "true","yes","y",
265 "false","no","n");
266 if(select < 3) io->quiet = YES;
267 }
268
269 s = XML_Get_Attribute_Value(p_elem,"memory.check");
270 if(s)
271 {
272 select = XML_Validate_Attr_Int(s,6,
273 "true","yes","y",
274 "false","no","n");
275 if(select >= 3) io->mem_question = NO;
276 }
277
278 /*! Read all partitionelem nodes and mixturelem nodes in each of them
279 */
280 do
281 {
282 p_elem = XML_Search_Node_Name("partitionelem",YES,p_elem);
283
284 if(p_elem == NULL) break;
285
286 buff = (option *)Make_Input();
287 Set_Defaults_Input(buff);
288 io->next = buff;
289 io->next->prev = io;
290
291 io = io->next;
292 if(first == YES)
293 {
294 io = io->prev;
295 io->next = NULL;
296 Free_Input(buff);
297 first = NO;
298 }
299
300
301 /*! Set the datatype (required when compressing data)
302 */
303 char *dt = NULL;
304 dt = XML_Get_Attribute_Value(p_elem,"data.type");
305 if(!dt)
306 {
307 PhyML_Fprintf(stderr,"\n. Please specify the type of data ('aa' or 'nt') for partition element '%s'",
308 XML_Get_Attribute_Value(p_elem,"id"));
309 PhyML_Fprintf(stderr,"\n. Syntax: 'data.type=\"aa\"' or 'data.type=\"nt\"'");
310 Exit("\n");
311 }
312
313 select = XML_Validate_Attr_Int(dt,2,"aa","nt");
314 switch(select)
315 {
316 case 0:
317 {
318 io->datatype = AA;
319 break;
320 }
321 case 1:
322 {
323 io->datatype = NT;
324 break;
325 }
326 default:
327 {
328 PhyML_Fprintf(stderr,"\n. Unknown data type. Must be either 'aa' or 'nt'.");
329 Exit("\n");
330 }
331 }
332
333
334 char *format = NULL;
335 format = XML_Get_Attribute_Value(p_elem,"format");
336 if(format)
337 {
338 if(!strcmp(format,"interleave") ||
339 !strcmp(format,"interleaved"))
340 {
341 io->interleaved = YES;
342 }
343 else if(!strcmp(format,"sequential"))
344 {
345 io->interleaved = NO;
346 }
347 }
348
349
350 /*! Attach a model to this io struct
351 */
352 io->mod = (t_mod *)Make_Model_Basic();
353 Set_Defaults_Model(io->mod);
354 io->mod->ras->n_catg = 1;
355 io->mod->io = io;
356 iomod = io->mod;
357
358 if(io->datatype == AA) io->mod->ns = 20;
359 else if(io->datatype == NT) io->mod->ns = 4;
360 else assert(FALSE);
361
362
363 /*! Attach an optimization structure to this model
364 */
365 iomod->s_opt = (t_opt *)Make_Optimiz();
366 Set_Defaults_Optimiz(iomod->s_opt);
367
368 iomod->s_opt->opt_kappa = NO;
369 iomod->s_opt->opt_lambda = NO;
370 iomod->s_opt->opt_rr = NO;
371 iomod->s_opt->opt_subst_param = NO;
372
373 /*! Input file
374 */
375 alignment = XML_Get_Attribute_Value(p_elem,"file.name");
376
377 if(!alignment)
378 {
379 PhyML_Fprintf(stderr,"\n. 'file.name' tag is mandatory. Please amend your");
380 PhyML_Fprintf(stderr,"\n. XML file accordingly.");
381 Exit("\n");
382 }
383
384 strcpy(io->in_align_file,alignment);
385
386 /*! Open pointer to alignment
387 */
388 io->fp_in_align = Openfile(io->in_align_file,READ);
389
390
391 /*! Load sequence file
392 */
393 io->data = Get_Seq(io);
394
395 /*! Close pointer to alignment
396 */
397 fclose(io->fp_in_align);
398
399 /*! Compress alignment
400 */
401 io->cdata = Compact_Data(io->data,io);
402
403 /*! Free uncompressed alignment
404 */
405 Free_Seq(io->data,io->n_otu);
406
407 /*! Create new mixture tree
408 */
409 buff = (t_tree *)Make_Tree_From_Scratch(io->cdata->n_otu,io->cdata);
410
411 if(mixt_tree)
412 {
413 mixt_tree->next_mixt = buff;
414 mixt_tree->next_mixt->prev_mixt = mixt_tree;
415 mixt_tree = mixt_tree->next_mixt;
416 mixt_tree->dp = mixt_tree->prev_mixt->dp+1;
417 }
418 else mixt_tree = buff;
419
420 /*! Connect mixt_tree to xml struct
421 */
422 mixt_tree->xml_root = root;
423
424 /*! Connect mixt_tree to io struct
425 */
426 mixt_tree->io = io;
427
428 /*! Connect mixt_tree to model struct
429 */
430 mixt_tree->mod = iomod;
431
432 /*! mixt_tree is a mixture tree
433 */
434 mixt_tree->is_mixt_tree = YES;
435
436 /*! mixt_tree is a mixture tree
437 */
438 mixt_tree->mod->is_mixt_mod = YES;
439
440 /*! Connect mixt_tree to compressed data
441 */
442 mixt_tree->data = io->cdata;
443
444 /*! Set total number of patterns
445 */
446 mixt_tree->n_pattern = io->cdata->crunch_len;
447
448 /*! Remove branch lengths from mixt_tree */
449 for(i=0;i<2*mixt_tree->n_otu-1;++i)
450 {
451 Free_Scalar_Dbl(mixt_tree->a_edges[i]->l);
452 Free_Scalar_Dbl(mixt_tree->a_edges[i]->l_old);
453 Free_Scalar_Dbl(mixt_tree->a_edges[i]->l_var);
454 Free_Scalar_Dbl(mixt_tree->a_edges[i]->l_var_old);
455 }
456
457 /*! Connect last tree of the mixture for the
458 previous partition element to the next mixture tree
459 */
460 if(tree)
461 {
462 tree->next = mixt_tree;
463 mixt_tree->prev = tree;
464 }
465
466 /*! Do the same for the model
467 */
468 if(mod)
469 {
470 mod->next = iomod;
471 iomod->prev = mod;
472 }
473
474 if(!root_tree) root_tree = mixt_tree;
475
476 /*! Tree size scaling factor */
477 char *scale_tree = NULL;
478 scale_tree = XML_Get_Attribute_Value(p_elem,"optimise.tree.scale");
479
480 if(scale_tree)
481 {
482 int select;
483
484 select = XML_Validate_Attr_Int(scale_tree,6,
485 "true","yes","y",
486 "false","no","n");
487
488 if(select < 3) mixt_tree->mod->s_opt->opt_br_len_mult = YES;
489 }
490
491 scale_tree = NULL;
492 scale_tree = XML_Get_Attribute_Value(p_elem,"tree.scale");
493
494 if(scale_tree)
495 {
496 char *scale_val;
497
498 scale_val = XML_Get_Attribute_Value(p_elem,"tree.scale");
499 if(scale_val)
500 {
501 mixt_tree->mod->br_len_mult->v = String_To_Dbl(scale_val);
502 mixt_tree->mod->br_len_mult_unscaled->v = String_To_Dbl(scale_val);
503 Free(scale_val);
504 }
505 }
506
507 /*! Process all the mixtureelem tags in this partition element
508 */
509 n_components = 0;
510 m_elem = p_elem;
511 first_m_elem = 0;
512 mod = NULL;
513 tree = NULL;
514 class_number = 0;
515 do
516 {
517 m_elem = XML_Search_Node_Name("mixtureelem",YES,m_elem);
518 if(m_elem == NULL) break;
519
520
521 if(!strcmp(m_elem->name,"mixtureelem"))
522 {
523 first_m_elem++;
524
525 /*! Rewind tree and model when processing a new mixtureelem node
526 */
527 if(first_m_elem > 1)
528 {
529 while(tree->prev && tree->prev->is_mixt_tree == NO) { tree = tree->prev; }
530 while(mod->prev && mod->prev->is_mixt_mod == NO) { mod = mod->prev; }
531 }
532
533 /*! Read and process model components
534 */
535 char *list;
536 list = XML_Get_Attribute_Value(m_elem,"list");
537
538 j = 0;
539 for(i=0;i<(int)strlen(list);++i) if(list[i] == ',') j++;
540
541 if(j != n_components && first_m_elem > 1)
542 {
543 PhyML_Fprintf(stderr,"\n. Discrepancy in the number of elements in nodes 'mixtureelem' partitionelem id '%s'",p_elem->id);
544 PhyML_Fprintf(stderr,"\n. Check 'mixturelem' node with list '%s'",list);
545 Exit("\n");
546 }
547 n_components = j;
548
549 i = j = 0;
550 component[0] = '\0';
551 while(1)
552 {
553 if(list[j] == ',' || j == (int)strlen(list))
554 {
555 /*! Reading a new component
556 */
557
558 if(first_m_elem == YES) /* Only true when processing the first mixtureelem node */
559 {
560 t_tree *this_tree;
561 t_mod *this_mod;
562
563 /*! Create new tree
564 */
565 this_tree = (t_tree *)Make_Tree_From_Scratch(io->cdata->n_otu,io->cdata);
566
567 /*! Update the number of mixtures */
568 iomod->n_mixt_classes++;
569
570 if(tree)
571 {
572 tree->next = this_tree;
573 tree->next->prev = tree;
574 }
575 else
576 {
577 mixt_tree->next = this_tree;
578 mixt_tree->next->prev = mixt_tree;
579 }
580
581 tree = this_tree;
582 tree->mixt_tree = mixt_tree;
583
584
585 /*! Create a new model
586 */
587 this_mod = (t_mod *)Make_Model_Basic();
588 Set_Defaults_Model(this_mod);
589 this_mod->ras->n_catg = 1;
590 this_mod->ns = iomod->ns;
591 /*! All br_len_multiplier point to the corresponding */
592 /*! parameter in the relevant mixt_tree */
593 Free_Scalar_Dbl(this_mod->br_len_mult);
594 this_mod->br_len_mult = iomod->br_len_mult;
595
596 Free_Scalar_Dbl(this_mod->br_len_mult_unscaled);
597 this_mod->br_len_mult_unscaled = iomod->br_len_mult_unscaled;
598
599 if(mod)
600 {
601 mod->next = this_mod;
602 mod->next->prev = mod;
603 }
604 else
605 {
606 this_mod->prev = iomod;
607 }
608
609 mod = this_mod;
610 if(!iomod->next) iomod->next = mod;
611 mod->io = io;
612
613 mod->s_opt = (t_opt *)Make_Optimiz();
614 Set_Defaults_Optimiz(mod->s_opt);
615
616 mod->s_opt->opt_alpha = NO;
617 mod->s_opt->opt_pinvar = NO;
618
619 tree->data = io->cdata;
620 tree->n_pattern = io->cdata->crunch_len;
621 tree->io = io;
622 tree->mod = mod;
623
624 if(tree->n_pattern != tree->prev->n_pattern) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
625 }
626
627 /*! Read a component */
628 component[i] = '\0';
629 if(j != (int)strlen(list)-1) i = 0;
630
631 /*! Find which node this ID corresponds to
632 */
633 instance = XML_Search_Node_ID(component,YES,root);
634
635 if(!instance)
636 {
637 PhyML_Fprintf(stderr,"\n. Could not find a node with id: '%s'.",component);
638 PhyML_Fprintf(stderr,"\n. Problem with 'mixtureelem' node, list '%s'.",list);
639 Exit("\n");
640 }
641
642 if(!instance->parent)
643 {
644 PhyML_Fprintf(stderr,"\n. Node '%s' with id:'%s' has no parent.",instance->name,component);
645 Exit("\n");
646 }
647
648 parent = instance->parent;
649
650 ////////////////////////////////////////
651 // SUBSTITUTION MODEL //
652 ////////////////////////////////////////
653
654 if(!strcmp(parent->name,"ratematrices"))
655 {
656 /* ! First time we process this 'instance' node which has this 'ratematrices' parent */
657 if(instance->ds->obj == NULL)
658 {
659 Make_Ratematrix_From_XML_Node(instance,io,mod);
660
661 ds = instance->ds;
662
663 /*! Connect the data structure n->ds to mod->r_mat */
664 ds->obj = (t_rmat *)(mod->r_mat);
665
666 /*! Create and connect the data structure n->ds->next to mod->kappa */
667 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
668 ds = ds->next;
669 ds->obj = (scalar_dbl *)(mod->kappa);
670
671 /*! Create and connect the data structure n->ds->next to mod->s_opt->opt_kappa */
672 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
673 ds = ds->next;
674 ds->obj = (int *)(&mod->s_opt->opt_kappa);
675
676 /*! Create and connect the data structure n->ds->next to mod->s_opt->opt_rr */
677 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
678 ds = ds->next;
679 ds->obj = (int *)(&mod->s_opt->opt_rr);
680
681 /*! Create and connect the data structure n->ds->next to mod->whichmodel */
682 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
683 ds = ds->next;
684 ds->obj = (int *)(&mod->whichmodel);
685
686 /*! Create and connect the data structure n->ds->next to mod->modelname */
687 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
688 ds = ds->next;
689 ds->obj = (t_string *)(mod->modelname);
690
691 /*! Create and connect the data structure n->ds->next to mod->ns */
692 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
693 ds = ds->next;
694 ds->obj = (int *)(&mod->ns);
695
696 /*! Create and connect the data structure n->ds->next to mod->modelname */
697 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
698 ds = ds->next;
699 ds->obj = (t_string *)(mod->custom_mod_string);
700
701
702 /*! Create and connect the data structure n->ds->next to mod->fp_aa_rate_mat */
703 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
704 ds = ds->next;
705 ds->obj = (FILE *)(mod->fp_aa_rate_mat);
706
707 /*! Create and connect the data structure n->ds->next to mod->aa_rate_mat_file */
708 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
709 ds = ds->next;
710 ds->obj = (t_string *)mod->aa_rate_mat_file;
711 }
712 else
713 {
714 /*! Connect to already extisting r_mat & kappa structs. */
715 t_ds *ds;
716
717
718 ds = instance->ds;
719 Free(mod->r_mat);
720 mod->r_mat = (t_rmat *)ds->obj;
721
722 ds = ds->next;
723 Free_Scalar_Dbl(mod->kappa);
724 mod->kappa = (scalar_dbl *)ds->obj;
725
726 ds = ds->next;
727 mod->s_opt->opt_kappa = *((int *)ds->obj);
728
729 ds = ds->next;
730 mod->s_opt->opt_rr = *((int *)ds->obj);
731
732 ds = ds->next;
733 mod->whichmodel = *((int *)ds->obj);
734
735 ds = ds->next;
736 Free_String(mod->modelname);
737 mod->modelname = (t_string *)ds->obj;
738
739 ds = ds->next;
740 mod->ns = *((int *)ds->obj);
741
742 ds = ds->next;
743 Free_String(mod->custom_mod_string);
744 mod->custom_mod_string = (t_string *)ds->obj;
745
746 ds = ds->next;
747 mod->fp_aa_rate_mat = (FILE *)ds->obj;
748
749 ds = ds->next;
750 Free_String(mod->aa_rate_mat_file);
751 mod->aa_rate_mat_file = (t_string *)ds->obj;
752
753 }
754 }
755
756 ////////////////////////////////////////
757 // STATE FREQS //
758 ////////////////////////////////////////
759
760 else if(!strcmp(parent->name,"equfreqs"))
761 {
762 /* If n->ds == NULL, the corrresponding node data structure, n->ds, has not */
763 /* been initialized. If not, do nothing. */
764 if(instance->ds->obj == NULL)
765 {
766 Make_Efrq_From_XML_Node(instance,io,mod);
767
768 t_ds *ds;
769
770 ds = instance->ds;
771 ds->obj = (t_efrq *)mod->e_frq;
772
773 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
774 ds = ds->next;
775 ds->obj = (int *)(&mod->s_opt->opt_state_freq);
776
777 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
778 ds = ds->next;
779 ds->obj = (int *)(&mod->e_frq->user_state_freq);
780
781 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
782 ds = ds->next;
783 ds->obj = (vect_dbl *)(mod->e_frq->user_b_freq);
784 }
785 else
786 {
787 /* Connect the data structure n->ds to mod->e_frq */
788
789 ds = instance->ds;
790 mod->e_frq = (t_efrq *)ds->obj;
791
792 ds = ds->next;
793 mod->s_opt->opt_state_freq = *((int *)ds->obj);
794
795 ds = ds->next;
796 mod->e_frq->user_state_freq = *((int *)ds->obj);
797
798 ds = ds->next;
799 mod->e_frq->user_b_freq = (vect_dbl *)ds->obj;
800 }
801 }
802
803 //////////////////////////////////////////
804 // TOPOLOGY //
805 //////////////////////////////////////////
806
807 else if(!strcmp(parent->name,"topology"))
808 {
809 if(parent->ds->obj == NULL) Make_Topology_From_XML_Node(instance,io,iomod);
810
811 ds = parent->ds;
812
813 int buff;
814 ds->obj = (int *)(& buff);
815 }
816
817 //////////////////////////////////////////
818 // RAS //
819 //////////////////////////////////////////
820
821 else if(!strcmp(parent->name,"siterates"))
822 {
823 char *rate_value = NULL;
824 /* scalar_dbl *r; */
825 phydbl val;
826
827
828 /*! First time we process this 'siterates' node, check that its format is valid.
829 and process it afterwards.
830 */
831 if(parent->ds->obj == NULL)
832 {
833 class_number = 0;
834
835 Make_RAS_From_XML_Node(parent,iomod);
836
837 ds = parent->ds;
838
839 ds->obj = (t_ras *)iomod->ras;
840
841 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
842 ds = ds->next;
843 ds->obj = (int *)(&iomod->s_opt->opt_alpha);
844
845 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
846 ds = ds->next;
847 ds->obj = (int *)(&iomod->s_opt->opt_free_mixt_rates);
848 }
849 else /*! Connect ras struct to an already defined one. Same for opt_alpha & opt_free_mixt_rates */
850 {
851 ds = parent->ds;
852
853 if(iomod->ras != (t_ras *)ds->obj) Free_RAS(iomod->ras);
854 iomod->ras = (t_ras *)ds->obj;
855
856 ds = ds->next;
857 iomod->s_opt->opt_alpha = *((int *)ds->obj);
858
859 ds = ds->next;
860 iomod->s_opt->opt_free_mixt_rates = *((int *)ds->obj);
861 }
862
863 rate_value = XML_Get_Attribute_Value(instance,"init.value");
864
865 val = 1.;
866 if(rate_value) val = String_To_Dbl(rate_value);
867
868 if(instance->ds->obj == NULL)
869 {
870 instance->ds->obj = (phydbl *)(&val);
871 instance->ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
872 instance->ds->next->obj = (int *)(class_num + class_number);
873
874 iomod->ras->gamma_rr->v[class_number] = val;
875 iomod->ras->init_rr = NO;
876
877 if(Are_Equal(val,0.0,1E-20) == NO) class_number++;
878 }
879
880
881 /*! Note: ras is already connected to the relevant t_ds stucture. No need
882 to connect ras->gamma_rr or ras->p_invar */
883
884 /*! Invariant */
885 if(Are_Equal(val,0.0,1E-20))
886 {
887 mod->ras->invar = YES;
888 }
889 else
890 {
891 mod->ras->parent_class_number = *((int *)instance->ds->next->obj);
892 }
893
894 xml_node *orig_w = NULL;
895 orig_w = XML_Search_Node_Attribute_Value("appliesto",instance->id,YES,instance->parent);
896
897
898 if(orig_w)
899 {
900 char *weight;
901 weight = XML_Get_Attribute_Value(orig_w,"value");
902 if(mod->ras->invar == YES)
903 {
904 iomod->ras->pinvar->v = String_To_Dbl(weight);
905 }
906 else
907 {
908 int class;
909 class = *((int *)instance->ds->next->obj);
910 iomod->ras->gamma_r_proba->v[class] = String_To_Dbl(weight);
911 iomod->ras->init_r_proba = NO;
912 }
913 }
914 }
915
916 //////////////////////////////////////////////
917 // BRANCH LENGTHS //
918 //////////////////////////////////////////////
919
920 else if(!strcmp(parent->name,"branchlengths"))
921 {
922 int i;
923 int n_otu;
924
925 n_otu = tree->n_otu;
926
927 if(instance->ds->obj == NULL)
928 {
929 if(!lens)
930 {
931 ori_lens = (scalar_dbl **)mCalloc(2*tree->n_otu-1,sizeof(scalar_dbl *));
932 ori_lens_old = (scalar_dbl **)mCalloc(2*tree->n_otu-1,sizeof(scalar_dbl *));
933
934 ori_lens_var = (scalar_dbl **)mCalloc(2*tree->n_otu-1,sizeof(scalar_dbl *));
935 ori_lens_var_old = (scalar_dbl **)mCalloc(2*tree->n_otu-1,sizeof(scalar_dbl *));
936
937 lens = ori_lens;
938 lens_old = ori_lens_old;
939
940 lens_var = ori_lens_var;
941 lens_var_old = ori_lens_var_old;
942
943 lens_size = 2*tree->n_otu-1;
944 }
945 else
946 {
947 ori_lens = (scalar_dbl **)mRealloc(ori_lens,2*tree->n_otu-1+lens_size,sizeof(scalar_dbl *));
948 ori_lens_old = (scalar_dbl **)mRealloc(ori_lens_old,2*tree->n_otu-1+lens_size,sizeof(scalar_dbl *));
949
950 ori_lens_var = (scalar_dbl **)mRealloc(ori_lens_var,2*tree->n_otu-1+lens_size,sizeof(scalar_dbl *));
951 ori_lens_var_old = (scalar_dbl **)mRealloc(ori_lens_var_old,2*tree->n_otu-1+lens_size,sizeof(scalar_dbl *));
952
953 lens = ori_lens + lens_size;;
954 lens_old = ori_lens_old + lens_size;;
955
956 lens_var = ori_lens_var + lens_size;;
957 lens_var_old = ori_lens_var_old + lens_size;;
958
959 lens_size += 2*tree->n_otu-1;
960 }
961
962 for(i=0;i<2*tree->n_otu-1;++i)
963 {
964 lens[i] = (scalar_dbl *)mCalloc(1,sizeof(scalar_dbl));
965 Init_Scalar_Dbl(lens[i]);
966
967 lens_old[i] = (scalar_dbl *)mCalloc(1,sizeof(scalar_dbl));
968 Init_Scalar_Dbl(lens_old[i]);
969
970 lens_var[i] = (scalar_dbl *)mCalloc(1,sizeof(scalar_dbl));
971 Init_Scalar_Dbl(lens_var[i]);
972
973 lens_var_old[i] = (scalar_dbl *)mCalloc(1,sizeof(scalar_dbl));
974 Init_Scalar_Dbl(lens_var_old[i]);
975
976 Free_Scalar_Dbl(tree->a_edges[i]->l);
977 Free_Scalar_Dbl(tree->a_edges[i]->l_old);
978
979 Free_Scalar_Dbl(tree->a_edges[i]->l_var);
980 Free_Scalar_Dbl(tree->a_edges[i]->l_var_old);
981
982 if(tree->prev &&
983 tree->prev->a_edges[i]->l == mixt_tree->a_edges[i]->l &&
984 tree->prev->is_mixt_tree == NO)
985 {
986 PhyML_Fprintf(stderr,"\n. %p %p",tree->a_edges[i]->l,mixt_tree->a_edges[i]->l);
987 PhyML_Fprintf(stderr,"\n. Only one set of edge lengths is allowed ");
988 PhyML_Fprintf(stderr,"\n. in a 'partitionelem'. Please fix your XML file.");
989 Exit("\n");
990 }
991 }
992
993 char *l_min = NULL;
994 l_min = XML_Get_Attribute_Value(instance,"l.min");
995
996 if(l_min)
997 {
998 iomod->l_min = String_To_Dbl(l_min);
999 }
1000
1001
1002 char *opt_bl = NULL;
1003 opt_bl = XML_Get_Attribute_Value(instance,"optimise.lens");
1004
1005 if(opt_bl)
1006 {
1007 if(!strcmp(opt_bl,"yes") || !strcmp(opt_bl,"true"))
1008 {
1009 iomod->s_opt->opt_bl = YES;
1010 }
1011 else
1012 {
1013 iomod->s_opt->opt_bl = NO;
1014 }
1015 }
1016
1017 ds = instance->ds;
1018
1019 ds->obj = (scalar_dbl **)lens;
1020
1021 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
1022 ds = ds->next;
1023 ds->obj = (scalar_dbl **)lens_old;
1024
1025 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
1026 ds = ds->next;
1027 ds->obj = (scalar_dbl **)lens_var;
1028
1029 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
1030 ds = ds->next;
1031 ds->obj = (scalar_dbl **)lens_var_old;
1032
1033 ds->next = (t_ds *)mCalloc(1,sizeof(t_ds));
1034 ds = ds->next;
1035 ds->obj = (int *)(&iomod->s_opt->opt_bl);
1036 }
1037 else
1038 {
1039 for(i=0;i<2*tree->n_otu-1;++i)
1040 {
1041 Free_Scalar_Dbl(tree->a_edges[i]->l);
1042 Free_Scalar_Dbl(tree->a_edges[i]->l_old);
1043 Free_Scalar_Dbl(tree->a_edges[i]->l_var);
1044 Free_Scalar_Dbl(tree->a_edges[i]->l_var_old);
1045 }
1046
1047 ds = instance->ds;
1048
1049 lens = (scalar_dbl **)ds->obj;
1050
1051 ds = ds->next;
1052 lens_old = (scalar_dbl **)ds->obj;
1053
1054 ds = ds->next;
1055 lens_var = (scalar_dbl **)ds->obj;
1056
1057 ds = ds->next;
1058 lens_var_old = (scalar_dbl **)ds->obj;
1059
1060 ds = ds->next;
1061 iomod->s_opt->opt_bl = *((int *)ds->obj);
1062 }
1063
1064 if(n_otu != tree->n_otu)
1065 {
1066 PhyML_Fprintf(stderr,"\n. All the data sets should display the same number of sequences.");
1067 PhyML_Fprintf(stderr,"\n. Found at least one data set with %d sequences and one with %d sequences.",n_otu,tree->n_otu);
1068 Exit("\n");
1069 }
1070
1071 For(i,2*tree->n_otu-1)
1072 {
1073 tree->a_edges[i]->l = lens[i];
1074 mixt_tree->a_edges[i]->l = lens[i];
1075 tree->a_edges[i]->l_old = lens_old[i];
1076 mixt_tree->a_edges[i]->l_old = lens_old[i];
1077
1078 tree->a_edges[i]->l_var = lens_var[i];
1079 mixt_tree->a_edges[i]->l_var = lens_var[i];
1080 tree->a_edges[i]->l_var_old = lens_var_old[i];
1081 mixt_tree->a_edges[i]->l_var_old = lens_var_old[i];
1082 }
1083 }
1084
1085 ///////////////////////////////////////////////
1086 ///////////////////////////////////////////////
1087 ///////////////////////////////////////////////
1088
1089 if(first_m_elem > 1) // Done with this component, move to the next tree and model
1090 {
1091 if(tree->next) tree = tree->next;
1092 if(mod->next) mod = mod->next;
1093 }
1094 }
1095 else if(list[j] != ' ')
1096 {
1097 component[i] = list[j];
1098 i++;
1099 }
1100 j++;
1101 if(j == (int)strlen(list)+1) break;
1102
1103 } // end of mixtureelem processing
1104 } // end of partitionelem processing
1105 }
1106 while(1);
1107 }
1108 while(1);
1109
1110
1111 if(ori_lens) Free(ori_lens);
1112 if(ori_lens_old) Free(ori_lens_old);
1113 if(ori_lens_var) Free(ori_lens_var);
1114 if(ori_lens_var_old) Free(ori_lens_var_old);
1115
1116 while(io->prev != NULL) io = io->prev;
1117 while(mixt_tree->prev != NULL) mixt_tree = mixt_tree->prev;
1118
1119
1120 /*! Finish making the models */
1121 mod = mixt_tree->mod;
1122 do
1123 {
1124 Make_Model_Complete(mod);
1125 mod = mod->next;
1126 }
1127 while(mod);
1128
1129 Check_Taxa_Sets(mixt_tree);
1130
1131 #if defined PHYML
1132 Check_Mandatory_XML_Node(root,"phyml");
1133 #elif defined PHYTIME
1134 Check_Mandatory_XML_Node(root,"phytime");
1135 #endif
1136
1137 Check_Mandatory_XML_Node(root,"topology");
1138 Check_Mandatory_XML_Node(root,"branchlengths");
1139 Check_Mandatory_XML_Node(root,"ratematrices");
1140 Check_Mandatory_XML_Node(root,"equfreqs");
1141 Check_Mandatory_XML_Node(root,"siterates");
1142 Check_Mandatory_XML_Node(root,"partitionelem");
1143 Check_Mandatory_XML_Node(root,"mixtureelem");
1144
1145 if(!io->mod->s_opt->random_input_tree) io->mod->s_opt->n_rand_starts = 1;
1146
1147
1148 /* int r_seed = (io->r_seed < 0)?(time(NULL)):(io->r_seed); */
1149 /* srand(r_seed); */
1150 /* io->r_seed = r_seed; */
1151
1152 Free(component);
1153 Free(class_num);
1154
1155 fclose(fp);
1156 return mixt_tree;
1157 }
1158
1159 //////////////////////////////////////////////////////////////
1160 //////////////////////////////////////////////////////////////
1161
XML_Load_File(FILE * fp)1162 xml_node *XML_Load_File(FILE *fp)
1163 {
1164 int c;
1165 char *buffer,*bufptr;
1166 int bufsize;
1167 xml_node *parent,*node;
1168
1169 buffer = (char *)mCalloc(T_MAX_XML_TAG,sizeof(char));
1170
1171 bufsize = T_MAX_XML_TAG;
1172 bufptr = buffer;
1173 parent = NULL;
1174 node = NULL;
1175
1176 while((c = fgetc(fp)) != EOF)
1177 {
1178 if(c == '<' && bufptr > buffer)
1179 {
1180 *bufptr = '\0';
1181
1182 /* PhyML_Printf("\n. Read value '%s' for node '%s'",buffer,node->name); */
1183 /* fflush(NULL); */
1184
1185 XML_Set_Node_Value(node,buffer);
1186 bufptr = buffer;
1187 }
1188
1189 if(c == '<')
1190 {
1191 bufptr = buffer;
1192
1193 while((c = fgetc(fp)) != EOF)
1194 {
1195 if(isspace(c) != NO || c == '>' || (c == '/' && bufptr > buffer)) break; // End of open or close tag
1196 else if(c == '<')
1197 {
1198 Exit("\n. Bare < in element!");
1199 }
1200 else if(XML_Add_Character(c,&bufptr,&buffer,&bufsize))
1201 {
1202 PhyML_Printf("\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1203 Exit("\n");
1204 }
1205 }
1206
1207 *bufptr = '\0';
1208
1209 if(!strcmp(buffer,"!--")) // Get the rest of the comment
1210 {
1211 while((c = fgetc(fp)) != EOF)
1212 {
1213
1214 if(c == '>' && bufptr > (buffer + 4) && bufptr[-3] != '-' &&
1215 bufptr[-2] == '-' && bufptr[-1] == '-') break;
1216 else if(XML_Add_Character(c,&bufptr,&buffer,&bufsize))
1217 {
1218 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1219 Exit("\n");
1220 }
1221 }
1222 *bufptr = '\0';
1223
1224 if(c != '>')
1225 {
1226 PhyML_Fprintf(stderr,"\n. Early EOF in comment node.");
1227 Exit("\n");
1228 }
1229 }
1230 else if(buffer[0] == '/') // Close tag
1231 {
1232 if(strcmp(buffer+1,parent->name))
1233 {
1234 PhyML_Fprintf(stderr,"\n. Opened tag with name '%s' and closed it with '%s'...",node->name,buffer+1);
1235 PhyML_Fprintf(stderr,"\n. Err. in file %s at line %d\n",__FILE__,__LINE__);
1236 Exit("\n");
1237 }
1238
1239 /* printf("\n. Closing node with name '%s'",node->name); */
1240
1241 if(node->parent)
1242 {
1243 parent = parent->parent;
1244 node = parent;
1245 }
1246 }
1247 else if(buffer[0] == '?')
1248 {
1249 while((c = fgetc(fp)) != EOF)
1250 {
1251 if (c == '>' && bufptr > buffer && bufptr[-1] == '?')
1252 break;
1253 else if (XML_Add_Character(c, &bufptr, &buffer, &bufsize))
1254 {
1255 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1256 Exit("\n");
1257 }
1258 }
1259
1260 if(c != '>')
1261 {
1262 PhyML_Fprintf(stderr,"\n. An error occurred when reading the processing instruction.");
1263 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1264 Exit("\n");
1265 }
1266
1267 *bufptr = '\0';
1268
1269 }
1270 else // Open tag
1271 {
1272 node = XML_Make_Node(buffer);
1273 XML_Init_Node(parent,node,buffer);
1274 if(!parent) parent = node;
1275
1276 if(isspace(c) != NO) c=XML_Parse_Element(fp,node);
1277 else if(c == '/')
1278 {
1279 if((c=fgetc(fp)) != '>')
1280 {
1281 PhyML_Fprintf(stderr,"\n. Expected '>' but read '%c' instead",c);
1282 Exit("\n");
1283 }
1284 c = '/';
1285 }
1286
1287 if(c != '/') parent = node;
1288
1289 buffer[0] = '\0';
1290 }
1291 bufptr = buffer;
1292 }
1293 else if(isspace(c) == NO)
1294 {
1295 if(XML_Add_Character(c,&bufptr,&buffer,&bufsize))
1296 {
1297 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1298 Exit("\n");
1299 }
1300 }
1301 }
1302 Free(buffer);
1303 return node;
1304 }
1305
1306 //////////////////////////////////////////////////////////////
1307 //////////////////////////////////////////////////////////////
1308
XML_Add_Character(int c,char ** bufptr,char ** buffer,int * bufsize)1309 int XML_Add_Character(int c, char **bufptr, char **buffer, int *bufsize)
1310 {
1311 char *newbuffer;
1312
1313 if(*bufptr >= (*buffer + *bufsize - 4))
1314 {
1315 // Increase the size of the buffer...
1316
1317 if (*bufsize < 1024)
1318 (*bufsize) *= 2;
1319 else
1320 (*bufsize) += 1024;
1321
1322 if((newbuffer = realloc(*buffer, *bufsize)) == NULL)
1323 {
1324 Free(*buffer);
1325 PhyML_Fprintf(stderr,"\n. Unable to expand string buffer to %d bytes!", *bufsize);
1326 Exit("\n");
1327 }
1328
1329 *bufptr = newbuffer + (*bufptr - *buffer);
1330 *buffer = newbuffer;
1331 }
1332
1333 /* *(*bufptr)++ = tolower(c); */
1334 *(*bufptr)++ = c;
1335 return 0;
1336 }
1337
1338 //////////////////////////////////////////////////////////////
1339 //////////////////////////////////////////////////////////////
1340
XML_Parse_Element(FILE * fp,xml_node * n)1341 int XML_Parse_Element(FILE *fp, xml_node *n)
1342 {
1343 int c;
1344 int quote;
1345 char *name, *value, *ptr;
1346 int namesize, valsize;
1347
1348 name = (char *)mCalloc(64,sizeof(char));
1349 value = (char *)mCalloc(64,sizeof(char));
1350
1351 namesize = 64;
1352 valsize = 64;
1353
1354 while((c = fgetc(fp)) != EOF)
1355 {
1356
1357 if(isspace(c) != NO) continue;
1358
1359 if(c == '/') // End of tag
1360 {
1361 /* printf("\n. Closing node '%s'.",n->name); */
1362
1363 quote = fgetc(fp);
1364 if(quote != '>')
1365 {
1366 PhyML_Fprintf(stderr,"\n. Expected '>' after '%c' but read '%c' instead",c,quote);
1367 Exit("\n");
1368 }
1369 break;
1370 }
1371 else if(c == '<')
1372 {
1373 Exit("\n. Bare < in element!");
1374 }
1375 else if(c == '>') // End of tag
1376 {
1377 break;
1378 }
1379
1380 name[0] = c;
1381 ptr = name + 1;
1382
1383 if(c == '\"' || c == '\'') // Name is in quotes
1384 {
1385 quote = c;
1386
1387 while((c = fgetc(fp)) != EOF)
1388 {
1389 if(XML_Add_Character(c,&ptr,&name,&namesize))
1390 {
1391 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1392 Exit("\n");
1393 }
1394 if(c == quote) break;
1395 }
1396 }
1397 else // Name not in quotes
1398 {
1399 while((c = fgetc(fp)) != EOF)
1400 {
1401 if(isspace(c) != NO || c == '=' || c == '/' || c == '>' || c == '?')
1402 break;
1403 else
1404 {
1405 if(XML_Add_Character(c,&ptr,&name,&namesize))
1406 {
1407 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1408 Exit("\n");
1409 }
1410 }
1411 }
1412 }
1413
1414 *ptr = '\0';
1415
1416 while(c != EOF && isspace(c) != NO) c = fgetc(fp);
1417
1418 if(c == '=') // Read the attribute value
1419 {
1420 while((c = fgetc(fp)) != EOF && isspace(c) != NO);
1421
1422 if(c == EOF)
1423 {
1424 PhyML_Fprintf(stderr,"\n. Missing value in attribute.");
1425 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1426 Exit("\n");
1427 }
1428
1429 if(c == '\'' || c == '\"')
1430 {
1431 quote = c;
1432 ptr = value;
1433
1434 while((c = fgetc(fp)) != EOF)
1435 {
1436 if(c == quote) break;
1437 else
1438 {
1439 if(XML_Add_Character(c,&ptr,&value,&valsize))
1440 {
1441 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1442 Exit("\n");
1443 }
1444 }
1445 }
1446 *ptr = '\0';
1447 }
1448 else
1449 {
1450 value[0] = c;
1451 ptr = value + 1;
1452
1453 while((c = fgetc(fp)) != EOF)
1454 {
1455 if(isspace(c) != NO || c == '=' || c == '/' || c == '>')
1456 break;
1457 else
1458 {
1459 if(XML_Add_Character(c,&ptr,&value,&valsize))
1460 {
1461 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1462 Exit("\n");
1463 }
1464 }
1465 }
1466 }
1467 }
1468
1469 /* printf("\n. Setting attribute '%s=%s' to node '%s'",name,value,n->name); */
1470 XML_Make_And_Set_Attribute(n,name,value);
1471
1472 if(c == '>') break;
1473
1474
1475 }
1476 Free(name);
1477 Free(value);
1478
1479 /* printf("\n. Return '%c'\n",c); */
1480 return(c);
1481 }
1482
1483 //////////////////////////////////////////////////////////////
1484 //////////////////////////////////////////////////////////////
1485
XML_Search_Attribute(xml_node * n,char * target_attr_name)1486 xml_attr *XML_Search_Attribute(xml_node *n, char *target_attr_name)
1487 {
1488 xml_attr *attr;
1489
1490 attr = n->attr;
1491 do
1492 {
1493 if(!strcmp(attr->name,target_attr_name)) return attr;
1494 attr = attr->next;
1495 }
1496 while(attr);
1497
1498 return(NULL);
1499 }
1500
1501 //////////////////////////////////////////////////////////////
1502 //////////////////////////////////////////////////////////////
1503
XML_Make_And_Set_Attribute(xml_node * n,char * attr_name,char * attr_value)1504 int XML_Make_And_Set_Attribute(xml_node *n, char *attr_name, char *attr_value)
1505 {
1506 xml_attr *prev;
1507 char *s;
1508
1509 prev = NULL;
1510 while(n->attr != NULL)
1511 {
1512 prev = n->attr;
1513 n->attr = n->attr->next;
1514 }
1515
1516 n->attr = XML_Make_Attribute(prev,attr_name,attr_value);
1517 XML_Init_Attribute(n->attr);
1518 n->n_attr++;
1519
1520 // rewind
1521 while(n->attr->prev != NULL) n->attr = n->attr->prev;
1522
1523 s = To_Lower_String(attr_name);
1524 if(!strcmp(s,"id"))
1525 {
1526 XML_Set_Node_Id(n,attr_value);
1527 /* printf("\n. Node '%s' id is '%s'",n->name,n->id); */
1528 }
1529 Free(s);
1530
1531 return(0);
1532 }
1533
1534 //////////////////////////////////////////////////////////////
1535 //////////////////////////////////////////////////////////////
1536
XML_Set_Attribute_Value(xml_node * n,char * attr_name,char * attr_value)1537 int XML_Set_Attribute_Value(xml_node *n, char *attr_name, char *attr_value)
1538 {
1539 xml_attr *attr;
1540 char *s;
1541
1542 attr = XML_Search_Attribute(n,attr_name);
1543 if(attr == NULL) return(-1);
1544
1545 s = To_Lower_String(attr_value);
1546 strcpy(attr->value,s);
1547
1548 Free(s);
1549
1550 return(0);
1551 }
1552
1553 //////////////////////////////////////////////////////////////
1554 //////////////////////////////////////////////////////////////
1555
XML_Add_Attribute(xml_node * n,char * attr_name,char * attr_value)1556 int XML_Add_Attribute(xml_node *n, char *attr_name, char *attr_value)
1557 {
1558 xml_attr *attr;
1559
1560 attr = n->attr;
1561 while(attr->next != NULL) attr = attr->next;
1562 assert(attr);
1563
1564 attr->next = XML_Make_Attribute(attr,attr_name,attr_value);
1565 XML_Init_Attribute(attr);
1566 n->n_attr++;
1567
1568 return(0);
1569 }
1570
1571 //////////////////////////////////////////////////////////////
1572 //////////////////////////////////////////////////////////////
1573
XML_Add_Node(xml_node * parent,char * nd_name)1574 xml_node *XML_Add_Node(xml_node *parent, char *nd_name)
1575 {
1576 xml_node *n;
1577 n = XML_Make_Node(nd_name);
1578 XML_Init_Node(parent,n,nd_name);
1579 return(n);
1580 }
1581
1582 //////////////////////////////////////////////////////////////
1583 //////////////////////////////////////////////////////////////
1584
XML_Set_Node_Id(xml_node * n,char * id)1585 int XML_Set_Node_Id(xml_node *n, char *id)
1586 {
1587 XML_Make_Node_Id(n,id);
1588 strcpy(n->id,id);
1589 return(0);
1590 }
1591
1592 //////////////////////////////////////////////////////////////
1593 //////////////////////////////////////////////////////////////
1594
XML_Set_Node_Value(xml_node * n,char * val)1595 int XML_Set_Node_Value(xml_node *n, char *val)
1596 {
1597 XML_Make_Node_Value(n,val);
1598 strcpy(n->value,val);
1599 return(0);
1600 }
1601
1602 //////////////////////////////////////////////////////////////
1603 //////////////////////////////////////////////////////////////
1604
XML_Search_Node_Generic(char * nd_name,char * attr_name,char * attr_val,int skip,xml_node * node)1605 xml_node *XML_Search_Node_Generic(char *nd_name, char *attr_name, char *attr_val, int skip, xml_node *node)
1606 {
1607
1608 xml_node *match;
1609
1610 /* if(nd_name) printf("\n. [1] nd_name:%s attr_name:%s attr_val:%s \n", nd_name, attr_name, attr_val); */
1611 /* else printf("\n. attr_name:%s attr_val:%s \n", attr_name, attr_val); */
1612
1613 /* printf("\n. name:%s child:%s next:%s ", */
1614 /* node?node->name:"xx", */
1615 /* node->child?node->child->name:"xx", */
1616 /* node->next?node->next->name:"xx"); fflush(NULL); */
1617
1618
1619 match = NULL;
1620 if(skip == NO && nd_name && attr_name && attr_val)
1621 {
1622 if(!strcmp(nd_name, node -> name))
1623 {
1624 xml_attr *attr = XML_Search_Attribute(node, attr_name);
1625 if(attr && !strcmp(attr -> value, attr_val)) match = node;
1626 }
1627 }
1628 else if(skip == NO && !nd_name && attr_name && attr_val)
1629 {
1630 xml_attr *attr = XML_Search_Attribute(node, attr_name);
1631 if(attr && !strcmp(attr -> value, attr_val)) match = node;
1632 }
1633 else if(skip == NO && nd_name && !attr_name && attr_val)
1634 {
1635 if(!strcmp(nd_name, node -> name))
1636 {
1637 do
1638 {
1639 if(!strcmp(node -> attr -> value, attr_val))
1640 {
1641 match = node;
1642 break;
1643 }
1644 node -> attr = node -> attr -> next;
1645 if(!node -> attr) break;
1646 }
1647 while(1);
1648 }
1649 }
1650 else if(skip == NO && nd_name && attr_name && !attr_val)
1651 {
1652 if(!strcmp(nd_name, node -> name))
1653 {
1654 do
1655 {
1656 if(!strcmp(node -> attr -> name, attr_name))
1657 {
1658 match = node;
1659 break;
1660 }
1661 node -> attr = node -> attr -> next;
1662 if(!node -> attr) break;
1663 }
1664 while(1);
1665 }
1666 }
1667 else if(skip == NO && nd_name && !attr_name && !attr_val)
1668 {
1669 if(!strcmp(nd_name, node -> name)) match = node;
1670 }
1671 else if(skip == NO && !nd_name && attr_name && !attr_val)
1672 {
1673 xml_attr *attr = XML_Search_Attribute(node, attr_name);
1674 if(attr) match = node;
1675 }
1676 else if(skip == NO && !nd_name && !attr_name && attr_val)
1677 {
1678 do
1679 {
1680 if(!strcmp(node -> attr -> value, attr_val))
1681 {
1682 match = node;
1683 break;
1684 }
1685 node -> attr = node -> attr -> next;
1686 if(!node -> attr) break;
1687 }
1688 while(1);
1689 }
1690
1691 // If node has a child, node = child, else if node has next, node = next, else if node
1692 // has parent, node = parent->next else node = NULL
1693
1694 if(match) return(match);
1695 if(node -> child)
1696 {
1697 match = XML_Search_Node_Generic(nd_name, attr_name, attr_val, NO, node -> child);
1698 }
1699 if(!match && node -> next)
1700 {
1701 match = XML_Search_Node_Generic(nd_name, attr_name, attr_val, NO, node -> next);
1702 }
1703 if(match == NULL && node -> parent)
1704 {
1705 if(node -> parent == NULL) // Reached the root
1706 {
1707 PhyML_Fprintf(stderr,"\n. Could not find a node with name '%s'.", attr_name);
1708 Exit("\n");
1709 }
1710 return NULL;
1711 }
1712 return match;
1713 }
1714
1715 //////////////////////////////////////////////////////////////
1716 //////////////////////////////////////////////////////////////
1717
XML_Search_Node_Name(char * name,int skip,xml_node * node)1718 xml_node *XML_Search_Node_Name(char *name, int skip, xml_node *node)
1719 {
1720 xml_node *match;
1721
1722 /* printf("\n. name:%s child:%s next:%s ", */
1723 /* node?node->name:"xx", */
1724 /* node->child?node->child->name:"xx", */
1725 /* node->next?node->next->name:"xx"); fflush(NULL); */
1726
1727 match = NULL;
1728 if(skip == NO && !strcmp(node->name,name)) match = node;
1729 else
1730 {
1731 // If node has a child, node = child, else if node has next, node = next, else if node
1732 // has parent, node = parent->next else node = NULL
1733 if(node->child)
1734 {
1735 match = XML_Search_Node_Name(name,NO,node->child);
1736 }
1737 if(match == NULL && node->next)
1738 {
1739 match = XML_Search_Node_Name(name,NO,node->next);
1740 }
1741 if(match == NULL && node->parent)
1742 {
1743 if(node->parent == NULL) // Reached the root
1744 {
1745 PhyML_Fprintf(stderr,"\n. Could not find a node with name '%s'.",name);
1746 Exit("\n");
1747 }
1748 return NULL;
1749 }
1750 }
1751 return match;
1752 }
1753
1754 //////////////////////////////////////////////////////////////
1755 //////////////////////////////////////////////////////////////
1756
XML_Search_Node_ID(char * id,int skip,xml_node * node)1757 xml_node *XML_Search_Node_ID(char *id, int skip, xml_node *node)
1758 {
1759 xml_node *match;
1760
1761 if(!node)
1762 {
1763 PhyML_Fprintf(stderr,"\n. Err. in file %s at line %d\n",__FILE__,__LINE__);
1764 Exit("\n");
1765 }
1766
1767
1768 match = NULL;
1769 if(skip == NO && node->id && !strcmp(node->id,id)) match = node;
1770 else
1771 {
1772 // If node has a child, node = child, else if node has next, node = next, else if node
1773 // has parent, node = parent->next else node = NULL
1774 if(node->child)
1775 {
1776 match = XML_Search_Node_ID(id,NO,node->child);
1777 }
1778 if(match == NULL && node->next)
1779 {
1780 match = XML_Search_Node_ID(id,NO,node->next);
1781 }
1782 if(match == NULL && node->parent)
1783 {
1784 if(node->parent == NULL) // Reached the root
1785 {
1786 PhyML_Fprintf(stderr,"\n. Could not find a node with id '%s'.",id);
1787 Exit("\n");
1788 }
1789
1790 return NULL;
1791 }
1792 }
1793 return match;
1794 }
1795
1796 //////////////////////////////////////////////////////////////
1797 //////////////////////////////////////////////////////////////
1798
XML_Search_Node_Attribute_Value(char * attr_name,char * value,int skip,xml_node * node)1799 xml_node *XML_Search_Node_Attribute_Value(char *attr_name, char *value, int skip, xml_node *node)
1800 {
1801 xml_node *match;
1802
1803
1804 if(!node)
1805 {
1806 PhyML_Fprintf(stderr,"\n. node: %p attr: %p",node,node?node->attr:NULL);
1807 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1808 Exit("\n");
1809 }
1810
1811 match = NULL;
1812
1813 if(skip)
1814 {
1815 match = XML_Search_Node_Attribute_Value(attr_name, value, NO, node->child);
1816 return match;
1817 }
1818
1819 if(skip == NO && node->attr)
1820 {
1821 xml_attr *attr;
1822 char *sname, *sval;
1823 char *qname, *qval;
1824
1825 qname = To_Lower_String(attr_name);
1826 qval = To_Lower_String(value);
1827
1828 attr = node->attr;
1829 do
1830 {
1831 sname = To_Lower_String(attr->name);
1832 sval = To_Lower_String(attr->value);
1833
1834 if(!strcmp(sname,qname) && !strcmp(sval,qval))
1835 {
1836 match = node;
1837 Free(sname);
1838 Free(sval);
1839 break;
1840 }
1841
1842 Free(sname);
1843 Free(sval);
1844
1845 attr = attr->next;
1846
1847 if(!attr) break;
1848 }
1849 while(1);
1850
1851 Free(qval);
1852 Free(qname);
1853 }
1854
1855 if(match) return(match);
1856
1857 if(node->child)
1858 {
1859 match = XML_Search_Node_Attribute_Value(attr_name,value,NO,node->child);
1860 return match;
1861 }
1862 if(node->next && !match)
1863 {
1864 match = XML_Search_Node_Attribute_Value(attr_name,value,NO,node->next);
1865 return match;
1866 }
1867 return NULL;
1868 }
1869
1870 //////////////////////////////////////////////////////////////
1871 //////////////////////////////////////////////////////////////
1872
XML_Get_Attribute_Value(xml_node * node,char * attr_name)1873 char *XML_Get_Attribute_Value(xml_node *node, char *attr_name)
1874 {
1875 xml_attr *attr;
1876
1877 attr = node->attr;
1878
1879 while(attr && strcmp(attr->name,attr_name)) attr = attr->next;
1880
1881 return(attr?attr->value:NULL);
1882 }
1883
1884 //////////////////////////////////////////////////////////////
1885 //////////////////////////////////////////////////////////////
1886
XML_Validate_Attr_Int(char * target,int num,...)1887 int XML_Validate_Attr_Int(char *target, int num, ...)
1888 {
1889 va_list args;
1890 int i;
1891 char *s,*sc_s;
1892 char *sc_target;
1893
1894 sc_target = To_Lower_String(target);
1895
1896 va_start(args,num);
1897 for(i=0;i<num;i++)
1898 {
1899 s = va_arg(args, char *);
1900 sc_s = To_Lower_String(s);
1901 if(!strcmp(sc_s,sc_target))
1902 {
1903 Free(sc_s);
1904 break;
1905 }
1906 Free(sc_s);
1907 }
1908 va_end(args);
1909
1910 if(i == num)
1911 {
1912 i = -1;
1913 PhyML_Fprintf(stderr,"\n. Attribute value '%s' is not valid",target);
1914 Exit("\n");
1915 }
1916
1917 Free(sc_target);
1918
1919 return(i);
1920 }
1921
1922 //////////////////////////////////////////////////////////////
1923 //////////////////////////////////////////////////////////////
1924
XML_Check_Siterates_Node(xml_node * parent)1925 void XML_Check_Siterates_Node(xml_node *parent)
1926 {
1927 xml_node *n;
1928 int n_weights_nodes;
1929 char *rate_value = NULL;
1930 phydbl buff;
1931 int n_zeros;
1932 char *endptr;
1933
1934 if(!parent)
1935 {
1936 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
1937 Exit("\n");
1938 }
1939
1940 if(strcmp(parent->name,"siterates"))
1941 {
1942 PhyML_Fprintf(stderr,"\n. Node name '%s' should be 'siterates'",parent->name);
1943 Exit("\n");
1944 }
1945
1946 // Check that only one 'weights' node is present
1947 n_weights_nodes = 0;
1948 n = parent->child;
1949 do
1950 {
1951 if(!strcmp(n->name,"weights")) n_weights_nodes++;
1952 if(n_weights_nodes > 1)
1953 {
1954 PhyML_Fprintf(stderr,"\n. Only one distribution is authorized for 'siterates' nodes.");
1955 Exit("\n");
1956 }
1957 n = n->next;
1958 if(!n) break;
1959 }
1960 while(1);
1961
1962 // Check that one rate value is set to zero if gamma+inv model is used
1963 n = XML_Search_Node_Attribute_Value("family","gamma+inv",YES,parent);
1964 if(!n) return;
1965 else
1966 {
1967 n_zeros = 0;
1968 n = parent->child;
1969 do
1970 {
1971 if(!strcmp(n->name,"instance"))
1972 {
1973 rate_value = NULL;
1974 rate_value = XML_Get_Attribute_Value(n,"init.value");
1975
1976 if(rate_value)
1977 {
1978 errno = 0;
1979 buff = strtod(rate_value,&endptr);
1980
1981 if(rate_value == endptr || errno == ERANGE)
1982 {
1983 PhyML_Fprintf(stderr,"\n. value: %s",rate_value);
1984 PhyML_Fprintf(stderr,"\n. Error in reading attribute 'init.value' in node 'instance'.");
1985 Exit("\n");
1986 }
1987
1988 if(buff < 1.E-20) n_zeros++;
1989 }
1990 }
1991 n = n->next;
1992 if(!n) break;
1993 }
1994 while(1);
1995
1996 if(n_zeros != 1)
1997 {
1998 PhyML_Fprintf(stderr,"\n. # of zero-rates: %d",n_zeros);
1999 PhyML_Fprintf(stderr,"\n. Exactly one rate value has to be set to zero when using the 'gamma+inv' model.");
2000 PhyML_Fprintf(stderr,"\n. Component id: %s",parent->id);
2001 Exit("\n");
2002 }
2003 }
2004 }
2005
2006 //////////////////////////////////////////////////////////////
2007 //////////////////////////////////////////////////////////////
2008
XML_Get_Number_Of_Classes_Siterates(xml_node * parent)2009 int XML_Get_Number_Of_Classes_Siterates(xml_node *parent)
2010 {
2011 xml_node *n;
2012 int n_classes;
2013
2014 if(!parent)
2015 {
2016 PhyML_Printf("\n. Err. in file %s at line %d\n",__FILE__,__LINE__);
2017 Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
2018 }
2019
2020 n_classes = 0;
2021 n = parent->child;
2022 do
2023 {
2024 if(!strcmp(n->name,"instance")) n_classes++;
2025 n = n->next;
2026 if(!n) break;
2027 }
2028 while(1);
2029
2030 n = NULL;
2031 n = XML_Search_Node_Attribute_Value("family","gamma+inv",YES,parent);
2032
2033 if(!n) return n_classes;
2034 else return n_classes-1;
2035 }
2036
2037 //////////////////////////////////////////////////////////////
2038
XML_Siterates_Has_Invariants(xml_node * parent)2039 int XML_Siterates_Has_Invariants(xml_node *parent)
2040 {
2041 xml_node *n;
2042
2043 if(!parent)
2044 {
2045 PhyML_Fprintf(stderr,"\n. Err in file %s at line %d\n",__FILE__,__LINE__);
2046 Exit("\n");
2047 }
2048
2049 n = NULL;
2050 n = XML_Search_Node_Attribute_Value("family","gamma+inv",YES,parent);
2051
2052 if(!n) return NO;
2053 else return YES;
2054 }
2055
2056 //////////////////////////////////////////////////////////////
2057
XML_Count_Number_Of_Node_With_ID(char * id,int * count,xml_node * n)2058 void XML_Count_Number_Of_Node_With_ID(char *id, int *count, xml_node *n)
2059 {
2060 if(!id) return;
2061 if(n->id && !strcmp(n->id,id)) (*count)++;
2062
2063 if(n->child) XML_Count_Number_Of_Node_With_ID(id,count,n->child);
2064 if(n->next) XML_Count_Number_Of_Node_With_ID(id,count,n->next);
2065
2066 }
2067
2068 //////////////////////////////////////////////////////////////
2069
XML_Count_Number_Of_Node_With_Name(char * name,int * count,xml_node * n)2070 void XML_Count_Number_Of_Node_With_Name(char *name, int *count, xml_node *n)
2071 {
2072 if(!name) return;
2073 if(n->name && !strcmp(n->name,name)) (*count)++;
2074
2075 if(n->child) XML_Count_Number_Of_Node_With_Name(name,count,n->child);
2076 if(n->next) XML_Count_Number_Of_Node_With_Name(name,count,n->next);
2077
2078 }
2079
2080 //////////////////////////////////////////////////////////////
2081
XML_Check_Duplicate_ID(xml_node * n)2082 void XML_Check_Duplicate_ID(xml_node *n)
2083 {
2084 int count;
2085
2086 count = 0;
2087 XML_Count_Number_Of_Node_With_ID(n->id,&count,n);
2088
2089 if(count > 1)
2090 {
2091 PhyML_Fprintf(stderr,"\n. Node ID '%s' was found more than once.",n->id);
2092 PhyML_Fprintf(stderr,"\n. Each ID must be unique. Please amend your XML");
2093 PhyML_Fprintf(stderr,"\n. file accordingly.");
2094 Exit("\n");
2095 }
2096
2097 if(n->child) XML_Check_Duplicate_ID(n->child);
2098 if(n->next) XML_Check_Duplicate_ID(n->next);
2099 }
2100
2101 //////////////////////////////////////////////////////////////
2102
XML_Copy_XML_Graph(xml_node * root)2103 xml_node *XML_Copy_XML_Graph(xml_node *root)
2104 {
2105 xml_node *cpy_root;
2106
2107 cpy_root = XML_Make_Node(root->name);
2108 XML_Copy_XML_Node(cpy_root,root);
2109
2110 return(cpy_root);
2111 }
2112
2113 //////////////////////////////////////////////////////////////
2114
XML_Copy_XML_Node(xml_node * cpy_root,xml_node * root)2115 void XML_Copy_XML_Node(xml_node *cpy_root, xml_node *root)
2116 {
2117 xml_attr *attr,*cpy_attr;
2118
2119 strcpy(cpy_root->name,root->name);
2120
2121 XML_Make_Node_Id(cpy_root,root->id);
2122 if(root->id) strcpy(cpy_root->id,root->id);
2123
2124 XML_Make_Node_Value(cpy_root,root->value);
2125 if(root->value) strcpy(cpy_root->value,root->value);
2126
2127 cpy_root->n_attr = root->n_attr;
2128
2129 if(root->attr)
2130 {
2131 cpy_root->attr = XML_Make_Attribute(NULL,root->attr->name,root->attr->value);
2132 XML_Init_Attribute(cpy_root->attr);
2133 attr = root->attr;
2134 cpy_attr = cpy_root->attr;
2135 while(attr->next)
2136 {
2137 fflush(NULL);
2138 cpy_attr->next = XML_Make_Attribute(cpy_attr,attr->next->name,attr->next->value);
2139 XML_Init_Attribute(cpy_attr->next);
2140 attr = attr->next;
2141 cpy_attr = cpy_attr->next;
2142 }
2143 }
2144
2145 if(root->child)
2146 {
2147 cpy_root->child = XML_Make_Node(root->child->name);
2148 cpy_root->child->parent = cpy_root;
2149 XML_Copy_XML_Node(cpy_root->child,root->child);
2150 }
2151
2152 if(root->next)
2153 {
2154 cpy_root->next = XML_Make_Node(root->next->name);
2155 cpy_root->next->prev = cpy_root;
2156 XML_Copy_XML_Node(cpy_root->next,root->next);
2157 }
2158 }
2159
2160 //////////////////////////////////////////////////////////////
2161
XML_Write_XML_Graph(FILE * fp,xml_node * root)2162 void XML_Write_XML_Graph(FILE *fp, xml_node *root)
2163 {
2164 int indent;
2165 indent = 0;
2166 XML_Write_XML_Node(fp,&indent,root);
2167 }
2168
2169 //////////////////////////////////////////////////////////////
2170
XML_Write_XML_Node(FILE * fp,int * indent,xml_node * root)2171 void XML_Write_XML_Node(FILE *fp, int *indent, xml_node *root)
2172 {
2173 xml_node *n;
2174 xml_attr *attr;
2175 char *s;
2176 int i;
2177
2178 s = (char *)mCalloc((*indent)+1,sizeof(char));
2179 for(i=0;i<(*indent);++i) s[i]='\t';
2180 s[i]='\0';
2181
2182
2183 PhyML_Fprintf(fp,"\n%s",s);
2184
2185 n = root;
2186
2187 PhyML_Fprintf(fp,"<%s",n->name);
2188
2189 attr = n->attr;
2190 while(attr)
2191 {
2192 PhyML_Fprintf(fp," %s=\"%s\"",attr->name,attr->value);
2193 attr = attr->next;
2194 }
2195 PhyML_Fprintf(fp,">");
2196
2197 if(n->value != NULL) XML_Write_Node_Value(fp,s,n);
2198
2199 if(n->child)
2200 {
2201 (*indent)++;
2202 XML_Write_XML_Node(fp,indent,n->child);
2203 (*indent)--;
2204 }
2205
2206 PhyML_Fprintf(fp,"\n%s</%s>",s,n->name);
2207
2208 PhyML_Fprintf(fp,"\n");
2209
2210 if(n->next) XML_Write_XML_Node(fp,indent,n->next);
2211
2212 Free(s);
2213 }
2214
2215 //////////////////////////////////////////////////////////////
2216 //////////////////////////////////////////////////////////////
2217
XML_Write_Node_Value(FILE * fp,char * indent,xml_node * n)2218 void XML_Write_Node_Value(FILE *fp, char *indent, xml_node *n)
2219 {
2220 PhyML_Fprintf(fp,"\n");
2221
2222 char *tk = strtok(n->value,"\n");
2223
2224 do
2225 {
2226 PhyML_Fprintf(fp,"%s%s",indent,tk);
2227 tk = strtok(NULL,"\n");
2228 if(tk != NULL) PhyML_Fprintf(fp,"\n");
2229 }
2230 while(tk != NULL);
2231
2232 }
2233 //////////////////////////////////////////////////////////////
2234 //////////////////////////////////////////////////////////////
2235
Check_Mandatory_XML_Node(xml_node * root,char * name)2236 void Check_Mandatory_XML_Node(xml_node *root, char *name)
2237 {
2238 if(!XML_Search_Node_Name(name,NO,root))
2239 {
2240 PhyML_Fprintf(stderr,"\n. Could not find mandatory XML node with name '%s'.",name);
2241 PhyML_Fprintf(stderr,"\n. Please amend your XML file.");
2242 Exit("\n");
2243 }
2244 }
2245
2246 //////////////////////////////////////////////////////////////
2247 //////////////////////////////////////////////////////////////
2248
XML_Number_Of_Taxa_In_Clade(xml_node * n_clade)2249 int XML_Number_Of_Taxa_In_Clade(xml_node *n_clade)
2250 {
2251 int clade_size = 0;
2252 if(n_clade)
2253 {
2254 do
2255 {
2256 clade_size++;
2257 if(n_clade->next) n_clade = n_clade -> next;
2258 else break;
2259 }
2260 while(n_clade);
2261 }
2262 else
2263 {
2264 PhyML_Fprintf(stderr,"\n. Clade is empty.");
2265 PhyML_Fprintf(stderr,"\n. Err. in file %s at line %d\n",__FILE__,__LINE__);
2266 Exit("\n");
2267 }
2268 return(clade_size);
2269 }
2270
2271 //////////////////////////////////////////////////////////////
2272 //////////////////////////////////////////////////////////////
2273
XML_Read_Clade(xml_node * xnd_clade,t_tree * tree)2274 char **XML_Read_Clade(xml_node *xnd_clade, t_tree *tree)
2275 {
2276 int i;
2277 char **clade;
2278
2279 clade = (char **)mCalloc(tree->n_otu, sizeof(char *));
2280
2281 if(xnd_clade)
2282 {
2283 i = 0;
2284 do
2285 {
2286 clade[i] = xnd_clade->attr->value;
2287 i++;
2288 if(xnd_clade->next) xnd_clade = xnd_clade->next;
2289 else break;
2290 }
2291 while(xnd_clade);
2292 }
2293 else
2294 {
2295 PhyML_Fprintf(stderr,"== Clade is empty. \n");
2296 PhyML_Fprintf(stderr,"\n. Err. in file %s at line %d\n",__FILE__,__LINE__);
2297 Exit("\n");
2298 }
2299
2300 return(clade);
2301 }
2302
2303 //////////////////////////////////////////////////////////////
2304 //////////////////////////////////////////////////////////////
2305
XML_Read_Calibration(xml_node * xroot,t_tree * tree)2306 void XML_Read_Calibration(xml_node *xroot, t_tree *tree)
2307 {
2308 xml_node *xnd_clade,*xnd_cal,*xnd,*xnd_dum;
2309 phydbl low,up,alpha_proba_dbl,t,t_ref;
2310 t_cal *cal;
2311 t_clad *clade;
2312 char *calib_id,*clade_id,*clade_name,*alpha_proba_string;
2313 int i,j;
2314
2315 clade = NULL;
2316 cal = NULL;
2317
2318 xnd = xroot->child;
2319 assert(xnd);
2320
2321 do
2322 {
2323 if(!strcmp(xnd->name,"calibration")) // Found a XML node <calibration>.
2324 {
2325 // TO DO: make sure calibs are shared across partition elements -> need to write chain function to
2326 // call once the calib struct on the first mixt_tree is initialized.
2327 /* mixt_tree->rates->tot_num_cal++; */
2328 /* if (mixt_tree->rates->calib == NULL) mixt_tree->rates->calib = Make_Calib(mixt_tree->n_otu); */
2329
2330 xnd_cal = xnd;
2331
2332 low = 0.0;
2333 up = BIG;
2334
2335 cal = Make_Calibration();
2336 Init_Calibration(cal);
2337
2338 xnd_dum = XML_Search_Node_Name("lower",YES,xnd_cal);
2339 if(xnd_dum != NULL) low = String_To_Dbl(xnd_dum->value);
2340
2341 xnd_dum = XML_Search_Node_Name("upper",YES,xnd_cal);
2342 if(xnd_dum != NULL) up = String_To_Dbl(xnd_dum->value);
2343
2344 calib_id = XML_Get_Attribute_Value(xnd_cal,"id");
2345 cal->id = (char *)mCalloc(strlen(calib_id)+1,sizeof(char));
2346 strcpy(cal->id,calib_id);
2347
2348 cal->clade_list_size = 0;
2349 cal->current_clade_idx = 0;
2350 cal->lower = low;
2351 cal->upper = up;
2352 cal->is_primary = YES;
2353
2354 tree->times->a_cal[tree->times->n_cal] = cal;
2355 tree->times->n_cal++;
2356
2357 if(tree->times->n_cal > 10 * tree->n_otu)
2358 {
2359 PhyML_Fprintf(stderr,"\n. There are too many clades defined that are not found in the tree...");
2360 PhyML_Fprintf(stderr,"\n. Please remove some of them.");
2361 Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
2362 }
2363
2364 xnd_dum = xnd_cal->child;
2365 do
2366 {
2367 clade_name = NULL;
2368 if(!strcmp("appliesto",xnd_dum->name))
2369 {
2370 clade_name = XML_Get_Attribute_Value(xnd_dum,"clade.id");
2371
2372 if(!clade_name)
2373 {
2374 PhyML_Fprintf(stderr,"\n. Attribute 'value=CLADE_NAME' is mandatory");
2375 PhyML_Fprintf(stderr,"\n. Please amend your XML file accordingly.");
2376 Exit("\n");
2377 }
2378
2379 alpha_proba_string = XML_Get_Attribute_Value(xnd_dum,"probability");
2380
2381 alpha_proba_dbl = -1;
2382 if(!alpha_proba_string) alpha_proba_dbl = 1.0;
2383 else alpha_proba_dbl = String_To_Dbl(alpha_proba_string);
2384 assert(!(alpha_proba_dbl < 0. || alpha_proba_dbl > 1.));
2385
2386 if(strcmp("root",clade_name))
2387 {
2388 xnd_clade = XML_Search_Node_Generic("clade","id",clade_name,YES,xroot);
2389
2390 if(xnd_clade != NULL) // found clade with a given name
2391 {
2392 char **xclade;
2393 int clade_size,nd_num;
2394 int i;
2395
2396 xclade = XML_Read_Clade(xnd_clade->child,tree);
2397 clade_size = XML_Number_Of_Taxa_In_Clade(xnd_clade->child);
2398
2399 clade = Make_Clade();
2400
2401 clade->n_tax = clade_size;
2402 clade->tax_list = (char **)mCalloc(clade_size,sizeof(char *));
2403 for(i=0;i<clade_size;++i) clade->tax_list[i] = (char *)mCalloc(strlen(xclade[i])+1,sizeof(char));
2404 for(i=0;i<clade_size;++i) strcpy(clade->tax_list[i],xclade[i]);
2405
2406 clade_id = XML_Get_Attribute_Value(xnd_dum,"clade.id");
2407 if(clade_id == NULL)
2408 {
2409 PhyML_Fprintf(stderr,"\n. Attribute \"clade.id\" is missing in \"appliesto\" tag.");
2410 Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
2411 }
2412 clade->id = (char *)mCalloc(strlen(clade_id)+1,sizeof(char));
2413 strcpy(clade->id,clade_id);
2414
2415 nd_num = Find_Clade(clade->tax_list,clade_size,tree);
2416
2417 if(nd_num > -1)
2418 {
2419 clade->target_nd = tree->a_nodes[nd_num];
2420 clade->tip_list = Make_Target_Tip(clade->n_tax);
2421 Init_Target_Tip(clade,tree);
2422
2423 if(cal->clade_list_size == 0) cal->clade_list = (t_clad **)mCalloc(1,sizeof(t_clad *));
2424 else cal->clade_list = (t_clad **)mRealloc(cal->clade_list,cal->clade_list_size+1,sizeof(t_clad *));
2425 if(cal->clade_list_size == 0) cal->alpha_proba_list = (phydbl *)mCalloc(1,sizeof(phydbl));
2426 else cal->alpha_proba_list = (phydbl *)mRealloc(cal->alpha_proba_list,cal->clade_list_size+1,sizeof(phydbl));
2427
2428 cal->clade_list[cal->clade_list_size] = clade;
2429 cal->alpha_proba_list[cal->clade_list_size] = alpha_proba_dbl;
2430 cal->clade_list_size++;
2431 }
2432 else
2433 {
2434 Free_Clade(clade);
2435 }
2436
2437 Free(xclade);
2438 }
2439 else
2440 {
2441 PhyML_Fprintf(stderr,"\n. Calibration information for clade [%s] was not found.", clade_name);
2442 PhyML_Fprintf(stderr,"\n. Err. in file %s at line %d\n",__FILE__,__LINE__);
2443 Exit("\n");
2444 }
2445 }
2446 }
2447 xnd_dum = xnd_dum->next;
2448 }
2449 while(xnd_dum != NULL);
2450 if(clade_name == NULL)
2451 {
2452 PhyML_Fprintf(stderr,"\n. Could not find calibration information for calibration `%s'",cal->id);
2453 PhyML_Fprintf(stderr,"\n. Please amend your XML file.");
2454 Exit("\n");
2455 }
2456 }
2457 xnd = xnd->next;
2458 }
2459 while(xnd != NULL);
2460
2461
2462 PhyML_Printf("\n\n.......................................................................");
2463 for(i=0;i<tree->times->n_cal;++i)
2464 {
2465 cal = tree->times->a_cal[i];
2466
2467 phydbl sum = 0.0;
2468 for(j=0;j<cal->clade_list_size;j++) sum += cal->alpha_proba_list[j];
2469 for(j=0;j<cal->clade_list_size;j++) cal->alpha_proba_list[j] /= sum;
2470
2471 for(j=0;j<cal->clade_list_size;j++)
2472 {
2473 clade = cal->clade_list[j];
2474 PhyML_Printf("\n Calibration id: %s.",cal->id);
2475 PhyML_Printf("\n Lower bound set to: %15f time units.",cal->lower);
2476 PhyML_Printf("\n Upper bound set to: %15f time units.",cal->upper);
2477 PhyML_Printf("\n This calibration applies to node %d with probability %G.",clade->target_nd->num,cal->alpha_proba_list[j]);
2478 PhyML_Printf("\n.......................................................................");
2479 }
2480 }
2481
2482
2483 t_ref = t = -1.;
2484
2485 for(i=0;i<tree->times->n_cal;++i)
2486 {
2487 cal = tree->times->a_cal[i];
2488 for(j=0;j<cal->clade_list_size;j++)
2489 {
2490 clade = cal->clade_list[j];
2491 if(clade->target_nd->tax == YES)
2492 {
2493 if(t_ref < -.5)
2494 {
2495 t_ref = Uni()*(cal->upper - cal->lower) + cal->lower;
2496 t_ref = fabs(t_ref);
2497 }
2498 else
2499 {
2500 t = Uni()*(cal->upper - cal->lower) + cal->lower;
2501 t = fabs(t);
2502
2503 if(Are_Equal(t,t_ref,1.E-10) == NO)
2504 {
2505 tree->times->is_asynchronous = YES;
2506 break;
2507 }
2508 }
2509 }
2510 }
2511 }
2512
2513 PhyML_Printf("\n\n");
2514 PhyML_Printf("\n. Is asynchronous: %s",tree->times->is_asynchronous ? "yes" : "no");
2515 }
2516
2517 //////////////////////////////////////////////////////////////
2518 //////////////////////////////////////////////////////////////
2519
XML_Update_XML_Struct_Given_Model_Params(t_tree * tree)2520 void XML_Update_XML_Struct_Given_Model_Params(t_tree *tree)
2521 {
2522 xml_node *parent,*child,*x;
2523 char *s,*val;
2524 t_ds *ds;
2525 phydbl v;
2526
2527 val = (char *)mCalloc(T_MAX_LINE,sizeof(char));
2528
2529 parent = XML_Search_Node_Name("ratematrices",YES,tree->xml_root);
2530 if(parent != NULL)
2531 {
2532 child = parent->child;
2533
2534 assert(child != NULL);
2535
2536 do
2537 {
2538 ds = child->ds;
2539
2540 assert(ds);
2541
2542 x = XML_Search_Node_Name("rr",YES,child);
2543
2544 if(x != NULL)
2545 {
2546 s = XML_Get_Attribute_Value(x,"AC");
2547 if(s != NULL)
2548 {
2549 sprintf(val,"%f",((t_rmat *)(ds->obj))->rr->v[0]);
2550 XML_Set_Attribute_Value(x,"AC",val);
2551 }
2552
2553 s = XML_Get_Attribute_Value(x,"AG");
2554 if(s != NULL)
2555 {
2556 sprintf(val,"%f",((t_rmat *)(ds->obj))->rr->v[1]);
2557 XML_Set_Attribute_Value(x,"AC",val);
2558 }
2559
2560 s = XML_Get_Attribute_Value(x,"AT");
2561 if(s != NULL)
2562 {
2563 sprintf(val,"%f",((t_rmat *)(ds->obj))->rr->v[2]);
2564 XML_Set_Attribute_Value(x,"AT",val);
2565 }
2566
2567 s = XML_Get_Attribute_Value(x,"CG");
2568 if(s != NULL)
2569 {
2570 sprintf(val,"%f",((t_rmat *)(ds->obj))->rr->v[3]);
2571 XML_Set_Attribute_Value(x,"CG",val);
2572 }
2573
2574 s = XML_Get_Attribute_Value(x,"CT");
2575 if(s != NULL)
2576 {
2577 sprintf(val,"%f",((t_rmat *)(ds->obj))->rr->v[4]);
2578 XML_Set_Attribute_Value(x,"CT",val);
2579 }
2580 }
2581
2582
2583 ds = ds->next;
2584 assert(ds);
2585
2586 s = XML_Get_Attribute_Value(child,"tstv");
2587 if(s != NULL)
2588 {
2589 sprintf(val,"%f",((scalar_dbl *)(ds->obj))->v);
2590 XML_Set_Attribute_Value(child,"tstv",val);
2591 }
2592
2593 child = child->next;
2594 }
2595 while(child);
2596 }
2597
2598
2599 parent = XML_Search_Node_Name("siterates",YES,tree->xml_root);
2600
2601 if(parent != NULL)
2602 {
2603 child = parent->child;
2604 assert(child != NULL);
2605
2606 int class = 0;
2607 do
2608 {
2609 s = XML_Get_Attribute_Value(child,"init.value");
2610 if(s != NULL)
2611 {
2612 v = atof(s);
2613 if(Are_Equal(v,0.0,1E-20) == NO)
2614 {
2615 sprintf(val,"%f",((t_ras *)(parent->ds->obj))->gamma_rr->v[class]);
2616 XML_Set_Attribute_Value(child,"init.value",val);
2617 class++;
2618 }
2619 }
2620
2621 child = child->next;
2622 }
2623 while(child);
2624
2625
2626 child = XML_Search_Node_Name("weights",YES,parent);
2627 assert(child != NULL);
2628 s = XML_Get_Attribute_Value(child,"family");
2629 assert(s != NULL);
2630 if(!strcmp(s,"gamma"))
2631 {
2632 s = XML_Get_Attribute_Value(child,"alpha");
2633 if(s != NULL)
2634 {
2635 sprintf(val,"%f",((t_ras *)(parent->ds->obj))->alpha->v);
2636 XML_Set_Attribute_Value(child,"alpha",val);
2637 }
2638 }
2639 else if(!strcmp(s,"gamma+inv"))
2640 {
2641 s = XML_Get_Attribute_Value(child,"alpha");
2642 if(s != NULL)
2643 {
2644 sprintf(val,"%f",((t_ras *)(parent->ds->obj))->alpha->v);
2645 XML_Set_Attribute_Value(child,"alpha",val);
2646
2647 sprintf(val,"%f",((t_ras *)(parent->ds->obj))->pinvar->v);
2648 XML_Set_Attribute_Value(child,"pinv",val);
2649 }
2650 }
2651
2652
2653 child = XML_Search_Node_Name("weights",YES,parent);
2654 child = child->child;
2655 assert(child);
2656
2657
2658 class = 0;
2659 do
2660 {
2661 // Find the rate class that this weight instance points to
2662 s = XML_Get_Attribute_Value(child,"appliesto");
2663 assert(s);
2664 x = XML_Search_Node_Attribute_Value("id",s,YES,parent);
2665 assert(x);
2666
2667 s = XML_Get_Attribute_Value(x,"init.value");
2668 if(s != NULL)
2669 {
2670 v = atof(s);
2671 if(Are_Equal(v,0.0,1E-20) == NO)
2672 {
2673 sprintf(val,"%f",((t_ras *)(parent->ds->obj))->gamma_r_proba->v[class]);
2674 XML_Set_Attribute_Value(child,"value",val);
2675 class++;
2676 }
2677 else
2678 {
2679 sprintf(val,"%f",0.0);
2680 XML_Set_Attribute_Value(child,"value",val);
2681 }
2682 }
2683 child = child->next;
2684 }
2685 while(child);
2686
2687
2688 }
2689
2690
2691
2692
2693
2694
2695 }
2696
2697
2698
2699
2700 //////////////////////////////////////////////////////////////
2701 //////////////////////////////////////////////////////////////
2702 //////////////////////////////////////////////////////////////
2703 //////////////////////////////////////////////////////////////
2704 //////////////////////////////////////////////////////////////
2705 //////////////////////////////////////////////////////////////
2706 //////////////////////////////////////////////////////////////
2707 //////////////////////////////////////////////////////////////
2708 //////////////////////////////////////////////////////////////
2709 //////////////////////////////////////////////////////////////
2710 //////////////////////////////////////////////////////////////
2711 //////////////////////////////////////////////////////////////
2712 //////////////////////////////////////////////////////////////
2713