...
| Code Block |
|---|
cd muTrkIso root -b -q muonIso.C++ # produce output of histograms root effMuon.C # produce comparison plot among original trk iso, PF iso, and trk iso with a veto of 2nd lepton in AN note |
Electron id
| Code Block |
|---|
for( int j = 0; j < nelec; j++ ) {
if( fabs(elec_etaSC[j]) < 1.47 ) isEB = true;
else isEE = true;
double _modTrkIso = elec_modIsoPtTrks[j];
double _modEMHADIso = elec_modIsoEMHADDepth1[j];
for( int k = 0; k < nTT; k++ ) {
if( track_pT[k] < _ptCut ) continue;
double tmpDr = deltaR(track_eta[k], track_phi[k], elec_gsfEta[j], elec_gsfPhi[j]);
if( tmpDr < 0.3 ) {
_modTrkIso -= track_pT[k];
if( _modTrkIso < 0 ) _modTrkIso = 0;
}
}
bool isElecSel = false;
if( fabs(elec_etaSC[j]) < _etaCut
&& elec_ecalDriven[j] == 1
&& fabs(elec_dPhiIn[j]) < 0.06
&& fabs(elec_dEtaIn[j]) < 0.02
&& elec_HoverE[j] > -99 && elec_HoverE[j] < 0.1 // default: 0.05
&& (elec_mHits[j] >= 0 && elec_mHits[j] <= 1)
&& _modTrkIso < 5 // iso
) {
if( isEB ) {
if( //fabs(elec_dEtaIn[j]) < 0.007
//(elec_25over55[j] > 0.90 || elec_15over55[j] > 0.83)
_modEMHADIso < (2.0+0.05*elec_et[j]) // iso
&& fabs(elec_dxyVTX[j]) < 0.02 ) {
if( elec_et[j] > _ptCutLeading ) {
elecEBLeading.push_back(tmpElec);
if( elec_charge[j] > 0 ) elPlus.push_back(tmpElec);
else elMinus.push_back(tmpElec);
isElecSel = true;
}
else if( elec_et[j] > _ptCut ) {
elecEB2ndLeading.push_back(tmpElec);
if( elec_charge[j] > 0 ) elPlus.push_back(tmpElec);
else elMinus.push_back(tmpElec);
isElecSel = true;
}
}
}
else if( isEE ) {
bool isPassEMH = false;
// elec iso
if( elec_et[j] < 50. ) {
if( _modEMHADIso < (2.5) ) isPassEMH = true;
}
else {
if( _modEMHADIso < (2.5+0.05*(elec_et[j]-50.)) ) isPassEMH = true;
}
if( isPassEMH
&& elec_sigmaIEtaIEta[j] > -99 && elec_sigmaIEtaIEta[j] < 0.03
&& fabs(elec_dxyVTX[j]) < 0.05 ) {
if( elec_et[j] > _ptCutLeading ) {
elecEELeading.push_back(tmpElec);
if( elec_charge[j] > 0 ) elPlus.push_back(tmpElec);
else elMinus.push_back(tmpElec);
isElecSel = true;
}
else if( elec_et[j] > _ptCut ) {
elecEE2ndLeading.push_back(tmpElec);
if( elec_charge[j] > 0 ) elPlus.push_back(tmpElec);
else elMinus.push_back(tmpElec);
isElecSel = true;
}
}
}
}
|
...
| Code Block |
|---|
cd preapproval/accEff root GetSyst.C # then print out the efficiency difference between signal and ZZ MC |
Limit calculation
Use HiggsCombined package in CMSSW: twiki
| Code Block |
|---|
setenv SCRAM_ARCH slc5_amd64_gcc472 cmsrel CMSSW_6_1_1 ### must be >= 6.1.1, as older versions have bugs (6.2.X and 7.0.X are NOT supported either) cd CMSSW_6_1_1/src cmsenv git clone https://github.com/cms-analysis/HiggsAnalysis-CombinedLimit.git HiggsAnalysis/CombinedLimit cd CMSSW_6_1_1/src/HiggsAnalysis/CombinedLimit/preapproval_final/ # working directory used for preapproval results |
Produce limit plots
Need to produce limit plots for both data and mc samples with fine bin size to produce datacard of HiggsCombined package
| Code Block |
|---|
cd analysisCodes/preapproval/Limits/
root -b limit.C++ # to produce all necessary histograms for datacard and limit calculation |
Background fit
In current analysis (presented in preapproval), we set the background rate to 0 in the limit calculation because the background is negligible. Therefore this background fit is not used in current analysis.
The main motivation for this background fit is to estimate backgrounds in M(4l) > 1 TeV region because we have almost no statistics from SM backgrounds in this region. So we extrapolate the distribution of backgrounds from low mass region (M(4l) < 500 GeV) where we observe good agreement between data and mc.
| Code Block |
|---|
cd analysisCodes/preapproval/Limits/
root fitPlotBg.C |
Produce datacard
We need large number of datacard as an input for limit calculation and therefore we need to automatize it.
| Code Block |
|---|
cd analysisCodes/preapproval/Limits/
root prodDataCard.C
# you can find the datacard in each channel in dataCard directory |
Run limit calculation
Need to carefully test the option of the method and then produce all channels
In currently analysis, we use "BayesianToyMC" method (same as MarkovChainMC, but simpler option), recommended by experts.
| Code Block |
|---|
# for mumumumu channel as an example (same for other channels)
cd CMSSW_6_1_1/src/HiggsAnalysis/CombinedLimit/preapproval_final/
cmsenv
./calc.csh
# high mass calculation is done quickly relatively and low mass calculation takes longer time
# need to adjust the option (for example, reduce the integration range rMax, if the signal strength is smaller)
# example output
-- BayesianToyMC --
Limit: r < 0.677481 @ 95% CL
Done in 1.86 min (cpu), 1.86 min (real)
mean expected limit: r < 0.676587 +/- 4.86033e-05 @ 95%CL (200 toyMC)
median expected limit: r < 0.676539 @ 95%CL (200 toyMC)
68% expected band : 0.675963 < r < 0.677278
95% expected band : 0.675415 < r < 0.678158
# there is signal strength as a limit and 1/2 sigma expected band
# model-independent calculation
combine -M BayesianToyMC ../../dataCard/zprime_mumumumu_ind.txt --rMax 3 -t 200 -i 100000 > & outputInd_mumumumu.txt & |
Limit plot
| Code Block |
|---|
cd analysisCodes/preapproval/Limits/
# When you finish to calculate the limits for all channels and mass bins, you need extract the limit results from the output and make the input files as the files in inputLimits/
root drawLimit.C # produce the figure in AN note
root calcGridLimit.C # produce the table in AN note (and paper) |