ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Utils/src/MuonTools.cc
(Generate patch)

Comparing UserCode/MitPhysics/Utils/src/MuonTools.cc (file contents):
Revision 1.1 by pharris, Mon Nov 3 13:19:19 2008 UTC vs.
Revision 1.6 by loizides, Fri Nov 28 09:13:35 2008 UTC

# Line 1 | Line 1
1 < #include "MitWlnu/MuonTools/interface/MuonTools.hh"
1 > // $Id$
2 >
3 > #include "MitPhysics/Utils/interface/MuonTools.h"
4 > #include <TFile.h>
5 >
6   using namespace mithep;
7  
8 < MuonTools::MuonTools() :
9 <  fpion_em_etaEmi(0),
10 <  fpion_had_etaEmi(0),
7 <  fpion_had_etaTmi(0),
8 <  fpion_em_etaB(0),
9 <  fpion_had_etaB(0),
10 <  fpion_ho_etaB(0),
11 <  fpion_had_etaTpl(0),
12 <  fpion_em_etaEpl(0),
13 <  fpion_had_etaEpl(0),
8 > //--------------------------------------------------------------------------------------------------
9 > MuonTools::MuonTools(const char *mutemp, const char *pitemp) :
10 >  fIsInit(kFALSE),
11    fmuon_em_etaEmi(0),
12    fmuon_had_etaEmi(0),
13    fmuon_had_etaTmi(0),
# Line 19 | Line 16 | MuonTools::MuonTools() :
16    fmuon_ho_etaB(0),
17    fmuon_had_etaTpl(0),
18    fmuon_em_etaEpl(0),
19 <  fmuon_had_etaEpl(0) {}
20 < MuonTools::~MuonTools() {
21 <  fPion_templates->Close();
22 <  fMuon_templates->Close();
23 <  delete fpion_em_etaEmi;
24 <  delete fpion_had_etaEmi;
25 <  delete fpion_had_etaTmi;
26 <  delete fpion_em_etaB;
27 <  delete fpion_had_etaB;
28 <  delete fpion_ho_etaB;
29 <  delete fpion_had_etaTpl;
30 <  delete fpion_em_etaEpl;
31 <  delete fpion_had_etaEpl;
32 <  delete fmuon_em_etaEmi;
33 <  delete fmuon_had_etaEmi;
34 <  delete fmuon_had_etaTmi;
35 <  delete fmuon_em_etaB;
36 <  delete fmuon_had_etaB;
37 <  delete fmuon_ho_etaB;
38 <  delete fmuon_had_etaTpl;
39 <  delete fmuon_em_etaEpl;
40 <  delete fmuon_had_etaEpl;
41 < }
42 <
43 < double MuonTools::getCaloCompatability(mithep::Muon* iMuon,bool iEMSpecial, bool iCorrectedHCAL) {
44 <  if(fpion_em_etaEmi == 0) {
45 <    TFile* fPion_templates = new TFile("../../data/PionCaloTemplate.root","READ");
46 <    TFile* fMuon_templates = new TFile("../../data/MuonCaloTemplate.root","READ");
47 <    fpion_em_etaEmi  = (TH2D*) fPion_templates->Get("em_etaEmi");
48 <    fpion_had_etaEmi = (TH2D*) fPion_templates->Get("had_etaEmi");
49 <    fpion_had_etaTmi = (TH2D*) fPion_templates->Get("had_etaTmi");
50 <    fpion_em_etaB    = (TH2D*) fPion_templates->Get("em_etaB")   ;
51 <    fpion_had_etaB   = (TH2D*) fPion_templates->Get("had_etaB")  ;
52 <    fpion_ho_etaB    = (TH2D*) fPion_templates->Get("ho_etaB")   ;
53 <    fpion_had_etaTpl = (TH2D*) fPion_templates->Get("had_etaTpl");
54 <    fpion_em_etaEpl  = (TH2D*) fPion_templates->Get("em_etaEpl") ;
55 <    fpion_had_etaEpl = (TH2D*) fPion_templates->Get("had_etaEpl");
56 <    fmuon_em_etaEmi  = (TH2D*) fMuon_templates->Get("em_etaEmi") ;
57 <    fmuon_had_etaEmi = (TH2D*) fMuon_templates->Get("had_etaEmi");
58 <    fmuon_had_etaTmi = (TH2D*) fMuon_templates->Get("had_etaTmi");
59 <    fmuon_em_etaB    = (TH2D*) fMuon_templates->Get("em_etaB")   ;
60 <    fmuon_had_etaB   = (TH2D*) fMuon_templates->Get("had_etaB")  ;
61 <    fmuon_ho_etaB    = (TH2D*) fMuon_templates->Get("ho_etaB")   ;
62 <    fmuon_had_etaTpl = (TH2D*) fMuon_templates->Get("had_etaTpl");
63 <    fmuon_em_etaEpl  = (TH2D*) fMuon_templates->Get("em_etaEpl");
64 <    fmuon_had_etaEpl = (TH2D*) fMuon_templates->Get("had_etaEpl");
65 <  }
66 <  double lEta = -1.; double lP = -1;
67 <  double lEM  = -5.;      double lHad = 0;      double lHO = 0;
68 <  lEta = iMuon->Eta();
69 <  lP   = iMuon->P();
70 <  if(lP >= 2000.) lP = 1999.9;
71 <  if(!iEMSpecial || iMuon->EmEnergy() != 0.) lEM  = iMuon->EmEnergy();
72 <  lHad = iMuon->HadEnergy();
73 <  lHO  = iMuon->HoEnergy();
74 <  if(lP < 0. )           return 0.5;
75 <  if(fabs(lEta) >  2.5 ) return 0.5;
76 <  TH2D* lTMuonHad = NULL;
77 <  TH2D* lTPionHad = NULL;
78 <  TH2D* lTMuonHo  = NULL;
79 <  TH2D* lTPionHo  = NULL;
80 <  TH2D* lTMuonEm  = NULL;
81 <  TH2D* lTPionEm  = NULL;
82 <  
83 <  if(fabs(lEta) >=  1.27) {
84 <    if(iCorrectedHCAL) lHad *= 1.8/2.2;
85 <    if(lEta > 0) {
19 >  fmuon_had_etaEpl(0),
20 >  fpion_em_etaEmi(0),
21 >  fpion_had_etaEmi(0),
22 >  fpion_had_etaTmi(0),
23 >  fpion_em_etaB(0),
24 >  fpion_had_etaB(0),
25 >  fpion_ho_etaB(0),
26 >  fpion_had_etaTpl(0),
27 >  fpion_em_etaEpl(0),
28 >  fpion_had_etaEpl(0)
29 > {
30 >  // Constructor.
31 >
32 >  if (mutemp && pitemp)
33 >    Init(mutemp, pitemp);
34 > }
35 >
36 > //--------------------------------------------------------------------------------------------------
37 > MuonTools::~MuonTools()
38 > {
39 >  // Destructor.
40 >
41 >  DeleteHistos();
42 > }
43 >
44 > //--------------------------------------------------------------------------------------------------
45 > void MuonTools::DeleteHistos()
46 > {
47 >  // Delete histograms.
48 >
49 >  if (fIsInit) {
50 >    delete fpion_em_etaEmi;
51 >    delete fpion_had_etaEmi;
52 >    delete fpion_had_etaTmi;
53 >    delete fpion_em_etaB;
54 >    delete fpion_had_etaB;
55 >    delete fpion_ho_etaB;
56 >    delete fpion_had_etaTpl;
57 >    delete fpion_em_etaEpl;
58 >    delete fpion_had_etaEpl;
59 >    delete fmuon_em_etaEmi;
60 >    delete fmuon_had_etaEmi;
61 >    delete fmuon_had_etaTmi;
62 >    delete fmuon_em_etaB;
63 >    delete fmuon_had_etaB;
64 >    delete fmuon_ho_etaB;
65 >    delete fmuon_had_etaTpl;
66 >    delete fmuon_em_etaEpl;
67 >    delete fmuon_had_etaEpl;
68 >    fpion_em_etaEmi  = 0;
69 >    fpion_had_etaEmi = 0;
70 >    fpion_had_etaTmi = 0;
71 >    fpion_em_etaB    = 0;
72 >    fpion_had_etaB   = 0;
73 >    fpion_ho_etaB    = 0;
74 >    fpion_had_etaTpl = 0;
75 >    fpion_em_etaEpl  = 0;
76 >    fpion_had_etaEpl = 0;
77 >    fmuon_em_etaEmi  = 0;
78 >    fmuon_had_etaEmi = 0;
79 >    fmuon_had_etaTmi = 0;
80 >    fmuon_em_etaB    = 0;
81 >    fmuon_had_etaB   = 0;
82 >    fmuon_ho_etaB    = 0;
83 >    fmuon_had_etaTpl = 0;
84 >    fmuon_em_etaEpl  = 0;
85 >    fmuon_had_etaEpl = 0;
86 >    fIsInit = kFALSE;
87 >  }
88 > }
89 >
90 > //--------------------------------------------------------------------------------------------------
91 > Double_t MuonTools::GetCaloCompatability(const Muon *iMuon,
92 >                                         Bool_t iEMSpecial, Bool_t iCorrectedHCAL) const
93 > {
94 >  // todo
95 >
96 >  Double_t lEta = iMuon->Eta();
97 >  Double_t aEta = TMath::Abs(lEta);
98 >  if (aEta > 2.5)
99 >    return 0.5;
100 >
101 >  Double_t lP = iMuon->P();
102 >  if (lP >= 2000.)
103 >    lP = 1999.9;
104 >  if(lP < 0. )          
105 >    return 0.5;
106 >
107 >  Double_t lEM  = -5.;      
108 >  if (!iEMSpecial || iMuon->EmEnergy() != 0.)
109 >    lEM  = iMuon->EmEnergy();
110 >
111 >  Double_t lHad = iMuon->HadEnergy();
112 >  Double_t lHO = iMuon->HoEnergy();;
113 >
114 >  TH2D *lTMuonHad = 0;
115 >  TH2D *lTPionHad = 0;
116 >  TH2D *lTMuonHo  = 0;
117 >  TH2D *lTPionHo  = 0;
118 >  TH2D *lTMuonEm  = 0;
119 >  TH2D *lTPionEm  = 0;
120 >    
121 >  if (aEta >= 1.27) {
122 >  if (iCorrectedHCAL)
123 >      lHad *= 1.8/2.2;
124 >    if (lEta > 0) {
125        lTPionHad = fpion_had_etaEpl;
126        lTMuonHad = fmuon_had_etaEpl;
127      } else {
# Line 93 | Line 129 | double MuonTools::getCaloCompatability(m
129        lTMuonHad = fmuon_had_etaEmi;
130      }
131    }
132 <  if(fabs(lEta) <  1.27  && fabs(lEta) >=  1.1 ) {
133 <    if(iCorrectedHCAL)    lHad *= (1.8/(-2.2*fabs(lEta)+5.5));
134 <    if(lEta > 0) {
132 >
133 >  if (aEta < 1.27 && aEta >= 1.1) {
134 >    if (iCorrectedHCAL)    
135 >      lHad *= (1.8/(-2.2*aEta+5.5));
136 >    if (lEta > 0) {
137        lTPionHad  = fpion_had_etaTpl;
138        lTMuonHad  = fmuon_had_etaTpl;
139      } else {
# Line 103 | Line 141 | double MuonTools::getCaloCompatability(m
141        lTMuonHad  = fmuon_had_etaTmi;
142      }
143    }
144 <  if(fabs(lEta) <  1.1) {
145 <    if(iCorrectedHCAL)    lHad *= sin(2*atan(exp(iMuon->Eta())));
144 >
145 >  if (aEta < 1.1) {
146 >    if(iCorrectedHCAL)    
147 >      lHad *= TMath::Sin(2*TMath::ATan(TMath::Exp(lEta))); //todo ask!
148      lTPionHad  = fpion_had_etaB;
149      lTMuonHad  = fmuon_had_etaB;
150    }
151 <  if(lEta >  1.479  ) {
152 <    lTPionEm  = fpion_em_etaEpl;
153 <    lTMuonEm  = fmuon_em_etaEpl;
151 >  if (lEta > 1.479) {
152 >    lTPionEm = fpion_em_etaEpl;
153 >    lTMuonEm = fmuon_em_etaEpl;
154    }
155 <  if(fabs(lEta) <=  1.479) {
155 >  if (aEta <= 1.479) {
156      lTPionEm  = fpion_em_etaB;
157      lTMuonEm  = fmuon_em_etaB;
158    }
159 <  if(lEta < -1.479 ) {
159 >  if (lEta < -1.479) {
160      lTPionEm  = fpion_em_etaEmi;
161      lTMuonEm  = fmuon_em_etaEmi;
162    }
163 <  if(fabs(lEta) < 1.28) {
163 >  if (aEta < 1.28) {
164      lTPionHo  = fpion_ho_etaB;
165      lTMuonHo  = fmuon_ho_etaB;
166    }
167    
168 <  double lPBX = 1.;     double lPSX = 1.;
169 <  double lPBY = 1.;     double lPSY = 1.;
170 <  double lPBZ = 1.;     double lPSZ = 1.;
171 <  if(!overflow(lTPionEm, lP,lEM))  lPBX =  lTPionEm ->GetBinContent(lTPionEm ->GetXaxis()->FindBin(lP),lTPionEm ->GetYaxis()->FindBin(lEM) );
172 <  if(!overflow(lTPionHad,lP,lHad)) lPBY =  lTPionHad->GetBinContent(lTPionHad->GetXaxis()->FindBin(lP),lTPionHad->GetYaxis()->FindBin(lHad));
173 <  if(!overflow(lTPionHo, lP,lHO))  lPBZ =  lTPionHo ->GetBinContent(lTPionHo ->GetXaxis()->FindBin(lP),lTPionHo ->GetYaxis()->FindBin(lHO) );
174 <  if(!overflow(lTMuonEm, lP,lEM )) lPSX =  lTMuonEm ->GetBinContent(lTMuonEm ->GetXaxis()->FindBin(lP),lTMuonEm ->GetYaxis()->FindBin(lEM) );
175 <  if(!overflow(lTMuonHad,lP,lHad)) lPSY =  lTMuonHad->GetBinContent(lTMuonHad->GetXaxis()->FindBin(lP),lTMuonHad->GetYaxis()->FindBin(lHad));
176 <  if(!overflow(lTMuonHo ,lP,lHO))  lPSZ =  lTMuonHo ->GetBinContent(lTMuonHo ->GetXaxis()->FindBin(lP),lTMuonHo ->GetYaxis()->FindBin(lHO) );
168 >  Double_t lPBX = 1.;    
169 >  Double_t lPSX = 1.;
170 >  Double_t lPBY = 1.;    
171 >  Double_t lPSY = 1.;
172 >  Double_t lPBZ = 1.;    
173 >  Double_t lPSZ = 1.;
174 >  if (!Overflow(lTPionEm, lP,lEM))  
175 >    lPBX = lTPionEm ->GetBinContent(lTPionEm ->GetXaxis()->FindBin(lP),
176 >                                     lTPionEm ->GetYaxis()->FindBin(lEM));
177 >  if (!Overflow(lTPionHad,lP,lHad))
178 >    lPBY = lTPionHad->GetBinContent(lTPionHad->GetXaxis()->FindBin(lP),
179 >                                     lTPionHad->GetYaxis()->FindBin(lHad));
180 >  if (!Overflow(lTPionHo, lP,lHO))  
181 >    lPBZ = lTPionHo ->GetBinContent(lTPionHo ->GetXaxis()->FindBin(lP),
182 >                                     lTPionHo ->GetYaxis()->FindBin(lHO));
183 >  if (!Overflow(lTMuonEm, lP,lEM ))
184 >    lPSX = lTMuonEm ->GetBinContent(lTMuonEm ->GetXaxis()->FindBin(lP),
185 >                                     lTMuonEm ->GetYaxis()->FindBin(lEM));
186 >  if (!Overflow(lTMuonHad,lP,lHad))
187 >    lPSY = lTMuonHad->GetBinContent(lTMuonHad->GetXaxis()->FindBin(lP),
188 >                                    lTMuonHad->GetYaxis()->FindBin(lHad));
189 >  if (!Overflow(lTMuonHo ,lP,lHO))  
190 >    lPSZ =  lTMuonHo ->GetBinContent(lTMuonHo ->GetXaxis()->FindBin(lP),
191 >                                     lTMuonHo ->GetYaxis()->FindBin(lHO));
192    
193 <  if(lPSX == 0. || lPBX == 0. || (lEM <= 0. && !iEMSpecial)) {lPSX = 1.; lPBX = 1.;}
194 <  if(lPSY == 0. || lPBY == 0. || lHad == 0.) {lPSY = 1.; lPBY = 1.;}
195 <  if(lPSZ == 0. || lPBZ == 0. || lHO  == 0.) {lPSZ = 1.; lPBZ = 1.;}
196 <  if((lPSX*lPSY*lPSZ+lPBX*lPBY*lPBZ) > 0.) return lPSX*lPSY*lPSZ/(lPSX*lPSY*lPSZ+lPBX*lPBY*lPBZ);
193 >  if (lPSX == 0. || lPBX == 0. || (lEM <= 0. && !iEMSpecial)) {
194 >    lPSX = 1.;
195 >    lPBX = 1.;
196 >  }
197 >  if (lPSY == 0. || lPBY == 0. || lHad == 0.) {
198 >    lPSY = 1.;
199 >    lPBY = 1.;
200 >  }
201 >  if (lPSZ == 0. || lPBZ == 0. || lHO  == 0.) {
202 >    lPSZ = 1.;
203 >    lPBZ = 1.;
204 >  }
205 >  if ((lPSX*lPSY*lPSZ+lPBX*lPBY*lPBZ) > 0.)
206 >    return lPSX*lPSY*lPSZ/(lPSX*lPSY*lPSZ+lPBX*lPBY*lPBZ);
207 >
208    return 0.5;
209   }
210  
211 < bool MuonTools::isGood(mithep::Muon *iMuon,selection iSelection) {
212 <  double lVal = 0;
213 <  switch(iSelection) {
214 <  case AllArbitrated:
215 <    if(iMuon->StandaloneTrk() != 0 || iMuon->GlobalTrk()!= 0)  return true;
216 <    if(iMuon->NSegments() > 0) return true;
217 <    return false;
218 <    break;
219 <  case PromptTight:
220 <    return promptTight(-1,iMuon);
221 <    break;
222 <   case TMOneStationLoose:
223 <     return TMOneStation(iMuon,99999,999999);
224 <     break;
225 <  case TMOneStationTight:
226 <    return TMOneStation(iMuon);
227 <    break;
228 <  case TMLastStationLoose:
229 <    return TMLastStation(iMuon,999999,999999);
230 <    break;
231 <  case TMLastStationTight:
232 <    return TMLastStation(iMuon);
233 <    break;
234 <  case TM2DCompatibilityLoose:
235 <    lVal             = 1.2*getSegmentCompatability(iMuon);
236 <    if(lVal/1.2 == 0.5) return false;
237 <    lVal += 0.8*getCaloCompatability(iMuon,true,true);
238 <    if(lVal > 0.7) return true;
239 <    return false;
240 <    break;
241 <  case TM2DCompatibilityTight:
242 <    lVal             = 1.2*getSegmentCompatability(iMuon);
243 <    if(lVal/1.2 == 0.5) return false;
244 <    lVal += 0.8*getCaloCompatability(iMuon,true,true);
245 <    if(lVal > 1.0) return true;
246 <    return false;
247 <    break;
248 <  default:
249 <    return false;
211 > //--------------------------------------------------------------------------------------------------
212 > Bool_t MuonTools::Init(const char *mutemp, const char *pitemp)
213 > {
214 >  // Read histograms from given files.
215 >
216 >  if (fIsInit) {
217 >    DeleteHistos();
218 >  }
219 >
220 >  TDirectory::TContext context(0);
221 >
222 >  TFile *muon_templates = TFile::Open(mutemp);
223 >  if (!muon_templates) {
224 >    Fatal("Init", "Could not open file %s", mutemp);
225 >    return kFALSE;
226 >  }
227 >  fmuon_em_etaEmi  = LoadHisto("em_etaEmi",  muon_templates);
228 >  fmuon_had_etaEmi = LoadHisto("had_etaEmi", muon_templates);
229 >  fmuon_had_etaTmi = LoadHisto("had_etaTmi", muon_templates);
230 >  fmuon_em_etaB    = LoadHisto("em_etaB",    muon_templates);
231 >  fmuon_had_etaB   = LoadHisto("had_etaB",   muon_templates);
232 >  fmuon_ho_etaB    = LoadHisto("ho_etaB",    muon_templates);
233 >  fmuon_had_etaTpl = LoadHisto("had_etaTpl", muon_templates);
234 >  fmuon_em_etaEpl  = LoadHisto("em_etaEpl",  muon_templates);
235 >  fmuon_had_etaEpl = LoadHisto("had_etaEpl", muon_templates);
236 >  muon_templates->Close();
237 >  delete muon_templates;
238 >
239 >  TFile *pion_templates = TFile::Open(pitemp);
240 >  if (!pion_templates) {
241 >    Fatal("Init", "Could not open file %s", pitemp);
242 >    return kFALSE;
243 >  }
244 >
245 >  fpion_em_etaEmi  = LoadHisto("em_etaEmi",  pion_templates);
246 >  fpion_had_etaEmi = LoadHisto("had_etaEmi", pion_templates);
247 >  fpion_had_etaTmi = LoadHisto("had_etaTmi", pion_templates);
248 >  fpion_em_etaB    = LoadHisto("em_etaB",    pion_templates);
249 >  fpion_had_etaB   = LoadHisto("had_etaB",   pion_templates);
250 >  fpion_ho_etaB    = LoadHisto("ho_etaB",    pion_templates);
251 >  fpion_had_etaTpl = LoadHisto("had_etaTpl", pion_templates);
252 >  fpion_em_etaEpl  = LoadHisto("em_etaEpl",  pion_templates);
253 >  fpion_had_etaEpl = LoadHisto("had_etaEpl", pion_templates);
254 >  pion_templates->Close();
255 >  delete pion_templates;
256 >
257 >  fIsInit = kTRUE;
258 >  return kTRUE;
259 > }
260 >
261 > //--------------------------------------------------------------------------------------------------
262 > Bool_t MuonTools::IsGood(const mithep::Muon *iMuon, ESelType iSel) const
263 > {
264 >  // Return true if given muon qualifies given selection criterium.
265 >
266 >  Double_t tm2dcut = 0.;
267 >
268 >  switch(iSel) {
269 >    case kAllArbitrated:
270 >      if (iMuon->StandaloneTrk() != 0 || iMuon->GlobalTrk()!= 0)  
271 >        return kTRUE;
272 >      if (iMuon->NSegments() > 0)
273 >        return kTRUE;
274 >      break;
275 >    case kPromptTight:
276 >      return iMuon->PromptTight(Muon::kAny);
277 >      break;
278 >    case kTMOneStationLoose:
279 >      return iMuon->TMOneStation(99999,999999);
280 >      break;
281 >    case kTMOneStationTight:
282 >      return iMuon->TMOneStation();
283 >      break;
284 >    case kTMLastStationLoose:
285 >      return iMuon->TMLastStation(999999,999999);
286 >      break;
287 >    case kTMLastStationTight:
288 >      return iMuon->TMLastStation();
289 >      break;
290 >    case kTM2DCompatibilityLoose:
291 >      tm2dcut = 0.7;
292 >      break;
293 >    case kTM2DCompatibilityTight:
294 >      tm2dcut = 1.0;
295 >      break;
296 >    default:
297 >      return kFALSE;
298 >      break;
299 >  }
300 >
301 >  Double_t lVal = 1.2*GetSegmentCompatability(iMuon);
302 >  if (lVal/1.2 == 0.5)
303 >    return kFALSE;
304 >
305 >  lVal += 0.8*GetCaloCompatability(iMuon,kTRUE,kTRUE);
306 >  if (lVal > tm2dcut)
307 >    return kTRUE;
308 >  return kFALSE;
309 > }
310 >
311 > //--------------------------------------------------------------------------------------------------
312 > Double_t MuonTools::GetSegmentCompatability(const mithep::Muon *iMuon) const
313 > {
314 >  // todo
315 >
316 >  Int_t lNStationsCrossed = 0;
317 >  Int_t lNStationsSegment = 0;
318 >  
319 >  Int_t lStSegmentmatch[8];
320 >  Int_t lStCrossed[8];
321 >  Double_t lStBoundary[8];
322 >
323 >  Double_t lWeight  = 0.;
324 >  Bool_t lAdjust    = kTRUE;
325 >  for (Int_t i0 = 0; i0 < 8; ++i0) {
326 >    lStBoundary[i0] = 0.;
327 >    if(iMuon->GetTrackDist(i0) < 999999. ) {
328 >      lNStationsCrossed++;
329 >      lStCrossed[i0]  = 1;
330 >      if (iMuon->GetTrackDist(i0) > -10. )
331 >        lStBoundary[i0] = iMuon->GetTrackDist(i0);
332 >    } else
333 >      lStCrossed[i0]  = 0;
334 >
335 >    if(iMuon->GetDX(i0) < 999999.) { //Use iMuon->GetSegmentX--> CHECK
336 >      lNStationsSegment++;
337 >      lStSegmentmatch[i0] = 1;
338 >    } else
339 >      lStSegmentmatch[i0] = 0;
340 >
341 >    if(iMuon->GetDY(i0) < 999999.)
342 >      lAdjust = kFALSE;
343 >  }
344 >
345 >  if (lNStationsCrossed == 0)
346 >    return 0.5;
347 >
348 >  Double_t lStWeight[8];
349 >  Int_t lPCross = -1;
350 >  const Double_t lAtWeight = 0.5;
351 >  for (Int_t i0 = 0; i0< 8; ++i0) {
352 >    lStWeight[i0] = 0;
353 >    if (lStCrossed[i0] > 0) {
354 >      lPCross++;
355 >
356 >      switch (lNStationsCrossed) {
357 >        case 1 :
358 >          lStWeight[i0] =  1.;
359 >          break;
360 >        case 2 :
361 >          if (lPCross == 0 )
362 >            lStWeight[i0] = 0.33;
363 >          else
364 >            lStWeight[i0] = 0.67;
365 >          break;
366 >        case 3 :
367 >          if (lPCross == 0)
368 >            lStWeight[i0] = 0.23;
369 >          else if (lPCross == 1)
370 >            lStWeight[i0] = 0.33;
371 >          else
372 >            lStWeight[i0] = 0.44;
373 >          break;
374 >        case 4 :
375 >          if (lPCross == 0)
376 >            lStWeight[i0] = 0.10;
377 >          else if (lPCross == 1)
378 >            lStWeight[i0] = 0.20;
379 >          else if (lPCross == 2)
380 >            lStWeight[i0] = 0.30;
381 >          else                    
382 >            lStWeight[i0] = 0.40;
383 >          break;
384 >        default :
385 >          lStWeight[i0] = 1./lNStationsCrossed;
386 >      }
387 >    
388 >      if (lStSegmentmatch[i0] <= 0 && lStBoundary[i0] != 0.)
389 >        lStWeight[i0] *= lAtWeight*0.5*(TMath::Erf(lStBoundary[i0]/6.)+1.);
390 >      else if (lStSegmentmatch[i0] <= 0 && lStBoundary[i0] == 0)
391 >        lStWeight[i0] = 0.;
392 >      
393 >      if (lStSegmentmatch[i0] > 0) {
394 >        Double_t lP2X = TMath::Power(iMuon->GetPullX(i0),2.);
395 >        Double_t lP2Y = TMath::Power(iMuon->GetPullY(i0),2.);
396 >        Double_t lD2X = TMath::Power(iMuon->GetDX(i0),2.);
397 >        Double_t lD2Y = TMath::Power(iMuon->GetDY(i0),2.);
398 >        if (iMuon->GetDY(i0) < 999999 && iMuon->GetDX(i0) < 999999)
399 >          lStWeight[i0] *= SigWeight(TMath::Sqrt(lD2X+lD2Y),TMath::Sqrt(lP2X+lP2Y));
400 >        else if (iMuon->GetDY(i0) >= 999999 && i0 < 4)
401 >          lStWeight[i0] *= SigWeight(iMuon->GetDX(i0),iMuon->GetPullX(i0));
402 >        else if(i0 < 4)
403 >          lStWeight[i0] *= SigWeight(iMuon->GetDY(i0),iMuon->GetPullY(i0));
404 >      }
405 >    }
406 >    lWeight += lStWeight[i0];
407 >  }
408 >
409 >  return lWeight;
410 > }
411 >
412 > //--------------------------------------------------------------------------------------------------
413 > TH2D *MuonTools::LoadHisto(const char *name, TFile *file) const
414 > {
415 >  // Load histogram with given name from given file and return it.
416 >
417 >  TH2D *ret = dynamic_cast<TH2D*>(file->Get(name));
418 >  if (!ret) {
419 >    Fatal("LoadHisto", "Could not load histogram %s from file %s", name, file->GetName());
420 >    return 0;
421    }
422 <  return false;
422 >  ret->SetDirectory(0);
423 >  return ret;
424   }
425 +

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines