X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=src%2Fme.c;h=8a00591ab680cb6846fc9ea7765c4a9b7e893c0f;hb=199144b0297c3fc76d76c29e561151235e39f0af;hp=5b383f33b2c685f0fd8887ada72514bb94e2b09c;hpb=c827059eeafc8cbe41c812b26979543ab287803e;p=ape.git diff --git a/src/me.c b/src/me.c index 5b383f3..8a00591 100644 --- a/src/me.c +++ b/src/me.c @@ -1,26 +1,39 @@ +/* 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; @@ -29,53 +42,37 @@ void me_b(double *X, int *N, char **labels, char **treeStr, int *nni) 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; @@ -84,11 +81,6 @@ void me_o(double *X, int *N, char **labels, char **treeStr, int *nni) 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); @@ -100,37 +92,24 @@ void me_o(double *X, int *N, char **labels, char **treeStr, int *nni) } 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) { @@ -146,9 +125,10 @@ 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; @@ -159,9 +139,10 @@ double **loadMatrix (double *X, char **labels, int n, set *S) for(i=0; iindex2 = i; S = addToSet(v,S); for (j=i; jlabel,label,NODE_LABEL_LENGTH); +// strncpy(newNode->label,label,NODE_LABEL_LENGTH); + newNode->label = label; newNode->index = index; newNode->index2 = -1; newNode->parentEdge = parentEdge; @@ -295,8 +279,7 @@ tree *detrifurcate(tree *T) 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 ) { @@ -325,7 +308,8 @@ void compareSets(tree *T, set *S) 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; @@ -338,7 +322,8 @@ void compareSets(tree *T, set *S) 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; @@ -347,8 +332,7 @@ void compareSets(tree *T, set *S) } 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) @@ -356,16 +340,14 @@ void compareSets(tree *T, set *S) 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; } @@ -516,4 +498,3 @@ void freeSubTree(edge *e) e->head = NULL; free(e); } -