diff --git a/v2.5/toolkit/fastnlotoolkit/fastNLOCoeffBase.cc b/v2.5/toolkit/fastnlotoolkit/fastNLOCoeffBase.cc index bf3c8ec237e06557bec11407d37262ac9a9fecf5..1b3156088f785e149b67d5f013ba8e48f617f35d 100644 --- a/v2.5/toolkit/fastnlotoolkit/fastNLOCoeffBase.cc +++ b/v2.5/toolkit/fastnlotoolkit/fastNLOCoeffBase.cc @@ -299,9 +299,12 @@ void fastNLOCoeffBase::AddCoeffInfoBlock(int fICoeffInfoBlockFlag1, int fICoeffI } void fastNLOCoeffBase::AddCoeffInfoBlock(int fICoeffInfoBlockFlag1, int fICoeffInfoBlockFlag2, std::vector<std::string> Description, - std::string filename, unsigned int icola, unsigned int icolb) { + std::string filename, unsigned int icola, unsigned int icolb, double relfac) { info["AddCoeffInfoBlocks"]<<"Adding additional InfoBlock reading data from file "<<filename<<endl; std::vector<double> Content = fastNLOTools::ReadContentFromFile(filename, icola, icolb); + for ( unsigned int i=0; i<Content.size(); i++ ) { + Content[i] = relfac*Content[i]; + } AddCoeffInfoBlock(fICoeffInfoBlockFlag1, fICoeffInfoBlockFlag2, Description, Content); } @@ -357,7 +360,7 @@ void fastNLOCoeffBase::WriteCoeffInfoBlocks(ostream& table, int ITabVersionWrite debug["WriteCoeffInfoBlocks"]<<"Writing additional InfoBlocks; NCoeffInfoBlocks = "<<NCoeffInfoBlocks<<endl; table << NCoeffInfoBlocks << sep; for (int i=0; i<NCoeffInfoBlocks; i++) { - debug["WriteCoeffInfoBlocks"]<<"ICoeffInfoBlockFlags1,2 = "<<ICoeffInfoBlockFlag1[i]<<", "<<ICoeffInfoBlockFlag1[i]<<endl; + debug["WriteCoeffInfoBlocks"]<<"ICoeffInfoBlockFlags1,2 = "<<ICoeffInfoBlockFlag1[i]<<", "<<ICoeffInfoBlockFlag2[i]<<endl; table << ICoeffInfoBlockFlag1[i] << sep; table << ICoeffInfoBlockFlag2[i] << sep; table << NCoeffInfoBlockDescr[i] << sep; diff --git a/v2.5/toolkit/fastnlotoolkit/fastNLOReader.cc b/v2.5/toolkit/fastnlotoolkit/fastNLOReader.cc index 4746d9787e727897fe9824695ddc690a31e0270d..0adbe39b32afcdd2ae6ebc629da0daa5e8b55820 100644 --- a/v2.5/toolkit/fastnlotoolkit/fastNLOReader.cc +++ b/v2.5/toolkit/fastnlotoolkit/fastNLOReader.cc @@ -1578,7 +1578,7 @@ void fastNLOReader::CalcCrossSectionv21(fastNLOCoeffAddFlex* c) { iCIBIndex = c->GetCoeffInfoBlockIndex(0); logger.debug["CalcCrossSectionv21"]<<"Found CoeffInfoBlock "<<iCIBIndex<<" with statistical/numerical uncertainties."<<endl; iCIBFlag2 = c->GetCoeffInfoBlockFlag2(iCIBIndex); - dCIBCont = c->GetCoeffInfoContent(iCIBIndex); + dCIBCont = c->GetCoeffInfoBlockContent(iCIBIndex); } else { logger.debug["CalcCrossSectionv21"]<<"No CoeffInfoBlock found; uncertainties are initialised to zero."<<endl; } @@ -1694,7 +1694,7 @@ void fastNLOReader::CalcCrossSectionv20(fastNLOCoeffAddFix* c) { iCIBIndex = c->GetCoeffInfoBlockIndex(0); logger.debug["CalcCrossSectionv20"]<<"Found CoeffInfoBlock "<<iCIBIndex<<" with statistical/numerical uncertainties."<<endl; iCIBFlag2 = c->GetCoeffInfoBlockFlag2(iCIBIndex); - dCIBCont = c->GetCoeffInfoContent(iCIBIndex); + dCIBCont = c->GetCoeffInfoBlockContent(iCIBIndex); } else { logger.debug["CalcCrossSectionv20"]<<"No CoeffInfoBlock found; uncertainties are initialised to zero."<<endl; } @@ -2211,7 +2211,7 @@ bool fastNLOReader::TestXFX() { // if ( sum== 0. ) printf("fastNLOReader. Error. All 13 pdf probabilities are 0. There might be sth. wrong in the pdf interface. Please check FastNLOUser::GetXFX().\n"); for (int i = 0 ; i<13 ; i++) { if (pdftest[i] > 1.e10 || (pdftest[i] < 1.e-10 && pdftest[i] > 1.e-15)) { - logger.warn["TestXFX"]<<"The pdf probability of the "<<i<<"'s flavor seeems to be unreasonably large/small (pdf="<<pdftest[i]<<") at x="<<xtest<<", mu="<<mutest<<".\n"; + logger.warn["TestXFX"]<<"The pdf probability of the "<<i<<"'s flavor seems to be unreasonably large/small (pdf="<<pdftest[i]<<") at x="<<xtest<<", mu="<<mutest<<".\n"; } } return true; @@ -3627,25 +3627,34 @@ XsUncertainty fastNLOReader::GetXsUncertainty(const ENumUncertaintyStyle eNumUnc vector < double > MydXSection; unsigned int NObsBin = GetNObsBin(); - //! For approximation bias get reference PDF and member + //! For approximation bias get reference PDF & member and cross section + std::vector < std::string > dCIBDescr; + std::vector < double > dCIBCont; + std::string PDFset; + std::string PDFmember; if (eNumUnc == kApproxBias) { - // Find perturbative (additive) contributions - for (unsigned int j = 0 ; j<BBlocksSMCalc.size() ; j++) { - for (unsigned int i = 0 ; i<BBlocksSMCalc[j].size() ; i++) { - if (BBlocksSMCalc[j][i] && BBlocksSMCalc[j][i]->IsEnabled()) { - if (fastNLOCoeffAddFlex::CheckCoeffConstants(BBlocksSMCalc[j][i],true)) - // CalcCrossSectionv21((fastNLOCoeffAddFlex*)BBlocksSMCalc[j][i]); - cout << "FLEX" << endl; - if ( BBlocksSMCalc[j][i].HasCoeffInfoBlock(1,0) ){ - cout << "FLEX YES" << endl; - } - else if (fastNLOCoeffAddFix::CheckCoeffConstants(BBlocksSMCalc[j][i],true)) - cout << "FIX" << endl; - } - } + for (unsigned int i=0; i<fCoeff.size() ; i++) { + fastNLOCoeffBase* c = GetCoeffTable(i); + if ( ! c->HasCoeffInfoBlock(1) ) { + logger.error["GetNumUncertainty"]<<"ERROR! No InfoBlock found for reference cross sections, exiting."<<endl; + exit(35); + } + cout << "EEEEE i = " << i << ", Has 0, 1 (0,1), (1,0) = " << c->HasCoeffInfoBlock(0) << c->HasCoeffInfoBlock(1) << c->HasCoeffInfoBlock(0,1) << c->HasCoeffInfoBlock(1,0) << endl; + int iCIBIndex = c->GetCoeffInfoBlockIndex(1); + logger.debug["GetNumUncertainty"]<<"Found CoeffInfoBlock "<<iCIBIndex<<" with reference cross sections."<<endl; + dCIBDescr = c->GetCoeffInfoBlockDescription(iCIBIndex); + dCIBCont = c->GetCoeffInfoBlockContent(iCIBIndex); + } + if ( dCIBDescr.size() < 3 ) { + logger.error["GetNumUncertainty"]<<"ERROR! InfoBlock description too short for reference cross sections, exiting."<<endl; + logger.error["GetNumUncertainty"]<<" Line two and three should contain the used PDF set and member."<<endl; + exit(36); + } else { + PDFset = dCIBDescr[1]; + PDFmember = dCIBDescr[2]; } + cout << "PDF set: " << PDFset << ", PDF member = " << PDFmember << endl; } - exit(2); //! Cross section and absolute uncertainties CalcCrossSection(); @@ -3667,6 +3676,14 @@ XsUncertainty fastNLOReader::GetXsUncertainty(const ENumUncertaintyStyle eNumUnc XsUnc.dxsu.push_back(MydXSection[iobs]); XsUnc.dxsl.push_back(-MydXSection[iobs]); } + } else if (eNumUnc == kApproxBias) { + logger.info["GetNumUncertainty"]<<"Interpolation bias selected."<<endl; + for (unsigned int iobs = 0; iobs < NObsBin; iobs++) { + XsUnc.xs.push_back(MyXSection[iobs]); + XsUnc.dxsu.push_back(dCIBCont[iobs]); + XsUnc.dxsl.push_back(MyXSection[iobs]-dCIBCont[iobs]); + cout << "XXX XS = " << MyXSection[iobs] << ", XSREF = " << dCIBCont[iobs] << ", delta = " << MyXSection[iobs]-dCIBCont[iobs] << endl; + } } else { logger.error["GetNumUncertainty"]<<"ERROR! No valid numerical uncertainty style selected, exiting."<<endl; logger.error["GetNumUncertainty"]<<"Style enum = "<<eNumUnc<<endl; @@ -3675,12 +3692,20 @@ XsUncertainty fastNLOReader::GetXsUncertainty(const ENumUncertaintyStyle eNumUnc //! Divide by cross section != 0 to give relative uncertainties for (unsigned int iobs = 0; iobs < NObsBin; iobs++) { - if (fabs(XsUnc.xs[iobs]) > DBL_MIN) { - XsUnc.dxsu[iobs] = +fabs(XsUnc.dxsu[iobs] / XsUnc.xs[iobs]); - XsUnc.dxsl[iobs] = -fabs(XsUnc.dxsl[iobs] / XsUnc.xs[iobs]); + if (eNumUnc == kApproxBias) { + if (fabs(dCIBCont[iobs]) > DBL_MIN) { + XsUnc.dxsl[iobs] = XsUnc.dxsl[iobs] / fabs(dCIBCont[iobs]); + } else { + XsUnc.dxsl[iobs] = 0.; + } } else { - XsUnc.dxsu[iobs] = 0.; - XsUnc.dxsl[iobs] = 0.; + if (fabs(XsUnc.xs[iobs]) > DBL_MIN) { + XsUnc.dxsu[iobs] = +fabs(XsUnc.dxsu[iobs] / XsUnc.xs[iobs]); + XsUnc.dxsl[iobs] = -fabs(XsUnc.dxsl[iobs] / XsUnc.xs[iobs]); + } else { + XsUnc.dxsu[iobs] = 0.; + XsUnc.dxsl[iobs] = 0.; + } } logger.debug["GetNumUncertainty"]<<"iobs = " << iobs << ", dxsl = " << XsUnc.dxsl[iobs] << ", dxsu = " << XsUnc.dxsu[iobs] <<endl; } diff --git a/v2.5/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOCoeffBase.h b/v2.5/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOCoeffBase.h index 162b8d52c4ab45d2b72d7f21f3a240498ac57349..b6d519476445f0852a34c614e566aa32156aa4a2 100644 --- a/v2.5/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOCoeffBase.h +++ b/v2.5/toolkit/fastnlotoolkit/include/fastnlotk/fastNLOCoeffBase.h @@ -90,14 +90,17 @@ public: int GetCoeffInfoBlockFlag2(int Index) const { return ICoeffInfoBlockFlag2[Index]; }; void SetCoeffInfoBlockFlag1(int Index, int iFlag1) { ICoeffInfoBlockFlag1[Index] = iFlag1; }; void SetCoeffInfoBlockFlag2(int Index, int iFlag2) { ICoeffInfoBlockFlag2[Index] = iFlag2; }; - std::vector < double > GetCoeffInfoContent(int Index) const { return CoeffInfoBlockContent[Index]; }; + std::vector< std::string > GetCoeffInfoBlockDescription(int Index) const { return CoeffInfoBlockDescript[Index]; } + std::vector < double > GetCoeffInfoBlockContent(int Index) const { return CoeffInfoBlockContent[Index]; }; int GetNCoeffInfoBlocks() const {return NCoeffInfoBlocks;} // Provide uncertainty via input vector void AddCoeffInfoBlock(int ICoeffInfoBlockFlag1, int ICoeffInfoBlockFlag2, std::vector<std::string> Description, std::vector<double> Content); // Provide uncertainty reading from filename void AddCoeffInfoBlock(int ICoeffInfoBlockFlag1, int ICoeffInfoBlockFlag2, std::vector<std::string> Description, - std::string filename, unsigned int icola = 0, unsigned int icolb = 0); + std::string filename, unsigned int icola = 0, unsigned int icolb = 0, double relfac = 1); + + protected: void ReadBase(std::istream& table, int ITabVersionRead); diff --git a/v2.5/toolkit/src/fnlo-tk-modify.cc b/v2.5/toolkit/src/fnlo-tk-modify.cc index b702a55478046e61ba02a3d7af9ceb30da4b16a5..e20c56a41ac008b5cc0722a37b2fe959af024297 100644 --- a/v2.5/toolkit/src/fnlo-tk-modify.cc +++ b/v2.5/toolkit/src/fnlo-tk-modify.cc @@ -367,6 +367,7 @@ int main(int argc, char** argv) { std::string KeyOrders = Key + "Orders"; std::string KeyFiles = Key + "Files"; std::string KeyColumns = Key + "Columns"; + std::string KeyFacs = Key + "Facs"; cout << "AAAAAAAAAA: Key = " << Key << endl; if ( !read_steer::getexist(KeyFlag1) || !read_steer::getexist(KeyFlag2) ){ info["fnlo-tk-modify"]<<"No InfoBlock information found for key " << Key << " Skipped."<<endl; @@ -493,6 +494,7 @@ int main(int argc, char** argv) { unsigned int NOrders = read_steer::getstringarray(KeyOrders).size(); unsigned int NFiles = read_steer::getstringarray(KeyFiles).size(); unsigned int NCols = read_steer::getintarray(KeyColumns).size(); + unsigned int NFacs = read_steer::getintarray(KeyFacs).size(); if ( NFiles != NOrders ) { error["fnlo-tk-modify"]<<"Need one order specification per file, aborted! Found NFiles = " << NFiles << ", and NOrders = " << NOrders <<endl; exit(37); @@ -510,6 +512,14 @@ int main(int argc, char** argv) { exit(38); } cout << "ZZZZZZZZZZZ icola, icolb = " << icola << ", " << icolb << endl; + double relfac = 1.; + if ( NFacs == 0 ) { + } else if ( NFacs == 1 ) { + relfac = read_steer::getdoublearray(KeyFacs)[0]; + } else { + error["fnlo-tk-modify"]<<"Only one multiplicative factor allowed, but found more. Aborted! NFacs = " << NFacs <<endl; + exit(39); + } for ( unsigned int i = 0; i < NFiles; i++ ){ info["fnlo-tk-modify"]<<"InfoBlock file no. " << i << " is: " << read_steer::getstringarray(KeyFiles)[i] << endl; } @@ -541,7 +551,7 @@ int main(int argc, char** argv) { } } } - c->AddCoeffInfoBlock(iFlag1,iFlag2,Description,read_steer::getstringarray(KeyFiles)[ilo],icola,icolb); + c->AddCoeffInfoBlock(iFlag1,iFlag2,Description,read_steer::getstringarray(KeyFiles)[ilo],icola,icolb,relfac); ic += 1; } else if ( c->IsNLO() ) { info["fnlo-tk-modify"]<<"Found NLO contribution " << i << endl; @@ -553,7 +563,7 @@ int main(int argc, char** argv) { } } } - c->AddCoeffInfoBlock(iFlag1,iFlag2,Description,read_steer::getstringarray(KeyFiles)[inlo],icola,icolb); + c->AddCoeffInfoBlock(iFlag1,iFlag2,Description,read_steer::getstringarray(KeyFiles)[inlo],icola,icolb,relfac); ic += 1; } else if ( c->IsNNLO() ) { info["fnlo-tk-modify"]<<"Found NNLO contribution " << i << endl; @@ -565,7 +575,7 @@ int main(int argc, char** argv) { } } } - c->AddCoeffInfoBlock(iFlag1,iFlag2,Description,read_steer::getstringarray(KeyFiles)[innlo],icola,icolb); + c->AddCoeffInfoBlock(iFlag1,iFlag2,Description,read_steer::getstringarray(KeyFiles)[innlo],icola,icolb,relfac); ic += 1; } else { info["fnlo-tk-modify"]<<"Unknown contribution " << i << endl;