+void neworder_cladewise(int *n, int *e1, int *e2, int *N, int *neworder)
+/* n: nb of tips
+ m: nb of nodes
+ N: nb of edges */
+{
+ int i, j, k, *L, *pos, m = *N - *n + 1, degrmax = *n - m + 1;
+
+/* degrmax is the largest value that a node degree can be */
+
+/* L is a 1-d array storing, for each node, the C indices of the rows of
+ the edge matrix where the node is ancestor (i.e., present in the 1st
+ column). It is used in the same way than a matrix (which is actually
+ a vector) is used in R as a 2-d structure. */
+
+ L = (int*)R_alloc(m * degrmax, sizeof(int));
+
+/* pos gives the position for each 'row' of L, that is the number of elements
+ which have already been stored for that 'row'. */
+
+ pos = (int*)R_alloc(m, sizeof(int));
+ memset(pos, 0, m * sizeof(int));
+
+/* we now go down along the edge matrix */
+
+ for (i = 0; i < *N; i++) {
+ k = e1[i] - *n - 1; /* k is the 'row' index in L corresponding to node e1[i] */
+ j = pos[k]; /* the current 'column' position corresping to k */
+ pos[k]++; /* increment in case the same node is found in another row of the edge matrix */
+ L[k + m * j] = i;
+ }
+
+/* L is now ready: we can start the recursive calls. */
+/* We start with the root 'n + 1': its index will be changed into
+ the corresponding C index inside the recursive function. */
+
+ iii = 0;
+ foo_reorder(*n + 1, *n, m, e1, e2, neworder, L, pos);
+}
+
+