fix small bug due to recent update

tothuhien · Oct 27, 2021 · 26ba127 · 26ba127
1 parent 590d761
commit 26ba127
Show file tree

Hide file tree

Showing 5 changed files with 32 additions and 6 deletions.
diff --git a/bin/lsd2_mac b/bin/lsd2_mac
diff --git a/bin/lsd2_unix b/bin/lsd2_unix
diff --git a/src/confidence_interval.cpp b/src/confidence_interval.cpp
@@ -819,12 +819,16 @@ void output(double br,int y, InputOutputStream *io, Pr* pr,Node** nodes,ostream&
 
         tree2<<"tree "<<y<<" = ";
         if (!diffTopology){
-            for (int j=0;j<=pr->nbBranches;j++){
-                nodes[j]->D = (nodes[j]->D * pr->rho)*(1 / pr->rhoLower + 1/ pr->rhoUpper)/2;
+            if (pr->rhoLower!=0 && pr->rhoUpper!=0){
+                for (int j=0;j<=pr->nbBranches;j++){
+                    nodes[j]->D = (nodes[j]->D * pr->rho)*(1 / pr->rhoLower + 1/ pr->rhoUpper)/2;
+                }    
             }
             tree2<<nexusICDate(0,pr,nodes,T_min,T_max,HD_min,HD_max).c_str();
         } else{
-            nodes[0]->D = (nodes[0]->D * pr->rho)*(1 / pr->rhoLower + 1 / pr->rhoUpper)/2;
+            if (pr->rhoLower!=0 && pr->rhoUpper!=0){
+                nodes[0]->D = (nodes[0]->D * pr->rho)*(1 / pr->rhoLower + 1 / pr->rhoUpper)/2;
+            }
             tree2<<nexusDate(0,pr,nodes).c_str();
         }
         int n=0;

diff --git a/src/options.cpp b/src/options.cpp
@@ -25,7 +25,7 @@ Pr* getOptions( int argc, char** argv )
 
 Pr* getCommandLine( int argc, char** argv)
 {
-    const string VERSION="v.2.2";
+    const string VERSION="v.2.3";
     Pr* opt = new Pr();
     int c;
     string s;
@@ -354,7 +354,7 @@ Pr* getInterface()
 
 void printInterface(ostream& in, Pr* opt)
 {
-    const string VERSION = "v.2.2";
+    const string VERSION = "v.2.3";
 
     in<<"\nLEAST-SQUARE METHODS TO ESTIMATE RATES AND DATES - "<<VERSION<<" \n\n";
     in<<"\nInput files:\n";
@@ -511,7 +511,7 @@ void printHelp( void )
     const string BOLD = "\033[00;01m";
     const string LINE = "\033[00;04m";
     const string FLAT = "\033[00;00m";
-    const string VERSION = "v.2.2";
+    const string VERSION = "v.2.3";
 
     cout<<BOLD<<"LSD: LEAST-SQUARES METHODS TO ESTIMATE RATES AND DATES - "<<VERSION<<"\n\n";
     cout<<BOLD<<"DESCRIPTION\n"

diff --git a/src/readData.cpp b/src/readData.cpp
@@ -136,6 +136,7 @@ void readInputDate(InputOutputStream* io, Pr* pr,Node** &nodes,bool& constraintC
     bool* tipHaveTime = new bool[pr->nbBranches+1];
     int nbUniqueDate = 0;
     double uniqueDate = 0;
+    double minINodeDate = 0;
     for (int i=pr->nbINodes; i<= pr->nbBranches; i++) tipHaveTime[i]=false;
     if (io->inDate){
         io->inDate->seekg(0);
@@ -216,6 +217,14 @@ void readInputDate(InputOutputStream* io, Pr* pr,Node** &nodes,bool& constraintC
                         pr->haveUnique = true;
                     }
                 }
+                if (k<pr->nbINodes){
+                    if (type=='p' || type=='l'){
+                        if (minINodeDate > v1) minINodeDate = v1;
+                    } else if (type == 'b'){
+                        if (minINodeDate > v1) minINodeDate = v1;
+                        if (minINodeDate > v2) minINodeDate = v2;
+                    }
+                }
             }
             else{
                 w1=w1+" "+s;
@@ -293,6 +302,19 @@ void readInputDate(InputOutputStream* io, Pr* pr,Node** &nodes,bool& constraintC
             }
         }
     }
+    if (nbUniqueDate==0){
+        std::ostringstream oss;
+        oss<<"- There's not any input tip dates, so all tip dates are automatically set to 0.\n";
+        for (int i=pr->nbINodes;i<=pr->nbBranches;i++){
+            nodes[i]->type = 'p';
+            nodes[i]->D = 0;
+        }
+        if (minINodeDate <0){
+            cerr<<"Error: ancestral dates have to be smaller than tips date (0). Hence, for example, set -20 for a calibration of 20 my ago."<<endl;
+            exit(EXIT_FAILURE);
+        }
+        pr->warningMessage.push_back(oss.str());
+    }
     if (pr->inDateFormat == 2 && dateFormat != 2){//the general day format is year-month-day but there're some only year or year-month
         int ii = 0;
         for (int i=0; i<kSave.size();i++){