Electroweak Utilities for MC corrections and systematics evaluation

Code and tags

Code is in ElectroWeakAnalysis/Utilities. Latest working tag as of today is V00-03-06 (38X - 39X series). For 34X versions or lower you will probably need an up to date version of the package MuonAnalysis/MomentumScaleCalibration (use V00-03-03 for instance). For 33X versions you will need older tags of ElectroWeakAnalysis/Utilities for technical reasons (use V00-01-14). The package ElectroWeakAnalysis/Utilities is not expected to work properly with versions below 33X.

How to compile, link, run

This code requires special preparation because full LHAPDF libraries are not linked in CMSSW by default. Possible recipes can be found in the first lines of the ElectroWeakAnalysis/Utilities/README file. The recommended way to run corresponds to:

   1) scram setup lhapdffull
   2) touch $CMSSW_BASE/src/ElectroWeakAnalysis/Utilities/BuildFile.xml (to ensure recompilation of the whole package)
   3) cmsenv
   4) scram b

Alternatively, you can try to plug your own LHAPDF libraries. Assuming that they are installed under "/usr/local":

    1) Edit the file $CMSSW_BASE/config/toolbox/$SCRAM_ARCH/tools/selected/lhapdf.xml
    2a) In that file, change the default for LHAPDF_BASE to "/usr/local"
    2b) In that file, change the default "version" by the new LHAPDF version
    3) scram setup lhapdf (and check that environment variables point to /usr/local/...)
    4) touch $CMSSW_BASE/src/ElectroWeakAnalysis/Utilities/CMS.BuildFile.xml
    5) cmsenv
    6) scram b

Keep in mind that full LHAPDF libraries are more expensive in terms of memory and disk space than the default light LHAPDF libraries. On lxplus, there are limitations on the virtual memory that you are using. You can bypass them by using on tcsh shell the command "limit vmem unlim" (on bash, use "ulimit -v hard").

The code only makes sense on MC (it produces nothing on data).

Jobs should run via CRAB if you are using the "lhapdffull" libraries. But as of today, CRAB does not pass properly the necessary environment. We have checked that a way to circumvent this is to add the following two commands to the CMSSW.sh executable in the local crab_.../job/ subdirectory before submitting the jobs:

scramv1 setup lhapdffull
scramv1 b
You can add them for instance after the 'eval `scramv1 runtime -sh ...`' command. It may also run with your own libraries if you get direct access to the necessary LHAPDF paths at your Tier2/Tier3 and set things accordingly (but we have not checked this possibility yet).

Documentation inside CMSSW

A short document describing all available utilities can be found in ElectroWeakAnalysis/Utilities/README.

PDF SYSTEMATICS (PDFWeightProducer)

(Important: an optimal usage of this utility requires user intervention in order to change the default set of PDF libraries. Please read the installation instructions in order to get an optimal performance out of it.)

Summary

The aim of the PDFWeightProducer is to write into the data structure "std::vector" products containing the LHAPDF weights "PDF_user(i)/PDF_in_generated_sample" for each event. "PDF_user" refers to the PDF values of the set chosen by the user in the configuration file, "PDF_in_generated_sample" is the PDF used when the sample was generated and "i" runs over the total number of members in the set. These weights are the basic ingredient to assign any PDF systematic uncertainty to any observable.

PDFSystematicAnalyzer is a complete analyzer that provides estimates on the rate and acceptance for any given selection filter.

General usage in the analysis

For PDFWeightProducer, the relevant piece of code in the configuration file looks like this:

 
      # Produce PDF weights (maximum is 3)
      process.pdfWeights = cms.EDProducer("PdfWeightProducer",
            # Fix POWHEG if buggy (this PDF set will also appear on output,
            # so only two more PDF sets can be added in PdfSetNames if not "")
            #FixPOWHEG = cms.untracked.string("cteq66.LHgrid"),
            #GenTag = cms.untracked.InputTag("genParticles"),
            PdfInfoTag = cms.untracked.InputTag("generator"),
            PdfSetNames = cms.untracked.vstring(
                    "cteq66.LHgrid"
                  , "MRST2006nnlo.LHgrid"
                  , "NNPDF10_100.LHgrid"
            )
      )
where "cteq66.LHgrid", "MRST2006nnlo.LHgrid", "NNPDF10_100.LHgrid" are some of the possible input PDF set files present in the directory $LHAPATH/. One can write a maximum of three different weight vectors in one job. The product tags written in the data structure will be called "pdfWeights:cteq66", "pdfWeights:MRST2006nnlo" and "pdfWeights:NNPDF10", respectively. (Note that the NNPDF product has no "_100" suffix: CMSSW does not allow underscores in product labels, because they are used to separte the different items in the full product tag.)

In order to determine systematic uncertainties on your own observable, you need to weight your MC events differently, according to the weights produced by PDFWeighttProducer. The following example shows how to pick up these weights in your analyzer for the product "pdfWeights:cteq66":

void YourAnalyzer::analyze(const edm::Event & ev, const edm::EventSetup&){
     edm::InputTag pdfWeightTag("pdfWeights:cteq66"); // or any other PDF set
     edm::Handle<std::vector<double> > weightHandle;
     ev.getByLabel(pdfWeightTag, weightHandle);

     std::vector<double> weights = (*weightHandle);
     std::cout << "Event weight for central PDF:" << weights[0] << std::endl;
     unsigned int nmembers = weights.size();
     for (unsigned int j=1; j<nmembers; j+=2) {
           std::cout << "Event weight for PDF variation +" << (j+1)/2 << ": " << weights[j] << std::endl;
           std::cout << "Event weight for PDF variation -" << (j+1)/2 << ": " << weights[j+1] << std::endl;
     }
} 

Assigning PDF systematics. Systematics on rates and acceptances (PDFSystematicAnalyzer)

Conventions to assign PDF systematics may still change in the future. The most widely used recipe for CTEQ and MSTW PDFs is based on the ``master'' equations' defined in: hep-ph/0605240. Variations of any observable ``X'' can be easily obtained via the event weights written by PDFWeightProducer. There are also specific recipes for each set (CTEQ6 PDFs correspond to 90% CL variations, NNPDFs use a different logic for the meaning of the sets). At present (Autumn 2010) we are using the PDF4LHC recommendations (http://www.hep.ucl.ac.uk/pdf4lhc/PDF4LHCrecom.pdf).

Let us note that the weights written by PDFWeightProducer consider PDF variations at the hard scattering level only. This is standard practice. Strictly speaking, there are more PDF dependences in the generator programs. Sudakov factors are employed to build up additional initial-state QCD radiation, and these factors are also governed by PDFs. From a quantitative point of view, these PDF dependences are second order effects compared with the dominant ones at the hard scattering level. Moreover, these Sudakov PDF dependences are absorbed in practice in the ISR systematic uncertainties. related with the implementation of the backward shower evolution.

Typical cases are the determination of PDF systematics on the expected rate or on the estimated acceptance for a given process after experimental cuts. These two cases have been considered in PDFSystematicAnalyzer. The aim is to automatize the evaluation of systematics for cross section measurements. The relevant piece of code in a configuration file that uses PDFSystematicAnalyzer.cc.

-

Full example

import FWCore.ParameterSet.Config as cms

# Process name
process = cms.Process("PDFANA")

# Max events
process.maxEvents = cms.untracked.PSet(
    input = cms.untracked.int32(-1)
)

# Printouts
process.MessageLogger = cms.Service("MessageLogger",
      cout = cms.untracked.PSet(
            default = cms.untracked.PSet(limit = cms.untracked.int32(100)),
            threshold = cms.untracked.string('INFO')
      ),
      destinations = cms.untracked.vstring('cout')
)

# Input files
process.source = cms.Source("PoolSource",
      debugVerbosity = cms.untracked.uint32(0),
      debugFlag = cms.untracked.bool(False),
      fileNames = cms.untracked.vstring("file:/data4/Wmunu_Summer09-MC_31X_V3_AODSIM-v1/0009/F82D4260-507F-DE11-B5D6-00093D128828.root")
)

# Produce PDF weights (maximum is 3)
process.pdfWeights = cms.EDProducer("PdfWeightProducer",
      FixPOWHEG = cms.untracked.bool(False), # fix POWHEG (it requires cteq66* PDFs in the list)
      PdfInfoTag = cms.untracked.InputTag("generator"),
      PdfSetNames = cms.untracked.vstring(
              "cteq65.LHgrid"
            , "MRST2006nnlo.LHgrid"
            , "MRST2007lomod.LHgrid"
      )
)

# Selector and parameters
# WMN fast selector, which requires
# the libraries and plugins fron the ElectroWeakAnalysis/WMuNu package
process.load("ElectroWeakAnalysis.WMuNu.WMuNuSelection_cff")

# Collect uncertainties for rate and acceptance
process.pdfSystematics = cms.EDFilter("PdfSystematicsAnalyzer",
      SelectorPath = cms.untracked.string('pdfana'),
      PdfWeightTags = cms.untracked.VInputTag(
              "pdfWeights:cteq65"
            , "pdfWeights:MRST2006nnlo"
            , "pdfWeights:MRST2007lomod"
      )
)

# Main path
process.pdfana = cms.Path(
       process.pdfWeights
      *process.selectCaloMetWMuNus
)

# Collect results
process.end = cms.EndPath(process.pdfSystematics)

The latest reference example can be found here: ElectroWeakAnalysis/Utilities/test/PdfSystematicsAnalyzer.py.

Example of output

>>>> Begin of PDF weight systematics summary >>>>
Total number of analyzed data: 2259 [events]
Total number of selected data: 1054 [events], corresponding to acceptance: [46.6578 +- 1.04964] %

>>>>> PDF UNCERTAINTIES ON RATE >>>>>>
RATE Results for PDF set cteq65 ---->
        Estimate for central PDF member: 1200 [events]
        i.e. [13.93 +- 2.566] % relative variation with respect to original PDF
        Number of eigenvectors for uncertainty estimation: 20
        Relative uncertainty with respect to central member: +4.28 / -4.381 [%]
RATE Results for PDF set MRST2006nnlo ---->
        Estimate for central PDF member: 1204 [events]
        i.e. [14.27 +- 2.574] % relative variation with respect to original PDF
        Number of eigenvectors for uncertainty estimation: 15
        Relative uncertainty with respect to central member: +1.743 / -2.073 [%]
RATE Results for PDF set MRST2007lomod ---->
        Estimate for central PDF member: 1209 [events]
        i.e. [14.77 +- 2.6] % relative variation with respect to original PDF
        NO eigenvectors for uncertainty estimation

>>>>> PDF UNCERTAINTIES ON ACCEPTANCE >>>>>>
ACCEPTANCE Results for PDF set cteq65 ---->
        Estimate for central PDF member acceptance: [47.3181 +- 1.06611] %
        i.e. [1.415 +- 1.106] % relative variation with respect to the original PDF
        Number of eigenvectors for uncertainty estimation: 20
        Relative uncertainty with respect to central member: +0.496 / -0.5614 [%]
ACCEPTANCE Results for PDF set MRST2006nnlo ---->
        Estimate for central PDF member acceptance: [47.3774 +- 1.06734] %
        i.e. [1.542 +- 1.115] % relative variation with respect to the original PDF
        Number of eigenvectors for uncertainty estimation: 15
        Relative uncertainty with respect to central member: +0.2287 / -0.3635 [%]
ACCEPTANCE Results for PDF set MRST2007lomod ---->
        Estimate for central PDF member acceptance: [47.0694 +- 1.06622] %
        i.e. [0.8822 +- 1.188] % relative variation with respect to the original PDF
        NO eigenvectors for uncertainty estimation
>>>> End of PDF weight systematics summary >>>>

INITIAL STATE RADIATION (ISR) SYSTEMATICS (ISRWeightProducer)

Summary

The aim of the ISRWeightProducer is to write into the data structure a product containing an ISR event weight of type "double". This weight can be used to study the variations due to a modification of the Z (or W) pt spectrum at the generator level. The input weights as a function of the boson Pt are proovided via vectors in the configuration file.

SimpleSystematicsAnalyzer is a complete analyzer that provides the systematics on the acceptance due to ISR uncertainties (and to other sources: this analyzer can process many input weights simultaneously).

Theoretical considerations

The input weights must be determined from generator-level studies, by comparing PT boson distributions with different generator settings. Concerning perturbative effects, a sensible approach is a reweight to the distributions for the PYTHIA setting PARP(64)=0.1, which is known to give a better agreement with NLO programs like RESBOS, POWHEG, ..., both for Z and W bosons. Today's PYTHIA reference in CMS MC production is PARP(64) = 0.2. The PYTHIA PARP(64) parameter is able to control the production of additional parton branchngs between the hard interaction and the iniitial proton-proton state. It is equivalent to a modification of the Lambda QCD parameter (PARP(61)) in ISR banchings. PARP(64) (or equivalently PARP(61) restricted to initial state) is the only parameter that seems to provide a significant distortion of the pt boson distribution able to match NLO perturbative predictions.

This should be completed with another set of weights modifying the non-perturbative part of the PT spectrum. Within the PYTHIA scheme, this can be accomplished by means of the PARP(91) and PARP(93) parameters, which control the primordial kT of partons in the proton. This is known to provide significant distortions only for pt<4 GeV (i.e.below the bulk of the event density as a function of PT, which happens above 4 GeV at the LHC). A suggested choice is to assign the variation between the PYTHIA default settings PARP(91)=1, PARP(93)=5 and the D6T settings that are reference for CMS: PARP(91)=2.1, PARP(93)=15.

A summary presentation addressing these issues can be found HERE.

General usage in the analysis

For ISRWeightProducer, the relevant piece of code in the configuration file looks like this:

      # Produce event weights according to generated boson Pt
      # Example corresponds to approximate weights to study
      # systematic effects due to ISR uncertainties (Z/W boson)
      # 
      # Reweighting from PARP(64)=0.2 to PARP(64)=0.1 ...
      process.isrWeight = cms.EDProducer("ISRWeightProducer",
            GenTag = cms.untracked.InputTag("genParticles"),
            ISRBinEdges = cms.untracked.vdouble(
                  0.,  1.,  2.,  3.,  4.,  5.,  6.,  7.,  8.,  9.
                  , 10., 11., 12., 13., 14., 15., 16., 17., 18., 19.
                  , 20., 21., 22., 23., 24., 25., 26., 27., 28., 29.
                  , 30., 31., 32., 33., 34., 35., 36., 37., 38., 39.
                  , 40., 41., 42., 43., 44., 45., 46., 47., 48., 49.
                  , 999999.
            ),
            PtWeights = cms.untracked.vdouble(
                  0.800665, 0.822121, 0.851249, 0.868285, 0.878733
                  , 0.953853, 0.928108, 0.982021, 1.00659 , 1.00648
                  , 1.03218 , 1.04924 , 1.03621 , 1.08743 , 1.01951
                  , 1.10519 , 0.984263, 1.04853 , 1.06724 , 1.10183
                  , 1.0503  , 1.13162 , 1.03837 , 1.12936 , 0.999173
                  , 1.01453 , 1.11435 , 1.10545 , 1.07199 , 1.04542
                  , 1.00828 , 1.0822  , 1.09667 , 1.16144 , 1.13906
                  , 1.27974 , 1.14936 , 1.23235 , 1.06667 , 1.06363
                  , 1.14225 , 1.22955 , 1.12674 , 1.03944 , 1.04639
                  , 1.13667 , 1.20493 , 1.09349 , 1.2107  , 1.21073
            )
      )

More precise PtWeights may be provided in the future.

Full example

The latest reference examples assigning ISR systematics (and more) can be found here: ElectroWeakAnalysis/Utilities/test/SimpleSystematicsAnalyzer.py.

Example of output

For SimpleSystematicsAnalyzer, the ISR resutls may look like this (Wmunu selection):

Total number of analyzed data: 2259 [events]
Total number of selected data: 1054 [events], corresponding to acceptance: [46.6578 +- 1.04964] %
...
Results for Weight Tag: isrWeight ---->
        Total Events after reweighting: 2267.93 [events]
        Events selected after reweighting: 1049.55 [events]
        Acceptance after reweighting: [46.2778 +- 1.05352] %
        i.e. [-0.8145 +- 0.3966] % relative variation with respect to the original acceptance
...
>>>> End of Weight systematics summary >>>>

QED INITIAL STATE RADIATION SYSTEMATICS (ISRGammaWeightProducer)

Summary

ISRGammaWeightProducer writes into the data structure a product containing an event weight of type "double". This weight can be used to study possible uncertainties in the description of QED radiation in the initial state.

SimpleSystematicsAnalyzer is a complete analyzer that provides the systematics on the acceptance due to QED ISR uncertainties (and to other sources: this analyzer can process many input weights simultaneously).

Theoretical considerations

This code accounts for missing initial-state photon radiation effects in Z and W production. Although initial-state radiation effects are dominated by QCD radiation, this effect must be really considered in V+gamma analyses, where photon kinematics is important. It assumes PYTHIA6 as input. The logic is taken from hep-ph/9812455, which is the one used for QCD matrix-element reweighting in W/Z inclusive production. PYTHIA8 should already include this QED ISR reweighting. Note that this code is not only useful to assign systematics, but also to ``correct'' the inclusive W/Z PYTHIA6 Monte Carlos. For instance, it should providie a rather accurate prediction of Z+gamma final sates, essentially in agreement with NLO QCD predictions for this process.

General usage in the analysis

For ISRGammaWeightProducer, the relevant piece of code in the configuration file is extremely simple:
# Produce event weights to estimate missing O(alpha) terms + NLO QED terms
process.isrGammaWeight = cms.EDProducer("ISRGammaWeightProducer",
      GenTag = cms.untracked.InputTag("genParticles"),
)

Full example

The latest reference example assigning QED ISR systematics (and more) can be found here: ElectroWeakAnalysis/Utilities/test/SimpleSystematicsAnalyzer.py.

Example of output

For SimpleSystematicsAnalyzer, the QED ISR resutls may look like this (Wmunu selection):

>>>> Begin of Weight systematics summary >>>>
Total number of analyzed data: 2259 [events]
Total number of selected data: 1054 [events], corresponding to acceptance: [46.6578 +- 1.04964] %
...
Results for Weight Tag: isrGammaWeight ---->
        Total Events after reweighting: 2258.03 [events]
        Events selected after reweighting: 1053.75 [events]
        Acceptance after reweighting: [46.6667 +- 1.04985] %
        i.e. [0.01895 +- 0] % relative variation with respect to the original acceptance
...
>>>> End of Weight systematics summary >>>>

QED FINAL STATE RADIATION (FSR) SYSTEMATICS (FSRWeightProducer)

Summary

FSRWeightProducer writes into the data structure a product containing an event weight of type "double". This weight can be used to study possible uncertainties in the description of QED radiation in the final state.

SimpleSystematicsAnalyzer is a complete analyzer that provides the systematics on the acceptance due to QED FSR uncertainties (and to other sources: this analyzer can process many input weights simultaneously).

This code may be improved in the future.

Theoretical considerations

The current code assumes that final state QED radiation is implemented within a parton-shower approach (i.e. a la PYTHIA). It is evident that different recipes will be required for other Monte Carlo approaches, like PHOTOS or HORACE. The PYTHIA approach leads to small QED FSR uncertainties ( < 1%) for inclusive cross-section measurements, but larger uncertainties for studies where a precise knowledge of the lepton kinematics is necessary (W mass measurement) or for analyses with high-PT photons, isolated from the leptons (V-gamma final states).

Currently we determine the FSR weight by correcting the following missing effects:

  • Go from a soft-collinear approach to the exact O(alpha) result. This is done by using the same weight used in PHOTOS (hep-ph/0303260). It is only relevant for W, and very small for low-PT or collinear radiation.
  • Use alpha(pt**2) coupling instead of just alpha(0) to absorb some higher-order corrections.

Note that this code does not account for the following effects: a) initial state photon radiation, b) initial-final state photon interference. Both effects are negligible in inclusive W/Z cross section analyses: a) is included in ISR uncertainties, dominated by QCD radiation; b) is negligible around W and Z resonances. Effect a) is relevant in V+gamma analyses, where photon kinematics is important, and it is discussed above. Effect b) may be important in Drell-Yan studies away from the Z mass, and must be studied with other full QED Monte Carlo approaches, like for instance the one of HORACE.

General usage in the analysis

For FSRWeightProducer, the relevant piece of code in the configuration file is extremely simple:
# Produce event weights to estimate missing O(alpha) terms + NLO QED terms
process.fsrWeight = cms.EDProducer("FSRWeightProducer",
      GenTag = cms.untracked.InputTag("genParticles"),
)

Full example

The latest reference example assigning Weak systematics (and more) can be found here: ElectroWeakAnalysis/Utilities/test/SimpleSystematicsAnalyzer.py.

Example of output

For SimpleSystematicsAnalyzer, the FSR resutls may look like this (Wmunu selection):

>>>> Begin of Weight systematics summary >>>>
Total number of analyzed data: 2259 [events]
Total number of selected data: 1054 [events], corresponding to acceptance: [46.6578 +- 1.04964] %
...
Results for Weight Tag: fsrWeight ---->
        Total Events after reweighting: 2261.98 [events]
        Events selected after reweighting: 1055.08 [events]
        Acceptance after reweighting: [46.644 +- 1.04934] %
        i.e. [-0.02967 +- 0.1123] % relative variation with respect to the original acceptance
...
>>>> End of Weight systematics summary >>>>

UNCERTAINTIES DUE TO WEAK EFFECTS (WeakEffectsWeightProducer)

Summary

WeakEffectsWeightProducer also writes into the data structure a product containing an event weight of type "double". This weight can be used to study the effect of missing weak contributions in the generation of Drell-Yan events. A future extension will also consider weak effects in W production.

Theoretical considerations

Usual Monte Carlo generators deal with Drell-Yan production in a Born-improved approximation. From the practical point of view, this implies working with Born-level expressions in which the QED coupling constant is determined at the scale of the Z mass and weak couplings are calculated for an effective value of sin2w of around 0.232. This provides an accurate ``weak description'' - at the percent level - of Drell-Yan production around the Z peak. In order to improve things, we can estimate the influence of two missing weak effects:

  • Non-universal corrections around the Z resonance: they can be accounted for by introducing an effective rho parameter in the Z couplings of 1.004 for all fermions (except for the b quark, but it has a minimal contribution at the parton level). This enchances the Drell-Yan cross section around the Z pole by 0.4% and leaves the low-mass Drell-Yan region unchanged. See for instance Baur et al., PRD 65, 033007 for further details on this missing effect.
  • Sudakov factors at high invariant mass: they essentially correspond to box-diagram corrections that grow as log(s/MZ**2), log(s/MW**2)), and are therefore relevant only for large Drell-Yan masses. We actually take just the leading order approximation in this logarithmic expansion from hep-ph/9809321.

General usage in the analysis

For WeightEffectsWeightProducer, the relevant piece of code in the configuration file looks like this:

# Produce event weights to estimate missing weak terms (=> include missing rho factor for Z diagrams)
process.weakWeight = cms.EDProducer("WeakEffectsWeightProducer",
      GenParticlesTag = cms.untracked.InputTag("genParticles"),
      RhoParameter = cms.untracked.double(1.004)
)

A previous step in practice is the selection of Drell-Yan events around the Z mass. This can be obtained with the configuration file in ElectroWeakUtilities/python/zmmSelection_cfi.py:

process.load("ElectroWeakAnalysis/Utilities/zmmSelection_cfi")
process.myPath(
        process.goldenZMMSelectionSequence
      *process.weakWeight
)

MUON CORRECTIONS TO MC (DistortedMuonProducer*)

Summary

The aim of the DistortedMuonProducer and DistortedMuonProducerFromDB is to write into the data structure a new reco::Muon collection distorted by smearing and efficiency effects. This new collection can be easily used as a realistic Monte Carlo prediction for the latest steps of the analysis. The input distortions are provided in the configuration file, either ad-hoc (DistortedMuonProducer) or from database (DistortedMuonProducerFromDB).

General usage in the analysis

For DistoretedMuonProducer, the relevant piece of code in the configuration file looks like this:

# GEN-REC muon matching
process.genMatchMap = cms.EDFilter("MCTruthDeltaRMatcherNew",
    src = cms.InputTag("muons"),
    matched = cms.InputTag("genParticles"),
    distMin = cms.double(0.15),
    matchPDGId = cms.vint32(13)
)   
    
# Create a new "distorted" Muon collection
process.distortedMuons = cms.EDFilter("DistortedMuonProducer",
      MuonTag = cms.untracked.InputTag("muons"),
      GenMatchTag = cms.untracked.InputTag("genMatchMap"),

      EtaBinEdges = cms.untracked.vdouble(-2.1,-1.2, 1.2, 2.1), # one more entry than next vectors

      ShiftOnOneOverPt = cms.untracked.vdouble(1.e-3), #in [1/GeV] units
      RelativeShiftOnPt = cms.untracked.vdouble(0.), # relative
      UncertaintyOnOneOverPt = cms.untracked.vdouble(3.e-3, 2.e-3, 3.e-3), #in [1/GeV] units
      RelativeUncertaintyOnPt = cms.untracked.vdouble(5.e-3, 5.e-3, 5.e-3), # relative

      EfficiencyRatioOverMC = cms.untracked.vdouble(0.90, 0.90, 0.90)
)

The different fields are almost self-explanatory. The idea is to use the collection "distortedMuons" instead of "muons" at the latest stages of the analysis.

The present version assumes that efficiencies and resolutions only change as a function of "eta". A possible extension in the future will allow to use distortions that also depend on "phi". Note that, in order to include distortions properly at the generator level, we perform a previous matching between reconstructed and generated muons using standard CMSSW matching utilities.

Inefficiencies are implemented by dropping muons in the event according to a hit-or-miss method.Therefore, all components of the vector field "EfficiencyRatioOverMC" must be less or equal than 1, i.e. we assume that the MC simulation provides a better muon efficiency than data everywhere.

For DistortedMuonProducerFromDB, the relevant piece of code in the configuration file looks like this:

# Create a new "distorted" Muon collection
process.distortedMuons = cms.EDFilter("DistortedMuonProducerFromDB",
      MuonTag = cms.untracked.InputTag("muons"),

      DBScaleLabel = cms.untracked.string(''),
      DBDataResolutionLabel = cms.untracked.string('MuScleFit_Resol_OctoberExercise_EWK_InnerTrack_WithLabel'),
      DBMCResolutionLabel = cms.untracked.string('MuScleFit_Resol_OctoberExercise_SherpaIdealMC_WithLabel'),
)

At present, the code only handles momentum resolution distortions (and not inefficiencies).

Since all the information is taken from database, the previous piece must be accompanied by the corresponding "PoolDBESSource"" statements, i.e.:

# Database source statements
process.load("CondCore.DBCommon.CondDBCommon_cfi")
process.poolDBESSource1 = cms.ESSource("PoolDBESSource",
      BlobStreamerName = cms.untracked.string('TBufferBlobStreamingService'),
      DBParameters = cms.PSet(
            messageLevel = cms.untracked.int32(2)
      ),
      timetype = cms.untracked.string('runnumber'),
      connect = cms.string('frontier://FrontierPrep/CMS_COND_PHYSICSTOOLS'),
      toGet = cms.VPSet(
            cms.PSet(
                  record = cms.string('MuScleFitDBobjectRcd'),
                  tag = cms.string('MuScleFit_Scale_OctoberExercise_EWK_InnerTrack'),
                  label = cms.untracked.string('')
            )
      )
)
process.poolDBESSource2 = cms.ESSource("PoolDBESSource",
      BlobStreamerName = cms.untracked.string('TBufferBlobStreamingService'),
      DBParameters = cms.PSet(
            messageLevel = cms.untracked.int32(2)
      ),
      timetype = cms.untracked.string('runnumber'),
      connect = cms.string('frontier://FrontierPrep/CMS_COND_PHYSICSTOOLS'),
      toGet = cms.VPSet(
            cms.PSet(
                  record = cms.string('MuScleFitDBobjectRcd'),
                  tag = cms.string('MuScleFit_Resol_OctoberExercise_EWK_InnerTrack_WithLabel'),
                  label = cms.untracked.string('MuScleFit_Resol_OctoberExercise_EWK_InnerTrack_WithLabel')
            )
      )
)
process.poolDBESSource3 = cms.ESSource("PoolDBESSource",
      BlobStreamerName = cms.untracked.string('TBufferBlobStreamingService'),
      DBParameters = cms.PSet(
            messageLevel = cms.untracked.int32(2)
      ),
      timetype = cms.untracked.string('runnumber'),
      connect = cms.string('frontier://FrontierPrep/CMS_COND_PHYSICSTOOLS'),
      toGet = cms.VPSet(
            cms.PSet(
                  record = cms.string('MuScleFitDBobjectRcd'),
                  tag = cms.string('MuScleFit_Resol_OctoberExercise_SherpaIdealMC_WithLabel'),
                  label = cms.untracked.string('MuScleFit_Resol_OctoberExercise_SherpaIdealMC_WithLabel')
            )
      )
)

Full examples

The latest reference examples can be found here: ElectroWeakAnalysis/Utilities/test/distortedMuons.py and ElectroWeakAnalysis/Utilities/test/distortedMuonsFromDB.py.

MET CORRECTIONS TO MC (DistortedMETProducer)

Summary

The aim of the DistortedMETProducer is to write into the data structure a new reco::MET collection after distortion. This new collection can be used as a realistic Monte Carlo prediction for (muon-uncorrected) MET at the latest steps of the analysis. At present, only the simplest possibility of scaling MET and sum(ET) by a global scale factor is implemented, and no matching with GenMET has been tried yet.

This utility must be considered as a rather naive approach. It expected that the CMS.JetMET group will provide more sophisticated tools to distort MET, and that most of the MET systematics will be extracted from comparisons between calorimetric MET, particle-flow MET and track-corrected MET.

General usage in the analysis

For DistortedMETProducer, the relevant piece of code in the configuration file looks like this:

      process.distortedMET = cms.EDFilter("DistortedMETProducer",
        MetTag = cms.untracked.InputTag("met"),
        MetScaleShift = cms.untracked.double(0.1)
      )

The paarmeter "MetScaleShift" must be interpreted as a relative shift. The idea is to use the collection "distortedMET" instead of "met" at the latest stages of the analysis. A possible extension in the future will be to distort MET by making use of the MET significance.

Full example

The latest reference example can be found here: ElectroWeakAnalysis/Utilities/test/distortedMET.py.

Links

Presentation at the EWK Muon Meeting

Contact

Please mail any comments, critics, requests to juan.alcaraz@cernNOSPAMPLEASE.ch

Responsible: Juan Alcaraz

Edit | Attach | Watch | Print version | History: r33 < r32 < r31 < r30 < r29 | Backlinks | Raw View | WYSIWYG | More topic actions
Topic revision: r33 - 2016-03-12 - IlyaGorbunov
 
    • Cern Search Icon Cern Search
    • TWiki Search Icon TWiki Search
    • Google Search Icon Google Search

    CMSPublic All webs login

This site is powered by the TWiki collaboration platform Powered by PerlCopyright & 2008-2019 by the contributing authors. All material on this collaboration platform is the property of the contributing authors.
Ideas, requests, problems regarding TWiki? Send feedback