Changes in / [70b9632:7bb13cd] in git

Files:

• CHANGELOG (modified) (1 diff)
• Makefile (modified) (17 diffs)
• README (modified) (1 diff)
• VERSION (modified) (1 diff)
• cards/CMS_PhaseII/CMS_PhaseII_Substructure_PIX4022_0PU.tcl (modified) (24 diffs)
• cards/CMS_PhaseII/CMS_PhaseII_Substructure_PIX4022_200PU.tcl (modified) (25 diffs)
• cards/converter_card.tcl (modified) (1 diff)
• cards/delphes_card_CMS.tcl (modified) (4 diffs)
• cards/delphes_card_CMS_simple.tcl (deleted)
• classes/DelphesClasses.cc (modified) (8 diffs)
• classes/DelphesClasses.h (modified) (10 diffs)
• classes/SortableObject.h (modified) (1 diff)
• doc/genMakefile.tcl (modified) (1 diff)
• examples/Pythia8/configParticleGun.cmnd (modified) (1 diff)
• examples/Validation.cpp (modified) (6 diffs)
• modules/BeamSpotFilter.cc (deleted)
• modules/BeamSpotFilter.h (deleted)
• modules/Calorimeter.cc (modified) (8 diffs)
• modules/Calorimeter.h (modified) (3 diffs)
• modules/ImpactParameterSmearing.cc (modified) (4 diffs)
• modules/ModulesLinkDef.h (modified) (6 diffs)
• modules/ParticlePropagator.cc (modified) (8 diffs)
• modules/ParticlePropagator.h (modified) (2 diffs)
• modules/PhotonConversions.cc (modified) (2 diffs)
• modules/PileUpMerger.cc (modified) (9 diffs)
• modules/RecoPuFilter.cc (deleted)
• modules/RecoPuFilter.h (deleted)
• modules/SimpleCalorimeter.cc (modified) (10 diffs)
• modules/SimpleCalorimeter.h (modified) (3 diffs)
• modules/TimeSmearing.cc (modified) (2 diffs)
• modules/TrackCountingBTagging.cc (modified) (2 diffs)
• modules/TrackSmearing.cc (deleted)
• modules/TrackSmearing.h (deleted)
• modules/TreeWriter.cc (modified) (7 diffs)
• modules/VertexFinder.cc (deleted)
• modules/VertexFinder.h (deleted)
• modules/VertexFinderDA4D.cc (deleted)
• modules/VertexFinderDA4D.h (deleted)
• modules/VertexSorter.cc (deleted)
• modules/VertexSorter.h (deleted)
• readers/DelphesCMSFWLite.cpp (modified) (7 diffs)
• readers/DelphesROOT.cpp (deleted)
• validation.sh (deleted)

CHANGELOG

@@ -1 @@
-3.3.3:
-- improved compatibility with ROOT >= 6.04
-- removed test of praticle stability based on ROOT pdgtable (#347, #753, #821, #855)
-- fixed bugs in DelphesCMSFWLite
-- fixed bugs in PUPPI
-- fixed compiler and linker flags for Pythia8
-- added CMS Phase II cards
-- added particle gun to DelphesPythia8
-- added jet merging to DelphesPythia8 based on Pythia8 main32.cc
-- added UseMiniCone option in Isolation module
-- added RecoPuFilter module
-- added back OldCalorimeter for CheckMate bwd compatibilty (#591)
-- updated tracking and calorimeter resolution in CMS cards
-
 3.3.2:
 - added pre-validated card for CMS upgrade Phase II studies

Makefile

@@ -250 +250 @@
         external/ExRootAnalysis/ExRootTreeBranch.h \
         external/ExRootAnalysis/ExRootProgressBar.h
-DelphesROOT$(ExeSuf): \
-        tmp/readers/DelphesROOT.$(ObjSuf)
-
-tmp/readers/DelphesROOT.$(ObjSuf): \
-        readers/DelphesROOT.cpp \
-        modules/Delphes.h \
-        classes/DelphesStream.h \
-        classes/DelphesClasses.h \
-        classes/DelphesFactory.h \
-        external/ExRootAnalysis/ExRootTreeWriter.h \
-        external/ExRootAnalysis/ExRootTreeReader.h \
-        external/ExRootAnalysis/ExRootTreeBranch.h \
-        external/ExRootAnalysis/ExRootProgressBar.h
 DelphesSTDHEP$(ExeSuf): \
         tmp/readers/DelphesSTDHEP.$(ObjSuf)
@@ -278 +265 @@
         DelphesHepMC$(ExeSuf) \
         DelphesLHEF$(ExeSuf) \
-        DelphesROOT$(ExeSuf) \
         DelphesSTDHEP$(ExeSuf)
 
@@ -284 +270 @@
         tmp/readers/DelphesHepMC.$(ObjSuf) \
         tmp/readers/DelphesLHEF.$(ObjSuf) \
-        tmp/readers/DelphesROOT.$(ObjSuf) \
         tmp/readers/DelphesSTDHEP.$(ObjSuf)
 
@@ -399 +384 @@
         modules/EnergySmearing.h \
         modules/MomentumSmearing.h \
-        modules/TrackSmearing.h \
         modules/ImpactParameterSmearing.h \
         modules/TimeSmearing.h \
@@ -423 +407 @@
         modules/StatusPidFilter.h \
         modules/PdgCodeFilter.h \
-        modules/BeamSpotFilter.h \
-        modules/RecoPuFilter.h \
         modules/Cloner.h \
         modules/Weighter.h \
@@ -430 +412 @@
         modules/JetFlavorAssociation.h \
         modules/JetFakeParticle.h \
-        modules/VertexSorter.h \
-        modules/VertexFinder.h \
-        modules/VertexFinderDA4D.h \
         modules/ExampleModule.h
 tmp/modules/ModulesDict$(PcmSuf): \
@@ -639 +618 @@
         classes/DelphesFactory.h \
         classes/DelphesFormula.h
-tmp/modules/BeamSpotFilter.$(ObjSuf): \
-        modules/BeamSpotFilter.$(SrcSuf) \
-        modules/BeamSpotFilter.h \
-        classes/DelphesClasses.h \
-        classes/DelphesFactory.h \
-        classes/DelphesFormula.h \
-        external/ExRootAnalysis/ExRootResult.h \
-        external/ExRootAnalysis/ExRootFilter.h \
-        external/ExRootAnalysis/ExRootClassifier.h
 tmp/modules/Calorimeter.$(ObjSuf): \
         modules/Calorimeter.$(SrcSuf) \
@@ -880 +850 @@
         external/ExRootAnalysis/ExRootFilter.h \
         external/ExRootAnalysis/ExRootClassifier.h
-tmp/modules/RecoPuFilter.$(ObjSuf): \
-        modules/RecoPuFilter.$(SrcSuf) \
-        modules/RecoPuFilter.h \
-        classes/DelphesClasses.h \
-        classes/DelphesFactory.h \
-        classes/DelphesFormula.h \
-        external/ExRootAnalysis/ExRootResult.h \
-        external/ExRootAnalysis/ExRootFilter.h \
-        external/ExRootAnalysis/ExRootClassifier.h
 tmp/modules/SimpleCalorimeter.$(ObjSuf): \
         modules/SimpleCalorimeter.$(SrcSuf) \
@@ -950 +911 @@
         modules/TrackPileUpSubtractor.$(SrcSuf) \
         modules/TrackPileUpSubtractor.h \
-        classes/DelphesClasses.h \
-        classes/DelphesFactory.h \
-        classes/DelphesFormula.h \
-        external/ExRootAnalysis/ExRootResult.h \
-        external/ExRootAnalysis/ExRootFilter.h \
-        external/ExRootAnalysis/ExRootClassifier.h
-tmp/modules/TrackSmearing.$(ObjSuf): \
-        modules/TrackSmearing.$(SrcSuf) \
-        modules/TrackSmearing.h \
         classes/DelphesClasses.h \
         classes/DelphesFactory.h \
@@ -981 +933 @@
         classes/DelphesFactory.h \
         classes/DelphesFormula.h \
-        external/ExRootAnalysis/ExRootResult.h \
-        external/ExRootAnalysis/ExRootFilter.h \
-        external/ExRootAnalysis/ExRootClassifier.h
-tmp/modules/VertexFinder.$(ObjSuf): \
-        modules/VertexFinder.$(SrcSuf) \
-        modules/VertexFinder.h \
-        classes/DelphesClasses.h \
-        classes/DelphesFactory.h \
-        classes/DelphesFormula.h \
-        classes/DelphesPileUpReader.h \
-        external/ExRootAnalysis/ExRootResult.h \
-        external/ExRootAnalysis/ExRootFilter.h \
-        external/ExRootAnalysis/ExRootClassifier.h
-tmp/modules/VertexFinderDA4D.$(ObjSuf): \
-        modules/VertexFinderDA4D.$(SrcSuf) \
-        modules/VertexFinderDA4D.h \
-        classes/DelphesClasses.h \
-        classes/DelphesFactory.h \
-        classes/DelphesFormula.h \
-        classes/DelphesPileUpReader.h \
-        external/ExRootAnalysis/ExRootResult.h \
-        external/ExRootAnalysis/ExRootFilter.h \
-        external/ExRootAnalysis/ExRootClassifier.h
-tmp/modules/VertexSorter.$(ObjSuf): \
-        modules/VertexSorter.$(SrcSuf) \
-        modules/VertexSorter.h \
-        classes/DelphesClasses.h \
-        classes/DelphesFactory.h \
-        classes/DelphesFormula.h \
-        classes/DelphesPileUpReader.h \
         external/ExRootAnalysis/ExRootResult.h \
         external/ExRootAnalysis/ExRootFilter.h \
@@ -1074 +996 @@
         tmp/modules/AngularSmearing.$(ObjSuf) \
         tmp/modules/BTagging.$(ObjSuf) \
-        tmp/modules/BeamSpotFilter.$(ObjSuf) \
         tmp/modules/Calorimeter.$(ObjSuf) \
         tmp/modules/Cloner.$(ObjSuf) \
@@ -1099 +1020 @@
         tmp/modules/PileUpJetID.$(ObjSuf) \
         tmp/modules/PileUpMerger.$(ObjSuf) \
-        tmp/modules/RecoPuFilter.$(ObjSuf) \
         tmp/modules/SimpleCalorimeter.$(ObjSuf) \
         tmp/modules/StatusPidFilter.$(ObjSuf) \
@@ -1108 +1028 @@
         tmp/modules/TrackCountingTauTagging.$(ObjSuf) \
         tmp/modules/TrackPileUpSubtractor.$(ObjSuf) \
-        tmp/modules/TrackSmearing.$(ObjSuf) \
         tmp/modules/TreeWriter.$(ObjSuf) \
         tmp/modules/UniqueObjectFinder.$(ObjSuf) \
-        tmp/modules/VertexFinder.$(ObjSuf) \
-        tmp/modules/VertexFinderDA4D.$(ObjSuf) \
-        tmp/modules/VertexSorter.$(ObjSuf) \
         tmp/modules/Weighter.$(ObjSuf)
 
@@ -1718 +1634 @@
         tmp/external/tcl/tclVar.$(ObjSuf)
 
-modules/VertexFinderDA4D.h: \
-        classes/DelphesModule.h \
-        classes/DelphesClasses.h
-        @touch $@
-
-modules/TrackSmearing.h: \
-        classes/DelphesModule.h
-        @touch $@
-
 external/fastjet/ClusterSequence.hh: \
         external/fastjet/PseudoJet.hh \
@@ -1987 +1894 @@
         @touch $@
 
-modules/VertexSorter.h: \
-        classes/DelphesModule.h \
-        classes/DelphesClasses.h
-        @touch $@
-
 modules/Delphes.h: \
         classes/DelphesModule.h
-        @touch $@
-
-modules/VertexFinder.h: \
-        classes/DelphesModule.h \
-        classes/DelphesClasses.h
         @touch $@
 
@@ -2070 +1967 @@
         @touch $@
 
-modules/RecoPuFilter.h: \
-        classes/DelphesModule.h
-        @touch $@
-
 modules/Hector.h: \
         classes/DelphesModule.h
@@ -2173 +2066 @@
 
 modules/FastJetFinder.h: \
-        classes/DelphesModule.h
-        @touch $@
-
-modules/BeamSpotFilter.h: \
         classes/DelphesModule.h
         @touch $@

README

@@ -4 +4 @@
 Commands to get the code:
 
-   wget http://cp3.irmp.ucl.ac.be/downloads/Delphes-3.3.3.tar.gz
+   wget http://cp3.irmp.ucl.ac.be/downloads/Delphes-3.3.2.tar.gz
 
-   tar -zxf Delphes-3.3.3.tar.gz
+   tar -zxf Delphes-3.3.2.tar.gz
 
 Commands to compile the code:
 
-   cd Delphes-3.3.3
+   cd Delphes-3.3.2
 
    make

VERSION

@@ -1 @@
-3.3.3
+3.3.2

cards/CMS_PhaseII/CMS_PhaseII_Substructure_PIX4022_0PU.tcl

@@ -33 +33 @@
   PhotonEnergySmearing
   ElectronFilter
-
   TrackPileUpSubtractor
-  RecoPuFilter
 
   TowerMerger
   NeutralEFlowMerger
-
+  EFlowMergerAllTracks
   EFlowMerger
-  EFlowMergerCHS
-  Rho
 
   LeptonFilterNoLep
@@ -52 +48 @@
 
   PhotonIsolation
-  PhotonIsolationCHS
   PhotonEfficiency
-  PhotonEfficiencyCHS
 
   ElectronIsolation
-  ElectronIsolationCHS
-
   ElectronEfficiency
-  ElectronEfficiencyCHS
 
   MuonIsolation
-  MuonIsolationCHS
-
   MuonLooseIdEfficiency
   MuonTightIdEfficiency
-
-  MuonLooseIdEfficiencyCHS
-  MuonTightIdEfficiencyCHS
 
   NeutrinoFilter
@@ -82 +68 @@
   FastJetFinder
   FastJetFinderAK8
-  JetPileUpSubtractor
-  JetPileUpSubtractorAK8
   FastJetFinderPUPPI
   FastJetFinderPUPPIAK8
@@ -137 +121 @@
 
   # maximum spread in the beam direction in m
-  set ZVertexSpread 0.25
+  set ZVertexSpread 0
 
   # maximum spread in time in s
-  set TVertexSpread 800E-12
+  set TVertexSpread 0
 
   # vertex smearing formula f(z,t) (z,t need to be respectively given in m,s) - {exp(-(t^2/160e-12^2/2))*exp(-(z^2/0.053^2/2))}
-  set VertexDistributionFormula {exp(-(t^2/160e-12^2/2))*exp(-(z^2/0.053^2/2))}
+  set VertexDistributionFormula 1
 
 }
@@ -289 +273 @@
   source muonMomentumResolution.tcl
 }
+
+
+
 
 ##############
@@ -540 +527 @@
 }
 
-########################
-# Reco PU filter
-########################
-
-module RecoPuFilter RecoPuFilter {
-  set InputArray HCal/eflowTracks
-  set OutputArray eflowTracks
-}
 
 ###################################################
@@ -560 +539 @@
 }
 
+
 ####################
 # Neutral eflow erger
@@ -571 +551 @@
 }
 
-#####################
+
+####################
 # Energy flow merger
-#####################
+####################
 
 module Merger EFlowMerger {
@@ -583 +564 @@
 }
 
-############################
-# Energy flow merger no PU
-############################
-
-module Merger EFlowMergerCHS {
+##################################
+# Energy flow merger (all tracks)
+##################################
+
+module Merger EFlowMergerAllTracks {
 # add InputArray InputArray
-  add InputArray RecoPuFilter/eflowTracks
+  add InputArray TrackMerger/tracks
   add InputArray PhotonEnergySmearing/eflowPhotons
   add InputArray HCal/eflowNeutralHadrons
@@ -714 +695 @@
 }
 
+
 #####################
 # MC truth jet finder
@@ -753 +735 @@
 }
 
-
-#############
-# Rho pile-up
-#############
-
-module FastJetFinder Rho {
-#  set InputArray Calorimeter/towers
-  set InputArray EFlowMergerCHS/eflow
-
-  set ComputeRho true
-  set RhoOutputArray rho
-
-  # area algorithm: 0 Do not compute area, 1 Active area explicit ghosts, 2 One ghost passive area, 3 Passive area, 4 Voronoi, 5 Active area
-  set AreaAlgorithm 5
-
-  # jet algorithm: 1 CDFJetClu, 2 MidPoint, 3 SIScone, 4 kt, 5 Cambridge/Aachen, 6 antikt
-  set JetAlgorithm 4
-  set ParameterR 0.4
-  set GhostEtaMax 5.0
-
-  add RhoEtaRange -5.0 -4.0
-  add RhoEtaRange -4.0 -1.5
-  add RhoEtaRange -1.5 1.5
-  add RhoEtaRange 1.5 4.0
-  add RhoEtaRange 4.0 5.0
-
-  set JetPTMin 0.0
-}
-
-
-##############
+############
 # Jet finder
-##############
+############
 
 module FastJetFinder FastJetFinder {
 #  set InputArray TowerMerger/towers
-  set InputArray EFlowMergerCHS/eflow
+  set InputArray EFlowMerger/eflow
 
   set OutputArray jets
-
-  set AreaAlgorithm 5
 
   # algorithm: 1 CDFJetClu, 2 MidPoint, 3 SIScone, 4 kt, 5 Cambridge/Aachen, 6 antikt
@@ -805 +755 @@
 module FastJetFinder FastJetFinderAK8 {
 #  set InputArray TowerMerger/towers
-  set InputArray EFlowMergerCHS/eflow
+  set InputArray EFlowMerger/eflow
 
   set OutputArray jets
-
-  set AreaAlgorithm 5
 
   # algorithm: 1 CDFJetClu, 2 MidPoint, 3 SIScone, 4 kt, 5 Cambridge/Aachen, 6 antikt
@@ -836 +784 @@
 }
 
-###########################
-# Jet Pile-Up Subtraction
-###########################
-
-module JetPileUpSubtractor JetPileUpSubtractor {
-  set JetInputArray FastJetFinder/jets
-  set RhoInputArray Rho/rho
-
-  set OutputArray jets
-
-  set JetPTMin 15.0
-}
-
-##############################
-# Jet Pile-Up Subtraction AK8
-##############################
-
-module JetPileUpSubtractor JetPileUpSubtractorAK8 {
-  set JetInputArray FastJetFinderAK8/jets
-  set RhoInputArray Rho/rho
-
-  set OutputArray jets
-
-  set JetPTMin 15.0
-}
-
 module FastJetFinder FastJetFinderPUPPI {
 #  set InputArray TowerMerger/towers
@@ -911 +833 @@
 
 module EnergyScale JetEnergyScale {
-  set InputArray JetPileUpSubtractor/jets
+  set InputArray FastJetFinder/jets
   set OutputArray jets
 
@@ -919 +841 @@
 
 module EnergyScale JetEnergyScaleAK8 {
-  set InputArray JetPileUpSubtractorAK8/jets
+  set InputArray FastJetFinderAK8/jets
   set OutputArray jets
 
@@ -985 +907 @@
 }
 
-
-########################
-# Photon isolation CHS #
-########################
-
-module Isolation PhotonIsolationCHS {
-
-  # particle for which calculate the isolation
-  set CandidateInputArray PhotonFilter/photons
-
-  # isolation collection
-  set IsolationInputArray EFlowMerger/eflow
-
-  # output array
-  set OutputArray photons
-
-  # isolation cone
-  set DeltaRMax 0.3
-
-  # minimum pT
-  set PTMin     1.0
-
-  # iso ratio to cut
-  set PTRatioMax 9999.
-
-}
 
 
@@ -1033 +929 @@
 
 
-#####################
-# Photon efficiency #
-#####################
-
-module Efficiency PhotonEfficiencyCHS {
-
-  ## input particles
-  set InputArray PhotonIsolationCHS/photons
-  ## output particles
-  set OutputArray photons
-  # set EfficiencyFormula {efficiency formula as a function of eta and pt}
-  # efficiency formula for photons
-  set EfficiencyFormula {                      (pt <= 10.0) * (0.00) + \
-                           (abs(eta) <= 1.5) * (pt > 10.0)  * (0.9635) + \
-         (abs(eta) > 1.5 && abs(eta) <= 4.0) * (pt > 10.0)  * (0.9624) + \
-         (abs(eta) > 4.0)                                   * (0.00)}
-
-}
-
 ######################
 # Electron isolation #
@@ -1072 +949 @@
 }
 
-
-##########################
-# Electron isolation CHS #
-##########################
-
-module Isolation ElectronIsolationCHS {
-
-  set CandidateInputArray ElectronFilter/electrons
-
-  # isolation collection
-  set IsolationInputArray EFlowMerger/eflow
-
-  set OutputArray electrons
-
-  set DeltaRMax 0.3
-  set PTMin 1.0
-  set PTRatioMax 9999.
-
-}
 
 
@@ -1137 +995 @@
 }
 
-###########################
-# Electron efficiency CHS #
-###########################
-
-module Efficiency ElectronEfficiencyCHS {
-
-  set InputArray ElectronIsolationCHS/electrons
-  set OutputArray electrons
-
-  # set EfficiencyFormula {efficiency formula as a function of eta and pt}
-  # efficiency formula for electrons
-  set EfficiencyFormula {
-                                      (pt <= 4.0)  * (0.00) + \
-                         (abs(eta) <= 1.45 ) * (pt >  4.0 && pt <= 6.0)   * (0.50) + \
-                         (abs(eta) <= 1.45 ) * (pt >  6.0 && pt <= 8.0)   * (0.70) + \
-                         (abs(eta) <= 1.45 ) * (pt >  8.0 && pt <= 10.0)  * (0.85) + \
-                         (abs(eta) <= 1.45 ) * (pt > 10.0 && pt <= 30.0)  * (0.94) + \
-                         (abs(eta) <= 1.45 ) * (pt > 30.0 && pt <= 50.0)  * (0.97) + \
-                         (abs(eta) <= 1.45 ) * (pt > 50.0 && pt <= 70.0)  * (0.98) + \
-                         (abs(eta) <= 1.45 ) * (pt > 70.0 )  * (1.0) + \
-                         (abs(eta) > 1.45  && abs(eta) <= 1.55) * (pt >  4.0 && pt <= 10.0)   * (0.35) + \
-                         (abs(eta) > 1.45  && abs(eta) <= 1.55) * (pt > 10.0 && pt <= 30.0)   * (0.40) + \
-                         (abs(eta) > 1.45  && abs(eta) <= 1.55) * (pt > 30.0 && pt <= 70.0)   * (0.45) + \
-                         (abs(eta) > 1.45  && abs(eta) <= 1.55) * (pt > 70.0 )  * (0.55) + \
-                         (abs(eta) >= 1.55 && abs(eta) <= 2.0 ) * (pt >  4.0 && pt <= 10.0)  * (0.75) + \
-                         (abs(eta) >= 1.55 && abs(eta) <= 2.0 ) * (pt > 10.0 && pt <= 30.0)  * (0.85) + \
-                         (abs(eta) >= 1.55 && abs(eta) <= 2.0 ) * (pt > 30.0 && pt <= 50.0)  * (0.95) + \
-                         (abs(eta) >= 1.55 && abs(eta) <= 2.0 ) * (pt > 50.0 && pt <= 70.0)  * (0.95) + \
-                         (abs(eta) >= 1.55 && abs(eta) <= 2.0 ) * (pt > 70.0 )  * (1.0) + \
-                         (abs(eta) >= 2.0 && abs(eta) <= 2.5 ) * (pt >  4.0 && pt <= 10.0)  * (0.65) + \
-                         (abs(eta) >= 2.0 && abs(eta) <= 2.5 ) * (pt > 10.0 && pt <= 30.0)  * (0.75) + \
-                         (abs(eta) >= 2.0 && abs(eta) <= 2.5 ) * (pt > 30.0 && pt <= 50.0)  * (0.90) + \
-                         (abs(eta) >= 2.0 && abs(eta) <= 2.5 ) * (pt > 50.0 && pt <= 70.0)  * (0.90) + \
-                         (abs(eta) >= 2.0 && abs(eta) <= 2.5 ) * (pt > 70.0 )  * (0.90) + \
-                         (abs(eta) > 2.5 && abs(eta) <= 4.0 ) * (pt > 4.0 && pt <= 10.0) * (0.65) + \
-                                          (abs(eta) > 2.5 && abs(eta) <= 4.0 ) * (pt > 10.0 && pt <= 30.0) * (0.75) + \
-                                          (abs(eta) > 2.5 && abs(eta) <= 4.0 ) * (pt > 30.0 && pt <= 50.0) * (0.90) + \
-                                          (abs(eta) > 2.5 && abs(eta) <= 4.0 ) * (pt > 50.0 && pt <= 70.0) * (0.90) + \
-                                          (abs(eta) > 2.5 && abs(eta) <= 4.0 ) * (pt > 70.0 ) * (0.90) + \
-                                          (abs(eta) > 4.0) * (0.00)
-
-  }
-}
-
-
 ##################
 # Muon isolation #
@@ -1200 +1013 @@
 }
 
-######################
-# Muon isolation CHS #
-######################
-
-module Isolation MuonIsolationCHS {
-  set CandidateInputArray MuonMomentumSmearing/muons
-
-  # isolation collection
-  set IsolationInputArray EFlowMerger/eflow
-
-  set OutputArray muons
-
-  set DeltaRMax 0.3
-  set PTMin 1.0
-  set PTRatioMax 9999.
-
-}
-
-
-#####################
-# Muon Loose Id     #
-#####################
+
+##################
+# Muon Loose Id  #
+##################
 
 module Efficiency MuonLooseIdEfficiency {
@@ -1243 +1038 @@
 }
 
-
-#####################
-# Muon Loose Id CHS #
-#####################
-
-module Efficiency MuonLooseIdEfficiencyCHS {
-    set InputArray MuonIsolationCHS/muons
-    set OutputArray muons
-    # tracking + LooseID efficiency formula for muons
-    source muonLooseId.tcl
-
-}
-
-
-######################
-# Muon Tight Id  CHS #
-######################
-
-module Efficiency MuonTightIdEfficiencyCHS {
-    set InputArray MuonIsolationCHS/muons
-    set OutputArray muons
-    # tracking + TightID efficiency formula for muons
-    source muonTightId.tcl
-}
 
 
@@ -1479 +1250 @@
 module StatusPidFilter GenParticleFilter {
 
-    set InputArray Delphes/allParticles
+    set InputArray  Delphes/allParticles
     set OutputArray filteredParticles
     set PTMin 5.0
@@ -1507 +1278 @@
   add Branch MuonLooseIdEfficiency/muons MuonLoose Muon
   add Branch MuonTightIdEfficiency/muons MuonTight Muon
-
-  add Branch PhotonEfficiencyCHS/photons PhotonCHS Photon
-  add Branch ElectronEfficiencyCHS/electrons ElectronCHS Electron
-  add Branch MuonLooseIdEfficiencyCHS/muons MuonLooseCHS Muon
-  add Branch MuonTightIdEfficiencyCHS/muons MuonTightCHS Muon
 
   add Branch JetEnergyScale/jets Jet Jet
@@ -1523 +1289 @@
   add Branch GenPileUpMissingET/momentum GenPileUpMissingET MissingET
   add Branch ScalarHT/energy ScalarHT ScalarHT
-
-}
+}

cards/CMS_PhaseII/CMS_PhaseII_Substructure_PIX4022_200PU.tcl

@@ -33 +33 @@
   PhotonEnergySmearing
   ElectronFilter
-
   TrackPileUpSubtractor
-  RecoPuFilter
 
   TowerMerger
   NeutralEFlowMerger
-
+  EFlowMergerAllTracks
   EFlowMerger
-  EFlowMergerCHS
-  Rho
 
   LeptonFilterNoLep
@@ -52 +48 @@
 
   PhotonIsolation
-  PhotonIsolationCHS
   PhotonEfficiency
-  PhotonEfficiencyCHS
 
   ElectronIsolation
-  ElectronIsolationCHS
-
   ElectronEfficiency
-  ElectronEfficiencyCHS
 
   MuonIsolation
-  MuonIsolationCHS
-
   MuonLooseIdEfficiency
   MuonTightIdEfficiency
-
-  MuonLooseIdEfficiencyCHS
-  MuonTightIdEfficiencyCHS
 
   NeutrinoFilter
@@ -82 +68 @@
   FastJetFinder
   FastJetFinderAK8
-  JetPileUpSubtractor
-  JetPileUpSubtractorAK8
   FastJetFinderPUPPI
   FastJetFinderPUPPIAK8
@@ -142 +126 @@
   set TVertexSpread 800E-12
 
-  # vertex smearing formula f(z,t) (z,t need to be respectively given in m,s) - {exp(-(t^2/160e-12^2/2))*exp(-(z^2/0.053^2/2))}
+  # vertex smearing formula f(z,t) (z,t need to be respectively given in m,s)
   set VertexDistributionFormula {exp(-(t^2/160e-12^2/2))*exp(-(z^2/0.053^2/2))}
 
@@ -154 +138 @@
 module ParticlePropagator ParticlePropagator {
   set InputArray PileUpMerger/stableParticles
-  #set InputArray Delphes/stableParticles
 
   set OutputArray stableParticles
@@ -289 +272 @@
   source muonMomentumResolution.tcl
 }
+
+
+
 
 ##############
@@ -540 +526 @@
 }
 
-########################
-# Reco PU filter
-########################
-
-module RecoPuFilter RecoPuFilter {
-  set InputArray HCal/eflowTracks
-  set OutputArray eflowTracks
-}
 
 ###################################################
@@ -560 +538 @@
 }
 
+
 ####################
 # Neutral eflow erger
@@ -571 +550 @@
 }
 
-#####################
+
+####################
 # Energy flow merger
-#####################
+####################
 
 module Merger EFlowMerger {
@@ -583 +563 @@
 }
 
-############################
-# Energy flow merger no PU
-############################
-
-module Merger EFlowMergerCHS {
+##################################
+# Energy flow merger (all tracks)
+##################################
+
+module Merger EFlowMergerAllTracks {
 # add InputArray InputArray
-  add InputArray RecoPuFilter/eflowTracks
+  add InputArray TrackMerger/tracks
   add InputArray PhotonEnergySmearing/eflowPhotons
   add InputArray HCal/eflowNeutralHadrons
@@ -714 +694 @@
 }
 
+
 #####################
 # MC truth jet finder
@@ -754 +735 @@
 
 
-#############
-# Rho pile-up
-#############
-
-module FastJetFinder Rho {
-#  set InputArray Calorimeter/towers
-  set InputArray EFlowMergerCHS/eflow
-
-  set ComputeRho true
-  set RhoOutputArray rho
-
-  # area algorithm: 0 Do not compute area, 1 Active area explicit ghosts, 2 One ghost passive area, 3 Passive area, 4 Voronoi, 5 Active area
-  set AreaAlgorithm 5
-
-  # jet algorithm: 1 CDFJetClu, 2 MidPoint, 3 SIScone, 4 kt, 5 Cambridge/Aachen, 6 antikt
-  set JetAlgorithm 4
-  set ParameterR 0.4
-  set GhostEtaMax 5.0
-
-  add RhoEtaRange -5.0 -4.0
-  add RhoEtaRange -4.0 -1.5
-  add RhoEtaRange -1.5 1.5
-  add RhoEtaRange 1.5 4.0
-  add RhoEtaRange 4.0 5.0
-
-  set JetPTMin 0.0
-}
-
-
-##############
+
+############
 # Jet finder
-##############
+############
 
 module FastJetFinder FastJetFinder {
 #  set InputArray TowerMerger/towers
-  set InputArray EFlowMergerCHS/eflow
+  set InputArray EFlowMerger/eflow
 
   set OutputArray jets
-
-  set AreaAlgorithm 5
 
   # algorithm: 1 CDFJetClu, 2 MidPoint, 3 SIScone, 4 kt, 5 Cambridge/Aachen, 6 antikt
@@ -805 +756 @@
 module FastJetFinder FastJetFinderAK8 {
 #  set InputArray TowerMerger/towers
-  set InputArray EFlowMergerCHS/eflow
+  set InputArray EFlowMerger/eflow
 
   set OutputArray jets
-
-  set AreaAlgorithm 5
 
   # algorithm: 1 CDFJetClu, 2 MidPoint, 3 SIScone, 4 kt, 5 Cambridge/Aachen, 6 antikt
@@ -836 +785 @@
 }
 
-###########################
-# Jet Pile-Up Subtraction
-###########################
-
-module JetPileUpSubtractor JetPileUpSubtractor {
-  set JetInputArray FastJetFinder/jets
-  set RhoInputArray Rho/rho
-
-  set OutputArray jets
-
-  set JetPTMin 15.0
-}
-
-##############################
-# Jet Pile-Up Subtraction AK8
-##############################
-
-module JetPileUpSubtractor JetPileUpSubtractorAK8 {
-  set JetInputArray FastJetFinderAK8/jets
-  set RhoInputArray Rho/rho
-
-  set OutputArray jets
-
-  set JetPTMin 15.0
-}
-
 module FastJetFinder FastJetFinderPUPPI {
 #  set InputArray TowerMerger/towers
@@ -911 +834 @@
 
 module EnergyScale JetEnergyScale {
-  set InputArray JetPileUpSubtractor/jets
+  set InputArray FastJetFinder/jets
   set OutputArray jets
 
@@ -919 +842 @@
 
 module EnergyScale JetEnergyScaleAK8 {
-  set InputArray JetPileUpSubtractorAK8/jets
+  set InputArray FastJetFinderAK8/jets
   set OutputArray jets
 
@@ -985 +908 @@
 }
 
-
-########################
-# Photon isolation CHS #
-########################
-
-module Isolation PhotonIsolationCHS {
-
-  # particle for which calculate the isolation
-  set CandidateInputArray PhotonFilter/photons
-
-  # isolation collection
-  set IsolationInputArray EFlowMerger/eflow
-
-  # output array
-  set OutputArray photons
-
-  # isolation cone
-  set DeltaRMax 0.3
-
-  # minimum pT
-  set PTMin     1.0
-
-  # iso ratio to cut
-  set PTRatioMax 9999.
-
-}
 
 
@@ -1033 +930 @@
 
 
-#####################
-# Photon efficiency #
-#####################
-
-module Efficiency PhotonEfficiencyCHS {
-
-  ## input particles
-  set InputArray PhotonIsolationCHS/photons
-  ## output particles
-  set OutputArray photons
-  # set EfficiencyFormula {efficiency formula as a function of eta and pt}
-  # efficiency formula for photons
-  set EfficiencyFormula {                      (pt <= 10.0) * (0.00) + \
-                           (abs(eta) <= 1.5) * (pt > 10.0)  * (0.9635) + \
-         (abs(eta) > 1.5 && abs(eta) <= 4.0) * (pt > 10.0)  * (0.9624) + \
-         (abs(eta) > 4.0)                                   * (0.00)}
-
-}
-
 ######################
 # Electron isolation #
@@ -1072 +950 @@
 }
 
-
-##########################
-# Electron isolation CHS #
-##########################
-
-module Isolation ElectronIsolationCHS {
-
-  set CandidateInputArray ElectronFilter/electrons
-
-  # isolation collection
-  set IsolationInputArray EFlowMerger/eflow
-
-  set OutputArray electrons
-
-  set DeltaRMax 0.3
-  set PTMin 1.0
-  set PTRatioMax 9999.
-
-}
 
 
@@ -1137 +996 @@
 }
 
-###########################
-# Electron efficiency CHS #
-###########################
-
-module Efficiency ElectronEfficiencyCHS {
-
-  set InputArray ElectronIsolationCHS/electrons
-  set OutputArray electrons
-
-  # set EfficiencyFormula {efficiency formula as a function of eta and pt}
-  # efficiency formula for electrons
-  set EfficiencyFormula {
-                                      (pt <= 4.0)  * (0.00) + \
-                         (abs(eta) <= 1.45 ) * (pt >  4.0 && pt <= 6.0)   * (0.50) + \
-                         (abs(eta) <= 1.45 ) * (pt >  6.0 && pt <= 8.0)   * (0.70) + \
-                         (abs(eta) <= 1.45 ) * (pt >  8.0 && pt <= 10.0)  * (0.85) + \
-                         (abs(eta) <= 1.45 ) * (pt > 10.0 && pt <= 30.0)  * (0.94) + \
-                         (abs(eta) <= 1.45 ) * (pt > 30.0 && pt <= 50.0)  * (0.97) + \
-                         (abs(eta) <= 1.45 ) * (pt > 50.0 && pt <= 70.0)  * (0.98) + \
-                         (abs(eta) <= 1.45 ) * (pt > 70.0 )  * (1.0) + \
-                         (abs(eta) > 1.45  && abs(eta) <= 1.55) * (pt >  4.0 && pt <= 10.0)   * (0.35) + \
-                         (abs(eta) > 1.45  && abs(eta) <= 1.55) * (pt > 10.0 && pt <= 30.0)   * (0.40) + \
-                         (abs(eta) > 1.45  && abs(eta) <= 1.55) * (pt > 30.0 && pt <= 70.0)   * (0.45) + \
-                         (abs(eta) > 1.45  && abs(eta) <= 1.55) * (pt > 70.0 )  * (0.55) + \
-                         (abs(eta) >= 1.55 && abs(eta) <= 2.0 ) * (pt >  4.0 && pt <= 10.0)  * (0.75) + \
-                         (abs(eta) >= 1.55 && abs(eta) <= 2.0 ) * (pt > 10.0 && pt <= 30.0)  * (0.85) + \
-                         (abs(eta) >= 1.55 && abs(eta) <= 2.0 ) * (pt > 30.0 && pt <= 50.0)  * (0.95) + \
-                         (abs(eta) >= 1.55 && abs(eta) <= 2.0 ) * (pt > 50.0 && pt <= 70.0)  * (0.95) + \
-                         (abs(eta) >= 1.55 && abs(eta) <= 2.0 ) * (pt > 70.0 )  * (1.0) + \
-                         (abs(eta) >= 2.0 && abs(eta) <= 2.5 ) * (pt >  4.0 && pt <= 10.0)  * (0.65) + \
-                         (abs(eta) >= 2.0 && abs(eta) <= 2.5 ) * (pt > 10.0 && pt <= 30.0)  * (0.75) + \
-                         (abs(eta) >= 2.0 && abs(eta) <= 2.5 ) * (pt > 30.0 && pt <= 50.0)  * (0.90) + \
-                         (abs(eta) >= 2.0 && abs(eta) <= 2.5 ) * (pt > 50.0 && pt <= 70.0)  * (0.90) + \
-                         (abs(eta) >= 2.0 && abs(eta) <= 2.5 ) * (pt > 70.0 )  * (0.90) + \
-                         (abs(eta) > 2.5 && abs(eta) <= 4.0 ) * (pt > 4.0 && pt <= 10.0) * (0.65) + \
-                                          (abs(eta) > 2.5 && abs(eta) <= 4.0 ) * (pt > 10.0 && pt <= 30.0) * (0.75) + \
-                                          (abs(eta) > 2.5 && abs(eta) <= 4.0 ) * (pt > 30.0 && pt <= 50.0) * (0.90) + \
-                                          (abs(eta) > 2.5 && abs(eta) <= 4.0 ) * (pt > 50.0 && pt <= 70.0) * (0.90) + \
-                                          (abs(eta) > 2.5 && abs(eta) <= 4.0 ) * (pt > 70.0 ) * (0.90) + \
-                                          (abs(eta) > 4.0) * (0.00)
-
-  }
-}
-
-
 ##################
 # Muon isolation #
@@ -1200 +1014 @@
 }
 
-######################
-# Muon isolation CHS #
-######################
-
-module Isolation MuonIsolationCHS {
-  set CandidateInputArray MuonMomentumSmearing/muons
-
-  # isolation collection
-  set IsolationInputArray EFlowMerger/eflow
-
-  set OutputArray muons
-
-  set DeltaRMax 0.3
-  set PTMin 1.0
-  set PTRatioMax 9999.
-
-}
-
-
-#####################
-# Muon Loose Id     #
-#####################
+
+##################
+# Muon Loose Id  #
+##################
 
 module Efficiency MuonLooseIdEfficiency {
@@ -1243 +1039 @@
 }
 
-
-#####################
-# Muon Loose Id CHS #
-#####################
-
-module Efficiency MuonLooseIdEfficiencyCHS {
-    set InputArray MuonIsolationCHS/muons
-    set OutputArray muons
-    # tracking + LooseID efficiency formula for muons
-    source muonLooseId.tcl
-
-}
-
-
-######################
-# Muon Tight Id  CHS #
-######################
-
-module Efficiency MuonTightIdEfficiencyCHS {
-    set InputArray MuonIsolationCHS/muons
-    set OutputArray muons
-    # tracking + TightID efficiency formula for muons
-    source muonTightId.tcl
-}
 
 
@@ -1479 +1251 @@
 module StatusPidFilter GenParticleFilter {
 
-    set InputArray Delphes/allParticles
+    set InputArray  Delphes/allParticles
     set OutputArray filteredParticles
     set PTMin 5.0
@@ -1507 +1279 @@
   add Branch MuonLooseIdEfficiency/muons MuonLoose Muon
   add Branch MuonTightIdEfficiency/muons MuonTight Muon
-
-  add Branch PhotonEfficiencyCHS/photons PhotonCHS Photon
-  add Branch ElectronEfficiencyCHS/electrons ElectronCHS Electron
-  add Branch MuonLooseIdEfficiencyCHS/muons MuonLooseCHS Muon
-  add Branch MuonTightIdEfficiencyCHS/muons MuonTightCHS Muon
 
   add Branch JetEnergyScale/jets Jet Jet
@@ -1523 +1290 @@
   add Branch GenPileUpMissingET/momentum GenPileUpMissingET MissingET
   add Branch ScalarHT/energy ScalarHT ScalarHT
-
-}
+}

cards/converter_card.tcl

@@ -13 +13 @@
 module TreeWriter TreeWriter {
 # add Branch InputArray BranchName BranchClass
-  add Branch Delphes/stableParticles Particle GenParticle
+  add Branch Delphes/allParticles Particle GenParticle
 }
 

cards/delphes_card_CMS.tcl

@@ -35 +35 @@
   
   FastJetFinder
-  FastCaloJetFinder
 
   JetEnergyScale
@@ -211 +210 @@
   set HCalEnergyMin 1.0
 
-  set ECalEnergySignificanceMin 2.0
-  set HCalEnergySignificanceMin 2.0
+  set ECalEnergySignificanceMin 1.0
+  set HCalEnergySignificanceMin 1.0
 
   set SmearTowerCenter true
@@ -498 +497 @@
 }
 
-################
-# caloJet finder
-################
-
-module FastJetFinder FastCaloJetFinder {
-  set InputArray Calorimeter/towers
-
-  set OutputArray jets
-
-  # algorithm: 1 CDFJetClu, 2 MidPoint, 3 SIScone, 4 kt, 5 Cambridge/Aachen, 6 antikt
-  set JetAlgorithm 6
-  set ParameterR 0.5
-
-  set JetPTMin 20.0
-}
-
 ##################
 # Jet Energy Scale
@@ -627 +610 @@
 
   add Branch UniqueObjectFinder/jets Jet Jet
-  add Branch FastCaloJetFinder/jets CaloJet Jet
   add Branch UniqueObjectFinder/electrons Electron Electron
   add Branch UniqueObjectFinder/photons Photon Photon

classes/DelphesClasses.cc

@@ -41 +41 @@
 CompBase *Tower::fgCompare = CompE<Tower>::Instance();
 CompBase *HectorHit::fgCompare = CompE<HectorHit>::Instance();
-CompBase *Vertex::fgCompare = CompSumPT2<Vertex>::Instance();
 CompBase *Candidate::fgCompare = CompMomentumPt<Candidate>::Instance();
 
@@ -122 +121 @@
   Charge(0), Mass(0.0),
   IsPU(0), IsRecoPU(0), IsConstituent(0), IsFromConversion(0),
-  ClusterIndex(-1), ClusterNDF(0), ClusterSigma(0), SumPT2(0), BTVSumPT2(0), GenDeltaZ(0), GenSumPT2(0),
   Flavor(0), FlavorAlgo(0), FlavorPhys(0),
   BTag(0), BTagAlgo(0), BTagPhys(0),
@@ -129 +127 @@
   Momentum(0.0, 0.0, 0.0, 0.0),
   Position(0.0, 0.0, 0.0, 0.0),
-  PositionError(0.0, 0.0, 0.0, 0.0),
-  InitialPosition(0.0, 0.0, 0.0, 0.0),
   Area(0.0, 0.0, 0.0, 0.0),
-  L(0),
-  D0(0), ErrorD0(0), 
-  DZ(0), ErrorDZ(0), 
-  P(0),  ErrorP(0), 
-  PT(0), ErrorPT(0), 
-  CtgTheta(0), ErrorCtgTheta(0), 
-  Phi(0), ErrorPhi(0),  
-  Xd(0), Yd(0), Zd(0), 
+  Dxy(0), SDxy(0), Xd(0), Yd(0), Zd(0),
   TrackResolution(0),
   NCharged(0),
@@ -256 +245 @@
   object.IsConstituent = IsConstituent;
   object.IsFromConversion = IsFromConversion;
-  object.ClusterIndex = ClusterIndex;
-  object.ClusterNDF = ClusterNDF;
-  object.ClusterSigma = ClusterSigma;
-  object.SumPT2 = SumPT2;
-  object.BTVSumPT2 = BTVSumPT2;
-  object.GenDeltaZ = GenDeltaZ;
-  object.GenSumPT2 = GenSumPT2;
   object.Flavor = Flavor;
   object.FlavorAlgo = FlavorAlgo;
@@ -280 +262 @@
   object.Momentum = Momentum;
   object.Position = Position;
-  object.InitialPosition = InitialPosition;
-  object.PositionError = PositionError;
   object.Area = Area;
-  object.L = L;
-  object.ErrorT = ErrorT;
-  object.D0 = D0;
-  object.ErrorD0 = ErrorD0;
-  object.DZ = DZ;
-  object.ErrorDZ = ErrorDZ;
-  object.P = P;
-  object.ErrorP = ErrorP;
-  object.PT = PT;
-  object.ErrorPT = ErrorPT;
-  object.CtgTheta = CtgTheta ;
-  object.ErrorCtgTheta = ErrorCtgTheta;
-  object.Phi = Phi;
-  object.ErrorPhi = ErrorPhi;  
+  object.Dxy = Dxy;
+  object.SDxy = SDxy;
   object.Xd = Xd;
   object.Yd = Yd;
@@ -314 +282 @@
   object.SumPtChargedPU = SumPtChargedPU;
   object.SumPt = SumPt;
-  object.ClusterIndex = ClusterIndex;
-  object.ClusterNDF = ClusterNDF;
-  object.ClusterSigma = ClusterSigma; 
-  object.SumPT2 = SumPT2;
-  
+
   object.FracPt[0] = FracPt[0];
   object.FracPt[1] = FracPt[1];
@@ -399 +363 @@
   Momentum.SetXYZT(0.0, 0.0, 0.0, 0.0);
   Position.SetXYZT(0.0, 0.0, 0.0, 0.0);
-  InitialPosition.SetXYZT(0.0, 0.0, 0.0, 0.0);
   Area.SetXYZT(0.0, 0.0, 0.0, 0.0);
-  L = 0.0;
-  ErrorT = 0.0;
-  D0 = 0.0;  
-  ErrorD0 = 0.0;
-  DZ = 0.0;
-  ErrorDZ = 0.0;
-  P =0.0;
-  ErrorP =0.0;
-  PT = 0.0;
-  ErrorPT = 0.0;
-  CtgTheta = 0.0;
-  ErrorCtgTheta = 0.0;
-  Phi = 0.0;
-  ErrorPhi = 0.0;
+  Dxy = 0.0;
+  SDxy = 0.0;
   Xd = 0.0;
   Yd = 0.0;
@@ -436 +387 @@
   SumPt = -999;
 
-  ClusterIndex = -1;
-  ClusterNDF = -99;
-  ClusterSigma = 0.0; 
-  SumPT2 = 0.0;
-  BTVSumPT2 = 0.0;
-  GenDeltaZ = 0.0;
-  GenSumPT2 = 0.0; 
-  
   FracPt[0] = 0.0;
   FracPt[1] = 0.0;

classes/DelphesClasses.h

@@ -147 +147 @@
   Float_t Pz; // particle momentum vector (z component) | hepevt.phep[number][2]
 
-  Float_t D0;
-  Float_t DZ;
-  Float_t P;
-  Float_t PT;
-  Float_t CtgTheta;
-  Float_t Phi;
+  Float_t PT; // particle transverse momentum
   Float_t Eta; // particle pseudorapidity
+  Float_t Phi; // particle azimuthal angle
+
   Float_t Rapidity; // particle rapidity
 
@@ -171 +168 @@
 //---------------------------------------------------------------------------
 
-class Vertex: public SortableObject
+class Vertex: public TObject
 {
 public:
@@ -179 +176 @@
   Float_t Z; // vertex position (z component)
 
-  Double_t ErrorX;
-  Double_t ErrorY;
-  Double_t ErrorZ;
-  Double_t ErrorT;
-
-  Int_t Index;
-  Int_t NDF;
-  Double_t Sigma;
-  Double_t SumPT2;
-  Double_t BTVSumPT2;
-  Double_t GenDeltaZ;
-  Double_t GenSumPT2;
-
-  TRefArray Constituents; // references to constituents
-
-  static CompBase *fgCompare; //!
-  const CompBase *GetCompare() const { return fgCompare; }
-
-  ClassDef(Vertex, 3)
+  ClassDef(Vertex, 1)
 };
 
@@ -418 +397 @@
   Int_t Charge; // track charge
 
+  Float_t PT; // track transverse momentum
+
   Float_t Eta; // track pseudorapidity
+  Float_t Phi; // track azimuthal angle
 
   Float_t EtaOuter; // track pseudorapidity at the tracker edge
@@ -433 +415 @@
   Float_t TOuter; // track position (z component) at the tracker edge
 
-  Float_t L; // track path length
-  Float_t ErrorT; // error on the time measurement
-
-  Float_t D0;     // track signed transverse impact parameter
-  Float_t ErrorD0;    // signed error on the track signed transverse impact parameter
-
-  Float_t DZ; // track transverse momentum
-  Float_t ErrorDZ; // track transverse momentum error
-
-  Float_t P; // track transverse momentum
-  Float_t ErrorP; // track transverse momentum error
-
-  Float_t PT; // track transverse momentum
-  Float_t ErrorPT; // track transverse momentum error
-
-  Float_t CtgTheta; // track transverse momentum
-  Float_t ErrorCtgTheta; // track transverse momentum error
-
-  Float_t Phi; // track azimuthal angle
-  Float_t ErrorPhi; // track azimuthal angle
-
+  Float_t Dxy;     // track signed transverse impact parameter
+  Float_t SDxy;    // signed error on the track signed transverse impact parameter
   Float_t Xd;      // X coordinate of point of closest approach to vertex
   Float_t Yd;      // Y coordinate of point of closest approach to vertex
@@ -460 +423 @@
   TRef Particle; // reference to generated particle
 
-  Int_t VertexIndex; // reference to vertex
-
   static CompBase *fgCompare; //!
   const CompBase *GetCompare() const { return fgCompare; }
@@ -565 +526 @@
   Float_t DeltaPhi;
 
-  TLorentzVector Momentum, Position, InitialPosition, PositionError, Area;
-
-  Float_t L; // path length
-  Float_t ErrorT; // path length
-  Float_t D0;
-  Float_t ErrorD0;
-  Float_t DZ;
-  Float_t ErrorDZ;
-  Float_t P;
-  Float_t ErrorP;
-  Float_t PT;
-  Float_t ErrorPT;
-  Float_t CtgTheta;
-  Float_t ErrorCtgTheta;
-  Float_t Phi;
-  Float_t ErrorPhi;
-
+  TLorentzVector Momentum, Position, Area;
+
+  Float_t Dxy;
+  Float_t SDxy;
   Float_t Xd;
   Float_t Yd;
@@ -587 +535 @@
 
   // tracking resolution
-
+  
   Float_t TrackResolution;
 
@@ -614 +562 @@
   Float_t SumPt;
 
-  // vertex variables
-
-  Int_t ClusterIndex;
-  Int_t ClusterNDF;
-  Double_t ClusterSigma;
-  Double_t SumPT2;
-  Double_t BTVSumPT2;
-  Double_t GenDeltaZ;
-  Double_t GenSumPT2;
-
   // N-subjettiness variables
 
@@ -657 +595 @@
   void SetFactory(DelphesFactory *factory) { fFactory = factory; }
 
-  ClassDef(Candidate, 5)
+  ClassDef(Candidate, 4)
 };
 

classes/SortableObject.h

@@ -156 +156 @@
       return -1;
     else if(t1->ET < t2->ET)
-      return 1;
-    else
-      return 0;
-  }
-};
-
-//---------------------------------------------------------------------------
-
-template <typename T>
-class CompSumPT2: public CompBase
-{
-  CompSumPT2() {}
-public:
-  static CompSumPT2 *Instance()
-  {
-    static CompSumPT2 single;
-    return &single;
-  }
-
-  Int_t Compare(const TObject *obj1, const TObject *obj2) const
-  {
-    const T *t1 = static_cast<const T*>(obj1);
-    const T *t2 = static_cast<const T*>(obj2);
-    if(t1->SumPT2 > t2->SumPT2)
-      return -1;
-    else if(t1->SumPT2 < t2->SumPT2)
       return 1;
     else

doc/genMakefile.tcl

@@ -263 +263 @@
 executableDeps {converters/*.cpp} {examples/*.cpp}
 
-executableDeps {readers/DelphesHepMC.cpp} {readers/DelphesLHEF.cpp} {readers/DelphesSTDHEP.cpp} {readers/DelphesROOT.cpp}
+executableDeps {readers/DelphesHepMC.cpp} {readers/DelphesLHEF.cpp} {readers/DelphesSTDHEP.cpp}
 
 puts {ifeq ($(HAS_CMSSW),true)}

examples/Pythia8/configParticleGun.cmnd

@@ -3 +3 @@
 ! 1) Settings used in the main program.
 
-Main:numberOfEvents = 100000         ! number of events to generate
+Main:numberOfEvents = 10000        ! number of events to generate
 Main:timesAllowErrors = 3          ! how many aborts before run stops
-Main:spareFlag1 = on                ! true means particle gun
-Main:spareMode1 = ID               ! 1-5 - di-quark, 21 - di-gluon, 11 - single electron, 13 - single muon, 22 - single photon
-Main:spareParm1 = 10000           ! max pt
+Main:spareFlag1 = on               ! true means particle gun
+Main:spareMode1 = ID               ! 1-5 - di-quark, 21 - di-gluon, 11 - single electron, 13 - single muon, 15 - single tau, 22 - single photon
+Main:spareParm1 = 10000            ! max pt
+Main:spareParm2 = 2.5              ! max eta
 
 ! 2) Settings related to output in init(), next() and stat().

examples/Validation.cpp

@@ -21 +21 @@
 #include <utility>
 #include <vector>
-#include <typeinfo>
 
 #include "TROOT.h"
@@ -29 +28 @@
 #include "TString.h"
 
-#include "TH1.h"
 #include "TH2.h"
-#include "TMath.h"
-#include "TStyle.h"
-#include "TGraph.h"
-#include "TCanvas.h"
 #include "THStack.h"
 #include "TLegend.h"
@@ -40 +34 @@
 #include "TClonesArray.h"
 #include "TLorentzVector.h"
-#include "TGraphErrors.h"
-#include "TMultiGraph.h"
 
 #include "classes/DelphesClasses.h"
@@ -55 +47 @@
 //------------------------------------------------------------------------------
 
-double ptrangemin = 10;
-double ptrangemax = 10000;
-static const int Nbins = 20;
+// Here you can put your analysis macro
 
-int objStyle = 1;
-int trackStyle = 7;
-int towerStyle = 3;
-
-Color_t objColor = kBlack;
-Color_t trackColor = kBlack;
-Color_t towerColor = kBlack;
-
-double effLegXmin = 0.22;
-double effLegXmax = 0.7;
-double effLegYmin = 0.22;
-double effLegYmax = 0.5;
-
-double resLegXmin = 0.62;
-double resLegXmax = 0.9;
-double resLegYmin = 0.52;
-double resLegYmax = 0.85;
-
-double topLeftLegXmin = 0.22;
-double topLeftLegXmax = 0.7;
-double topLeftLegYmin = 0.52;
-double topLeftLegYmax = 0.85;
-
-
-struct resolPlot
-{
-  TH1 *cenResolHist;
-  TH1 *fwdResolHist;
-  double ptmin;
-  double ptmax;
-  double xmin;
-  double xmax;
-  TString obj;
-
-  resolPlot();
-  resolPlot(double ptdown, double ptup, TString object);
-  void set(double ptdown, double ptup, TString object, double xmin = 0, double xmax = 2);
-  void print(){std::cout << ptmin << std::endl;}
-};
-
-
-resolPlot::resolPlot()
-{
-}
-
-resolPlot::resolPlot(double ptdown, double ptup, TString object)
-{
-  this->set(ptdown,ptup,object);
-}
-
-void resolPlot::set(double ptdown, double ptup, TString object, double xmin, double xmax)
-{
-  ptmin = ptdown;
-  ptmax = ptup;
-  obj = object;
-
-  cenResolHist = new TH1D(obj+"_delta_pt_"+Form("%4.2f",ptmin)+"_"+Form("%4.2f",ptmax)+"_cen", obj+"_delta_pt_"+Form("%4.2f",ptmin)+"_"+Form("%4.2f",ptmax)+"_cen", 200,  xmin, xmax);
-  fwdResolHist = new TH1D(obj+"_delta_pt_"+Form("%4.2f",ptmin)+"_"+Form("%4.2f",ptmax)+"_fwd", obj+"_delta_pt_"+Form("%4.2f",ptmin)+"_"+Form("%4.2f",ptmax)+"_fwd", 200,  xmin, xmax);
-}
-
-void HistogramsCollection(std::vector<resolPlot> *histos, double ptmin, double ptmax, TString obj, double xmin = 0, double xmax = 2)
-{
-  double width;
-  double ptdown;
-  double ptup;
-  resolPlot ptemp;
-
-  for (int i = 0; i < Nbins; i++)
-  {
-    width = (ptmax - ptmin) / Nbins;
-    ptdown = TMath::Power(10,ptmin + i * width );
-    ptup = TMath::Power(10,ptmin + (i+1) * width );
-    ptemp.set(ptdown, ptup, obj, xmin, xmax);
-    histos->push_back(ptemp);
-  }
-}
-
-//------------------------------------------------------------------------------
-
-class ExRootResult;
-class ExRootTreeReader;
-
-//------------------------------------------------------------------------------
-
-void BinLogX(TH1*h)
-{
-  TAxis *axis = h->GetXaxis();
-  int bins = axis->GetNbins();
-
-  Axis_t from = axis->GetXmin();
-  Axis_t to = axis->GetXmax();
-  Axis_t width = (to - from) / bins;
-  Axis_t *new_bins = new Axis_t[bins + 1];
-
-  for (int i = 0; i <= bins; i++)
-  {
-    new_bins[i] = TMath::Power(10, from + i * width);
-  }
-  axis->Set(bins, new_bins);
-  delete new_bins;
-}
-
-
-//------------------------------------------------------------------------------
-
-template<typename T>
-std::pair<TH1D*, TH1D*> GetEff(TClonesArray *branchReco, TClonesArray *branchParticle, TString name, int pdgID, ExRootTreeReader *treeReader)
-{
-
-  cout << "** Computing Efficiency of reconstructing "<< branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
-
-  Long64_t allEntries = treeReader->GetEntries();
-
-  GenParticle *particle;
-  T *recoObj;
-
-  TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
-
-  Float_t deltaR;
-  Float_t pt, eta;
-  Long64_t entry;
-
-  Int_t i, j;
-
-  TH1D *histGenPtcen = new TH1D(name+" gen spectra Pt",name+" gen spectra cen", Nbins, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax));
-  TH1D *histRecoPtcen = new TH1D(name+" reco spectra Pt",name+" reco spectra cen", Nbins, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax));
-  TH1D *histGenPtfwd  = new TH1D(name+" gen spectra Eta",name+" gen spectra fwd", Nbins, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax));
-  TH1D *histRecoPtfwd = new TH1D(name+" reco spectra Eta",name+" reco spectra fwd", Nbins, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax));
-
-
-  BinLogX(histGenPtcen);
-  BinLogX(histRecoPtcen);
-  BinLogX(histGenPtfwd);
-  BinLogX(histRecoPtfwd);
-
-  // Loop over all events
-  for(entry = 0; entry < allEntries; ++entry)
-  {
-    // Load selected branches with data from specified event
-    treeReader->ReadEntry(entry);
-
-    // Loop over all generated particle in event
-    for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
-    {
-
-      particle = (GenParticle*) branchParticle->At(i);
-      genMomentum = particle->P4();
-
-      deltaR = 999;
-
-      if (particle->PID == pdgID && genMomentum.Pt() > ptrangemin && genMomentum.Pt() < ptrangemax )
-      {
-
-        // Loop over all reco object in event
-        for(j = 0; j < branchReco->GetEntriesFast(); ++j)
-        {
-          recoObj = (T*)branchReco->At(j);
-          recoMomentum = recoObj->P4();
-          // this is simply to avoid warnings from initial state particle
-          // having infite rapidity ...
-        //if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
-
-          // take the closest parton candidate
-          if(TMath::Abs(pdgID) == 5)
-          {
-            Jet *jet = (Jet *)recoObj;
-            if(jet->BTag != 1) continue;
-          }
-          if(TMath::Abs(pdgID) == 15)
-          {
-            Jet *jet = (Jet *)recoObj;
-            if(jet->TauTag != 1) continue;
-          }
-          if(genMomentum.DeltaR(recoMomentum) < deltaR)
-          {
-            deltaR = genMomentum.DeltaR(recoMomentum);
-            bestRecoMomentum = recoMomentum;
-          }
-        }
-
-        pt  = genMomentum.Pt();
-        eta = genMomentum.Eta();
-
-        if (TMath::Abs(eta) < 1.5)
-        {
-          histGenPtcen->Fill(pt);
-          if(deltaR < 0.3) { histRecoPtcen->Fill(pt); }
-        }
-        else if (TMath::Abs(eta) < 2.5)
-        {
-          histGenPtfwd->Fill(pt);
-          if(deltaR < 0.3) { histRecoPtfwd->Fill(pt); }
-
-        }
-      }
-    }
-  }
-
-
-  std::pair<TH1D*,TH1D*> histos;
-
-  histRecoPtcen->Divide(histGenPtcen);
-  histRecoPtfwd->Divide(histGenPtfwd);
-
-  histos.first = histRecoPtcen;
-  histos.second = histRecoPtfwd;
-
-  return histos;
-}
-
-template<typename T>
-void GetEres(std::vector<resolPlot> *histos, TClonesArray *branchReco, TClonesArray *branchParticle, int pdgID, ExRootTreeReader *treeReader)
-{
-  Long64_t allEntries = treeReader->GetEntries();
-
-  cout << "** Computing resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
-
-  GenParticle *particle;
-  T* recoObj;
-
-  TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
-
-  Float_t deltaR;
-  Float_t pt, eta;
-  Long64_t entry;
-
-  Int_t i, j, bin;
-
-  // Loop over all events
-  for(entry = 0; entry < allEntries; ++entry)
-  {
-    // Load selected branches with data from specified event
-    treeReader->ReadEntry(entry);
-
-    // Loop over all reconstructed jets in event
-    for(i = 0; i < branchReco->GetEntriesFast(); ++i)
-    {
-      recoObj = (T*) branchReco->At(i);
-      recoMomentum = recoObj->P4();
-
-      deltaR = 999;
-
-     // Loop over all hard partons in event
-     for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
-     {
-        particle = (GenParticle*) branchParticle->At(j);
-        if (particle->PID == pdgID && particle->Status == 1)
-        {
-          genMomentum = particle->P4();
-
-          // this is simply to avoid warnings from initial state particle
-          // having infite rapidity ...
-          if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
-
-          // take the closest parton candidate
-          if(genMomentum.DeltaR(recoMomentum) < deltaR)
-          {
-             deltaR = genMomentum.DeltaR(recoMomentum);
-             bestGenMomentum = genMomentum;
-          }
-        }
-      }
-
-      if(deltaR < 0.3)
-      {
-        pt  = bestGenMomentum.Pt();
-        eta = TMath::Abs(bestGenMomentum.Eta());
-
-        for (bin = 0; bin < Nbins; bin++)
-        {
-          if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta > 0.0 && eta < 2.5)
-          {
-            if (eta < 1.5) {histos->at(bin).cenResolHist->Fill(recoMomentum.Pt()/bestGenMomentum.Pt());}
-            else if (eta < 2.5) {histos->at(bin).fwdResolHist->Fill(recoMomentum.Pt()/bestGenMomentum.Pt());}
-          }
-        }
-      }
-    }
-  }
-}
-void GetJetsEres(std::vector<resolPlot> *histos, TClonesArray *branchJet, TClonesArray *branchGenJet, ExRootTreeReader *treeReader)
-{
-
-  Long64_t allEntries = treeReader->GetEntries();
-
-  cout << "** Computing resolution of " << branchJet->GetName() << " induced by " << branchGenJet->GetName() << endl;
-
-  Jet *jet, *genjet;
-
-  TLorentzVector jetMomentum, genJetMomentum, bestGenJetMomentum;
-
-  Float_t deltaR;
-  Float_t pt, eta;
-  Long64_t entry;
-
-  Int_t i, j, bin;
-
-  // Loop over all events
-  for(entry = 0; entry < allEntries; ++entry)
-  {
-    // Load selected branches with data from specified event
-    treeReader->ReadEntry(entry);
-
-    if(entry%10000 == 0) cout << "Event number: "<< entry <<endl;
-
-    // Loop over all reconstructed jets in event
-    for(i = 0; i < TMath::Min(2,branchJet->GetEntriesFast()); ++i) //branchJet->GetEntriesFast(); ++i)
-    {
-
-      jet = (Jet*) branchJet->At(i);
-      jetMomentum = jet->P4();
-
-      deltaR = 999;
-
-     // Loop over all hard partons in event
-     for(j = 0; j < TMath::Min(2,branchGenJet->GetEntriesFast()); ++j)
-     {
-        genjet = (Jet*) branchGenJet->At(j);
-
-        genJetMomentum = genjet->P4();
-
-        // this is simply to avoid warnings from initial state particle
-        // having infite rapidity ...
-        if(genJetMomentum.Px() == 0 && genJetMomentum.Py() == 0) continue;
-
-        // take the closest parton candidate
-        if(genJetMomentum.DeltaR(jetMomentum) < deltaR)
-        {
-           deltaR = genJetMomentum.DeltaR(jetMomentum);
-           bestGenJetMomentum = genJetMomentum;
-        }
-      }
-
-      if(deltaR < 0.25)
-      {
-        pt  = genJetMomentum.Pt();
-        eta = TMath::Abs(genJetMomentum.Eta());
-
-        for (bin = 0; bin < Nbins; bin++)
-        {
-            if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta < 1.5)
-            {
-                histos->at(bin).cenResolHist->Fill(jetMomentum.Pt()/bestGenJetMomentum.Pt());
-            }
-            else if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta < 2.5)
-            {
-                histos->at(bin).fwdResolHist->Fill(jetMomentum.Pt()/bestGenJetMomentum.Pt());
-            }
-
-        }
-      }
-    }
-  }
-}
-
-void GetMetres(std::vector<resolPlot> *histos, TClonesArray *branchScalarHT, TClonesArray *branchMet, TClonesArray *branchJet, ExRootTreeReader *treeReader)
-{
-
-  Long64_t allEntries = treeReader->GetEntries();
-
-  cout << "** Computing resolution of " << branchMet->GetName() << " vs " << branchScalarHT->GetName() << endl;
-
-  MissingET *met;
-  ScalarHT *scalarHT;
-
-  Long64_t entry;
-
-  Int_t bin;
-  Double_t ht;
-
-  Jet *jet;
-  TLorentzVector p1, p2;
-
-  // Loop over all events
-  for(entry = 0; entry < allEntries; ++entry)
-  {
-    // Load selected branches with data from specified event
-    treeReader->ReadEntry(entry);
-
-    if(entry%10000 == 0) cout << "Event number: "<< entry <<endl;
-
-    if (branchJet->GetEntriesFast() > 1)
-    {
-
-      jet = (Jet*) branchJet->At(0);
-      p1 = jet->P4();
-      jet = (Jet*) branchJet->At(1);
-      p2 = jet->P4();
-
-      met = (MissingET*) branchMet->At(0);
-      scalarHT = (ScalarHT*) branchScalarHT->At(0);
-      ht = scalarHT->HT;
-
-      if(p1.Pt() < 0.75*ht/2) continue;
-      if(p2.Pt() < 0.75*ht/2) continue;
-
-      for (bin = 0; bin < Nbins; bin++)
-      {
-        if(ht > histos->at(bin).ptmin && ht < histos->at(bin).ptmax )
-        {
-          histos->at(bin).cenResolHist->Fill(met->P4().Px());
-          histos->at(bin).fwdResolHist->Fill(met->P4().Py());
-        }
-      }
-    }
-  }
-}
-
-
-std::pair<Double_t, Double_t> GausFit(TH1* hist)
-{
-  TF1 *f1 = new TF1("f1", "gaus", hist->GetMean()-2*hist->GetRMS(), hist->GetMean()+2*hist->GetRMS());
-  hist->Fit("f1","RQ");
-  Double_t sig = f1->GetParameter(2);
-  Double_t sigErr = f1->GetParError(2);
-  delete f1;
-  return make_pair (sig, sigErr);
-}
-
-
-TGraphErrors EresGraph(std::vector<resolPlot> *histos, bool central, bool rms = false)
-{
-  Int_t bin;
-  Int_t count = 0;
-  TGraphErrors gr = TGraphErrors(Nbins/2);
-  Double_t sig = 0;
-  Double_t sigErr = 0;
-  for (bin = 0; bin < Nbins; bin++)
-  {
-    if (central == true && histos->at(bin).cenResolHist->GetEntries() > 100)
-    {
-      std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).cenResolHist);
-      if (rms == true)
-      {
-        gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, sigvalues.second);
-        gr.SetPointError(count,0, sigvalues.second); // to correct
-      }
-      else
-      {
-        gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, sigvalues.first);
-        gr.SetPointError(count,0, sigvalues.second);
-      }
-      count++;
-    }
-
-    else if (central == false && histos->at(bin).fwdResolHist->GetEntries() > 10)
-    {
-      std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).fwdResolHist);
-      if (rms == true)
-      {
-        gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, sigvalues.second);
-        gr.SetPointError(count,0, sigvalues.second); // to correct
-      }
-      else
-      {
-        gr.SetPoint(count,(histos->at(bin).ptmin+histos->at(bin).ptmax)/2.0, sigvalues.first);
-        gr.SetPointError(count,0, sigvalues.second);
-      }
-      count++;
-    }
-
-  }
-  return gr;
-}
-
-
-//------------------------------------------------------------------------------
-
-
-// type 1 : object, 2 : track, 3 : tower
-
-void addGraph(TMultiGraph *mg, TGraphErrors *gr, TLegend *leg, int type)
-{
-
-  gr->SetLineWidth(2);
-
-  switch ( type )
-  {
-    case 1:
-      gr->SetLineColor(objColor);
-      gr->SetLineStyle(objStyle);
-      std::cout << "Adding " << gr->GetName() << std::endl;
-      mg->Add(gr);
-      leg->AddEntry(gr,"Reco","l");
-      break;
-
-    case 2:
-      gr->SetLineColor(trackColor);
-      gr->SetLineStyle(trackStyle);
-      mg->Add(gr);
-      leg->AddEntry(gr,"Track","l");
-      break;
-
-    case 3:
-      gr->SetLineColor(towerColor);
-      gr->SetLineStyle(towerStyle);
-      mg->Add(gr);
-      leg->AddEntry(gr,"Tower","l");
-      break;
-
-    case 0:
-      gr->SetLineColor(objColor);
-      gr->SetLineStyle(objStyle);
-      mg->Add(gr);
-      break;
-
-    default:
-      std::cout << "wrong type, possibles choices are Object, Track and Tower" << std::endl;
-      break;
-  }
-}
-
-void addHist(TH1D *h, TLegend *leg, int type)
-{
-  h->SetLineWidth(2);
-
-  switch ( type )
-  {
-    case 1:
-      h->SetLineColor(objColor);
-      h->SetLineStyle(objStyle);
-      leg->AddEntry(h,"Reco","l");
-      break;
-
-    case 2:
-      h->SetLineColor(trackColor);
-      h->SetLineStyle(trackStyle);
-      leg->AddEntry(h,"Track","l");
-      break;
-
-    case 3:
-      h->SetLineColor(towerColor);
-      h->SetLineStyle(towerStyle);
-      leg->AddEntry(h,"Tower","l");
-      break;
-
-    case 0:
-      h->SetLineColor(objColor);
-      h->SetLineStyle(objStyle);
-      break;
-
-    default:
-      std::cout << "wrong type, possibles choices are Object, Track and Tower" << std::endl;
-      break;
-  }
-}
-
-void DrawAxis(TMultiGraph *mg, TLegend *leg, double max)
-{
-  mg->SetMinimum(0.);
-  mg->SetMaximum(max);
-  mg->GetXaxis()->SetLimits(ptrangemin,ptrangemax);
-  mg->GetYaxis()->SetTitle("resolution");
-  mg->GetXaxis()->SetTitle("p_{T} [GeV]");
-  mg->GetYaxis()->SetTitleSize(0.07);
-  mg->GetXaxis()->SetTitleSize(0.07);
-  mg->GetYaxis()->SetLabelSize(0.06);
-  mg->GetXaxis()->SetLabelSize(0.06);
-  mg->GetYaxis()->SetLabelOffset(0.03);
-  mg->GetYaxis()->SetTitleOffset(1.4);
-  mg->GetXaxis()->SetTitleOffset(1.4);
-
-  mg->GetYaxis()->SetNdivisions(505);
-
-  leg->SetBorderSize(0);
-  leg->SetShadowColor(0);
-  leg->SetFillColor(0);
-  leg->SetFillStyle(0);
-
-  gStyle->SetOptTitle(0);
-  gPad->SetLogx();
-  gPad->SetBottomMargin(0.2);
-  gPad->SetLeftMargin(0.2);
-  gPad->Modified();
-  gPad->Update();
-
-}
-
-void DrawAxis(TH1D *h, TLegend *leg, int type)
-{
-
-  h->GetYaxis()->SetRangeUser(0,1.0);
-  if (type == 0) h->GetXaxis()->SetTitle("p_{T} [GeV]");
-  else h->GetXaxis()->SetTitle("#eta");
-  h->GetYaxis()->SetTitle("efficiency");
-  h->GetYaxis()->SetTitleSize(0.07);
-  h->GetXaxis()->SetTitleSize(0.07);
-  h->GetYaxis()->SetLabelSize(0.06);
-  h->GetXaxis()->SetLabelSize(0.06);
-  h->GetYaxis()->SetLabelOffset(0.03);
-  h->GetYaxis()->SetTitleOffset(1.3);
-  h->GetXaxis()->SetTitleOffset(1.4);
-
-  h->GetYaxis()->SetNdivisions(505);
-
-  leg->SetBorderSize(0);
-  leg->SetShadowColor(0);
-  leg->SetFillColor(0);
-  leg->SetFillStyle(0);
-
-  gStyle->SetOptTitle(0);
-  gStyle->SetOptStat(0);
-  gPad->SetBottomMargin(0.2);
-  gPad->SetLeftMargin(0.2);
-
-  gPad->Modified();
-  gPad->Update();
-
-}
-
-
-void Validation(const char *inputFileElectron, const char *inputFileMuon, const char *inputFilePhoton, const char *inputFileJet, const char *inputFileBJet, const char *inputFileTauJet, const char *outputFile)
-{
-  TChain *chainElectron = new TChain("Delphes");
-  chainElectron->Add(inputFileElectron);
-  ExRootTreeReader *treeReaderElectron = new ExRootTreeReader(chainElectron);
-
-  TChain *chainMuon = new TChain("Delphes");
-  chainMuon->Add(inputFileMuon);
-  ExRootTreeReader *treeReaderMuon = new ExRootTreeReader(chainMuon);
-
-  TChain *chainPhoton = new TChain("Delphes");
-  chainPhoton->Add(inputFilePhoton);
-  ExRootTreeReader *treeReaderPhoton = new ExRootTreeReader(chainPhoton);
-
-  TChain *chainJet = new TChain("Delphes");
-  chainJet->Add(inputFileJet);
-  ExRootTreeReader *treeReaderJet = new ExRootTreeReader(chainJet);
-
-  TChain *chainBJet = new TChain("Delphes");
-  chainBJet->Add(inputFileBJet);
-  ExRootTreeReader *treeReaderBJet = new ExRootTreeReader(chainBJet);
-
-  TChain *chainTauJet = new TChain("Delphes");
-  chainTauJet->Add(inputFileTauJet);
-  ExRootTreeReader *treeReaderTauJet = new ExRootTreeReader(chainTauJet);
-
-  TClonesArray *branchParticleElectron = treeReaderElectron->UseBranch("Particle");
-  TClonesArray *branchTrackElectron = treeReaderElectron->UseBranch("Track");
-  TClonesArray *branchTowerElectron = treeReaderElectron->UseBranch("Tower");
-  TClonesArray *branchElectron = treeReaderElectron->UseBranch("Electron");
-
-  TClonesArray *branchParticleMuon = treeReaderMuon->UseBranch("Particle");
-  TClonesArray *branchTrackMuon = treeReaderMuon->UseBranch("Track");
-  TClonesArray *branchMuon = treeReaderMuon->UseBranch("Muon");
-
-  TClonesArray *branchParticlePhoton = treeReaderPhoton->UseBranch("Particle");
-  TClonesArray *branchTowerPhoton = treeReaderPhoton->UseBranch("Tower");
-  TClonesArray *branchPhoton = treeReaderPhoton->UseBranch("Photon");
-
-  TClonesArray *branchGenJet = treeReaderJet->UseBranch("GenJet");
-  TClonesArray *branchPFJet = treeReaderJet->UseBranch("Jet");
-  TClonesArray *branchCaloJet = treeReaderJet->UseBranch("CaloJet");
-
-  TClonesArray *branchParticleBJet = treeReaderBJet->UseBranch("Particle");
-  TClonesArray *branchPFBJet = treeReaderBJet->UseBranch("Jet");
-
-  TClonesArray *branchParticleTauJet = treeReaderTauJet->UseBranch("Particle");
-  TClonesArray *branchPFTauJet = treeReaderTauJet->UseBranch("Jet");
-
-  TClonesArray *branchScalarHT = treeReaderJet->UseBranch("ScalarHT");
-  TClonesArray *branchMet = treeReaderJet->UseBranch("MissingET");
-
-  ///////////////
-  // Electrons //
-  ///////////////
-
-  // Reconstruction efficiency
-  TString elecs = "Electron";
-  int elID = 11;
-  std::pair<TH1D*,TH1D*> histos_el = GetEff<Electron>(branchElectron, branchParticleElectron, "Electron", elID, treeReaderElectron);
-
-  // tracking reconstruction efficiency
-  std::pair <TH1D*,TH1D*> histos_eltrack = GetEff<Track>(branchTrackElectron, branchParticleElectron, "electronTrack", elID, treeReaderElectron);
-
-  // Tower reconstruction efficiency
-  std::pair <TH1D*,TH1D*> histos_eltower = GetEff<Tower>(branchTowerElectron, branchParticleElectron, "electronTower", elID, treeReaderElectron);
-
-  // Electron Energy Resolution
-  std::vector<resolPlot> plots_el;
-  HistogramsCollection(&plots_el, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "electrons");
-  GetEres<Electron>(&plots_el, branchElectron, branchParticleElectron, elID, treeReaderElectron);
-  TGraphErrors gr_el = EresGraph(&plots_el, true);
-  TGraphErrors gr_elFwd = EresGraph(&plots_el, false);
-  gr_el.SetName("Electron");
-  gr_elFwd.SetName("ElectronFwd");
-
-  // Electron Track Energy Resolution
-  std::vector<resolPlot> plots_eltrack;
-  HistogramsCollection(&plots_eltrack, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "electronsTracks");
-  GetEres<Track>(&plots_eltrack, branchTrackElectron, branchParticleElectron, elID, treeReaderElectron);
-  TGraphErrors gr_eltrack = EresGraph(&plots_eltrack, true);
-  TGraphErrors gr_eltrackFwd = EresGraph(&plots_eltrack, false);
-  gr_eltrack.SetName("ElectronTracks");
-  gr_eltrackFwd.SetName("ElectronTracksFwd");
-
-  // Electron Tower Energy Resolution
-  std::vector<resolPlot> plots_eltower;
-  HistogramsCollection(&plots_eltower, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "electronsTower");
-  GetEres<Tower>(&plots_eltower, branchTowerElectron, branchParticleElectron, elID, treeReaderElectron);
-  TGraphErrors gr_eltower = EresGraph(&plots_eltower, true);
-  TGraphErrors gr_eltowerFwd = EresGraph(&plots_eltower, false);
-  gr_eltower.SetName("ElectronTower");
-  gr_eltrackFwd.SetName("ElectronTracksFwd");
-
-  // Canvases
-  TString elEff = "electronEff";
-  TCanvas *C_el1 = new TCanvas(elEff,elEff, 1600, 600);
-  C_el1->Divide(2);
-  C_el1->cd(1);
-  TLegend *leg_el1 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
-  leg_el1->SetHeader("#splitline{electrons}{|#eta| < 1.5}");
-  leg_el1->AddEntry("","","");
-
-  gPad->SetLogx();
-  histos_eltrack.first->Draw("][");
-  addHist(histos_eltrack.first, leg_el1, 2);
-  histos_el.first->Draw("same ][");
-  addHist(histos_el.first, leg_el1, 1);
-  DrawAxis(histos_eltrack.first, leg_el1, 0);
-
-  leg_el1->Draw();
-
-  C_el1->cd(2);
-  TLegend *leg_el2 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
-  leg_el2->SetHeader("#splitline{electrons}{1.5 < |#eta| < 2.5}");
-  leg_el2->AddEntry("","","");
-
-  gPad->SetLogx();
-  histos_eltrack.second->Draw("][");
-  addHist(histos_eltrack.second, leg_el2, 2);
-  histos_el.second->Draw("same ][");
-  addHist(histos_el.second, leg_el2, 1);
-
-  DrawAxis(histos_eltrack.second, leg_el2, 0);
-  leg_el2->Draw();
-
-  C_el1->cd(0);
-
-  TString elRes = "electronERes";
-  TString elResFwd = "electronEResForward";
-  TCanvas *C_el2 = new TCanvas(elRes,elRes, 1600, 600);
-  C_el2->Divide(2);
-  C_el2->cd(1);
-  TMultiGraph *mg_el = new TMultiGraph(elRes,elRes);
-  TLegend *leg_el = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
-  leg_el->SetHeader("#splitline{electrons}{|#eta| < 1.5}");
-  leg_el->AddEntry("","","");
-
-  addGraph(mg_el, &gr_eltower, leg_el, 3);
-  addGraph(mg_el, &gr_eltrack, leg_el, 2);
-  addGraph(mg_el, &gr_el, leg_el, 1);
-
-  mg_el->Draw("ACX");
-  leg_el->Draw();
-
-  DrawAxis(mg_el, leg_el, 0.1);
-
-  C_el2->cd(2);
-  TMultiGraph *mg_elFwd = new TMultiGraph(elResFwd,elResFwd);
-  TLegend *leg_elFwd = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
-  leg_elFwd->SetHeader("#splitline{electrons}{1.5 < |#eta| < 2.5}");
-  leg_elFwd->AddEntry("","","");
-
-  addGraph(mg_elFwd, &gr_eltowerFwd, leg_elFwd, 3);
-  addGraph(mg_elFwd, &gr_eltrackFwd, leg_elFwd, 2);
-  addGraph(mg_elFwd, &gr_elFwd, leg_elFwd, 1);
-
-  mg_elFwd->Draw("ACX");
-  leg_elFwd->Draw();
-
-  DrawAxis(mg_elFwd, leg_elFwd, 0.2);
-
-  C_el1->Print("validation.pdf(","pdf");
-  C_el2->Print("validation.pdf","pdf");
-
-  gDirectory->cd(0);
-
-  ///////////
-  // Muons //
-  ///////////
-
-  // Reconstruction efficiency
-  int muID = 13;
-  std::pair<TH1D*,TH1D*> histos_mu = GetEff<Muon>(branchMuon, branchParticleMuon,"Muon", muID, treeReaderMuon);
-
-  // muon tracking reconstruction efficiency
-  std::pair <TH1D*,TH1D*> histos_mutrack = GetEff<Track>(branchTrackMuon, branchParticleMuon, "muonTrack", muID, treeReaderMuon);
-
-  // Muon Energy Resolution
-  std::vector<resolPlot> plots_mu;
-  HistogramsCollection(&plots_mu, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "muons");
-  GetEres<Muon>(&plots_mu, branchMuon, branchParticleMuon, muID, treeReaderMuon);
-  TGraphErrors gr_mu = EresGraph(&plots_mu, true);
-  TGraphErrors gr_muFwd = EresGraph(&plots_mu, false);
-  gr_mu.SetName("Muon");
-  gr_muFwd.SetName("MuonFwd");
-
-  // Muon Track Energy Resolution
-  std::vector<resolPlot> plots_mutrack;
-  HistogramsCollection(&plots_mutrack, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "muonsTracks");
-  GetEres<Track>(&plots_mutrack, branchTrackMuon, branchParticleMuon, muID, treeReaderMuon);
-  TGraphErrors gr_mutrack = EresGraph(&plots_mutrack, true);
-  TGraphErrors gr_mutrackFwd = EresGraph(&plots_mutrack, false);
-  gr_mutrackFwd.SetName("MuonTracksFwd");
-
-  // Canvas
-
-  TString muEff = "muonEff";
-  TCanvas *C_mu1 = new TCanvas(muEff,muEff, 1600, 600);
-  C_mu1->Divide(2);
-  C_mu1->cd(1);
-  TLegend *leg_mu1 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
-  leg_mu1->SetHeader("#splitline{muons}{|#eta| < 1.5}");
-  leg_mu1->AddEntry("","","");
-
-
-  gPad->SetLogx();
-  histos_mutrack.first->Draw("][");
-  addHist(histos_mutrack.first, leg_mu1, 2);
-  histos_mu.first->Draw("same ][");
-  addHist(histos_mu.first, leg_mu1, 1);
-
-  DrawAxis(histos_mutrack.first, leg_mu1, 0);
-
-  leg_mu1->Draw();
-
-  C_mu1->cd(2);
-  TLegend *leg_mu2 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
-  leg_mu2->SetHeader("#splitline{muons}{1.5 < |#eta| < 2.5}");
-  leg_mu2->AddEntry("","","");
-
-  gPad->SetLogx();
-  histos_mutrack.second->Draw("][");
-  addHist(histos_mutrack.second, leg_mu2, 2);
-  histos_mu.second->Draw("same ][");
-  addHist(histos_mu.second, leg_mu2, 1);
-
-  DrawAxis(histos_mutrack.second, leg_mu2, 1);
-  leg_mu2->Draw();
-
-  TString muRes = "muonERes";
-  TString muResFwd = "muonEResFwd";
-
-  TCanvas *C_mu = new TCanvas(muRes,muRes, 1600, 600);
-  C_mu->Divide(2);
-  C_mu->cd(1);
-  TMultiGraph *mg_mu = new TMultiGraph(muRes,muRes);
-  TLegend *leg_mu = new TLegend(topLeftLegXmin,topLeftLegYmin,topLeftLegXmax,topLeftLegYmax);
-  leg_mu->SetHeader("#splitline{muons}{|#eta| < 1.5}");
-  leg_mu->AddEntry("","","");
-
-  addGraph(mg_mu, &gr_mutrack, leg_mu, 2);
-  addGraph(mg_mu, &gr_mu, leg_mu, 1);
-
-  mg_mu->Draw("ACX");
-  leg_mu->Draw();
-
-  DrawAxis(mg_mu, leg_mu, 0.3);
-
-  C_mu->cd(2);
-  TMultiGraph *mg_muFwd = new TMultiGraph(muResFwd,muResFwd);
-  TLegend *leg_muFwd = new TLegend(topLeftLegXmin,topLeftLegYmin,topLeftLegXmax,topLeftLegYmax);
-  leg_muFwd->SetHeader("#splitline{muons}{1.5 < |#eta| < 2.5}");
-  leg_muFwd->AddEntry("","","");
-
-  addGraph(mg_muFwd, &gr_mutrackFwd, leg_muFwd, 2);
-  addGraph(mg_muFwd, &gr_muFwd, leg_muFwd, 1);
-
-  mg_muFwd->Draw("ACX");
-  leg_muFwd->Draw();
-
-  DrawAxis(mg_muFwd, leg_muFwd, 0.3);
-
-  //C_mu1->SaveAs(muEff+".eps");
-  //C_mu->SaveAs(muRes+".eps");
-
-  C_mu1->Print("validation.pdf","pdf");
-  C_mu->Print("validation.pdf","pdf");
-
-  gDirectory->cd(0);
-
-  /////////////
-  // Photons //
-  /////////////
-
-  // Reconstruction efficiency
-  int phID = 22;
-  std::pair<TH1D*,TH1D*> histos_ph = GetEff<Electron>(branchPhoton, branchParticlePhoton, "Photon", phID, treeReaderPhoton);
-  std::pair<TH1D*,TH1D*> histos_phtower = GetEff<Electron>(branchTowerPhoton, branchParticlePhoton, "Photon", phID, treeReaderPhoton);
-
-  // Photon Energy Resolution
-  std::vector<resolPlot> plots_ph;
-  HistogramsCollection(&plots_ph, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "photons");
-  GetEres<Photon>(&plots_ph, branchPhoton, branchParticlePhoton, phID, treeReaderPhoton);
-  TGraphErrors gr_ph = EresGraph(&plots_ph, true);
-  TGraphErrors gr_phFwd = EresGraph(&plots_ph, false);
-  gr_ph.SetName("Photon");
-  gr_phFwd.SetName("PhotonFwd");
-
-
-  // Photon Tower Energy Resolution
-  std::vector<resolPlot> plots_phtower;
-  HistogramsCollection(&plots_phtower, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "photonsTower");
-  GetEres<Tower>(&plots_phtower, branchTowerPhoton, branchParticlePhoton, phID, treeReaderPhoton);
-  TGraphErrors gr_phtower = EresGraph(&plots_phtower, true);
-  TGraphErrors gr_phtowerFwd = EresGraph(&plots_phtower, false);
-  gr_phtower.SetName("PhotonTower");
-  gr_phtowerFwd.SetName("PhotonTowerFwd");
-
-  // Canvas
-
-  TString phEff = "photonEff";
-  TCanvas *C_ph1 = new TCanvas(phEff,phEff, 1600, 600);
-  C_ph1->Divide(2);
-  C_ph1->cd(1);
-  TLegend *leg_ph1 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
-  leg_ph1->SetHeader("#splitline{photons}{|#eta| < 1.5}");
-  leg_ph1->AddEntry("","","");
-
-
-  gPad->SetLogx();
-  histos_phtower.first->Draw("][");
-  addHist(histos_phtower.first, leg_ph1, 3);
-  histos_ph.first->Draw("same ][");
-  addHist(histos_ph.first, leg_ph1, 1);
-
-  DrawAxis(histos_phtower.first, leg_ph1, 0);
-  leg_ph1->Draw();
-
-  C_ph1->cd(2);
-  TLegend *leg_ph2 = new TLegend(effLegXmin,effLegYmin,effLegXmax,effLegYmax);
-  leg_ph2->SetHeader("#splitline{photons}{1.5 < |#eta| < 2.5}");
-  leg_ph2->AddEntry("","","");
-
-
-  gPad->SetLogx();
-  histos_phtower.second->Draw("][");
-  addHist(histos_phtower.second, leg_ph2, 3);
-  histos_ph.second->Draw("same ][");
-  addHist(histos_ph.second, leg_ph2, 1);
-
-  DrawAxis(histos_phtower.second, leg_ph2, 1);
-  leg_ph2->Draw();
-
-  C_ph1->SaveAs(phEff+".eps");
-
-  TString phRes = "phERes";
-  TString phResFwd = "phEResFwd";
-
-  TCanvas *C_ph = new TCanvas(phRes,phRes, 1600, 600);
-  C_ph->Divide(2);
-  C_ph->cd(1);
-  TMultiGraph *mg_ph = new TMultiGraph(phRes,phRes);
-  TLegend *leg_ph = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
-  leg_ph->SetHeader("#splitline{photons}{|#eta| < 1.5}");
-  leg_ph->AddEntry("","","");
-
-  addGraph(mg_ph, &gr_phtower, leg_ph, 3);
-  addGraph(mg_ph, &gr_ph, leg_ph, 1);
-
-  mg_ph->Draw("ACX");
-  leg_ph->Draw();
-
-  DrawAxis(mg_ph, leg_ph, 0.3);
-
-  C_ph->cd(2);
-  TMultiGraph *mg_phFwd = new TMultiGraph(phResFwd,phResFwd);
-  TLegend *leg_phFwd = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
-  leg_phFwd->SetHeader("#splitline{photons}{1.5 < |#eta| < 2.5}");
-  leg_phFwd->AddEntry("","","");
-
-  addGraph(mg_phFwd, &gr_phtowerFwd, leg_phFwd, 3);
-  addGraph(mg_phFwd, &gr_phFwd, leg_phFwd, 1);
-
-  mg_phFwd->Draw("ACX");
-  leg_phFwd->Draw();
-
-  DrawAxis(mg_phFwd, leg_phFwd, 0.3);
-
-  C_ph->SaveAs(phRes+".eps");
-
-  C_ph1->Print("validation.pdf","pdf");
-  C_ph->Print("validation.pdf","pdf");
-
-  gDirectory->cd(0);
-
-  //////////
-  // Jets //
-  //////////
-
-  // BJets Reconstruction efficiency
-  int bID = 5;
-  std::pair<TH1D*,TH1D*> histos_btag = GetEff<Jet>(branchPFBJet, branchParticleBJet,"BTag", bID, treeReaderBJet);
-
-  // TauJets Reconstruction efficiency
-  int tauID = 15;
-  std::pair<TH1D*,TH1D*> histos_tautag = GetEff<Jet>(branchPFTauJet, branchParticleTauJet,"TauTag", tauID, treeReaderTauJet);
-
-  // PFJets Energy Resolution
-  std::vector<resolPlot> plots_pfjets;
-  HistogramsCollection(&plots_pfjets, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "PFJet");
-  GetJetsEres(&plots_pfjets, branchPFJet, branchGenJet, treeReaderJet);
-  TGraphErrors gr_pfjets = EresGraph(&plots_pfjets, true);
-  TGraphErrors gr_pfjetsFwd = EresGraph(&plots_pfjets, false);
-  gr_pfjets.SetName("pfJet");
-  gr_pfjetsFwd.SetName("pfJetFwd");
-
-  // CaloJets Energy Resolution
-  std::vector<resolPlot> plots_calojets;
-  HistogramsCollection(&plots_calojets, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "CaloJet");
-  GetJetsEres(&plots_calojets, branchCaloJet, branchGenJet, treeReaderJet);
-  TGraphErrors gr_calojets = EresGraph(&plots_calojets, true);
-  TGraphErrors gr_calojetsFwd = EresGraph(&plots_calojets, false);
-  gr_calojets.SetName("caloJet");
-  gr_calojetsFwd.SetName("caloJetFwd");
-
-  // MET Resolution vs HT
-  std::vector<resolPlot> plots_met;
-  HistogramsCollection(&plots_met, TMath::Log10(ptrangemin), TMath::Log10(ptrangemax), "MET", -500, 500);
-  GetMetres(&plots_met, branchScalarHT, branchMet, branchPFJet, treeReaderJet);
-  TGraphErrors gr_met = EresGraph(&plots_met, true);
-  gr_calojets.SetName("MET");
-
-  // Canvas
-  TString btagEff = "btagEff";
-  TCanvas *C_btag1 = new TCanvas(btagEff,btagEff, 1600, 600);
-  C_btag1->Divide(2);
-  C_btag1->cd(1);
-  TLegend *leg_btag1 = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
-  leg_btag1->SetHeader("#splitline{B-tagging}{|#eta| < 1.5}");
-  leg_btag1->AddEntry("","","");
-
-  gPad->SetLogx();
-  histos_btag.first->Draw();
-  addHist(histos_btag.first, leg_btag1, 0);
-
-  DrawAxis(histos_btag.first, leg_btag1, 0);
-  leg_btag1->Draw();
-
-  C_btag1->cd(2);
-  TLegend *leg_btag2 = new TLegend(resLegXmin,resLegYmin,resLegXmax+0.05,resLegYmax);
-  leg_btag2->SetHeader("#splitline{B-tagging}{1.5 < |#eta| < 2.5}");
-  leg_btag2->AddEntry("","","");
-
-  gPad->SetLogx();
-  histos_btag.second->Draw();
-  addHist(histos_btag.second, leg_btag2, 0);
-
-  DrawAxis(histos_btag.second, leg_btag2, 0);
-  leg_btag2->Draw();
-  C_btag1->cd(0);
-
-  TString tautagEff = "tautagEff";
-  TCanvas *C_tautag1 = new TCanvas(tautagEff,tautagEff, 1600, 600);
-  C_tautag1->Divide(2);
-  C_tautag1->cd(1);
-  TLegend *leg_tautag1 = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
-  leg_tautag1->SetHeader("#splitline{#tau-tagging}{|#eta| < 1.5}");
-  leg_tautag1->AddEntry("","","");
-
-  gPad->SetLogx();
-  histos_tautag.first->Draw();
-  addHist(histos_tautag.first, leg_tautag1, 0);
-
-  DrawAxis(histos_tautag.first, leg_tautag1, 0);
-  leg_tautag1->Draw();
-
-  C_tautag1->cd(2);
-  TLegend *leg_tautag2 = new TLegend(resLegXmin,resLegYmin,resLegXmax+0.05,resLegYmax);
-  leg_tautag2->SetHeader("#splitline{#tau-tagging}{1.5 < |#eta| < 2.5}");
-  leg_tautag2->AddEntry("","","");
-
-  gPad->SetLogx();
-  histos_tautag.second->Draw();
-  addHist(histos_tautag.second, leg_tautag2, 0);
-
-  DrawAxis(histos_tautag.second, leg_tautag2, 0);
-  leg_tautag2->Draw();
-  C_tautag1->cd(0);
-
-  TString jetRes = "jetERes";
-  TString jetResFwd = "jetEResFwd";
-  TCanvas *C_jet = new TCanvas(jetRes,jetRes, 1600, 600);
-  C_jet->Divide(2);
-
-  C_jet->cd(1);
-  TMultiGraph *mg_jet = new TMultiGraph(jetRes,jetRes);
-  TLegend *leg_jet = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
-  leg_jet->SetHeader("#splitline{jets}{|#eta| < 1.5}");
-  leg_jet->AddEntry("","","");
-
-  addGraph(mg_jet, &gr_calojets, leg_jet, 3);
-  addGraph(mg_jet, &gr_pfjets, leg_jet, 1);
-
-  mg_jet->Draw("ACX");
-  leg_jet->Draw();
-
-  DrawAxis(mg_jet, leg_jet, 0.25);
-
-  C_jet->cd(2);
-  TMultiGraph *mg_jetFwd = new TMultiGraph(jetResFwd,jetResFwd);
-  TLegend *leg_jetFwd = new TLegend(resLegXmin,resLegYmin,resLegXmax,resLegYmax);
-  leg_jetFwd->SetHeader("#splitline{jets}{|#eta| < 1.5}");
-  leg_jetFwd->AddEntry("","","");
-
-  addGraph(mg_jetFwd, &gr_calojetsFwd, leg_jetFwd, 3);
-  addGraph(mg_jetFwd, &gr_pfjetsFwd, leg_jetFwd, 1);
-
-  mg_jetFwd->Draw("ACX");
-  leg_jetFwd->Draw();
-
-  DrawAxis(mg_jetFwd, leg_jetFwd, 0.25);
-
-  C_btag1->SaveAs(btagEff+".eps");
-  C_jet->SaveAs(jetRes+".eps");
-
-  TString metRes = "MetRes";
-  TCanvas *C_met = new TCanvas(metRes,metRes, 800, 600);
-
-  TMultiGraph *mg_met = new TMultiGraph(metRes,metRes);
-  TLegend *leg_met = new TLegend(topLeftLegXmin,topLeftLegYmin+0.2,topLeftLegXmax-0.2,topLeftLegYmax);
-  leg_met->SetHeader("E_{T}^{miss}");
-  leg_met->AddEntry("","","");
-
-
-  addGraph(mg_met, &gr_met, leg_met, 0);
-
-  mg_met->Draw("ACX");
-  leg_met->Draw();
-
-  DrawAxis(mg_met, leg_met, 300);
-
-  mg_met->GetXaxis()->SetTitle("H_{T} [GeV]");
-  mg_met->GetYaxis()->SetTitle("#sigma(ME_{x}) [GeV]");
-
-  C_met->SaveAs(metRes+".eps");
-
-  C_jet->Print("validation.pdf","pdf");
-  C_btag1->Print("validation.pdf","pdf");
-  C_tautag1->Print("validation.pdf","pdf");
-  C_met->Print("validation.pdf)","pdf");
-
-  TFile *fout = new TFile(outputFile,"recreate");
-
-  for (int bin = 0; bin < Nbins; bin++)
-  {
-    plots_el.at(bin).cenResolHist->Write();
-    plots_eltrack.at(bin).cenResolHist->Write();
-    plots_eltower.at(bin).cenResolHist->Write();
-    plots_el.at(bin).fwdResolHist->Write();
-    plots_eltrack.at(bin).fwdResolHist->Write();
-    plots_eltower.at(bin).fwdResolHist->Write();
-  }
-
-  histos_el.first->Write();
-  histos_el.second->Write();
-  histos_eltrack.first->Write();
-  histos_eltrack.second->Write();
-  histos_eltower.first->Write();
-  histos_eltower.second->Write();
-  C_el1->Write();
-  C_el2->Write();
-
-  histos_mu.first->Write();
-  histos_mu.second->Write();
-  histos_mutrack.first->Write();
-  histos_mutrack.second->Write();
-  C_mu1->Write();
-  C_mu->Write();
-
-  histos_ph.first->Write();
-  histos_ph.second->Write();
-  C_ph1->Write();
-  C_ph->Write();
-
-  for (int bin = 0; bin < Nbins; bin++)
-  {
-    plots_pfjets.at(bin).cenResolHist->Write();
-    plots_pfjets.at(bin).fwdResolHist->Write();
-    plots_calojets.at(bin).cenResolHist->Write();
-    plots_calojets.at(bin).fwdResolHist->Write();
-    plots_met.at(bin).cenResolHist->Write();
-  }
-  histos_btag.first->Write();
-  histos_btag.second->Write();
-  histos_tautag.first->Write();
-  histos_tautag.second->Write();
-  C_btag1->Write();
-  C_tautag1->Write();
-  C_jet->Write();
-  C_met->Write();
-
-  fout->Write();
-
-  cout << "** Exiting..." << endl;
-
-  delete treeReaderElectron;
-  delete treeReaderMuon;
-  delete treeReaderPhoton;
-  delete treeReaderJet;
-  delete treeReaderBJet;
-  delete treeReaderTauJet;
-  delete chainElectron;
-  delete chainMuon;
-  delete chainPhoton;
-  delete chainJet;
-  delete chainBJet;
-  delete chainTauJet;
-}
+#include "Validation.C"
 
 //------------------------------------------------------------------------------
@@ -1278 +59 @@
   if(argc != 3)
   {
-    cout << " Usage: " << appName << " input_file_electron input_file_muon input_file_photon input_file_jet input_file_bjet input_file_taujet output_file" << endl;
-    cout << " input_file_electron  - input file in ROOT format ('Delphes' tree)," << endl;
-    cout << " input_file_muon - input file in ROOT format ('Delphes' tree)," << endl;
-    cout << " input_file_photon - input file in ROOT format ('Delphes' tree)," << endl;
-    cout << " input_file_jet - input file in ROOT format ('Delphes' tree)," << endl;
-    cout << " input_file_bjet - input file in ROOT format ('Delphes' tree)," << endl;
-    cout << " input_file_taujet - input file in ROOT format ('Delphes' tree)," << endl;
+    cout << " Usage: " << appName << " input_file output_file" << endl;
+    cout << " input_file  - input file in ROOT format ('Delphes' tree)," << endl;
     cout << " output_file - output file in ROOT format" << endl;
     return 1;
@@ -1295 +71 @@
   TApplication app(appName, &appargc, appargv);
 
-  Validation(argv[1], argv[2], argv[3], argv[4], argv[5], argv[6], argv[7]);
+  TString inputFile(argv[1]);
+  TString outputFile(argv[2]);
+
+
+//------------------------------------------------------------------------------
+// Here you call your macro's main function
+
+  Validation(inputFile, outputFile);
+
+//------------------------------------------------------------------------------
+
 }
 

modules/Calorimeter.cc

@@ -58 +58 @@
   fItParticleInputArray(0), fItTrackInputArray(0)
 {
-  
+  Int_t i;
+
   fECalResolutionFormula = new DelphesFormula;
   fHCalResolutionFormula = new DelphesFormula;
 
-  fECalTowerTrackArray = new TObjArray;
-  fItECalTowerTrackArray = fECalTowerTrackArray->MakeIterator();
-
-  fHCalTowerTrackArray = new TObjArray;
-  fItHCalTowerTrackArray = fHCalTowerTrackArray->MakeIterator();
-  
+  for(i = 0; i < 2; ++i)
+  {
+    fECalTowerTrackArray[i] = new TObjArray;
+    fItECalTowerTrackArray[i] = fECalTowerTrackArray[i]->MakeIterator();
+
+    fHCalTowerTrackArray[i] = new TObjArray;
+    fItHCalTowerTrackArray[i] = fHCalTowerTrackArray[i]->MakeIterator();
+  }
 }
 
@@ -74 +77 @@
 Calorimeter::~Calorimeter()
 {
-  
+  Int_t i;
+
   if(fECalResolutionFormula) delete fECalResolutionFormula;
   if(fHCalResolutionFormula) delete fHCalResolutionFormula;
 
-  if(fECalTowerTrackArray) delete fECalTowerTrackArray;
-  if(fItECalTowerTrackArray) delete fItECalTowerTrackArray;
-
-  if(fHCalTowerTrackArray) delete fHCalTowerTrackArray;
-  if(fItHCalTowerTrackArray) delete fItHCalTowerTrackArray;
- 
+  for(i = 0; i < 2; ++i)
+  {
+    if(fECalTowerTrackArray[i]) delete fECalTowerTrackArray[i];
+    if(fItECalTowerTrackArray[i]) delete fItECalTowerTrackArray[i];
+
+    if(fHCalTowerTrackArray[i]) delete fHCalTowerTrackArray[i];
+    if(fItHCalTowerTrackArray[i]) delete fItHCalTowerTrackArray[i];
+  }
 }
 
@@ -212 +218 @@
   Double_t ecalFraction, hcalFraction;
   Double_t ecalEnergy, hcalEnergy;
+  Double_t ecalSigma, hcalSigma;
   Int_t pdgCode;
 
@@ -361 +368 @@
       fHCalTowerEnergy = 0.0;
 
-      fECalTrackEnergy = 0.0;
-      fHCalTrackEnergy = 0.0;
-      
-      fECalTrackSigma = 0.0;
-      fHCalTrackSigma = 0.0;
-      
+      fECalTrackEnergy[0] = 0.0;
+      fECalTrackEnergy[1] = 0.0;
+
+      fHCalTrackEnergy[0] = 0.0;
+      fHCalTrackEnergy[1] = 0.0;
+
       fTowerTrackHits = 0;
       fTowerPhotonHits = 0;
 
-      fECalTowerTrackArray->Clear();
-      fHCalTowerTrackArray->Clear();
-    
+      fECalTowerTrackArray[0]->Clear();
+      fECalTowerTrackArray[1]->Clear();
+
+      fHCalTowerTrackArray[0]->Clear();
+      fHCalTowerTrackArray[1]->Clear();
     }
 
@@ -397 +406 @@
       if(fECalTrackFractions[number] > 1.0E-9 && fHCalTrackFractions[number] < 1.0E-9)
       {
-        fECalTrackEnergy += ecalEnergy;
-        fECalTrackSigma += (track->TrackResolution)*momentum.E()*(track->TrackResolution)*momentum.E();
-        fECalTowerTrackArray->Add(track);
+        ecalSigma = fECalResolutionFormula->Eval(0.0, fTowerEta, 0.0, momentum.E());
+        if(ecalSigma/momentum.E() < track->TrackResolution)
+        {
+          fECalTrackEnergy[0] += ecalEnergy;
+          fECalTowerTrackArray[0]->Add(track);
+        }
+        else
+        {
+          fECalTrackEnergy[1] += ecalEnergy;
+          fECalTowerTrackArray[1]->Add(track);
+        }
       }
-     
       else if(fECalTrackFractions[number] < 1.0E-9 && fHCalTrackFractions[number] > 1.0E-9)
       {
-        fHCalTrackEnergy += hcalEnergy;
-        fHCalTrackSigma += (track->TrackResolution)*momentum.E()*(track->TrackResolution)*momentum.E();
-        fHCalTowerTrackArray->Add(track);
+        hcalSigma = fHCalResolutionFormula->Eval(0.0, fTowerEta, 0.0, momentum.E());
+        if(hcalSigma/momentum.E() < track->TrackResolution)
+        {
+          fHCalTrackEnergy[0] += hcalEnergy;
+          fHCalTowerTrackArray[0]->Add(track);
+        }
+        else
+        {
+          fHCalTrackEnergy[1] += hcalEnergy;
+          fHCalTowerTrackArray[1]->Add(track);
+        }
       }
-      
       else if(fECalTrackFractions[number] < 1.0E-9 && fHCalTrackFractions[number] < 1.0E-9)
       {
@@ -453 +476 @@
   Double_t energy, pt, eta, phi;
   Double_t ecalEnergy, hcalEnergy;
-  Double_t ecalNeutralEnergy, hcalNeutralEnergy;
-  
   Double_t ecalSigma, hcalSigma;
-  Double_t ecalNeutralSigma, hcalNeutralSigma;
-
-  Double_t weightTrack, weightCalo, bestEnergyEstimate, rescaleFactor;
-  
+
   TLorentzVector momentum;
   TFractionMap::iterator itFractionMap;
@@ -466 +484 @@
 
   if(!fTower) return;
+
 
   ecalSigma = fECalResolutionFormula->Eval(0.0, fTowerEta, 0.0, fECalTowerEnergy);
@@ -535 +554 @@
     fTowerOutputArray->Add(fTower);
   }
-  
+
   // fill energy flow candidates
-  fECalTrackSigma = TMath::Sqrt(fECalTrackSigma);
-  fHCalTrackSigma = TMath::Sqrt(fHCalTrackSigma);
-
-  //compute neutral excesses
-  ecalNeutralEnergy = max( (ecalEnergy - fECalTrackEnergy) , 0.0);
-  hcalNeutralEnergy = max( (hcalEnergy - fHCalTrackEnergy) , 0.0);
-  
-  ecalNeutralSigma = ecalNeutralEnergy / TMath::Sqrt(fECalTrackSigma*fECalTrackSigma + ecalSigma*ecalSigma);
-  hcalNeutralSigma = hcalNeutralEnergy / TMath::Sqrt(fHCalTrackSigma*fHCalTrackSigma + hcalSigma*hcalSigma);
-  
-   // if ecal neutral excess is significant, simply create neutral EflowPhoton tower and clone each track into eflowtrack
-  if(ecalNeutralEnergy > fECalEnergyMin && ecalNeutralSigma > fECalEnergySignificanceMin)
+
+  ecalEnergy -= fECalTrackEnergy[1];
+  hcalEnergy -= fHCalTrackEnergy[1];
+
+  fItECalTowerTrackArray[0]->Reset();
+  while((track = static_cast<Candidate*>(fItECalTowerTrackArray[0]->Next())))
+  {
+    mother = track;
+    track = static_cast<Candidate*>(track->Clone());
+    track->AddCandidate(mother);
+
+    track->Momentum *= ecalEnergy/fECalTrackEnergy[0];
+
+    fEFlowTrackOutputArray->Add(track);
+  }
+
+  fItECalTowerTrackArray[1]->Reset();
+  while((track = static_cast<Candidate*>(fItECalTowerTrackArray[1]->Next())))
+  {
+    mother = track;
+    track = static_cast<Candidate*>(track->Clone());
+    track->AddCandidate(mother);
+
+    fEFlowTrackOutputArray->Add(track);
+  }
+
+  fItHCalTowerTrackArray[0]->Reset();
+  while((track = static_cast<Candidate*>(fItHCalTowerTrackArray[0]->Next())))
+  {
+    mother = track;
+    track = static_cast<Candidate*>(track->Clone());
+    track->AddCandidate(mother);
+
+    track->Momentum *= hcalEnergy/fHCalTrackEnergy[0];
+
+    fEFlowTrackOutputArray->Add(track);
+  }
+
+  fItHCalTowerTrackArray[1]->Reset();
+  while((track = static_cast<Candidate*>(fItHCalTowerTrackArray[1]->Next())))
+  {
+    mother = track;
+    track = static_cast<Candidate*>(track->Clone());
+    track->AddCandidate(mother);
+
+    fEFlowTrackOutputArray->Add(track);
+  }
+
+  if(fECalTowerTrackArray[0]->GetEntriesFast() > 0) ecalEnergy = 0.0;
+  if(fHCalTowerTrackArray[0]->GetEntriesFast() > 0) hcalEnergy = 0.0;
+
+  ecalSigma = fECalResolutionFormula->Eval(0.0, fTowerEta, 0.0, ecalEnergy);
+  hcalSigma = fHCalResolutionFormula->Eval(0.0, fTowerEta, 0.0, hcalEnergy);
+
+  if(ecalEnergy < fECalEnergyMin || ecalEnergy < fECalEnergySignificanceMin*ecalSigma) ecalEnergy = 0.0;
+  if(hcalEnergy < fHCalEnergyMin || hcalEnergy < fHCalEnergySignificanceMin*hcalSigma) hcalEnergy = 0.0;
+
+  energy = ecalEnergy + hcalEnergy;
+
+  if(ecalEnergy > 0.0)
   {
     // create new photon tower
     tower = static_cast<Candidate*>(fTower->Clone());
-    pt =  ecalNeutralEnergy / TMath::CosH(eta);
-    
-    tower->Momentum.SetPtEtaPhiE(pt, eta, phi, ecalNeutralEnergy);
-    tower->Eem = ecalNeutralEnergy;
+
+    pt = ecalEnergy / TMath::CosH(eta);
+
+    tower->Momentum.SetPtEtaPhiE(pt, eta, phi, ecalEnergy);
+    tower->Eem = ecalEnergy;
     tower->Ehad = 0.0;
     tower->PID = 22;
-    
+
     fEFlowPhotonOutputArray->Add(tower);
-   
-    //clone tracks
-    fItECalTowerTrackArray->Reset();
-    while((track = static_cast<Candidate*>(fItECalTowerTrackArray->Next())))
-    {
-      mother = track;
-      track = static_cast<Candidate*>(track->Clone());
-      track->AddCandidate(mother);
-
-      fEFlowTrackOutputArray->Add(track);
-    }
-  
-  }
- 
-  // if neutral excess is not significant, rescale eflow tracks, such that the total charged equals the best measurement given by the calorimeter and tracking
-  else if(fECalTrackEnergy > 0.0)
-  {
-    weightTrack = (fECalTrackSigma > 0.0) ? 1 / (fECalTrackSigma*fECalTrackSigma) : 0.0;
-    weightCalo  = (ecalSigma > 0.0) ? 1 / (ecalSigma*ecalSigma) : 0.0;
-  
-    bestEnergyEstimate = (weightTrack*fECalTrackEnergy + weightCalo*ecalEnergy) / (weightTrack + weightCalo); 
-    rescaleFactor = bestEnergyEstimate/fECalTrackEnergy;
-
-    //rescale tracks
-    fItECalTowerTrackArray->Reset();
-    while((track = static_cast<Candidate*>(fItECalTowerTrackArray->Next())))
-    {  
-      mother = track;
-      track = static_cast<Candidate*>(track->Clone());
-      track->AddCandidate(mother);
-
-      track->Momentum *= rescaleFactor;
-
-      fEFlowTrackOutputArray->Add(track);
-    }
-  }
-
-
-  // if hcal neutral excess is significant, simply create neutral EflowNeutralHadron tower and clone each track into eflowtrack
-  if(hcalNeutralEnergy > fHCalEnergyMin && hcalNeutralSigma > fHCalEnergySignificanceMin)
-  {
-    // create new photon tower
+  }
+  if(hcalEnergy > 0.0)
+  {
+    // create new neutral hadron tower
     tower = static_cast<Candidate*>(fTower->Clone());
-    pt =  hcalNeutralEnergy / TMath::CosH(eta);
-    
-    tower->Momentum.SetPtEtaPhiE(pt, eta, phi, hcalNeutralEnergy);
-    tower->Ehad = hcalNeutralEnergy;
+
+    pt = hcalEnergy / TMath::CosH(eta);
+
+    tower->Momentum.SetPtEtaPhiE(pt, eta, phi, hcalEnergy);
     tower->Eem = 0.0;
-    
+    tower->Ehad = hcalEnergy;
+
     fEFlowNeutralHadronOutputArray->Add(tower);
-   
-    //clone tracks
-    fItHCalTowerTrackArray->Reset();
-    while((track = static_cast<Candidate*>(fItHCalTowerTrackArray->Next())))
-    {
-      mother = track;
-      track = static_cast<Candidate*>(track->Clone());
-      track->AddCandidate(mother);
-
-      fEFlowTrackOutputArray->Add(track);
-    }
-  
-  }
- 
-  // if neutral excess is not significant, rescale eflow tracks, such that the total charged equals the best measurement given by the calorimeter and tracking
-  else if(fHCalTrackEnergy > 0.0)
-  {
-    weightTrack = (fHCalTrackSigma > 0.0) ? 1 / (fHCalTrackSigma*fHCalTrackSigma) : 0.0;
-    weightCalo  = (hcalSigma > 0.0) ? 1 / (hcalSigma*hcalSigma) : 0.0;
-  
-    bestEnergyEstimate = (weightTrack*fHCalTrackEnergy + weightCalo*hcalEnergy) / (weightTrack + weightCalo); 
-    rescaleFactor = bestEnergyEstimate / fHCalTrackEnergy;
-
-    //rescale tracks
-    fItHCalTowerTrackArray->Reset();
-    while((track = static_cast<Candidate*>(fItHCalTowerTrackArray->Next())))
-    {  
-      mother = track;
-      track = static_cast<Candidate*>(track->Clone());
-      track->AddCandidate(mother);
-
-      track->Momentum *= rescaleFactor;
-
-      fEFlowTrackOutputArray->Add(track);
-    }
-  }
-  
-  
+  }
 }
 

modules/Calorimeter.h

@@ -58 +58 @@
   Double_t fTowerEta, fTowerPhi, fTowerEdges[4];
   Double_t fECalTowerEnergy, fHCalTowerEnergy;
-  Double_t fECalTrackEnergy, fHCalTrackEnergy;
+  Double_t fECalTrackEnergy[2], fHCalTrackEnergy[2];
 
   Double_t fTimingEnergyMin;
@@ -70 +70 @@
   Double_t fECalEnergySignificanceMin;
   Double_t fHCalEnergySignificanceMin;
-
-  Double_t fECalTrackSigma;
-  Double_t fHCalTrackSigma;
 
   Bool_t fSmearTowerCenter;
@@ -106 +103 @@
   TObjArray *fEFlowNeutralHadronOutputArray; //!
 
-  TObjArray *fECalTowerTrackArray; //!
-  TIterator *fItECalTowerTrackArray; //!
+  TObjArray *fECalTowerTrackArray[2]; //!
+  TIterator *fItECalTowerTrackArray[2]; //!
 
-  TObjArray *fHCalTowerTrackArray; //!
-  TIterator *fItHCalTowerTrackArray; //!
+  TObjArray *fHCalTowerTrackArray[2]; //!
+  TIterator *fItHCalTowerTrackArray[2]; //!
 
   void FinalizeTower();

modules/ImpactParameterSmearing.cc

@@ -96 +96 @@
 {
   Candidate *candidate, *particle, *mother;
-  Double_t xd, yd, zd, d0, sx, sy, sz, dd0;
+  Double_t xd, yd, zd, dxy, sx, sy, sz, ddxy;
   Double_t pt, eta, px, py, phi, e;
 
@@ -103 +103 @@
   {
 
-    // take momentum before smearing (otherwise apply double smearing on d0)
+    // take momentum before smearing (otherwise apply double smearing on dxy)
     particle = static_cast<Candidate*>(candidate->GetCandidates()->At(0));
 
@@ -131 +131 @@
 
     // calculate impact parameter (after-smearing)
-    d0 = (xd*py - yd*px)/pt;
+    dxy = (xd*py - yd*px)/pt;
 
-    dd0 = gRandom->Gaus(0.0, fFormula->Eval(pt, eta, phi, e));
+    ddxy = gRandom->Gaus(0.0, fFormula->Eval(pt, eta, phi, e));
 
     // fill smeared values in candidate
@@ -143 +143 @@
     candidate->Zd = zd;
 
-    candidate->D0 = d0;
-    candidate->ErrorD0 = dd0;
+    candidate->Dxy = dxy;
+    candidate->SDxy = ddxy;
 
     candidate->AddCandidate(mother);

modules/ModulesLinkDef.h

@@ -35 +35 @@
 #include "modules/EnergySmearing.h"
 #include "modules/MomentumSmearing.h"
-#include "modules/TrackSmearing.h"
 #include "modules/ImpactParameterSmearing.h"
 #include "modules/TimeSmearing.h"
@@ -59 +58 @@
 #include "modules/StatusPidFilter.h"
 #include "modules/PdgCodeFilter.h"
-#include "modules/BeamSpotFilter.h"
-#include "modules/RecoPuFilter.h"
 #include "modules/Cloner.h"
 #include "modules/Weighter.h"
@@ -66 +63 @@
 #include "modules/JetFlavorAssociation.h"
 #include "modules/JetFakeParticle.h"
-#include "modules/VertexSorter.h"
-#include "modules/VertexFinder.h"
-#include "modules/VertexFinderDA4D.h"
 #include "modules/ExampleModule.h"
 
@@ -86 +80 @@
 #pragma link C++ class EnergySmearing+;
 #pragma link C++ class MomentumSmearing+;
-#pragma link C++ class TrackSmearing+;
 #pragma link C++ class ImpactParameterSmearing+;
 #pragma link C++ class TimeSmearing+;
@@ -110 +103 @@
 #pragma link C++ class StatusPidFilter+;
 #pragma link C++ class PdgCodeFilter+;
-#pragma link C++ class BeamSpotFilter+;
-#pragma link C++ class RecoPuFilter+;
 #pragma link C++ class Cloner+;
 #pragma link C++ class Weighter+;
@@ -117 +108 @@
 #pragma link C++ class JetFlavorAssociation+;
 #pragma link C++ class JetFakeParticle+;
-#pragma link C++ class VertexSorter+;
-#pragma link C++ class VertexFinder+;
-#pragma link C++ class VertexFinderDA4D+;
 #pragma link C++ class ExampleModule+;
 

modules/ParticlePropagator.cc

@@ -66 +66 @@
 {
 }
-
 
 //------------------------------------------------------------------------------
@@ -92 +91 @@
   fItInputArray = fInputArray->MakeIterator();
 
-  // import beamspot
-  try
-  {
-    fBeamSpotInputArray = ImportArray(GetString("BeamSpotInputArray", "BeamSpotFilter/beamSpotParticle"));
-  }
-  catch(runtime_error &e)
-  {
-    fBeamSpotInputArray = 0;
-  }
   // create output arrays
 
@@ -121 +111 @@
 {
   Candidate *candidate, *mother;
-  TLorentzVector candidatePosition, candidateMomentum, beamSpotPosition;
+  TLorentzVector candidatePosition, candidateMomentum;
   Double_t px, py, pz, pt, pt2, e, q;
   Double_t x, y, z, t, r, phi;
@@ -130 +120 @@
   Double_t tmp, discr, discr2;
   Double_t delta, gammam, omega, asinrho;
-  Double_t rcu, rc2, xd, yd, zd;
-  Double_t l, d0, dz, p, ctgTheta, phip, etap, alpha;
-  Double_t bsx, bsy, bsz;
+  Double_t rcu, rc2, dxy, xd, yd, zd;
 
   const Double_t c_light = 2.99792458E8;
-
-  if (!fBeamSpotInputArray || fBeamSpotInputArray->GetSize () == 0)
-    beamSpotPosition.SetXYZT(0.0, 0.0, 0.0, 0.0);
-  else
-  {
-    Candidate &beamSpotCandidate = *((Candidate *) fBeamSpotInputArray->At(0));
-    beamSpotPosition = beamSpotCandidate.Position;
-  }
 
   fItInputArray->Reset();
@@ -152 +132 @@
     y = candidatePosition.Y()*1.0E-3;
     z = candidatePosition.Z()*1.0E-3;
-
-    bsx = beamSpotPosition.X()*1.0E-3;
-    bsy = beamSpotPosition.Y()*1.0E-3;
-    bsz = beamSpotPosition.Z()*1.0E-3;
-
     q = candidate->Charge;
 
@@ -207 +182 @@
       z_t = z + pz*t;
 
-      l = TMath::Sqrt( (x_t - x)*(x_t - x) + (y_t - y)*(y_t - y) + (z_t - z)*(z_t - z));
-
       mother = candidate;
       candidate = static_cast<Candidate*>(candidate->Clone());
 
-      candidate->InitialPosition = candidatePosition;
       candidate->Position.SetXYZT(x_t*1.0E3, y_t*1.0E3, z_t*1.0E3, candidatePosition.T() + t*e*1.0E3);
-      candidate->L = l*1.0E3;
 
       candidate->Momentum = candidateMomentum;
@@ -268 +239 @@
       zd = z + (TMath::Sqrt(xd*xd + yd*yd) - TMath::Sqrt(x*x + y*y))*pz/pt;
 
-      // use perigee momentum rather than original particle
-      // momentum, since the orignal particle momentum isn't known
-
-      px = TMath::Sign(1.0,r) * pt * (-y_c / r_c);
-      py = TMath::Sign(1.0,r) * pt * (x_c / r_c);
-      etap = candidateMomentum.Eta();
-      phip = TMath::ATan2(py, px);
-
-      candidateMomentum.SetPtEtaPhiE(pt, etap, phip, candidateMomentum.E());
-
-      // calculate additional track parameters (correct for beamspot position)
-
-      d0        = (  (x - bsx) * py - (y - bsy) * px) / pt;
-      dz        = z - ((x - bsx) * px + (y - bsy) * py) / pt * (pz / pt);
-      p         = candidateMomentum.P();
-      ctgTheta  = 1.0 / TMath::Tan (candidateMomentum.Theta ());
-
+      // calculate impact paramater
+      dxy = (xd*py - yd*px)/pt;
 
       // 3. time evaluation t = TMath::Min(t_r, t_z)
@@ -331 +287 @@
       r_t = TMath::Hypot(x_t, y_t);
 
-
-      // compute path length for an helix
-
-      alpha = pz*1.0E9 / c_light / gammam;
-      l = t * TMath::Sqrt(alpha*alpha + r*r*omega*omega);
-
       if(r_t > 0.0)
       {
-
-        // store these variables before cloning
-        candidate->D0 = d0*1.0E3;
-        candidate->DZ = dz*1.0E3;
-        candidate->P  = p;
-        candidate->PT = pt;
-        candidate->CtgTheta = ctgTheta;
-        candidate->Phi = phip;
-
         mother = candidate;
         candidate = static_cast<Candidate*>(candidate->Clone());
 
-        candidate->InitialPosition = candidatePosition;
         candidate->Position.SetXYZT(x_t*1.0E3, y_t*1.0E3, z_t*1.0E3, candidatePosition.T() + t*c_light*1.0E3);
 
         candidate->Momentum = candidateMomentum;
-
-            candidate->L  =  l*1.0E3;
-
-            candidate->Xd = xd*1.0E3;
+        candidate->Dxy = dxy*1.0E3;
+        candidate->Xd = xd*1.0E3;
         candidate->Yd = yd*1.0E3;
         candidate->Zd = zd*1.0E3;

modules/ParticlePropagator.h

@@ -35 +35 @@
 class TClonesArray;
 class TIterator;
-class TLorentzVector;
 
 class ParticlePropagator: public DelphesModule
@@ -56 +55 @@
 
   const TObjArray *fInputArray; //!
-  const TObjArray *fBeamSpotInputArray; //!
 
   TObjArray *fOutputArray; //!

modules/PhotonConversions.cc

@@ -59 +59 @@
   fItInputArray(0), fConversionMap(0), fDecayXsec(0)
 {
-  fDecayXsec = new TF1("decayXsec","1.0 - 4.0/3.0 * x * (1.0 - x)", 0.0, 1.0);
+  fDecayXsec = new TF1;
   fConversionMap = new DelphesCylindricalFormula;
 }
@@ -84 +84 @@
 
   fConversionMap->Compile(GetString("ConversionMap", "0.0"));
+
+#if  ROOT_VERSION_CODE < ROOT_VERSION(6,04,00)
+  fDecayXsec->Compile("1.0 - 4.0/3.0 * x * (1.0 - x)");
+#else
+  fDecayXsec->GetFormula()->Compile("1.0 - 4.0/3.0 * x * (1.0 - x)");
+#endif
+  fDecayXsec->SetRange(0.0, 1.0);
 
   // import array with output from filter/classifier module

modules/PileUpMerger.cc

@@ -115 +115 @@
   TDatabasePDG *pdg = TDatabasePDG::Instance();
   TParticlePDG *pdgParticle;
-  Int_t pid, nch, nvtx = -1;
+  Int_t pid;
   Float_t x, y, z, t, vx, vy;
-  Float_t px, py, pz, e, pt;
-  Double_t dz, dphi, dt, sumpt2, dz0, dt0;
+  Float_t px, py, pz, e;
+  Double_t dz, dphi, dt;
   Int_t numberOfEvents, event, numberOfParticles;
   Long64_t allEntries, entry;
@@ -132 +132 @@
   fFunction->GetRandom2(dz, dt);
 
-  dz0 = -1.0e6;
-  dt0 = -1.0e6;
-
   dt *= c_light*1.0E3; // necessary in order to make t in mm/c
   dz *= 1.0E3; // necessary in order to make z in mm
-
-  //cout<<dz<<","<<dt<<endl;
-
   vx = 0.0;
   vy = 0.0;
-
   numberOfParticles = fInputArray->GetEntriesFast();
-  nch = 0;
-  sumpt2 = 0.0;
-
-  factory = GetFactory();
-  vertex = factory->NewCandidate();
-
   while((candidate = static_cast<Candidate*>(fItInputArray->Next())))
   {
@@ -156 +143 @@
     z = candidate->Position.Z();
     t = candidate->Position.T();
-    pt = candidate->Momentum.Pt();
-
-    // take postion and time from first stable particle
-    if (dz0 < -999999.0)
-      dz0 = z;
-    if (dt0 < -999999.0)
-      dt0 = t;
-
-    // cancel any possible offset in position and time the input file
-    candidate->Position.SetZ(z - dz0 + dz);
-    candidate->Position.SetT(t - dt0 + dt);
-
-    candidate->IsPU = 0;
-
+    candidate->Position.SetZ(z + dz);
+    candidate->Position.SetT(t + dt);
     fParticleOutputArray->Add(candidate);
-
-    if(TMath::Abs(candidate->Charge) >  1.0E-9)
-    {
-      nch++;
-      sumpt2 += pt*pt;
-      vertex->AddCandidate(candidate);
-    }
   }
 
   if(numberOfParticles > 0)
   {
-    vx /= sumpt2;
-    vy /= sumpt2;
-  }
-
-  nvtx++;
+    vx /= numberOfParticles;
+    vy /= numberOfParticles;
+  }
+
+  factory = GetFactory();
+
+  vertex = factory->NewCandidate();
   vertex->Position.SetXYZT(vx, vy, dz, dt);
-  vertex->ClusterIndex = nvtx;
-  vertex->ClusterNDF = nch;
-  vertex->SumPT2 = sumpt2;
-  vertex->GenSumPT2 = sumpt2;
   fVertexOutputArray->Add(vertex);
 
@@ -204 +170 @@
       numberOfEvents = gRandom->Integer(2*fMeanPileUp + 1);
       break;
-    case 2:
-      numberOfEvents = fMeanPileUp;
-      break;
     default:
       numberOfEvents = gRandom->Poisson(fMeanPileUp);
@@ -213 +176 @@
 
   allEntries = fReader->GetEntries();
-
 
   for(event = 0; event < numberOfEvents; ++event)
@@ -236 +198 @@
     vx = 0.0;
     vy = 0.0;
-
     numberOfParticles = 0;
-    sumpt2 = 0.0;
-
-    //factory = GetFactory();
-    vertex = factory->NewCandidate();
-
     while(fReader->ReadParticle(pid, x, y, z, t, px, py, pz, e))
     {
@@ -259 +215 @@
       candidate->Momentum.SetPxPyPzE(px, py, pz, e);
       candidate->Momentum.RotateZ(dphi);
-      pt = candidate->Momentum.Pt();
 
       x -= fInputBeamSpotX;
@@ -269 +224 @@
       vx += candidate->Position.X();
       vy += candidate->Position.Y();
-
       ++numberOfParticles;
-      if(TMath::Abs(candidate->Charge) >  1.0E-9)
-      {
-        nch++;
-        sumpt2 += pt*pt;
-        vertex->AddCandidate(candidate);
-      }
 
       fParticleOutputArray->Add(candidate);
@@ -287 +235 @@
     }
 
-    nvtx++;
-
+    vertex = factory->NewCandidate();
     vertex->Position.SetXYZT(vx, vy, dz, dt);
-
-    vertex->ClusterIndex = nvtx;
-    vertex->ClusterNDF = nch;
-    vertex->SumPT2 = sumpt2;
-    vertex->GenSumPT2 = sumpt2;
-
     vertex->IsPU = 1;
 
     fVertexOutputArray->Add(vertex);
-
-  }
-}
-
-//------------------------------------------------------------------------------
+  }
+}
+
+//------------------------------------------------------------------------------

modules/SimpleCalorimeter.cc

@@ -58 +58 @@
   fItParticleInputArray(0), fItTrackInputArray(0)
 {
-  
+  Int_t i;
+
   fResolutionFormula = new DelphesFormula;
-  fTowerTrackArray = new TObjArray;
-  fItTowerTrackArray = fTowerTrackArray->MakeIterator();
- 
+
+  for(i = 0; i < 2; ++i)
+  {
+    fTowerTrackArray[i] = new TObjArray;
+    fItTowerTrackArray[i] = fTowerTrackArray[i]->MakeIterator();
+  }
 }
 
@@ -69 +73 @@
 SimpleCalorimeter::~SimpleCalorimeter()
 {
- 
+  Int_t i;
+
   if(fResolutionFormula) delete fResolutionFormula;
-  if(fTowerTrackArray) delete fTowerTrackArray;
-  if(fItTowerTrackArray) delete fItTowerTrackArray;
- 
+
+  for(i = 0; i < 2; ++i)
+  {
+    if(fTowerTrackArray[i]) delete fTowerTrackArray[i];
+    if(fItTowerTrackArray[i]) delete fItTowerTrackArray[i];
+  }
 }
 
@@ -190 +198 @@
   Double_t fraction;
   Double_t energy;
+  Double_t sigma;
   Int_t pdgCode;
 
@@ -331 +340 @@
       fTowerEnergy = 0.0;
 
-      fTrackEnergy = 0.0;
-      fTrackSigma = 0.0;
+      fTrackEnergy[0] = 0.0;
+      fTrackEnergy[1] = 0.0;
 
       fTowerTime = 0.0;
@@ -342 +351 @@
       fTowerPhotonHits = 0;
 
-      fTowerTrackArray->Clear();
-     }
+      fTowerTrackArray[0]->Clear();
+      fTowerTrackArray[1]->Clear();
+    }
 
     // check for track hits
@@ -361 +371 @@
       if(fTrackFractions[number] > 1.0E-9)
       {
-             
-       // compute total charged energy   
-       fTrackEnergy += energy;
-       fTrackSigma += ((track->TrackResolution)*momentum.E())*((track->TrackResolution)*momentum.E());
-       
-       fTowerTrackArray->Add(track);
-      
+        sigma = fResolutionFormula->Eval(0.0, fTowerEta, 0.0, momentum.E());
+        if(sigma/momentum.E() < track->TrackResolution)
+        {
+          fTrackEnergy[0] += energy;
+          fTowerTrackArray[0]->Add(track);
+        }
+        else
+        {
+          fTrackEnergy[1] += energy;
+          fTowerTrackArray[1]->Add(track);
+        }
       }
-       
       else
       {
@@ -390 +403 @@
     fTowerEnergy += energy;
 
-    fTowerTime += energy*position.T();
-    fTowerTimeWeight += energy;
+    fTowerTime += TMath::Sqrt(energy)*position.T();
+    fTowerTimeWeight += TMath::Sqrt(energy);
 
     fTower->AddCandidate(particle);
@@ -405 +418 @@
 {
   Candidate *tower, *track, *mother;
-  Double_t energy,neutralEnergy, pt, eta, phi;
-  Double_t sigma, neutralSigma;
+  Double_t energy, pt, eta, phi;
+  Double_t sigma;
   Double_t time;
-   
-  Double_t weightTrack, weightCalo, bestEnergyEstimate, rescaleFactor;
 
   TLorentzVector momentum;
@@ -425 +436 @@
 
   if(energy < fEnergyMin || energy < fEnergySignificanceMin*sigma) energy = 0.0;
-
 
   if(fSmearTowerCenter)
@@ -454 +464 @@
   if(energy > 0.0) fTowerOutputArray->Add(fTower);
 
- 
-  // e-flow candidates
-
-  //compute neutral excess
-  
-  fTrackSigma = TMath::Sqrt(fTrackSigma);
-  neutralEnergy = max( (energy - fTrackEnergy) , 0.0);
-  
-  //compute sigma_trk total
-  neutralSigma = neutralEnergy / TMath::Sqrt(fTrackSigma*fTrackSigma+ sigma*sigma);
-   
-  // if neutral excess is significant, simply create neutral Eflow tower and clone each track into eflowtrack
-  if(neutralEnergy > fEnergyMin && neutralSigma > fEnergySignificanceMin)
+  // fill e-flow candidates
+
+  energy -= fTrackEnergy[1];
+
+  fItTowerTrackArray[0]->Reset();
+  while((track = static_cast<Candidate*>(fItTowerTrackArray[0]->Next())))
+  {
+    mother = track;
+    track = static_cast<Candidate*>(track->Clone());
+    track->AddCandidate(mother);
+
+    track->Momentum *= energy/fTrackEnergy[0];
+
+    fEFlowTrackOutputArray->Add(track);
+  }
+
+  fItTowerTrackArray[1]->Reset();
+  while((track = static_cast<Candidate*>(fItTowerTrackArray[1]->Next())))
+  {
+    mother = track;
+    track = static_cast<Candidate*>(track->Clone());
+    track->AddCandidate(mother);
+
+    fEFlowTrackOutputArray->Add(track);
+  }
+
+  if(fTowerTrackArray[0]->GetEntriesFast() > 0) energy = 0.0;
+
+  sigma = fResolutionFormula->Eval(0.0, fTowerEta, 0.0, energy);
+  if(energy < fEnergyMin || energy < fEnergySignificanceMin*sigma) energy = 0.0;
+
+  // save energy excess as an energy flow tower
+  if(energy > 0.0)
   {
     // create new photon tower
     tower = static_cast<Candidate*>(fTower->Clone());
-    pt = neutralEnergy / TMath::CosH(eta);
-
-    tower->Eem = (!fIsEcal) ? 0 : neutralEnergy;
-    tower->Ehad = (fIsEcal) ? 0 : neutralEnergy;
+    pt = energy / TMath::CosH(eta);
+
+    tower->Eem = (!fIsEcal) ? 0 : energy;
+    tower->Ehad = (fIsEcal) ? 0 : energy;
+    
+    tower->Momentum.SetPtEtaPhiE(pt, eta, phi, energy);
+    
     tower->PID = (fIsEcal) ? 22 : 0;
- 
-    tower->Momentum.SetPtEtaPhiE(pt, eta, phi, neutralEnergy);
+    
     fEFlowTowerOutputArray->Add(tower);
-    
-    fItTowerTrackArray->Reset();
-    while((track = static_cast<Candidate*>(fItTowerTrackArray->Next())))
-    {
-      mother = track;
-      track = static_cast<Candidate*>(track->Clone());
-      track->AddCandidate(mother);
-
-      fEFlowTrackOutputArray->Add(track);
-    }
-  }
-    
-  // if neutral excess is not significant, rescale eflow tracks, such that the total charged equals the best measurement given by the calorimeter and tracking
-  else if (fTrackEnergy > 0.0)
-  {
-    weightTrack = (fTrackSigma > 0.0) ? 1 / (fTrackSigma*fTrackSigma) : 0.0;
-    weightCalo  = (sigma > 0.0) ? 1 / (sigma*sigma) : 0.0;
-  
-    bestEnergyEstimate = (weightTrack*fTrackEnergy + weightCalo*energy) / (weightTrack + weightCalo); 
-    rescaleFactor = bestEnergyEstimate/fTrackEnergy;
-   
-    fItTowerTrackArray->Reset();
-    while((track = static_cast<Candidate*>(fItTowerTrackArray->Next())))
-    {
-      mother = track;
-      track = static_cast<Candidate*>(track->Clone());
-      track->AddCandidate(mother);
-
-      track->Momentum *= rescaleFactor;
-
-      fEFlowTrackOutputArray->Add(track);
-    }
-  }
-   
- } 
+  }
+}
 
 //------------------------------------------------------------------------------

modules/SimpleCalorimeter.h

@@ -59 +59 @@
   Double_t fTowerEta, fTowerPhi, fTowerEdges[4];
   Double_t fTowerEnergy;
-  Double_t fTrackEnergy;
+  Double_t fTrackEnergy[2];
 
   Double_t fTowerTime;
@@ -72 +72 @@
 
   Double_t fEnergySignificanceMin;
-
-  Double_t fTrackSigma;
 
   Bool_t fSmearTowerCenter;
@@ -104 +102 @@
   TObjArray *fEFlowTowerOutputArray; //!
 
-  TObjArray *fTowerTrackArray; //!
-  TIterator *fItTowerTrackArray; //!
+  TObjArray *fTowerTrackArray[2]; //!
+  TIterator *fItTowerTrackArray[2]; //!
 
   void FinalizeTower();

modules/TimeSmearing.cc

@@ -93 +93 @@
 {
   Candidate *candidate, *mother;
-  Double_t ti, tf_smeared, tf;
+  Double_t t;
   const Double_t c_light = 2.99792458E8;
 
@@ -99 +99 @@
   while((candidate = static_cast<Candidate*>(fItInputArray->Next())))
   {
-    const TLorentzVector &candidateInitialPosition = candidate->InitialPosition;
-    const TLorentzVector &candidateFinalPosition = candidate->Position;
-
-    ti = candidateInitialPosition.T()*1.0E-3/c_light;
-    tf = candidateFinalPosition.T()*1.0E-3/c_light;
+    const TLorentzVector &candidatePosition = candidate->Position;
+    t = candidatePosition.T()*1.0E-3/c_light;
 
     // apply smearing formula
-    tf_smeared = gRandom->Gaus(tf, fTimeResolution);
-    ti = ti + tf_smeared - tf;
+    t = gRandom->Gaus(t, fTimeResolution);
 
     mother = candidate;
     candidate = static_cast<Candidate*>(candidate->Clone());
-    candidate->InitialPosition.SetT(ti*1.0E3*c_light);
-    candidate->Position.SetT(tf*1.0E3*c_light);
-
-    candidate->ErrorT = fTimeResolution*1.0E3*c_light;
+    candidate->Position.SetT(t*1.0E3*c_light);
 
     candidate->AddCandidate(mother);

modules/TrackCountingBTagging.cc

@@ -96 +96 @@
 
   Double_t jpx, jpy;
-  Double_t dr, tpt;
-  Double_t xd, yd, d0, dd0, ip, sip;
+  Double_t dr, tpx, tpy, tpt;
+  Double_t xd, yd, dxy, ddxy, ip, sip;
 
   Int_t sign;
@@ -117 +117 @@
     {
       const TLorentzVector &trkMomentum = track->Momentum;
-      
+
       dr = jetMomentum.DeltaR(trkMomentum);
+
       tpt = trkMomentum.Pt();
+      tpx = trkMomentum.Px();
+      tpy = trkMomentum.Py();
+
       xd = track->Xd;
       yd = track->Yd;
-      d0 = TMath::Hypot(xd, yd);
-      dd0 = track->ErrorD0;
+      dxy = TMath::Hypot(xd, yd);
+      ddxy = track->SDxy;
 
       if(tpt < fPtMin) continue;
       if(dr > fDeltaR) continue;
-      if(d0 > fIPmax) continue;
+      if(dxy > fIPmax) continue;
 
       sign = (jpx*xd + jpy*yd > 0.0) ? 1 : -1;
 
-      ip = sign*d0;
-      sip = ip / TMath::Abs(dd0);
+      ip = sign*dxy;
+      sip = ip / TMath::Abs(ddxy);
 
       if(sip > fSigMin) count++;

modules/TreeWriter.cc

@@ -215 +215 @@
     entry->Pz = momentum.Pz();
 
-    entry->D0            = candidate->D0;
-    entry->DZ            = candidate->DZ;
-    entry->P             = candidate->P;
-    entry->PT            = candidate->PT;
-    entry->CtgTheta      = candidate->CtgTheta;
-    entry->Phi           = candidate->Phi;
-
     entry->Eta = eta;
     entry->Phi = momentum.Phi();
@@ -240 +233 @@
 {
   TIter iterator(array);
-  Candidate *candidate = 0, *constituent = 0;
+  Candidate *candidate = 0;
   Vertex *entry = 0;
 
   const Double_t c_light = 2.99792458E8;
 
-  Double_t x, y, z, t, xError, yError, zError, tError, sigma, sumPT2, btvSumPT2, genDeltaZ, genSumPT2;
-  UInt_t index, ndf;
-
-  CompBase *compare = Candidate::fgCompare;
-  Candidate::fgCompare = CompSumPT2<Candidate>::Instance();
-  array->Sort();
-  Candidate::fgCompare = compare;
-
   // loop over all vertices
   iterator.Reset();
   while((candidate = static_cast<Candidate*>(iterator.Next())))
   {
-
-    index = candidate->ClusterIndex;
-    ndf = candidate->ClusterNDF;
-    sigma = candidate->ClusterSigma;
-    sumPT2 = candidate->SumPT2;
-    btvSumPT2 = candidate->BTVSumPT2;
-    genDeltaZ = candidate->GenDeltaZ;
-    genSumPT2 = candidate->GenSumPT2;
-
-    x = candidate->Position.X();
-    y = candidate->Position.Y();
-    z = candidate->Position.Z();
-    t = candidate->Position.T()*1.0E-3/c_light;
-
-    xError = candidate->PositionError.X ();
-    yError = candidate->PositionError.Y ();
-    zError = candidate->PositionError.Z ();
-    tError = candidate->PositionError.T ()*1.0E-3/c_light;
+    const TLorentzVector &position = candidate->Position;
 
     entry = static_cast<Vertex*>(branch->NewEntry());
 
-    entry->Index = index;
-    entry->NDF = ndf;
-    entry->Sigma = sigma;
-    entry->SumPT2 = sumPT2;
-    entry->BTVSumPT2 = btvSumPT2;
-    entry->GenDeltaZ = genDeltaZ;
-    entry->GenSumPT2 = genSumPT2;
-
-    entry->X = x;
-    entry->Y = y;
-    entry->Z = z;
-    entry->T = t;
-
-    entry->ErrorX = xError;
-    entry->ErrorY = yError;
-    entry->ErrorZ = zError;
-    entry->ErrorT = tError;
-
-
-    TIter itConstituents(candidate->GetCandidates());
-    itConstituents.Reset();
-    entry->Constituents.Clear();
-    while((constituent = static_cast<Candidate*>(itConstituents.Next())))
-    {
-      entry->Constituents.Add(constituent);
-    }
-
-  }
-}
-
+    entry->X = position.X();
+    entry->Y = position.Y();
+    entry->Z = position.Z();
+    entry->T = position.T()*1.0E-3/c_light;
+  }
+}
 
 //------------------------------------------------------------------------------
@@ -317 +261 @@
   Candidate *particle = 0;
   Track *entry = 0;
-  Double_t pt, signz, cosTheta, eta, rapidity, p, ctgTheta, phi;
+  Double_t pt, signz, cosTheta, eta, rapidity;
   const Double_t c_light = 2.99792458E8;
 
@@ -348 +292 @@
     entry->TOuter = position.T()*1.0E-3/c_light;
 
-    entry->L = candidate->L;
-
-    entry->D0            = candidate->D0;
-    entry->ErrorD0       = candidate->ErrorD0;
-    entry->DZ            = candidate->DZ;
-    entry->ErrorDZ       = candidate->ErrorDZ;
-
-    entry->ErrorP        = candidate->ErrorP;
-    entry->ErrorPT       = candidate->ErrorPT;
-    entry->ErrorCtgTheta = candidate->ErrorCtgTheta;
-    entry->ErrorPhi      = candidate->ErrorPhi;
-
+    entry->Dxy = candidate->Dxy;
+    entry->SDxy = candidate->SDxy ;
     entry->Xd = candidate->Xd;
     entry->Yd = candidate->Yd;
@@ -367 +301 @@
 
     pt = momentum.Pt();
-    p = momentum.P();
-    phi = momentum.Phi();
-    ctgTheta = (TMath::Tan(momentum.Theta()) != 0) ? 1/TMath::Tan(momentum.Theta()) : 1e10;
-
     cosTheta = TMath::Abs(momentum.CosTheta());
     signz = (momentum.Pz() >= 0.0) ? 1.0 : -1.0;
@@ -376 +306 @@
     rapidity = (cosTheta == 1.0 ? signz*999.9 : momentum.Rapidity());
 
-    entry->PT  = pt;
     entry->Eta = eta;
-    entry->Phi = phi;
-    entry->CtgTheta = ctgTheta;
+    entry->Phi = momentum.Phi();
+    entry->PT = pt;
 
     particle = static_cast<Candidate*>(candidate->GetCandidates()->At(0));
@@ -390 +319 @@
 
     entry->Particle = particle;
-
-    entry->VertexIndex = candidate->ClusterIndex;
-
   }
 }

readers/DelphesCMSFWLite.cpp

@@ -65 +65 @@
   ExRootTreeBranch *branchEvent, ExRootTreeBranch *branchRwgt,
   DelphesFactory *factory, TObjArray *allParticleOutputArray,
-  TObjArray *stableParticleOutputArray, TObjArray *partonOutputArray, Bool_t firstEvent)
+  TObjArray *stableParticleOutputArray, TObjArray *partonOutputArray)
 {
-
   fwlite::Handle< GenEventInfoProduct > handleGenEventInfo;
+
   fwlite::Handle< LHEEventProduct > handleLHEEvent;
+
   fwlite::Handle< vector< reco::GenParticle > > handleParticle;
-
   vector< reco::GenParticle >::const_iterator itParticle;
 
@@ -78 +78 @@
 
   handleGenEventInfo.getByLabel(event, "generator");
-
-  if (!((handleLHEEvent.getBranchNameFor(event, "source")).empty()))
-  { 
-    handleLHEEvent.getByLabel(event, "source");
-  }
-  else if (!((handleLHEEvent.getBranchNameFor(event, "externalLHEProducer")).empty()))
-  {
-    handleLHEEvent.getByLabel(event, "externalLHEProducer");
-  }
-  else if (firstEvent)
-  {
-    std::cout<<"Wrong LHEEvent Label! Please, check the input file."<<std::endl;
-  }
-
-  if (!((handleParticle.getBranchNameFor(event, "genParticles")).empty()))
-  {
-    handleParticle.getByLabel(event, "genParticles");
-  }
-  else if (!((handleParticle.getBranchNameFor(event, "prunedGenParticles")).empty()))
-  {
-    handleParticle.getByLabel(event, "prunedGenParticles");
-  }
-  else
-  {
-    std::cout<<"Wrong GenParticle Label! Please, check the input file."<<std::endl;
-    exit(-1);
-  }
-
-  Bool_t foundLHE = !((handleLHEEvent.getBranchNameFor(event, "source")).empty()) || !((handleLHEEvent.getBranchNameFor(event, "externalLHEProducer")).empty());
+  handleLHEEvent.getByLabel(event, "externalLHEProducer");
+  handleParticle.getByLabel(event, "genParticles");
 
   HepMCEvent *element;
@@ -118 +91 @@
   Double_t px, py, pz, e, mass;
   Double_t x, y, z;
+
+  const vector< gen::WeightsInfo > &vectorWeightsInfo = handleLHEEvent->weights();
+  vector< gen::WeightsInfo >::const_iterator itWeightsInfo;
 
   element = static_cast<HepMCEvent *>(branchEvent->NewEntry());
@@ -141 +117 @@
   element->ProcTime = 0.0;
 
-
-  if(foundLHE)
-  {
-    const vector< gen::WeightsInfo > &vectorWeightsInfo = handleLHEEvent->weights();
-    vector< gen::WeightsInfo >::const_iterator itWeightsInfo;
-    
-    for(itWeightsInfo = vectorWeightsInfo.begin(); itWeightsInfo != vectorWeightsInfo.end(); ++itWeightsInfo)
-    {
-      weight = static_cast<Weight *>(branchRwgt->NewEntry());
-      weight->Weight = itWeightsInfo->wgt;
-    }  
+  for(itWeightsInfo = vectorWeightsInfo.begin(); itWeightsInfo != vectorWeightsInfo.end(); ++itWeightsInfo)
+  {
+    weight = static_cast<Weight *>(branchRwgt->NewEntry());
+    weight->Weight = itWeightsInfo->wgt;
   }
 
@@ -238 +207 @@
   Int_t i;
   Long64_t eventCounter, numberOfEvents;
-  Bool_t firstEvent = kTRUE;
 
   if(argc < 4)
@@ -272 +240 @@
 
     branchEvent = treeWriter->NewBranch("Event", HepMCEvent::Class());
-    branchRwgt = treeWriter->NewBranch("Weight", Weight::Class());
+    branchRwgt = treeWriter->NewBranch("Rwgt", Weight::Class());
 
     confReader = new ExRootConfReader;
@@ -313 +281 @@
       modularDelphes->Clear();
       treeWriter->Clear();
-
       for(event.toBegin(); !event.atEnd() && !interrupted; ++event)
       {
         ConvertInput(event, eventCounter, branchEvent, branchRwgt, factory,
-          allParticleOutputArray, stableParticleOutputArray, partonOutputArray, firstEvent);
+          allParticleOutputArray, stableParticleOutputArray, partonOutputArray);
         modularDelphes->ProcessTask();
-          
-        firstEvent = kFALSE;
 
         treeWriter->Fill();