+/* me.c 2012-04-30 */
+
+/* Copyright 2007-2008 Olivier Gascuel, Rick Desper,
+ R port by Vincent Lefort and Emmanuel Paradis */
+
+/* This file is part of the R-package `ape'. */
+/* See the file ../COPYING for licensing issues. */
+
#include "me.h"
//functions from me_balanced.c
-tree *BMEaddSpecies(tree *T,node *v, double **D, double **A);
+tree *BMEaddSpecies(tree *T, node *v, double **D, double **A);
void assignBMEWeights(tree *T, double **A);
void makeBMEAveragesTable(tree *T, double **D, double **A);
//functions from me_ols.c
-tree *GMEaddSpecies(tree *T,node *v, double **D, double **A);
+tree *GMEaddSpecies(tree *T, node *v, double **D, double **A);
void assignOLSWeights(tree *T, double **A);
void makeOLSAveragesTable(tree *T, double **D, double **A);
//functions from bNNI.c
void bNNI(tree *T, double **avgDistArray, int *count, double **D, int numSpecies);
//functions from NNI.c
void NNI(tree *T, double **avgDistArray, int *count, double **D, int numSpecies);
+//functions from SPR.c
+void SPR(tree *T, double **D, double **A, int *count);
+//functions from TBR.c
+void TBR(tree *T, double **D, double **A);
-void me_b(double *X, int *N, char **labels, char **treeStr, int *nni)
+void me_b(double *X, int *N, int *labels,
+ int *nni, int *spr, int *tbr,
+ int *edge1, int *edge2, double *el)
{
double **D, **A;
set *species, *slooper;
node *addNode;
tree *T;
- char *str;
int n, nniCount;
n = *N;
species = (set *) malloc(sizeof(set));
species->firstNode = NULL;
species->secondNode = NULL;
- str = (char *)R_alloc(MAX_INPUT_SIZE, sizeof(char));
- /* added by EP */
- if (strlen(str))
- strncpy(str, "", strlen(str));
- /* end */
- D = loadMatrix (X, labels, n, species);
- A = initDoubleMatrix(2 * n - 2);
+ D = loadMatrix(X, labels, n, species);
+ A = initDoubleMatrix(2*n - 2);
for(slooper = species; NULL != slooper; slooper = slooper->secondNode)
{
addNode = copyNode(slooper->firstNode);
- T = BMEaddSpecies(T,addNode,D,A);
+ T = BMEaddSpecies(T, addNode, D, A);
}
// Compute bNNI
- if (*nni == 1)
- bNNI(T,A,&nniCount,D,n);
+ if (*nni) bNNI(T, A, &nniCount, D, n);
assignBMEWeights(T,A);
- NewickPrintTreeStr(T,str);
+ if (*spr) SPR(T, D, A, &nniCount);
+ if (*tbr) TBR(T, D, A);
- if (strlen (str) < MAX_INPUT_SIZE -1)
- {
- *treeStr = (char *)R_alloc(MAX_INPUT_SIZE, sizeof(char));
- /* added by EP */
- if (strlen(*treeStr))
- strncpy(*treeStr, "", strlen(*treeStr));
- /* end */
- strncpy (*treeStr, str, strlen(str));
- }
+ tree2phylo(T, edge1, edge2, el, labels, n);
-/* free (str); */
freeMatrix(D,n);
freeMatrix(A,2*n - 2);
freeSet(species);
freeTree(T);
T = NULL;
-
- /* return; */
}
-void me_o(double *X, int *N, char **labels, char **treeStr, int *nni)
+void me_o(double *X, int *N, int *labels, int *nni,
+ int *edge1, int *edge2, double *el)
{
double **D, **A;
set *species, *slooper;
node *addNode;
tree *T;
- char *str;
int n, nniCount;
n = *N;
species = (set *) malloc(sizeof(set));
species->firstNode = NULL;
species->secondNode = NULL;
- str = (char *)R_alloc(MAX_INPUT_SIZE, sizeof(char));
- /* added by EP */
- if (strlen(str))
- strncpy(str, "", strlen(str));
- /* end */
D = loadMatrix (X, labels, n, species);
A = initDoubleMatrix(2 * n - 2);
}
makeOLSAveragesTable(T,D,A);
// Compute NNI
- if (*nni == 1)
+ if (*nni)
NNI(T,A,&nniCount,D,n);
assignOLSWeights(T,A);
- NewickPrintTreeStr(T,str);
-
- if (strlen (str) < MAX_INPUT_SIZE -1)
- {
- *treeStr = (char *)R_alloc(MAX_INPUT_SIZE, sizeof(char));
- /* added by EP */
- if (strlen(*treeStr))
- strncpy(*treeStr, "", strlen(*treeStr));
- /* end */
- strncpy (*treeStr, str, strlen (str));
- }
+ tree2phylo(T, edge1, edge2, el, labels, n);
- /* free (str); */
freeMatrix(D,n);
freeMatrix(A,2*n - 2);
freeSet(species);
freeTree(T);
T = NULL;
-
- return;
}
-/*************************************************************************
+/*
- MATRIX FUNCTIONS
+ -- MATRIX FUNCTIONS --
-*************************************************************************/
+*/
double **initDoubleMatrix(int d)
{
return(A);
}
-double **loadMatrix (double *X, char **labels, int n, set *S)
+//double **loadMatrix (double *X, char **labels, int n, set *S)
+double **loadMatrix (double *X, int *labels, int n, set *S)
{
- char nextString[MAX_LABEL_LENGTH];
+// char nextString[MAX_LABEL_LENGTH];
node *v;
double **table;
int i, j, a, b;
for(i=0; i<n; i++)
{
- strncpy (nextString, labels[i], MAX_LABEL_LENGTH);
+// strncpy (nextString, labels[i], MAX_LABEL_LENGTH);
// ReplaceForbiddenChars (nextString, '_');
- v = makeNewNode(nextString,-1);
+// v = makeNewNode(nextString,-1);
+ v = makeNewNode(labels[i], -1);
v->index2 = i;
S = addToSet(v,S);
for (j=i; j<n; j++) {
return (table);
}
-/*************************************************************************
+/*
- GRAPH FUNCTIONS
+ -- GRAPH FUNCTIONS --
-*************************************************************************/
+*/
set *addToSet(node *v, set *X)
{
return(X);
}
-node *makeNewNode(char *label, int i)
+//node *makeNewNode(char *label, int i)
+node *makeNewNode(int label, int i)
{
return(makeNode(label,NULL,i));
}
-node *makeNode(char *label, edge *parentEdge, int index)
+//node *makeNode(char *label, edge *parentEdge, int index)
+node *makeNode(int label, edge *parentEdge, int index)
{
node *newNode; /*points to new node added to the graph*/
newNode = (node *) malloc(sizeof(node));
- strncpy(newNode->label,label,NODE_LABEL_LENGTH);
+// strncpy(newNode->label,label,NODE_LABEL_LENGTH);
+ newNode->label = label;
newNode->index = index;
newNode->index2 = -1;
newNode->parentEdge = parentEdge;
return(T);
if (NULL != v->parentEdge)
{
- Rprintf ("Error: root %s is poorly rooted.\n",v->label);
- exit(0);
+ error("root %d is poorly rooted.", v->label);
}
for(e = v->middleEdge, v->middleEdge = NULL; NULL != e; e = f )
{
for(X = S; NULL != X; X = X->secondNode)
{
w = X->firstNode;
- if (0 == strcmp(v->label,w->label))
+// if (0 == strcmp(v->label,w->label))
+ if (v->label == w->label)
{
v->index2 = w->index2;
w->index2 = -1;
for(X = S; NULL != X; X = X->secondNode)
{
w = X->firstNode;
- if (0 == strcmp(v->label,w->label))
+// if (0 == strcmp(v->label,w->label))
+ if (v->label == w->label)
{
v->index2 = w->index2;
w->index2 = -1;
}
if (-1 == v->index2)
{
- Rprintf("Error leaf %s in tree not in distance matrix.\n",v->label);
- exit(0);
+ error("leaf %d in tree not in distance matrix.", v->label);
}
e = depthFirstTraverse(T,NULL);
while (NULL != e)
v = e->head;
if ((leaf(v)) && (-1 == v->index2))
{
- Rprintf("Error leaf %s in tree not in distance matrix.\n",v->label);
- exit(0);
+ error("leaf %d in tree not in distance matrix.", v->label);
}
e = depthFirstTraverse(T,e);
}
for(X = S; NULL != X; X = X->secondNode)
if (X->firstNode->index2 > -1)
{
- Rprintf("Error node %s in matrix but not a leaf in tree.\n",X->firstNode->label);
- exit(0);
+ error("node %d in matrix but not a leaf in tree.", X->firstNode->label);
}
return;
}
e->head = NULL;
free(e);
}
-