From b1705525aeb8e7d285bbfdba2d8fd258d973bc54 Mon Sep 17 00:00:00 2001
From: Klaus Rabbertz <klaus.rabbertz@cern.ch>
Date: Wed, 19 May 2021 22:12:16 +0200
Subject: [PATCH] Read also stat. uncertainties from fnlo-stat-unc (NLOJet++)

---
 v2.3/toolkit/README                           |  4 +-
 v2.3/toolkit/TODO                             |  7 +-
 v2.3/toolkit/fastnlotoolkit/fastNLOTools.cc   | 65 ++++++++++++++-----
 .../include/fastnlotk/fastNLOConstants.h.in   |  8 +++
 v2.3/toolkit/src/fnlo-tk-statunc.cc           | 39 ++---------
 5 files changed, 68 insertions(+), 55 deletions(-)

diff --git a/v2.3/toolkit/README b/v2.3/toolkit/README
index 2f4287d8..607317cd 100644
--- a/v2.3/toolkit/README
+++ b/v2.3/toolkit/README
@@ -1,12 +1,12 @@
 # -*-sh-*-
-# fastNLO_toolkit version 2.3.1:
+# fastNLO_toolkit version 2.5.0:
 # ==============================
 # Please contact the authors for any questions or problems.
 
 
 # Content:
 # --------
-# Generic C++ linkable library and code to create, fill, read, and evaluate fastNLO v2 tables
+# Generic C++ linkable library and code to create, fill, modify, read, and evaluate fastNLO v2 tables
 
 
 # Requirements:
diff --git a/v2.3/toolkit/TODO b/v2.3/toolkit/TODO
index 6c2a2ec6..cd423536 100644
--- a/v2.3/toolkit/TODO
+++ b/v2.3/toolkit/TODO
@@ -2,20 +2,18 @@
 TODO list for fastNLO_toolkit (ordered by priority)
 -----------------------------------------------------------------------
 
+- Creator for Referenztabellen
 - Add checks for HOPPET, CRunDec, QCDNUM
 - Change everything to use LHAPDF6 features instead of deprecated LHAGLUE;
   partially done in HOPPET interface
 - Drop LHAPDF5 support
 - Use binary read/write table to save storage space and read in times (Mark's suggestion)
+  (Not clear that binary storage actually saves anything ...! Test first.)
 - Provide features to allow -> APPLgrid conversion and the inverse -> fastNLO
 - Test installation with lgcenv environment
    - Problem with zlib?
    - Problem with HOPPET?
- - Include stat. uncertainties into table
- - Creator for Referenztabellen
  - Streamline usage of logger.messagelevel, shout|cout, error|warn|silent etc. in speaker
-
- - KR: Check, which points are actually done ...
  - Getter for PhaseSpace-check: a la: GetIsWithinPhaseSpace(double obs1, double obs2, etc.) (cache probably observables for later usage)
  - Write documentation (e.g. twiki)
  - converter aktualisieren ??
@@ -34,6 +32,7 @@ TODO list for fastNLO_toolkit (ordered by priority)
 =======================================================================
 DONE
 =======================================================================
+ - Include stat. uncertainties into table
  - Modifier & Concatenator for tables/bins
  - fast Warmup-run (just some events and then 'guess' PS boundaries)
  - 11x11 = 121 PDF-LiCos
diff --git a/v2.3/toolkit/fastnlotoolkit/fastNLOTools.cc b/v2.3/toolkit/fastnlotoolkit/fastNLOTools.cc
index 667b74f7..fea7ffad 100644
--- a/v2.3/toolkit/fastnlotoolkit/fastNLOTools.cc
+++ b/v2.3/toolkit/fastnlotoolkit/fastNLOTools.cc
@@ -393,35 +393,70 @@ namespace fastNLOTools {
 
    //______________________________________________________________________________
    std::vector <double> ReadInfoBlockContent(std::string filename) {
+      std::string extension = "";
       std::ifstream infile;
       std::string line;
       std::vector <double> Uncertainty;
+
+      //! Determine extension to differentiate for parsing among
+      //  - NLOJet++ --> fnlo-tk-statunc .log files
+      //  - NNLOJET  --> .dat files
+      if ( filename.find_last_of(".") != std::string::npos ) {
+         extension = filename.substr(filename.find_last_of(".")+1);
+      }
       info["ReadInfoBlockContent"]<<"Reading additional InfoBlock content from file: " << filename <<endl;
       infile.open(filename);
       if (infile.is_open()) {
-         int iline = 0;
+         int  iline = 0;
+         bool lline = false;
          // Read line-by-line
          while(std::getline(infile, line)) {
+            cout << "line = " << line << endl;
             // Put line into stringstream and read word-by-word
             std::istringstream iss(line);
-            std::string word;
+            std::string word, word1, word2;
             iss >> word;
-            // For now assume NNLOJET dat file format:
+            // For 'dat' extension assume NNLOJET dat file format:
             // - Skip all lines starting with comment symbol '#'
             // - Read cross section and absolute statistical uncertainty from 4th and 5th columns
-            if ( word.at(0) != '#' ) {
-               // Skip first three words of each line
-               iss >> word;
-               iss >> word;
-               double xs, dxs;
-               iss >> xs;
-               iss >> dxs;
-               if ( fabs(xs) > DBL_MIN ) {
-                  Uncertainty.push_back(dxs/xs);
-               } else {
-                  Uncertainty.push_back(0.);
+            if ( extension == "dat" ) {
+               if ( word.at(0) != '#' ) {
+                  // Skip first three words of each line
+                  iss >> word;
+                  iss >> word;
+                  double xs, dxs;
+                  iss >> xs;
+                  iss >> dxs;
+                  if ( fabs(xs) > DBL_MIN ) {
+                     Uncertainty.push_back(dxs/xs);
+                  } else {
+                     Uncertainty.push_back(0.);
+                  }
+                  iline += 1;
+               }
+            }
+            // For 'log' extension assume fnlo-tk-stat v2.5 log file format:
+            // (New v2.5 separator lines starting with #- - - - - - -)
+            // - Start at first line with "#-" as 1st word and
+            // - stop again at next line starting with "#-" in 1st word
+            else if ( extension == "log" ) {
+               if ( word == "#-" ) {
+                  lline = ! lline;
+               } else if ( lline ) {
+                  cout << "AAAAA word = " << word << endl;
+                  // Skip second & third word of each uncertainty line (x section; lower uncertainty)
+                  iss >> word;
+                  iss >> word;
+                  double dxsrel;
+                  iss >> dxsrel;
+                  cout << "dxsrel = " << dxsrel << endl;
+                  Uncertainty.push_back(dxsrel);
+                  iline += 1;
+                  cout << "BBBBB iline = " << iline << endl;
                }
-               iline += 1;
+            } else {
+               error["ReadInfoBlockContent"]<<"Unkown file extension, aborted! Filename is: " << filename <<endl;
+               exit(32);
             }
          }
       } else {
diff --git a/v2.3/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOConstants.h.in b/v2.3/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOConstants.h.in
index 92f51300..3493edfb 100644
--- a/v2.3/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOConstants.h.in
+++ b/v2.3/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOConstants.h.in
@@ -167,27 +167,35 @@ namespace fastNLO {
    const std::string _CSEP20("####################");
    const std::string _DSEP20("====================");
    const std::string _SSEP20("--------------------");
+   const std::string _TSEP20(" - - - - - - - - - -");
    const std::string _CSEP20C(" ######################");
    const std::string _DSEP20C(" #=====================");
    const std::string _SSEP20C(" #---------------------");
+   const std::string _TSEP20C(" #- - - - - - - - - - -");
    const std::string _CSEP40  = _CSEP20  + _CSEP20;
    const std::string _DSEP40  = _DSEP20  + _DSEP20;
    const std::string _SSEP40  = _SSEP20  + _SSEP20;
+   const std::string _TSEP40  = _TSEP20  + _TSEP20;
    const std::string _CSEP40C = _CSEP20C + _CSEP20;
    const std::string _DSEP40C = _DSEP20C + _DSEP20;
    const std::string _SSEP40C = _SSEP20C + _SSEP20;
+   const std::string _TSEP40C = _TSEP20C + _TSEP20;
    const std::string _CSEPS  = _CSEP40  + _CSEP40;
    const std::string _DSEPS  = _DSEP40  + _DSEP40;
    const std::string _SSEPS  = _SSEP40  + _SSEP40;
+   const std::string _TSEPS  = _TSEP40  + _TSEP40;
    const std::string _CSEPSC = _CSEP40C + _CSEP40;
    const std::string _DSEPSC = _DSEP40C + _DSEP40;
    const std::string _SSEPSC = _SSEP40C + _SSEP40;
+   const std::string _TSEPSC = _TSEP40C + _TSEP40;
    const std::string _CSEPL  = _CSEPS  + _CSEPS ;
    const std::string _DSEPL  = _DSEPS  + _DSEPS ;
    const std::string _SSEPL  = _SSEPS  + _SSEPS ;
+   const std::string _TSEPL  = _TSEPS  + _TSEPS ;
    const std::string _CSEPLC = _CSEPSC + _CSEPS ;
    const std::string _DSEPLC = _DSEPSC + _DSEPS ;
    const std::string _SSEPLC = _SSEPSC + _SSEPS ;
+   const std::string _TSEPLC = _TSEPSC + _TSEPS ;
 #endif
 }
 
diff --git a/v2.3/toolkit/src/fnlo-tk-statunc.cc b/v2.3/toolkit/src/fnlo-tk-statunc.cc
index 77ea0523..5ec74b88 100644
--- a/v2.3/toolkit/src/fnlo-tk-statunc.cc
+++ b/v2.3/toolkit/src/fnlo-tk-statunc.cc
@@ -108,7 +108,7 @@ int main(int argc, char** argv) {
          man << "   - Specify the PDF set including the absolute path." << endl;
          man << "   - Download the desired PDF set from the LHAPDF web site." << endl;
          man << "[order]: Fixed-order precision to use, def. = NLO" << endl;
-         man << "   Alternatives: LO, NNLO, NLO-ONLY, NNLO-ONLY (if available)" << endl;
+         man << "   Alternatives: LO, NNLO, NLO_only, NNLO_only (if available)" << endl;
          man << "[nmin]: Smallest table number nnnn to start with, def. = 0000." << endl;
          man << "[nmax]: Largest  table number nnnn to end with, def. = 1000." << endl;
          man << "[Verbosity]: Set verbosity level of table evaluation [DEBUG,INFO,WARNING,ERROR], def. = WARNING" << endl;
@@ -158,11 +158,11 @@ int main(int argc, char** argv) {
       } else if ( FixedOrderChoice == "NNLO" ) {
          eOrder = kNextToNextToLeading;
          shout["fnlo-tk-statunc"] << "Deriving NNLO cross sections for comparison." << endl;
-      } else if ( FixedOrderChoice == "NLO-ONLY" ) {
+      } else if ( FixedOrderChoice == "NLO_only" ) {
          eOrder = kNextToLeading;
          lexclusive = true;
          shout["fnlo-tk-statunc"] << "Deriving NLO contributions for comparison." << endl;
-      } else if ( FixedOrderChoice == "NNLO-ONLY" ) {
+      } else if ( FixedOrderChoice == "NNLO_only" ) {
          eOrder = kNextToNextToLeading;
          lexclusive = true;
          shout["fnlo-tk-statunc"] << "Deriving NNLO contributions for comparison." << endl;
@@ -365,35 +365,6 @@ int main(int argc, char** argv) {
       exit(1);
    }
 
-   // KR TEST MERGING
-
-   fastNLOTable* resultTable = NULL;
-   int nmerge = 0;
-   EMerge moption = kMerge;
-   for ( const string& path : alltables ) {
-      nmerge++;
-      if (nmerge == 1) {
-         resultTable = new fastNLOTable(path);
-      } else {
-         fastNLOTable tab(path);
-         resultTable->MergeTable(tab, moption); // merge
-      }
-   }
-   // --- Write result
-   resultTable->SetFilename("mergetest.tab.gz");
-   resultTable->SetITabVersionWrite(25000);
-   // Loop over contributions in table
-   for ( int ic=0; ic<resultTable->GetNcontrib()+resultTable->GetNdata(); ic++ ) {
-      fastNLOCoeffAddBase* coeff = (fastNLOCoeffAddBase*)resultTable->GetCoeffTable(ic);
-      int ib = coeff->GetNCoeffInfoBlocks();
-      cout<< "IIIIIIIII = " << ic << ", ib = " << ib <<endl;
-   }
-   //   AddCoeffInfoBlock
-   info["Blubber"]<<"Write merged results to file "<<resultTable->GetFilename()<<"."<<endl;
-   resultTable->WriteTable();
-
-   // KR TEST MERGING
-
    //! Write out statistical fluctuation info
    yell << _CSEPSC << endl;
    shout << "Special info on statistical fluctuations:" << endl;
@@ -475,7 +446,7 @@ int main(int argc, char** argv) {
    shout << "Relative Statistical Uncertainties" << endl;
    yell  << _SSEPSC << endl;
    shout << "bin      cross section           lower uncertainty       upper uncertainty" << endl;
-   yell  << _SSEPSC << endl;
+   yell  << _TSEPSC << endl;
 
    vector < double > dxsu;
    vector < double > dxsl;
@@ -484,7 +455,7 @@ int main(int argc, char** argv) {
       dxsu.push_back(dxs[iobs]);
       dxsl.push_back(dxs[iobs]);
    }
-   yell << _SSEPSC << endl;
+   yell << _TSEPSC << endl;
 
    //! Without YODA we can stop here
 #ifndef WITH_YODA
-- 
GitLab