Mods/src/ElectronIDMod.cc

// $Id: ElectronIDMod.cc,v 1.18 2009/04/12 13:56:30 ceballos Exp $

#include "MitPhysics/Mods/interface/ElectronIDMod.h"
#include "MitPhysics/Init/interface/ModNames.h"

using namespace mithep;

ClassImp(mithep::ElectronIDMod)

//--------------------------------------------------------------------------------------------------
ElectronIDMod::ElectronIDMod(const char *name, const char *title) : 
  BaseMod(name,title),
  fElectronBranchName(Names::gkElectronBrn),
  fConversionBranchName(Names::gkMvfConversionBrn),
  fGoodElectronsName(ModNames::gkGoodElectronsName),  
  fVertexName(string("PrimaryVertexesBeamSpot").c_str()),
  fElectronIDType("Tight"),
  fElectronIsoType("TrackJuraSliding"),
  fElectronPtMin(10),
  fIDLikelihoodCut(0.9),
  fTrackIsolationCut(5.0),
  fCaloIsolationCut(5.0),
  fEcalJuraIsoCut(5.0),
  fHcalIsolationCut(5.0),
  fApplyConvFilter(kTRUE),
  fApplyD0Cut(kTRUE),
  fD0Cut(0.025),
  fElIdType(kIdUndef),
  fElIsoType(kIsoUndef),
  fElectrons(0),
  fConversions(0),
  fVertices(0)
{
  // Constructor.
}

//--------------------------------------------------------------------------------------------------
void ElectronIDMod::Process()
{
  // Process entries of the tree. 

  LoadBranch(fElectronBranchName);

  ElectronOArr *GoodElectrons = new ElectronOArr;
  GoodElectrons->SetName(fGoodElectronsName);

  for (UInt_t i=0; i<fElectrons->GetEntries(); ++i) {
    const Electron *e = fElectrons->At(i);        

    if (e->Pt() <= fElectronPtMin) 
      continue;
    
    Bool_t idcut = kFALSE;
    switch (fElIdType) {
      case kTight:
        idcut = e->PassTightID();
        break;
      case kLoose:
        idcut = e->PassLooseID();
       break;
      case kLikelihood:
        idcut = (e->IDLikelihood() > fIDLikelihoodCut);
        break;
      case kNoId:
        idcut = kTRUE;
        break;
      case kCustomId:
      default:
        break;
    }

    if (!idcut) 
      continue;

    Bool_t isocut = kFALSE;
    switch (fElIsoType) {
      case kTrackCalo:
        isocut = (e->TrackIsolation() < fTrackIsolationCut) &&
                 (e->CaloIsolation() < fCaloIsolationCut);
        break;
      case kTrackJura:
        isocut = (e->TrackIsolation() < fTrackIsolationCut) &&
                 (e->EcalJurassicIsolation() < fEcalJuraIsoCut) &&
                 (e->HcalIsolation() < fHcalIsolationCut);
        break;
      case kTrackJuraSliding:
        {
          Double_t totalIso = e->TrackIsolation() + e->EcalJurassicIsolation() - 1.5;
          if ((totalIso < (e->Pt()-10.0)*5.0/15.0 && e->Pt() <= 25) ||
              (totalIso < 5.0 && e->Pt() > 25) ||
               totalIso <= 0)
            isocut = kTRUE;
        }
        break;
      case kNoIso:
        isocut = kTRUE;
        break;
      case kCustomIso:
      default:
        break;
    }

    if (!isocut) 
      continue;

    // apply conversion filter
    Bool_t isGoodConversion = kFALSE;
    if (fApplyConvFilter) {
      LoadBranch(fConversionBranchName);
      for (UInt_t ifc=0; ifc<fConversions->GetEntries(); ifc++) {
        
        Bool_t ConversionMatchFound = kFALSE;
        for (UInt_t d=0; d<fConversions->At(ifc)->NDaughters(); d++) {
          const Track *trk = dynamic_cast<const ChargedParticle*>
                       (fConversions->At(ifc)->Daughter(d))->Trk();
          if (e->GsfTrk() == trk) {
            ConversionMatchFound = kTRUE;
            break;
          }
        }

        // if match between the e-track and one of the conversion legs
        if (ConversionMatchFound == kTRUE){
          isGoodConversion =  (fConversions->At(ifc)->Prob() > 0.0005) &&
                              (fConversions->At(ifc)->Lxy() > 0) &&
                              (fConversions->At(ifc)->Lz() > 0);

          if (isGoodConversion == kTRUE) {
            for (UInt_t d=0; d<fConversions->At(ifc)->NDaughters(); d++) {
              const Track *trk = dynamic_cast<const ChargedParticle*>
                           (fConversions->At(ifc)->Daughter(d))->Trk();
                      
              if (trk) {
                // These requirements are not used for the GSF track (i == 1)
                if (!(trk->NHits() > 8 && trk->Prob() > 0.005 && i == 0))
                  isGoodConversion = kFALSE;
              
                const StableData *sd = dynamic_cast<const StableData*>
                                  (fConversions->At(ifc)->DaughterDat(d));
                if (sd->NWrongHits() != 0)
                  isGoodConversion = kFALSE;
              
              } else {
                isGoodConversion = kFALSE;
              }
            }
          }
        }

        if (isGoodConversion == kTRUE) break;

      } // loop over all conversions 
      
    }
    if (isGoodConversion == kTRUE) continue;

    if (fApplyD0Cut) {
      LoadBranch(fVertexName);
      // d0 cut
      double d0_real = 99999;
      for(uint i0 = 0; i0 < fVertices->GetEntries(); i0++) {
        double pD0 = e->GsfTrk()->D0Corrected(*fVertices->At(i0));
        if(TMath::Abs(pD0) < TMath::Abs(d0_real)) d0_real = TMath::Abs(pD0);
      }
      if(d0_real >= fD0Cut) continue;
    }

    // add good electron
    GoodElectrons->Add(e);
  }

  // sort according to pt
  GoodElectrons->Sort();

  // add to event for other modules to use
  AddObjThisEvt(GoodElectrons);  
}

//--------------------------------------------------------------------------------------------------
void ElectronIDMod::SlaveBegin()
{
  // Run startup code on the computer (slave) doing the actual analysis. Here,
  // we just request the electron collection branch.

  ReqBranch(fElectronBranchName, fElectrons);

  if (fApplyConvFilter)
    ReqBranch(fConversionBranchName, fConversions);

  if (fApplyD0Cut)
    ReqBranch(fVertexName, fVertices);

  if (fElectronIDType.CompareTo("Tight") == 0) 
    fElIdType = kTight;
  else if (fElectronIDType.CompareTo("Loose") == 0) 
    fElIdType = kLoose;
  else if (fElectronIDType.CompareTo("Likelihood") == 0) 
    fElIdType = kLikelihood;
  else if (fElectronIDType.CompareTo("NoId") == 0) 
    fElIdType = kNoId;
  else if (fElectronIDType.CompareTo("Custom") == 0) {
    fElIdType = kCustomId;
    SendError(kWarning, "SlaveBegin",
              "Custom electron identification is not yet implemented.");
  } else {
    SendError(kAbortAnalysis, "SlaveBegin",
              "The specified electron identification %s is not defined.",
              fElectronIDType.Data());
    return;
  }

  if (fElectronIsoType.CompareTo("TrackCalo") == 0 )
    fElIsoType = kTrackCalo;
  else if (fElectronIsoType.CompareTo("TrackJura") == 0) 
    fElIsoType = kTrackJura;
  else if(fElectronIsoType.CompareTo("TrackJuraSliding") == 0)
    fElIsoType = kTrackJuraSliding;
  else if (fElectronIsoType.CompareTo("NoIso") == 0 )
    fElIsoType = kNoIso;
  else if (fElectronIsoType.CompareTo("Custom") == 0 ) {
    fElIsoType = kCustomIso;
    SendError(kWarning, "SlaveBegin",
              "Custom electron isolation is not yet implemented.");
  } else {
    SendError(kAbortAnalysis, "SlaveBegin",
              "The specified electron isolation %s is not defined.",
              fElectronIsoType.Data());
    return;
  }
}
Revision:	1.19
Committed:	Wed Apr 15 18:21:33 2009 UTC (16 years ago) by ceballos
Content type:	text/plain
Branch:	MAIN
Changes since 1.18:	+2 -2 lines
Log Message:	improving id
#	User	Rev	Content
1	ceballos	1.19	// $Id: ElectronIDMod.cc,v 1.18 2009/04/12 13:56:30 ceballos Exp $
2	loizides	1.1
3			#include "MitPhysics/Mods/interface/ElectronIDMod.h"
4	loizides	1.5	#include "MitPhysics/Init/interface/ModNames.h"
5	loizides	1.1
6			using namespace mithep;
7
8			ClassImp(mithep::ElectronIDMod)
9
10			//--------------------------------------------------------------------------------------------------
11	loizides	1.5	ElectronIDMod::ElectronIDMod(const char name, const char title) :
12	loizides	1.1	BaseMod(name,title),
13	loizides	1.5	fElectronBranchName(Names::gkElectronBrn),
14	ceballos	1.12	fConversionBranchName(Names::gkMvfConversionBrn),
15	loizides	1.5	fGoodElectronsName(ModNames::gkGoodElectronsName),
16	ceballos	1.12	fVertexName(string("PrimaryVertexesBeamSpot").c_str()),
17	loizides	1.1	fElectronIDType("Tight"),
18	ceballos	1.8	fElectronIsoType("TrackJuraSliding"),
19	loizides	1.1	fElectronPtMin(10),
20			fIDLikelihoodCut(0.9),
21			fTrackIsolationCut(5.0),
22			fCaloIsolationCut(5.0),
23	loizides	1.5	fEcalJuraIsoCut(5.0),
24			fHcalIsolationCut(5.0),
25	loizides	1.14	fApplyConvFilter(kTRUE),
26	ceballos	1.18	fApplyD0Cut(kTRUE),
27	loizides	1.14	fD0Cut(0.025),
28	loizides	1.5	fElIdType(kIdUndef),
29	ceballos	1.12	fElIsoType(kIsoUndef),
30	loizides	1.14	fElectrons(0),
31			fConversions(0),
32			fVertices(0)
33	loizides	1.1	{
34			// Constructor.
35			}
36
37			//--------------------------------------------------------------------------------------------------
38			void ElectronIDMod::Process()
39			{
40			// Process entries of the tree.
41
42	loizides	1.5	LoadBranch(fElectronBranchName);
43	loizides	1.1
44	loizides	1.6	ElectronOArr *GoodElectrons = new ElectronOArr;
45			GoodElectrons->SetName(fGoodElectronsName);
46	loizides	1.1
47	ceballos	1.18	for (UInt_t i=0; i<fElectrons->GetEntries(); ++i) {
48	loizides	1.5	const Electron *e = fElectrons->At(i);
49	loizides	1.1
50	loizides	1.5	if (e->Pt() <= fElectronPtMin)
51			continue;
52	loizides	1.1
53	loizides	1.5	Bool_t idcut = kFALSE;
54			switch (fElIdType) {
55			case kTight:
56			idcut = e->PassTightID();
57			break;
58			case kLoose:
59			idcut = e->PassLooseID();
60			break;
61			case kLikelihood:
62			idcut = (e->IDLikelihood() > fIDLikelihoodCut);
63			break;
64	loizides	1.11	case kNoId:
65			idcut = kTRUE;
66			break;
67	loizides	1.5	case kCustomId:
68			default:
69			break;
70	loizides	1.1	}
71
72	loizides	1.5	if (!idcut)
73			continue;
74
75			Bool_t isocut = kFALSE;
76			switch (fElIsoType) {
77			case kTrackCalo:
78			isocut = (e->TrackIsolation() < fTrackIsolationCut) &&
79			(e->CaloIsolation() < fCaloIsolationCut);
80			break;
81			case kTrackJura:
82			isocut = (e->TrackIsolation() < fTrackIsolationCut) &&
83			(e->EcalJurassicIsolation() < fEcalJuraIsoCut) &&
84			(e->HcalIsolation() < fHcalIsolationCut);
85			break;
86			case kTrackJuraSliding:
87	ceballos	1.18	{
88	loizides	1.5	Double_t totalIso = e->TrackIsolation() + e->EcalJurassicIsolation() - 1.5;
89	ceballos	1.12	if ((totalIso < (e->Pt()-10.0)*5.0/15.0 && e->Pt() <= 25) \|\|
90			(totalIso < 5.0 && e->Pt() > 25) \|\|
91	ceballos	1.7	totalIso <= 0)
92	loizides	1.5	isocut = kTRUE;
93			}
94			break;
95			case kNoIso:
96			isocut = kTRUE;
97			break;
98			case kCustomIso:
99			default:
100			break;
101	loizides	1.1	}
102
103	loizides	1.5	if (!isocut)
104			continue;
105
106	loizides	1.14	// apply conversion filter
107	ceballos	1.12	Bool_t isGoodConversion = kFALSE;
108	loizides	1.14	if (fApplyConvFilter) {
109	ceballos	1.12	LoadBranch(fConversionBranchName);
110			for (UInt_t ifc=0; ifc<fConversions->GetEntries(); ifc++) {
111
112			Bool_t ConversionMatchFound = kFALSE;
113			for (UInt_t d=0; d<fConversions->At(ifc)->NDaughters(); d++) {
114			const Track trk = dynamic_cast<const ChargedParticle>
115			(fConversions->At(ifc)->Daughter(d))->Trk();
116	ceballos	1.13	if (e->GsfTrk() == trk) {
117	ceballos	1.12	ConversionMatchFound = kTRUE;
118			break;
119			}
120			}
121
122			// if match between the e-track and one of the conversion legs
123			if (ConversionMatchFound == kTRUE){
124	ceballos	1.19	isGoodConversion = (fConversions->At(ifc)->Prob() > 0.0005) &&
125	ceballos	1.12	(fConversions->At(ifc)->Lxy() > 0) &&
126			(fConversions->At(ifc)->Lz() > 0);
127
128			if (isGoodConversion == kTRUE) {
129			for (UInt_t d=0; d<fConversions->At(ifc)->NDaughters(); d++) {
130			const Track trk = dynamic_cast<const ChargedParticle>
131			(fConversions->At(ifc)->Daughter(d))->Trk();
132
133			if (trk) {
134			// These requirements are not used for the GSF track (i == 1)
135			if (!(trk->NHits() > 8 && trk->Prob() > 0.005 && i == 0))
136			isGoodConversion = kFALSE;
137
138			const StableData sd = dynamic_cast<const StableData>
139			(fConversions->At(ifc)->DaughterDat(d));
140			if (sd->NWrongHits() != 0)
141			isGoodConversion = kFALSE;
142
143			} else {
144			isGoodConversion = kFALSE;
145			}
146			}
147			}
148			}
149
150			if (isGoodConversion == kTRUE) break;
151
152			} // loop over all conversions
153
154			}
155			if (isGoodConversion == kTRUE) continue;
156
157	ceballos	1.15	if (fApplyD0Cut) {
158			LoadBranch(fVertexName);
159			// d0 cut
160			double d0_real = 99999;
161			for(uint i0 = 0; i0 < fVertices->GetEntries(); i0++) {
162			double pD0 = e->GsfTrk()->D0Corrected(*fVertices->At(i0));
163			if(TMath::Abs(pD0) < TMath::Abs(d0_real)) d0_real = TMath::Abs(pD0);
164			}
165			if(d0_real >= fD0Cut) continue;
166	ceballos	1.12	}
167
168	loizides	1.5	// add good electron
169	ceballos	1.12	GoodElectrons->Add(e);
170	loizides	1.5	}
171	loizides	1.1
172	loizides	1.9	// sort according to pt
173			GoodElectrons->Sort();
174
175	loizides	1.5	// add to event for other modules to use
176	loizides	1.6	AddObjThisEvt(GoodElectrons);
177	loizides	1.1	}
178
179			//--------------------------------------------------------------------------------------------------
180			void ElectronIDMod::SlaveBegin()
181			{
182			// Run startup code on the computer (slave) doing the actual analysis. Here,
183	loizides	1.5	// we just request the electron collection branch.
184	loizides	1.1
185	loizides	1.5	ReqBranch(fElectronBranchName, fElectrons);
186	loizides	1.17
187	ceballos	1.16	if (fApplyConvFilter)
188	ceballos	1.15	ReqBranch(fConversionBranchName, fConversions);
189	loizides	1.17
190	ceballos	1.15	if (fApplyD0Cut)
191			ReqBranch(fVertexName, fVertices);
192	loizides	1.1
193	loizides	1.5	if (fElectronIDType.CompareTo("Tight") == 0)
194			fElIdType = kTight;
195			else if (fElectronIDType.CompareTo("Loose") == 0)
196			fElIdType = kLoose;
197			else if (fElectronIDType.CompareTo("Likelihood") == 0)
198			fElIdType = kLikelihood;
199	loizides	1.10	else if (fElectronIDType.CompareTo("NoId") == 0)
200			fElIdType = kNoId;
201	loizides	1.5	else if (fElectronIDType.CompareTo("Custom") == 0) {
202			fElIdType = kCustomId;
203			SendError(kWarning, "SlaveBegin",
204			"Custom electron identification is not yet implemented.");
205			} else {
206			SendError(kAbortAnalysis, "SlaveBegin",
207			"The specified electron identification %s is not defined.",
208			fElectronIDType.Data());
209			return;
210			}
211
212			if (fElectronIsoType.CompareTo("TrackCalo") == 0 )
213			fElIsoType = kTrackCalo;
214			else if (fElectronIsoType.CompareTo("TrackJura") == 0)
215			fElIsoType = kTrackJura;
216			else if(fElectronIsoType.CompareTo("TrackJuraSliding") == 0)
217			fElIsoType = kTrackJuraSliding;
218			else if (fElectronIsoType.CompareTo("NoIso") == 0 )
219			fElIsoType = kNoIso;
220			else if (fElectronIsoType.CompareTo("Custom") == 0 ) {
221			fElIsoType = kCustomIso;
222			SendError(kWarning, "SlaveBegin",
223			"Custom electron isolation is not yet implemented.");
224			} else {
225			SendError(kAbortAnalysis, "SlaveBegin",
226			"The specified electron isolation %s is not defined.",
227			fElectronIsoType.Data());
228			return;
229			}
230	loizides	1.1	}