shiba
 All Data Structures Files Functions Variables Macros Pages
io.c
Go to the documentation of this file.
1 
5 #define _GNU_SOURCE
6 #include <mxml.h>
7 #include "shiba.h"
8 
13 void readXML()
14 {
15 
16  FILE *fp;
17  mxml_node_t *tree;
18  mxml_node_t *node;
19 
20  char **TimeLabel; // The ID label in the XML file for the time period
21  char **SpaceLabel; // The ID label in the XML file for the spaces
22  char **PhyloLabel; // The ID label in the XML file for the phylogenies
23  char **TaxonLabel; // The ID label in the XML file for the taxa
24 
25  // Counters for occurrences of data elements
26  int t = 0; int s = 0; int p = 0; int x = 0;
27  int a = 0; int d = 0; int e = 0; int ss = 0;
28 
29  fp = fopen(DataFile, "r");
30  tree = mxmlLoadFile(NULL, fp, MXML_OPAQUE_CALLBACK);
31  fclose(fp);
32 
33  if (!tree) error("Unable to parse input file");
34 
35 
36  // ----------- Get core dimensioning elements first using Xpath-ish --------
37 
38  // Note the quirk that mxmlFindPath returns the child not the requested node:
39  node = mxmlGetParent(mxmlFindPath(tree, "shibaInput/time"));
40  Times = atoi(mxmlElementGetAttr(node, "n"));
41  if (!Times) error("time/@n not found in input");
42 
43  node = mxmlGetParent(mxmlFindPath(tree, "shibaInput/space"));
44  Spaces = atoi(mxmlElementGetAttr(node, "n"));
45  if (!Spaces) error("space/@n not found in input");
46 
47  node = mxmlGetParent(mxmlFindPath(tree, "shibaInput/phylo"));
48  Phylos = atoi(mxmlElementGetAttr(node, "n"));
49  if (!Phylos) error("phylo/@n not found in input");
50 
51  node = mxmlGetParent(mxmlFindPath(tree, "shibaInput/taxa"));
52  Taxa = atoi(mxmlElementGetAttr(node, "n"));
53  if (!Taxa) error("taxa/@n not found in input");
54 
55  // --- single config options
56 
57  node = mxmlGetParent(mxmlFindPath(tree, "*/config/nStartSpaces"));
58  Cfg.nStartSpaces = atoi(mxmlGetOpaque(mxmlGetFirstChild(node)));
59  if (!Cfg.nStartSpaces) error("//config/nStartSpaces not found in input");
60 
61  node = mxmlGetParent(mxmlFindPath(tree, "*/config/maxRuns"));
62  Cfg.maxRuns = atoi(mxmlGetOpaque(mxmlGetFirstChild(node)));
63  if (!Cfg.maxRuns) error("//config/maxRuns not found in input");
64  if (Cfg.maxRuns < RUN_BATCH) error("maxRuns must be more than specified");
65 
66  node = mxmlGetParent(mxmlFindPath(tree, "*/config/stopAfterSuccess"));
67  Cfg.stopAfterSuccess = atoi(mxmlGetOpaque(mxmlGetFirstChild(node)));
68  if (!Cfg.stopAfterSuccess) error("//stopAfterSuccess not found in input");
69 
70  if (Cfg.probSurvA == -1.0) {
71  node = mxmlGetParent(mxmlFindPath(tree, "*/config/probSurvA"));
72  Cfg.probSurvA = atof(mxmlGetOpaque(mxmlGetFirstChild(node)));
73  if (!Cfg.probSurvA) error("//config/probSurvA not found in input");
74  }
75 
76  if (Cfg.probDispA == -1.0) {
77  node = mxmlGetParent(mxmlFindPath(tree, "*/config/probDispA"));
78  Cfg.probDispA = atof(mxmlGetOpaque(mxmlGetFirstChild(node)));
79  if (!Cfg.probDispA) error("//config/probDispA not found in input");
80  }
81 
82  node = mxmlGetParent(mxmlFindPath(tree, "*/config/probSurvB"));
83  Cfg.probSurvB = atof(mxmlGetOpaque(mxmlGetFirstChild(node)));
84  if (!Cfg.probSurvB) error("//config/probSurvB not found in input");
85 
86  node = mxmlGetParent(mxmlFindPath(tree, "*/config/probDispB"));
87  Cfg.probDispB = atof(mxmlGetOpaque(mxmlGetFirstChild(node)));
88  if (!Cfg.probDispB) error("//config/probDispB not found in input");
89 
90  // ------------- Dimension arrays in the heap -----------------------------
91 
92  Area = mem2d_d(Times, Spaces);
93  Dist = mem3d_d(Times, Spaces, Spaces);
94  RealTime = mem1d_d(Times);
95  TimeLabel = mem2d1_c(Times);
96  SpaceName = mem2d1_c(Spaces);
97  SpaceLabel = mem2d1_c(Spaces);
98  Phylo = mem2d1_c(Phylos);
99  PhyloLabel = mem2d1_c(Phylos);
100  Taxon = mem2d1_c(Taxa);
101  TaxonLabel = mem2d1_c(Taxa);
102  Extant = mem3d_i(Taxa, Times, Spaces);
103  Cfg.startSpace = mem1d_i(Spaces);
104 
105  // --------------------- Crawl the XML tree to fill in the data ---------
106 
107  node = tree;
108  while(node)
109  {
110  if (mxmlGetElement(node))
111  {
112  // periodStart elements
113  if (!strcmp(mxmlGetElement(node), "periodStart"))
114  {
115  if (t >= Times) error("Too many periodStart elements");
116  if (!mxmlGetOpaque(mxmlGetFirstChild(node))) \
117  error("A periodStart element is missing");
118  else RealTime[t] = atof(mxmlGetOpaque(mxmlGetFirstChild(node)));
119  // test for declining times
120  if (t > 0)
121  {
126  if (RealTime[t] > RealTime[t-1])
127  error("Time sequence is not sorted in declining age");
128  }
129  if (!mxmlElementGetAttr(node, "id"))
130  error("An id attr is missing in periodStart");
131  else asprintf(&TimeLabel[t], "%s",
132  mxmlElementGetAttr(node, "id"));
133  t++;
134  }
135 
136  // spaceName elements
137  else if (!strcmp(mxmlGetElement(node), "spaceName"))
138  {
139  if (s >= Spaces) error("Too many spaceName elements");
140  if (!mxmlGetOpaque(mxmlGetFirstChild(node)))
141  error("A spaceName datum is missing");
142  else asprintf(&SpaceName[s], "%s",
143  mxmlGetOpaque(mxmlGetFirstChild(node)));
144  if (!mxmlElementGetAttr(node, "id"))
145  error("An id attr is missing in spaceName");
146  else asprintf(&SpaceLabel[s], "%s",
147  mxmlElementGetAttr(node, "id"));
148  s++;
149  }
150 
151  // phylo elements
152  else if (!strcmp(mxmlGetElement(node), "newick"))
153  {
154  if (p >= Phylos) error("Too many newick elements");
155  if (!mxmlGetOpaque(mxmlGetFirstChild(node)))
156  error("A newick datum is missing");
157  else
158  {
159  // scrub whitespace
160  char *tmp; int tmpn = 0;
161  asprintf(&tmp, "%s",
162  mxmlGetOpaque(mxmlGetFirstChild(node)));
163  for (int i = 0; i < strlen(tmp); i++)
164  if ((tmp[i] != 32) &&
165  (tmp[i] != 10) &&
166  (tmp[i] != 12) &&
167  (tmp[i] != 13)) {
168  if (!tmpn) {
169  asprintf(&Phylo[p], "%c", tmp[i]); tmpn++;
170  } else {
171  Sasprintf(Phylo[p], "%s%c", Phylo[p], tmp[i]); tmpn++;
172  }
173  }
174  free(tmp);
175  }
176  if (!mxmlElementGetAttr(node, "id"))
177  error("An id attr is missing in newick");
178  else asprintf(&PhyloLabel[p], "%s",
179  mxmlElementGetAttr(node, "id"));
180  p++;
181  }
182 
183  // taxa elements
184  else if (!strcmp(mxmlGetElement(node), "taxon"))
185  {
186  if (x >= Taxa) error("Too many taxon elements");
187  if (!mxmlGetOpaque(mxmlGetFirstChild(node))) \
188  error("A taxon datum is missing");
189  else asprintf(&Taxon[x], "%s", \
190  mxmlGetOpaque(mxmlGetFirstChild(node)));
191  if (!mxmlElementGetAttr(node, "id")) \
192  error("An id attr is missing in taxon");
193  else asprintf(&TaxonLabel[x], "%s", \
194  mxmlElementGetAttr(node, "id"));
195  x++;
196  }
197 
198  // area elements
199  else if (!strcmp(mxmlGetElement(node), "area"))
200  {
201  int a_t; int a_s;
202  if (a >= Spaces * Times) error("Too many area elements");
203  // Test for text child element
204  if (!mxmlGetOpaque(mxmlGetFirstChild(node)))
205  error("An area datum is missing");
206  // Read the time attr
207  if (!mxmlElementGetAttr(node, "time")) \
208  error("A time attr is missing in area");
209  else
210  {
211  a_t = -1;
212  for (int i = 0; i < Times; i++)
213  {
214  if (!strcmp(mxmlElementGetAttr(node, "time"),
215  TimeLabel[i])) a_t = i;
216  }
217  if (a_t == -1)
218  {
219  fprintf(stderr,
220  "XML parse error: //area/@time '%s' not IDREF\n",
221  mxmlElementGetAttr(node, "time"));
222  exit(1);
223  }
224  }
225  // Read the space attr
226  if (!mxmlElementGetAttr(node, "space")) {
227  error("A space attr is missing in area");
228  } else {
229  a_s = -1;
230  for (int i = 0; i < Spaces; i++)
231  if (!strcmp(mxmlElementGetAttr(node, "space"),
232  SpaceLabel[i])) a_s = i;
233  if (a_s == -1) {
234  fprintf(stderr,
235  "XML parse error: //area/@space '%s' not IDREF\n",
236  mxmlElementGetAttr(node, "space"));
237  exit(1);
238  }
239  }
240  // Set the datum
241  Area[a_t][a_s] = atof(mxmlGetOpaque(mxmlGetFirstChild(node)));
242  a++;
243  }
244 
245  // dist elements
246  else if (!strcmp(mxmlGetElement(node), "dist"))
247  {
248  int d_t; int d_s1; int d_s2;
249  if (a > Spaces * Spaces * Times)
250  error("Too many dist elements");
251  // Test for text child element
252  if (!mxmlGetOpaque(mxmlGetFirstChild(node)))
253  error("An dist datum is missing");
254  // Read the time attr
255  if (!mxmlElementGetAttr(node, "time")) \
256  error("A time attr is missing in dist");
257  else
258  {
259  d_t = -1;
260  for (int i = 0; i < Times; i++)
261  {
262  if (!strcmp(mxmlElementGetAttr(node, "time"),
263  TimeLabel[i])) d_t = i;
264  }
265  if (d_t == -1)
266  {
267  fprintf(stderr,
268  "XML parse error: //dist/@time '%s' not IDREF\n",
269  mxmlElementGetAttr(node, "time"));
270  exit(1);
271  }
272  }
273  // Read the space from attr
274  if (!mxmlElementGetAttr(node, "from"))
275  error("A from attr is missing in dist");
276  else
277  {
278  d_s1 = -1;
279  for (int i = 0; i < Spaces; i++)
280  {
281  if (!strcmp(mxmlElementGetAttr(node, "from"),
282  SpaceLabel[i])) d_s1 = i;
283  }
284  if (d_s1 == -1)
285  {
286  fprintf(stderr,
287  "XML parse error: //dist/@from '%s' not IDREF\n",
288  mxmlElementGetAttr(node, "from"));
289  exit(1);
290  }
291  }
292  // Read the space to attr
293  if (!mxmlElementGetAttr(node, "to"))
294  error("A to attr is missing in dist");
295  else
296  {
297  d_s2 = -1;
298  for (int i = 0; i < Spaces; i++)
299  {
300  if (!strcmp(mxmlElementGetAttr(node, "to"),
301  SpaceLabel[i])) d_s2 = i;
302  }
303  if (d_s2 == -1)
304  {
305  fprintf(stderr,
306  "XML parse error: //dist/@to '%s' not IDREF\n",
307  mxmlElementGetAttr(node, "to"));
308  exit(1);
309  }
310  }
311  // Set the datum (symetrically)
312  Dist[d_t][d_s1][d_s2] =
313  atof(mxmlGetOpaque(mxmlGetFirstChild(node)));
314  Dist[d_t][d_s2][d_s1] =
315  atof(mxmlGetOpaque(mxmlGetFirstChild(node)));
316  d++;
317  }
318 
319  // extant elements
320  else if (!strcmp(mxmlGetElement(node), "extant"))
321  {
322  int e_x; int e_t; int e_s;
323  // NB, No Test for text child element
324 
325  // Read the taxon from attr
326  if (!mxmlElementGetAttr(node, "taxon"))
327  error("A taxon attr is missing in extant");
328  else
329  {
330  e_x = -1;
331  for (int i = 0; i < Taxa; i++)
332  {
333  if (!strcmp(mxmlElementGetAttr(node, "taxon"),
334  TaxonLabel[i])) e_x = i;
335  }
336  if (e_x == -1)
337  {
338  fprintf(stderr,
339  "XML parse error: //extant/@taxon '%s' not IDREF\n",
340  mxmlElementGetAttr(node, "taxon"));
341  exit(1);
342  }
343  }
344  // Read the time attr
345  if (!mxmlElementGetAttr(node, "time")) \
346  error("A time attr is missing in extant");
347  else
348  {
349  e_t = -1;
350  for (int i = 0; i < Times; i++)
351  {
352  if (!strcmp(mxmlElementGetAttr(node, "time"),
353  TimeLabel[i])) e_t = i;
354  }
355  if (e_t == -1)
356  {
357  fprintf(stderr,
358  "XML parse error: //extant/@time '%s' not IDREF\n",
359  mxmlElementGetAttr(node, "time"));
360  exit(1);
361  }
362  }
363  // Read the space to attr
364  if (!mxmlElementGetAttr(node, "space"))
365  error("A space attr is missing in extant");
366  else
367  {
368  e_s = -1;
369  for (int i = 0; i < Spaces; i++)
370  {
371  if (!strcmp(mxmlElementGetAttr(node, "space"),
372  SpaceLabel[i])) e_s = i;
373  }
374  if (e_s == -1)
375  {
376  fprintf(stderr,
377  "XML parse error: //extant/@space '%s' not IDREF\n",
378  mxmlElementGetAttr(node, "space"));
379  exit(1);
380  }
381  }
382  // Set the datum (symetrically)
383  Extant[e_x][e_t][e_s] = 1;
384  e++;
385  }
386 
387  // ------------- Config --------------------------------------
388 
389  // startSpaces/allowed
390  else if (!strcmp(mxmlGetElement(node), "allowed"))
391  {
392  int ss_s;
393  if (ss >= Spaces) error("Too many startSpaces/allowed elements");
394  // Read the attr
395  if (!mxmlElementGetAttr(node, "space"))
396  error("A space attr is missing in allowed");
397  else
398  {
399  ss_s = -1;
400  for (int i = 0; i < Spaces; i++)
401  {
402  if (!strcmp(mxmlElementGetAttr(node, "space"),
403  SpaceLabel[i])) ss_s = i;
404  }
405  if (ss_s == -1)
406  {
407  fprintf(stderr,
408  "XML parse error: //allowed/@space '%s' not IDREF\n",
409  mxmlElementGetAttr(node, "space"));
410  exit(1);
411  }
412  }
413  Cfg.startSpace[ss_s] = 1;
414  ss++;
415  }
416  }
417 
418  node = mxmlWalkNext(node, tree, MXML_DESCEND);
419 
420  }
421 
422  // Check for correct dimensioning
423  if (t != Times) error("Wrong number of periodStart elements");
424  if (s != Spaces) error("Wrong number of spaceName elements");
425  if (p != Phylos) error("Wrong number of newick elements");
426  if (x != Taxa) error("Wrong number of taxon elements");
427  if (a != (Spaces * Times)) error("Wrong number of area elements");
428  if ((d != ((Spaces * (Spaces-1))/2) * Times) &&
429  (d != Spaces * Spaces * Times ))
430  error("Wrong number of dist elements");
431  if (e < 1) error("Need at least one extant taxon");
432 
433 
434  // free unused mem
435  mxmlDelete(tree);
436 
437  free2d1_c(PhyloLabel, Phylos);
438  free2d1_c(SpaceLabel, Spaces);
439  free2d1_c(TimeLabel, Times);
440  free2d1_c(TaxonLabel, Taxa);
441 
442 }
443 
446 void printIndata()
447 {
448  printf("# Data as read from %s\n", DataFile);
449 
450  printf("\n## Spaces:\n");
451  for (int i = 0; i < Spaces; i++) printf(" %3d: %s\n", i, SpaceName[i]);
452 
453  printf("\n## Times:\n");
454  for (int i = 0; i < Times; i++) printf(" %3d: %7.2f\n", i, RealTime[i]);
455 
456  printf("\n## Areas:\n");
457  printf(" Spaces :");
458  for (int i = 0; i < Spaces; i++) printf(" %6d", i);
459  printf("\n");
460  for (int i = 0; i < Times; i++)
461  {
462  printf(" Time %3d:", i);
463  for (int j = 0; j < Spaces; j++) printf(" %6.2f", Area[i][j]);
464  printf("\n");
465  }
466 
467  printf("\n## Distances:\n");
468  for (int t = 0; t < Times; t++)
469  {
470  printf(" Time: %3d\n", t);
471  printf(" From :");
472  for (int i = 0; i < Spaces; i++) printf(" %6d", i);
473  printf("\n");
474  for (int i = 0; i < Spaces; i++)
475  {
476  printf(" To %3d:", i);
477  for (int j = 0; j < Spaces; j++) printf(" %6.1f", Dist[t][i][j]);
478  printf("\n");
479  }
480  }
481 
482  printf("\n## Taxa:\n");
483  for (int i = 0; i < Taxa; i++) printf(" %d: %s\n", i, Taxon[i]);
484 
485  printf("\n## Phylos:\n");
486  for (int i = 0; i < Phylos; i++) printf(" %d: %s\n", i, Phylo[i]);
487 
488  printf("\n## LineagePeriods:\n");
489  printf(" Lineage : ");
490  for (int i = 0; i < Lineages; i++)
491  {
492  if ((i % 5) == 0) printf("|");
493  else printf(" ");
494  }
495  printf("\n");
496  for (int t = 0; t < Times; t++)
497  {
498  printf(" Period %3d: ", t);
499  for (int i = 0; i < Lineages; i++)
500  {
501  printf("%1d", LineagePeriod[i][t]);
502  }
503  printf("\n");
504  }
505 
506  printf("\n## LineageDaughters:\n");
507  for (int i = 0; i < Lineages; i++) {
508  printf(" Lineage %3d: ", i);
509  //printf(" %3d", LineageNDaughters[i]);
510  for (int j = 0; j < LineageDaughtersN[i]; j++)
511  printf(" %3d", LineageDaughters[i][j]);
512  printf("\n");
513  }
514 
515  printf("\n## ExtantTaxa:\n");
516  printf(" Spaces :");
517  for (int i = 0; i < Spaces; i++) printf(" %3d", i);
518  printf("\n");
519  for (int t = 0; t < Times; t++)
520  {
521  for (int i = 0; i < Taxa; i++)
522  {
523  int isOne = 0;
524  for (int j = 0; j < Spaces; j++)
525  if (Extant[i][t][j]) isOne++;
526  if (isOne)
527  {
528  printf(" Time %3d, Taxon %3d:", t, i);
529  for (int j = 0; j < Spaces; j++) printf(" %3d", Extant[i][t][j]);
530  printf("\n");
531  }
532  }
533  }
534 
535  printf("\n## LineageExtant:\n");
536  printf(" Spaces :");
537  for (int i = 0; i < Spaces; i++) printf(" %3d", i);
538  printf("\n");
539  for (int t = 0; t < Times; t++)
540  {
541  for (int i = 0; i < Lineages; i++)
542  {
543  int isOne = 0;
544  for (int j = 0; j < Spaces; j++)
545  if (LineageExtant[i][t][j]) isOne++;
546  if (isOne)
547  {
548  printf(" Time %3d, Lineage %3d:", t, i);
549  for (int j = 0; j < Spaces; j++) printf(" %3d", LineageExtant[i][t][j]);
550  printf("\n");
551  }
552  }
553  }
554 
555  printf("\n## StartSpaces:\n");
556  printf(" Spaces :");
557  for (int i = 0; i < Spaces; i++) printf(" %3d", i);
558  printf("\n");
559  printf(" Allowed :");
560  for (int j = 0; j < Spaces; j++) printf(" %3d", Cfg.startSpace[j]);
561  printf("\n");
562 
563  printf("\n## Dispersal Function:\n");
564  printf(" 0.0 0.5 1.0\n");
565  printf(
566  " +----+----+----+----+----+----+----+----+----+----+-> y (pDisp)\n");
567  for (int i = 0; i <= 10; i++) {
568  double d = pDispersal((double) i / 10.0);
569  if (!(i % 2)) printf(" %3.1f +", (double) i / 10.0);
570  else printf(" |");
571  for (int j = 0; j < ((int) (d * 50.0)); j++) printf("*");
572  printf("\n");
573  }
574  printf(" |\n");
575  printf(" v x (proportion of max distance)\n");
576 
577 
578  printf("\n## Survival Function:\n");
579  printf(" 0.0 0.5 1.0\n");
580  printf(
581  " +----+----+----+----+----+----+----+----+----+----+-> y (pSurv)\n");
582  for (int i = 0; i <= 10; i++) {
583  double s = pSurvival((double) i / 10.0 );
584  if (!(i % 2)) printf(" %3.1f +", (double) i / 10.0);
585  else printf(" |");
586  for (int j = 0; j < ((int) (s * 50.0)); j++) printf("*");
587  printf("\n");
588  }
589  printf(" |\n");
590  printf(" v x (proportion of total land area)\n");
591 
592 
593  printf("\n## Config:\n");
594  printf(" nStartSpaces : %10d\n", Cfg.nStartSpaces);
595  printf(" maxRuns : %10d\n", Cfg.maxRuns);
596  printf(" stopAfterSuccess : %10d\n", Cfg.stopAfterSuccess);
597  printf(" probSurvA : %10.2f\n", Cfg.probSurvA);
598  printf(" probSurvB : %10.2f\n", Cfg.probSurvB);
599  printf(" probDispA : %10.2f\n", Cfg.probDispA);
600  printf(" probDispB : %10.2f\n", Cfg.probDispB);
601 
602  printf("\n");
603 
604 }
605 
612 void error(char *error_msg)
613 {
614  printf("ERROR:\n %s\n", error_msg);
615  exit(1);
616  // NB: In the tcsh, `echo $status` returns the exit status of
617  // the preceding command.
618 }
619 
620 
628 double* mem1d_d(int dimx)
629 {
630  double *ptr;
631  ptr = calloc(dimx, sizeof(double));
632  if (!ptr) error("allocation failure in mem1d_d()");
633  return ptr;
634 }
635 
643 int* mem1d_i(int dimx)
644 {
645  int *ptr;
646  ptr = calloc(dimx, sizeof(int));
647  if (!ptr) error("allocation failure in mem1d_i()");
648  return ptr;
649 }
650 
658 int** mem2d1_i(int dimx)
659 {
660  int **ptr;
661  ptr = calloc(dimx, sizeof(int *));
662  if (!ptr) error("allocation failure in mem2d1_i()");
663  return ptr;
664 }
665 
670 void free2d1_i(int **ptr, int dimx)
671 {
672  for(int i = 0; i < dimx; i++)
673  free(ptr[i]);
674  free(ptr);
675 }
676 
685 char** mem2d1_c(int dimx)
686 {
687  char **ptr;
688  ptr = calloc(dimx, sizeof(char *));
689  if (!ptr) error("allocation failure in mem2d1_c()");
690  return ptr;
691 
692  // Direct, non-function alternative, e.g.:
693  // SpaceName = (char**) malloc( Spaces * sizeof(char*));
694 }
695 
700 void free2d1_c(char **ptr, int dimx)
701 {
702  for(int i = 0; i < dimx; i++)
703  free(ptr[i]);
704  free(ptr);
705 }
706 
707 
712 int** mem2d_i(int dimx, int dimy)
713 {
714  int **ptr;
715  ptr = calloc(dimx, sizeof(int *));
716  if (!ptr) error("allocation failure in mem2d_i() pass 1");
717  for(int i = 0; i < dimx; i++)
718  {
719  ptr[i] = calloc(dimy, sizeof(int));
720  if (!ptr[i]) error("allocation failure in mem2d_i() pass 2");
721  }
722  return ptr;
723 }
724 
725 void free2d_i(int **ptr, int dimx)
726 {
727  for(int i = 0; i < dimx; i++)
728  free(ptr[i]);
729  free(ptr);
730 }
731 
736 double** mem2d_d(int dimx, int dimy)
737 {
738  double **ptr;
739  ptr = calloc(dimx, sizeof(double *));
740  if (!ptr) error("allocation failure in mem2d_d() pass 1");
741  for(int i = 0; i < dimx; i++)
742  {
743  ptr[i] = calloc(dimy, sizeof(double));
744  if (!ptr[i]) error("allocation failure in mem2d_d() pass 2");
745  }
746  return ptr;
747 }
748 
749 void free2d_d(double **ptr, int dimx)
750 {
751  for(int i = 0; i < dimx; i++)
752  free(ptr[i]);
753  free(ptr);
754 }
755 
760 double*** mem3d_d(int dimx, int dimy, int dimz)
761 {
762  double ***ptr;
763  ptr = calloc(dimx, sizeof(double **));
764  if (!ptr) error("allocation failure in mem3d_d() pass 1");
765  for(int i = 0; i < dimx; i++)
766  {
767  ptr[i] = calloc(dimy, sizeof(double *));
768  if (!ptr[i]) error("allocation failure in mem3d_d() pass 2");
769  for(int j = 0; j < dimy; j++)
770  {
771  ptr[i][j] = calloc(dimz, sizeof(double));
772  if (!ptr[i][j]) error("allocation failure in mem3d_d() pass 3");
773  }
774  }
775  return ptr;
776 }
777 
778 void free3d_d(double ***ptr, int dimx, int dimy)
779 {
780  for(int i = 0; i < dimx; i++)
781  {
782  for(int j = 0; j < dimy; j++)
783  free(ptr[i][j]);
784  free(ptr[i]);
785  }
786  free(ptr);
787 }
788 
793 int*** mem3d_i(int dimx, int dimy, int dimz)
794 {
795  int ***ptr;
796  ptr = calloc(dimx, sizeof(int **));
797  if (!ptr) error("allocation failure in mem3d_d() pass 1");
798  for(int i = 0; i < dimx; i++)
799  {
800  ptr[i] = calloc(dimy, sizeof(int *));
801  if (!ptr[i]) error("allocation failure in mem3d_d() pass 2");
802  for(int j = 0; j < dimy; j++)
803  {
804  ptr[i][j] = calloc(dimz, sizeof(int));
805  if (!ptr[i][j]) error("allocation failure in mem3d_d() pass 3");
806  }
807  }
808  return ptr;
809 }
810 
811 void free3d_i(int ***ptr, int dimx, int dimy)
812 {
813  for(int i = 0; i < dimx; i++)
814  {
815  for(int j = 0; j < dimy; j++)
816  free(ptr[i][j]);
817  free(ptr[i]);
818  }
819  free(ptr);
820 }
821 
823 
824 void help()
825 {
826  char *version = VERSION;
827  printf("\n SHIBA (Simulated Historical Island Biogeography Analysis)\n");
828  printf(" (c) Campbell Webb 2013; Version %s\n", version);
829  printf(" Usage: shiba [ -dfhlpsv ]\n");
830  printf(" Options:\n");
831  printf(" -d NUM Use this value as the prob. of dispersal (probDispA)\n");
832  printf(" -f FILE Use this file as input. Default file: shibaInput.xml\n");
833  printf(" -h Print this help list\n");
834  printf(" -l Print input data summary\n");
835  printf(" -p INT Use this phylogeny (0...n). Default: 0\n");
836  printf(" -s NUM Use this value as the prob. of survival (probSurvA)\n");
837  printf(" -v Print highly verbose event descriptions\n\n");
838  exit(0);
839 }
840 
Generated on Wed Jan 16 2013 16:19:13 for shiba by   doxygen 1.8.2