ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/UfCode/UserArea/StatTools/cuts.C
Revision: 1.1
Committed: Wed Dec 8 10:57:12 2010 UTC (14 years, 4 months ago) by kkotov
Content type: text/plain
Branch: MAIN
Log Message:
*** empty log message ***

File Contents

# User Rev Content
1 kkotov 1.1 #include <TChain.h>
2     #include <iostream>
3    
4     // muon info
5     typedef struct {
6     int charge;
7     float pt;
8     float eta;
9     float phi;
10    
11     float normChiSquare;
12     float d0;
13    
14     int numValidMuonHits;
15     int numValidPixelHits;
16     int numValidTrackerHits;
17     int numValidStripHits;
18     int numSegmentMatches;
19    
20     float trackIsoSumPt;
21    
22     int isHltMatched;
23     char hltPath[128];
24    
25     } _MuonInfo;
26    
27     // true info
28     typedef struct {
29     int charge;
30     float pt;
31     float eta;
32     float phi;
33    
34     } _TrueInfo;
35    
36     bool isLoose(_MuonInfo& muon) {
37     bool isLoose=false;
38    
39     // kinematic cuts
40     if (muon.pt < 20. ) return isLoose; // pt cut
41     if (fabs(muon.eta) > 2.4) return isLoose; // eta cut
42     if (muon.numValidTrackerHits < 10) return isLoose; // # hits in tracker
43    
44     // iso cut
45     if (muon.trackIsoSumPt>10) return isLoose;
46    
47     // beam spot cut
48     if (fabs(muon.d0) > 0.2) return isLoose;
49    
50     // if all the cuts are passed
51     isLoose=true;
52     return isLoose;
53     }
54    
55     bool isKinTight(_MuonInfo& muon) {
56     bool isKinTight=false;
57    
58     // minimum requirement to start investigating: to be loose
59     if ( !isLoose(muon) ) return isKinTight;
60    
61     if ( muon.eta > 2.1 ) return isKinTight;
62     if ( muon.numValidMuonHits < 1 ) return isKinTight;
63     if ( muon.numValidPixelHits < 1 ) return isKinTight;
64     if ( muon.numSegmentMatches < 2 ) return isKinTight;
65     if ( muon.normChiSquare > 10) return isKinTight;
66    
67     isKinTight=true;
68     return isKinTight;
69     }
70    
71     bool cuts(void){
72     TChain *tree = new TChain("tree");
73     // tree->AddFile("tup_ci_M20.root");
74     // tree->AddFile("tup_qgf_M15.root");
75     // tree->AddFile("Zprime_1.5.root");
76     tree->AddFile("dimuons_DYToMuMu_M-20_CT10_TuneZ2_7TeV-Fall10.root");
77     // tree->AddFile("dimuons_run132440to147195.root");
78     // tree->AddFile("dimuons_run147196to148058.root");
79     // tree->AddFile("dimuons_run148059to149442.root");
80    
81     _MuonInfo reco1, reco2;
82    
83     tree->SetBranchAddress("reco1", &reco1);
84     tree->SetBranchAddress("reco2", &reco2);
85    
86     float recoMass, recoPt;
87     tree->SetBranchAddress("recoCandMassCktl", &recoMass);
88     tree->SetBranchAddress("recoCandPtCktl", &recoPt );
89     // tree->SetBranchAddress("recoCandMass", &recoMass);
90     // tree->SetBranchAddress("recoCandPt", &recoPt );
91    
92     _TrueInfo true1, true2;
93    
94     tree->SetBranchAddress("true1", &true1);
95     tree->SetBranchAddress("true2", &true2);
96    
97     float trueMass, truePt;
98     tree->SetBranchAddress("trueMass", &trueMass);
99     tree->SetBranchAddress("truePt", &truePt );
100    
101     TH1F *qT = new TH1F("qT","", 10000,0,1700);
102     TH1F *mass = new TH1F("mass","",100,60,120);
103    
104     const double mLow = 60., mHigh = 120.;
105    
106     int pass=0, total=0, event=0;
107     for(event=0; tree->GetEvent(event); event++){
108     // geometrical reconstruction acceptance
109    
110     if (reco1.charge*reco2.charge>0) continue;
111     if (recoMass < mLow) continue;
112     if (recoMass > mHigh) continue;
113     mass->Fill(recoMass);
114    
115     // if( recoPt<400 ) continue;
116     //if( fabs(reco1.eta) > 2.1 || fabs(reco2.eta) > 2.1) continue;
117     /*
118     bool match=false;
119     if( reco1.charge*reco2.charge<0 && reco1.charge>-10 && reco2.charge>-10 ){
120     if( fabs(reco1.eta-true1.eta)<0.1 && reco1.charge*true1.charge>0 &&
121     fabs(reco2.eta-true2.eta)<0.1 && reco2.charge*true2.charge>0 ){
122     reco1.pt -= (reco1.pt-true1.pt)*0.25;
123     reco2.pt -= (reco2.pt-true2.pt)*0.25;
124     match=true;
125     } else {
126     if( fabs(reco1.eta-true2.eta)<0.1 && reco1.charge*true2.charge>0 &&
127     fabs(reco2.eta-true1.eta)<0.1 && reco2.charge*true1.charge>0 ){
128     reco1.pt -= (reco1.pt-true2.pt)*0.25;
129     reco2.pt -= (reco2.pt-true1.pt)*0.25;
130     match=true;
131     }
132     }
133     }
134     if(!match)continue; else
135     */
136     total++;
137     //kinematical selection efficiency (no trigger)
138     if( !isLoose(reco1) || !isLoose(reco2) ) continue; // both loose
139     if( !isKinTight(reco1) && !isKinTight(reco2) ) continue; // at least one tight
140    
141     qT->Fill(recoPt);
142    
143     pass++;
144     }
145    
146     double eff = double(pass)/double(total), err = sqrt(eff*(1-eff)/total);
147     cout<<"Pass="<<pass<<" total="<<total<<" eff="<<eff<<"+-"<<err<<endl;
148    
149     }