File indexing completed on 2025-01-18 10:06:36
0001
0002
0003
0004
0005
0006
0007
0008 #ifndef Pythia8_Weights_H
0009 #define Pythia8_Weights_H
0010
0011 #include "Pythia8/Basics.h"
0012 #include "Pythia8/LHEF3.h"
0013 #include "Pythia8/PythiaStdlib.h"
0014 #include "Pythia8/SharedPointers.h"
0015
0016 namespace Pythia8 {
0017
0018
0019 class Info;
0020 class History;
0021 class PartonLevel;
0022 class Merging;
0023 class WeightContainer;
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035 class WeightsBase {
0036
0037 public:
0038
0039
0040 friend class History;
0041 friend class PartonLevel;
0042
0043
0044 virtual void init() {weightValues.resize(0); weightNames.resize(0);
0045 bookWeight("Baseline");}
0046 virtual void init(bool) {}
0047
0048
0049 virtual void clear() {fill(weightValues.begin(), weightValues.end(), 1.);}
0050
0051
0052
0053 virtual void bookVectors(vector<double> weights, vector<string> names) {
0054 for (int i = 0; i < (int)weights.size(); ++i) {
0055 replace(names[i].begin(), names[i].end(), ' ', '_');
0056 bookWeight(names[i], weights[i]);
0057 }
0058 }
0059
0060
0061
0062 virtual void collectWeightValues(vector<double>& outputWeights,
0063 double norm = 1.);
0064
0065
0066 virtual void collectWeightNames(vector<string>& outputNames);
0067
0068
0069
0070
0071
0072
0073
0074
0075 string getWeightsName(int iPos) const {
0076 string name = iPos >= 0
0077 && iPos < (int)weightNames.size() ? weightNames[iPos] : "";
0078 if (name.find(":") != string::npos)
0079 replace(name.begin(), name.end(), ':', '.');
0080 return name == "" ? to_string(iPos) : name;}
0081 virtual double getWeightsValue(int iPos) const { return weightValues[iPos]; }
0082 int getWeightsSize() const { return weightValues.size(); }
0083
0084
0085 void bookWeight(string name, double defaultValue = 1.) {
0086
0087 if (findIndexOfName(name) != -1) setValueByName(name, defaultValue);
0088 else {
0089 weightNames.push_back(name);
0090 weightValues.push_back(defaultValue);
0091 }
0092 }
0093
0094
0095 void setValueByIndex(int iPos, double val) {
0096 if (iPos < 0 || iPos >= (int)weightValues.size()) return;
0097 weightValues[iPos] = val;
0098 }
0099 void setValueByName(string name, double val) {
0100 setValueByIndex(findIndexOfName(name), val);}
0101
0102
0103 virtual void reweightValueByIndex(int iPos, double val) {
0104 if (iPos < 0 || iPos >= (int)weightValues.size()) return;
0105 weightValues[iPos] *= val;
0106 };
0107 virtual void reweightValueByName(string name, double val) {
0108
0109 int iPos = findIndexOfName(name);
0110 reweightValueByIndex(iPos, val);
0111 }
0112
0113
0114
0115 int findIndexOfName(string name) {
0116 vector<string>::iterator it
0117 = find(weightNames.begin(), weightNames.end(), name);
0118 unsigned long int index = distance(weightNames.begin(), it);
0119 if (index == weightNames.size()) return -1;
0120 return distance(weightNames.begin(), it);
0121 }
0122
0123
0124 void setPtrs(Info* infoPtrIn) { infoPtr = infoPtrIn; }
0125
0126 protected:
0127
0128
0129 void parse(string wvecKey, map<string, map<string, double> > &dct);
0130
0131
0132 vector<double> weightValues;
0133 vector<string> weightNames;
0134
0135
0136 Info* infoPtr{};
0137
0138 };
0139
0140
0141
0142
0143
0144 class WeightsShower : public WeightsBase {
0145
0146 public:
0147
0148
0149 void init(bool) override {}
0150
0151
0152 virtual void initWeightGroups(bool = false) {}
0153 virtual int nWeightGroups() const {return 0;}
0154 virtual string getGroupName( int ) const {return "";}
0155 virtual double getGroupWeight( int ) const {return 1.;}
0156
0157 };
0158
0159
0160
0161
0162
0163 class WeightsSimpleShower : public WeightsShower {
0164
0165 public:
0166
0167
0168 void init( bool doMerging) override;
0169
0170
0171
0172 void collectWeightNames(vector<string>& outputNames) override;
0173 void collectWeightValues(vector<double>& outputWeights,
0174 double norm = 1.) override;
0175
0176
0177 void initWeightGroups(bool = false) override;
0178
0179
0180 string getGroupName(int iGN) const override;
0181
0182
0183 double getGroupWeight(int iGW) const override;
0184 int nWeightGroups() const override { return externalVariations.size(); }
0185
0186
0187 bool initUniqueShowerVars();
0188
0189
0190 void setEnhancedTrial( double pTIn, double wtIn) { pTEnhanced = pTIn;
0191 wtEnhanced = wtIn; }
0192 double getEnhancedTrialPT() { return pTEnhanced;}
0193 double getEnhancedTrialWeight() { return wtEnhanced;}
0194
0195
0196 vector<string> getUniqueShowerVars() { return uniqueShowerVars; }
0197 vector<string> getUniqueShowerVars(vector<string> keys);
0198
0199
0200 bool initEnhanceFactors();
0201 unordered_map<string,double> getEnhanceFactors() { return enhanceFactors; }
0202
0203 string getInitialName(int iG) const { return initialNameSave[iG]; }
0204
0205
0206 map<int,double> varPDFplus, varPDFminus, varPDFmember;
0207
0208
0209 vector<string> uniqueShowerVars;
0210 unordered_map<string,double> enhanceFactors;
0211
0212
0213 vector<string> externalVariations;
0214 vector<vector<string> > externalVarNames;
0215 vector<string> externalGroupNames;
0216 vector<string> initialNameSave;
0217 vector<vector<int> > externalMap;
0218 int externalVariationsSize{};
0219
0220
0221 vector<double> getMuRWeightVector();
0222 vector<vector<string> > mergingVarNames;
0223 double pTEnhanced{}, wtEnhanced{};
0224
0225 };
0226
0227
0228
0229
0230
0231
0232
0233
0234
0235 class WeightsMerging : public WeightsBase {
0236
0237 public:
0238
0239
0240 friend class Merging;
0241 friend class WeightContainer;
0242
0243
0244 void init() override;
0245
0246
0247 void clear() override;
0248
0249
0250 void bookWeight(string name, double value, double valueFirst);
0251
0252
0253 void bookVectors(vector<double> weights, vector<string> names) override;
0254 void bookVectors(vector<double> weights, vector<double> weightsFirst,
0255 vector<string> names);
0256
0257 double getWeightsValue(int iPos) const override {
0258 return weightValues[iPos] - weightValuesFirst[iPos]; }
0259
0260 double getWeightsValueP(int iPos) const {
0261 return weightValuesP[iPos] - weightValuesFirstP[iPos]; }
0262 double getWeightsValuePC(int iPos) const {
0263 return weightValuesPC[iPos] - weightValuesFirstPC[iPos]; }
0264
0265
0266 void reweightValueByIndex(int iPos, double val) override;
0267 void reweightValueByName(string name, double val) override;
0268
0269
0270 void setValueFirstByIndex(int iPos, double val);
0271 void setValueFirstByName(string name, double val);
0272
0273
0274 void setValueVector(vector<double> ValueVector);
0275 void setValueFirstVector(vector<double> ValueFirstVector);
0276
0277
0278 vector<double> getMuRVarFactors();
0279
0280
0281 void setLHEFvariationMapping();
0282
0283
0284 void collectWeightNames(vector<string>& outputNames) override;
0285
0286
0287 void collectWeightValues(vector<double>& outputWeights,
0288 double norm = 1.) override;
0289
0290 protected:
0291
0292
0293 map<int,int> muRVarLHEFindex;
0294
0295
0296 vector<double> weightValuesFirst;
0297
0298
0299 vector<double> weightValuesP, weightValuesPC,
0300 weightValuesFirstP, weightValuesFirstPC;
0301
0302
0303 bool isNLO;
0304
0305 };
0306
0307
0308
0309
0310
0311
0312
0313 class WeightsLHEF : public WeightsBase {
0314
0315 public:
0316
0317
0318 friend class WeightsMerging;
0319
0320
0321 double centralWeight;
0322
0323
0324 void clear() override;
0325
0326
0327 void bookVectors(vector<double> weights, vector<string> names) override;
0328
0329
0330
0331 void collectWeightValues(vector<double>& outputWeights,
0332 double norm = 1.) override;
0333 void collectWeightNames(vector<string>& outputNames) override;
0334
0335
0336
0337 vector<string> convertNames(vector<string> names);
0338
0339 void identifyVariationsFromLHAinit(map<string,LHAweight>* weights);
0340
0341 protected:
0342
0343
0344 map<int,double> muRvars;
0345
0346 };
0347
0348
0349
0350
0351
0352
0353 class WeightsFragmentation : public WeightsBase {
0354
0355 public:
0356
0357
0358 void init() override;
0359
0360 int nWeightGroups() const {return externalGroupNames.size();}
0361
0362 string getGroupName(int iGN) const {return iGN < 0 || iGN >= nWeightGroups()
0363 ? "Null" : externalGroupNames[iGN];}
0364
0365
0366 double getGroupWeight(int iGW) const {
0367 if (iGW < 0 || iGW >= nWeightGroups()) return 1.0;
0368 double wgt(1.);
0369 for (const int &iWgt : externalMap[iGW]) wgt *= getWeightsValue(iWgt);
0370 return wgt;}
0371
0372
0373
0374 void collectWeightNames(vector<string>& outputNames) override;
0375 void collectWeightValues(vector<double>& outputWeights,
0376 double norm = 1.) override;
0377
0378
0379 vector<map<vector<double>, int> > weightParms{};
0380 vector<string> externalGroupNames{};
0381 vector<vector<int> > externalMap{};
0382
0383
0384 enum FactIndex{Z, Flav, PT};
0385
0386 };
0387
0388
0389
0390
0391
0392
0393
0394 class WeightContainer {
0395
0396 public:
0397
0398
0399
0400 WeightContainer() : weightNominal(1.0), xsecIsInit(false) {}
0401
0402
0403 double weightNominal;
0404 void setWeightNominal( double weightNow );
0405 double collectWeightNominal();
0406
0407
0408 WeightsLHEF weightsLHEF{};
0409
0410
0411
0412
0413
0414 WeightsShower* weightsShowerPtr{};
0415 WeightsSimpleShower weightsSimpleShower{};
0416
0417
0418 WeightsMerging weightsMerging{};
0419
0420
0421 WeightsFragmentation weightsFragmentation{};
0422
0423
0424 WeightsBase weightsUserHooks{};
0425
0426
0427 int numberOfWeights();
0428 double weightValueByIndex(int key = 0);
0429 string weightNameByIndex(int key = 0);
0430
0431
0432
0433
0434
0435 vector<double> weightValueVector();
0436
0437
0438
0439 vector<string> weightNameVector();
0440
0441
0442 void clear();
0443
0444
0445 void clearTotal();
0446
0447
0448 void init(bool doMerging);
0449
0450
0451 void initPtrs(Info* infoPtrIn);
0452
0453
0454 void initXsecVec();
0455
0456
0457 vector<double> getSampleXsec();
0458 vector<double> getTotalXsec();
0459 vector<double> getSampleXsecErr();
0460 vector<double> getTotalXsecErr();
0461
0462
0463 void accumulateXsec(double norm = 1.);
0464
0465 private:
0466
0467
0468 Info* infoPtr{};
0469
0470
0471 bool doSuppressAUXweights{};
0472
0473
0474 vector<double> sigmaTotal, sigmaSample, errorTotal, errorSample;
0475 bool xsecIsInit;
0476
0477 };
0478
0479
0480
0481 }
0482
0483 #endif