ts/wf/teseco [FEATURE]: provide new tolerance option -trim

src/ts/wf/teseco.cc

@@ -4,7 +4,7 @@
  * ----------------------------------------------------------------------------
  * 
  * \author Thomas Forbriger
- * \date 28/11/2016
+ * \date 29/11/2016
  * 
  * time series corrections
  * 
@@ -32,11 +32,12 @@
  *  - 30/07/2015   V1.2   provide option -datetolerance
  *  - 28/11/2016   V1.3   allow scaling of signal
  *  - 29/11/2016   V1.4   fix: make use of datetolerance
- * 
+ *                 V1.5   provide new tolerance option -trim
+ *
  * ============================================================================
  */
 #define TESECO_VERSION \
-  "TESECO   V1.4   time series corrections/ linear combination of signals"
+  "TESECO   V1.5   time series corrections/ linear combination of signals"
 
 #include <fstream>
 #include <sstream>
@@ -59,7 +60,7 @@ using std::endl;
 typedef aff::Series<double> Tseries;
 
 struct Options {
-  bool verbose, overwrite, add, debug;
+  bool verbose, overwrite, add, debug, trimseries;
   double datetolerance;
   std::string inputformat, outputformat;
 }; // struct Options
@@ -72,7 +73,7 @@ int main(int iargc, char* argv[])
   {
     TESECO_VERSION "\n"
     "usage: teseco [-v] [-o] [-a] [-type type] [-Type type]" "\n"
-    "              [-datetolerance t]" "\n"
+    "              [-datetolerance t] [-trim]" "\n"
     "              outfile signal [t:T] [f:F]\n"
     "              infile [t:T] [f:F] [infile [t:T] [f:F] ... ]" "\n"
     "   or: teseco --help|-h" "\n"
@@ -107,6 +108,9 @@ int main(int iargc, char* argv[])
     "             input time series; dates will be considered as different\n"
     "             if they differ by more than a fraction f of the average\n"
     "             sampling interval\n"
+    "-trim        trim all series to the shortest series provided\n"
+    "             if not set, providing series of different number of samples\n"
+    "             is considered to be an error\n"
   };
 
   // define commandline options
@@ -131,6 +135,8 @@ int main(int iargc, char* argv[])
     {"xhelp",arg_no,"-"},
     // 8: report details regarding data types
     {"datetolerance",arg_yes,"0."},
+    // 9: trim size of series to the least number of samples
+    {"trim",arg_no,"-"},
     {NULL}
   };
 
@@ -181,6 +187,7 @@ int main(int iargc, char* argv[])
   opt.debug=cmdline.optset(5);
   opt.outputformat=cmdline.string_arg(6);
   opt.datetolerance=cmdline.double_arg(8);
+  opt.trimseries=cmdline.optset(9);
 
   if (opt.verbose)
   { cout << TESECO_VERSION << endl; }
@@ -211,6 +218,10 @@ int main(int iargc, char* argv[])
   std::ofstream ofs(outfile.c_str(), oopenmode);
   datrw::oanystream os(ofs, opt.outputformat);
 
+  // prepare the largest number of samples found in one of the input series so
+  // far
+  unsigned int maxsamples=0;
+
   // prepare file FREE block
   sff::FREE filefree;
   filefree.append(TESECO_VERSION);
@@ -237,6 +248,9 @@ int main(int iargc, char* argv[])
     else
     { signal.read(is.idatstream(), opt.verbose); }
 
+    TFXX_assert(signal.size()>0,
+                "set of input signals is empty");
+
     // apply factor to each of the traces, if selected
     if (sigfile.haskey(factorkey))
     {
@@ -248,14 +262,29 @@ int main(int iargc, char* argv[])
       {
         *isig *= factor;
       }
-    }
+    } // if (sigfile.haskey(factorkey))
+
+    // set max samples if series are to be trimmed
+    if (opt.trimseries)
+    {
+      for (Tfile::Ttracevector::iterator isig=signal.begin();
+           isig != signal.end(); ++isig)
+      {
+        TFXX_assert(isig->size()>0, "missing samples in input signal");
+        if ((maxsamples == 0) || (maxsamples > isig->size()))
+        {
+          maxsamples = isig->size();
+        }
+      }
+    } // if (opt.trimseries)
+
   }
   TFXX_assert((signal.size()>0), "missing signal");
   signal.fileheader.append(filefree);
 
 
   // SFF WID2 compare
-  sff::WID2compare compare(sff::Fnsamples | sff::Fdt | sff::Fdate);
+  sff::WID2compare compare(sff::Fdt | sff::Fdate);
   compare.setdatetolerance(opt.datetolerance);
 
   // cycle through all input files
@@ -306,6 +335,7 @@ int main(int iargc, char* argv[])
         for (Tfile::Ttracevector::iterator isig=signal.begin();
              isig != signal.end(); ++isig)
         {
+          // check matching headers
           if (!compare (isig->header.wid2(),wid2))
           {
             cerr << "ERROR: header signature mismatch:" << endl;
@@ -313,6 +343,8 @@ int main(int iargc, char* argv[])
             cerr << isig->header.wid2().line();
             cerr << "correction signal:" << endl;
             cerr << wid2.line();
+
+            // check start date explicitly and provide a hint
             sff::WID2compare cmpdate(sff::Fdate);
             cmpdate.setdatetolerance(opt.datetolerance);
             if (!cmpdate(isig->header.wid2(),wid2))
@@ -325,7 +357,32 @@ int main(int iargc, char* argv[])
               cerr << wid2.date.timestring() << endl;
             }
             TFXX_abort("bailing out...");
-          }
+          } // if (!compare (isig->header.wid2(),wid2))
+
+          // check matching size
+          if (isig->size() != series.size())
+          {
+            if (opt.trimseries)
+            {
+              if (maxsamples > isig->size()) { maxsamples = isig->size(); }
+              if (maxsamples > series.size()) { maxsamples = series.size(); }
+              isig->setlastindex(isig->first()-1+maxsamples);
+              series.setlastindex(series.first()-1+maxsamples);
+            }
+            else
+            {
+              cerr << "ERROR: inconsistent number of samples:" << endl;
+              cerr << isig->size() << " samples for signal" << endl;
+              cerr << isig->header.wid2().line();
+              cerr << series.size() << " samples for correction signal" 
+                << endl;
+              cerr << wid2.line();
+              cerr << "consider to use option -trim" << endl;
+              TFXX_abort("bailing out...");
+            } // else, if (opt.trimseries)
+          } // (isig->size() != series.size())
+
+          // apply correction
           *isig -= series;
         }