some big fixes for ape 2.7-1

[ape.git] / src / rTrait.c

diff --git a/src/rTrait.c b/src/rTrait.c
index 6441e161b9c07413bb97e2cb4cdba48edce20788..11cc63e2a0bdb2a66105873e690ab1e03d115074 100644 (file)
--- a/src/rTrait.c
+++ b/src/rTrait.c
@@ -1,4 +1,4 @@
-/* rTrait.c       2010-01-11 */
+/* rTrait.c       2010-07-26 */
 
 /* Copyright 2010 Emmanuel Paradis */
 
@@ -16,13 +16,19 @@ void rTraitCont(int *model, int *Nedge, int *edge1, int *edge2, double *el,
        switch(*model) {
        case 1 : for (i = *Nedge - 1; i >= 0; i--) {
                        GetRNGstate();
-                       x[edge2[i]] = x[edge1[i]] + el[i] * sigma[i] * norm_rand();
+                       x[edge2[i]] = x[edge1[i]] + sqrt(el[i]) * sigma[i] * norm_rand();
                        PutRNGstate();
                }
                break;
        case 2 : for (i = *Nedge - 1; i >= 0; i--) {
                        GetRNGstate();
-                       x[edge2[i]] = x[edge1[i]] + (sigma[i]*norm_rand() - alpha[i]*(x[edge1[i]] - theta[i])) * el[i];
+                       x[edge2[i]] = x[edge1[i]] - alpha[i]*el[i]*(x[edge1[i]] - theta[i]) + sqrt(el[i])*sigma[i]*norm_rand();
+                       PutRNGstate();
+               }
+               break;
+       case 3 : for (i = *Nedge - 1; i >= 0; i--) {
+                       GetRNGstate();
+                       x[edge2[i]] = x[edge1[i]] - (1 - exp(alpha[i]*el[i]))*(x[edge1[i]] - theta[i]) + sqrt(el[i])*sigma[i]*norm_rand();
                        PutRNGstate();
                }
                break;