-
Notifications
You must be signed in to change notification settings - Fork 12
Usage of external codes
This page contains descriptions of external codes in nano/external
The codes are src/lumiTool.c and interface/lumiTool.h.
The following code is for initialization. The file Cert_271036-284044_13TeV_PromptReco_Collisions16_JSON.txt can be changed if you want to use your own lumi JSON file.
// lumiTool m_lumi; is already given above
string env = getenv("CMSSW_BASE");
string lumi = env+"/src/nano/analysis/data/Cert_271036-284044_13TeV_PromptReco_Collisions16_JSON.txt";
m_lumi = new lumiTool(lumi);
This code is from nano/analysis/src/nanoBase.cc. We already have this in the constructor of nanoBase class (into m_lumi). Once you run the above and keep m_lumi, in your event loop you can reject events not in the JSON file as following:
for ( [EVENT LOOP] ) {
// run and luminosityBlock are given in nanoBase class
if ( !( m_lumi->LumiCheck(run, luminosityBlock) ) ) continue;
...
}
Calculator of pile-up weight. The codes are src/pileUpTool.cc and interface/pileUpTool.h. This calculator is adjusted to Moriond17 MC/RD.
Initialization is done by just one line;
pileUpTool *m_pileUp;
m_pileUp = new pileUpTool();
It is already called in the constructor of nanoBase. It is also easy to get the pile-up weight. The following is enough:
// If you want to evaluate uncertainty, change the second parameter to +1 or -1
b_puweight = m_pileUp->getWeight(Pileup_nTrueInt, 0); // Pileup_nTrueInt is provided in the nanoBase class
Although the usage is quite easy, the numbers in src/pileUpTool.cc, contained in four arraries pileUpTool::Moriond17MC, pileUpTool::Moriond17RD, pileUpTool::Moriond17RD_up, pileUpTool::Moriond17RD_dn, come from somewhat complicated paths. The origin of the MC array is somewhat simple; see mix.input.nbPileupEvents.probValue in https://github.com/cms-sw/cmssw/blob/master/SimGeneral/MixingModule/python/mix_2016_25ns_Moriond17MC_PoissonOOTPU_cfi.py, as following https://twiki.cern.ch/twiki/bin/view/CMS/PdmVPileUpDescription#Startup2015 (although this page is for 2015 MC/data, the prescription still works).
The RD tables are somewhat complicated. The way for this is described in here. According to this page, first, we need the JSON file containing informations for pileup. The file is /afs/cern.ch/cms/CAF/CMSCOMM/COMM_DQM/certification/Collisions16/13TeV/PileUp/pileup_latest.txt (you can find it in lxplus). Then, the following command will give a root file:
pileupCalc.py -i [Cert_..._Collision16_JSON.txt] --inputLumiJSON pileup_latest.txt --calcMode true --minBiasXsec 69200 --maxPileupBin 75 --numPileupBins 75 [A_NANE_YOU_WANT].root
The lumi JSON file [Cert_..._Collision16_JSON.txt] can be found in nano/analysis/data if you ran getFiles. Since the 'pileup' number is an integer, maxPileupBin and numPileupBins have to be same. (In the given table the number of bins is 75, so I gave same). This command will create [A_NANE_YOU_WANT].root, which contains a histogram. Each of bin contents is the number in the RD table.
On the other hand, minBiasXsec is the measured cross section of MinBias events in 2016, which is 69.200 mb +- 5%. (In the CATTools, 69.000 mb is used; see sec. Dataset selections in here.) The uncertainty from pileup is estimated by giving the cross section which got up (69.200 mb + 5% = 72.660 mb) and down (69.200 mb - 5% = 65.740 mb). The RD_up and RD_dn tables are filled by the bin contents of them.
In the interface
There is computeTrigSF(lepton1, lepton2, direction = central:0 or up:1 or down:-1)
we select an lepton1 and an lepton 2 to the computeTrigSF,
sumId will be checked like this abs(lepton1.PdgCode()) + abs(lepton2.PdgCode())
11+11 = ee channel, 13+13 = mumu channel, else = emu channel.
aeta1 = abs(lepton1.Eta()), aeta2 = abs(lepton2.Eta())
and then, return lepton trigger weight for each Eta.
if ( channel == CH_ELEL ) {
if ( aeta1 < 0.3 ) {
if ( aeta2 < 0.3 ) return 0.965 + direction*0.010;
else if ( aeta2 < 0.6 ) return 0.989 + 0.5*direction*0.007;
else if ( aeta2 < 1.2 ) return 1.003 + 0.5*direction*0.004;
else if ( aeta2 < 1.7 ) return 0.990 + 0.5*direction*0.010;
else if ( aeta2 < 2.4 ) return 0.997 + 0.5*direction*0.011;
}
else if ( aeta1 < 0.6 ) {
if ( aeta2 < 0.3 ) return 0.981 + direction*0.007;
else if ( aeta2 < 0.6 ) return 0.972 + 0.5*direction*0.009;
else if ( aeta2 < 1.2 ) return 0.990 + 0.5*direction*0.005;
else if ( aeta2 < 1.7 ) return 1.000 + 0.5*direction*0.007;
else if ( aeta2 < 2.4 ) return 0.997 + 0.5*direction*0.011;
}
...
else if ( aeta1 < 2.4 ) {
if ( aeta2 < 0.3 ) return 1.005 + direction*0.013;
else if ( aeta2 < 0.6 ) return 0.995 + 0.5*direction*0.016;
else if ( aeta2 < 1.2 ) return 0.989 + 0.5*direction*0.011;
else if ( aeta2 < 1.7 ) return 0.985 + 0.5*direction*0.018;
else if ( aeta2 < 2.4 ) return 0.996 + 0.5*direction*0.014;
}
}
CH_MUMU and CH_MUEL are calculated, also.
CSVv2_Moriond_B_H.csv is used into BTagCalibration,
There are BTagCalibrationReader(OperatingPoint = LOOSE,MEDIUM,TIGHT , sysType = "central", "up", "down")
and eval_auto_bounds(sysType, JetFlavor, eta, pt, dicrimination=CSVv2).
declare, in the nanoBase.cc
string csvFileName = "CSVv2_Moriond17_B_H.csv";
std::string csvFile = env+"/src/nano/analysis/data/btagSF/"+csvFileName;
BTagCalibration calib("csvv2", csvFile);
m_btagSF = BTagCalibrationReader(BTagEntry::OP_MEDIUM, "central");
m_btagSF.load(calib, BTagEntry::FLAV_B, "comb");
calculate, in the topObjectSelection.cc
btagweight will be calculated after considering b, c and usdg.
double csvWgtHF = 1.0, csvWgtLF = 1.0, csvWgtC = 1.0, csvWgtTotal = 1.0;
#in the jet selection loop,
if (abs(Jet_hadronFlavour[i]) == 5) {
JF = BTagEntry::FLAV_B;
double iCSVWgtHF = m_btagSF.eval_auto_bounds("central", JF, abs(Jet_eta[i]), Jet_pt[i], Jet_btagCSVV2[i]);
if (iCSVWgtHF != 0) csvWgtHF *= iCSVWgtHF;
}
else if (abs(Jet_hadronFlavour[i]) == 4) {
JF = BTagEntry::FLAV_C;
double iCSVWgtC = m_btagSF.eval_auto_bounds("central", JF, abs(Jet_eta[i]), Jet_pt[i], Jet_btagCSVV2[i]);
if (iCSVWgtC != 0) csvWgtC *= iCSVWgtC;
}
else {
JF = BTagEntry::FLAV_USDG;
double iCSVWgtLF = m_btagSF.eval_auto_bounds("central", JF, abs(Jet_eta[i]), Jet_pt[i], Jet_btagCSVV2[i]);
if (iCSVWgtLF != 0) csvWgtLF *= iCSVWgtLF;
}
csvWgtTotal = csvWgtHF * csvWgtC * csvWgtLF;
#this is btagweight which we use
b_btagweight = csvWgtTotal;