Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:29:00

0001 {
0002 gSystem->Load("libFWCoreFWLite.so");
0003 FWLiteEnabler::enable();
0004 /////////SET THESE VALUES///////////////
0005 TFile file("TracksParTest2.root");
0006 TTree * tree = (TTree *) gROOT->FindObject("Events");
0007 TFile outFile("graphsParTest2.root","recreate");
0008 float intervals[] = {0,0.5,1.0,1.5,2.0,2.5};
0009 int bins = 5;
0010 const char * ptinterval  = "4.99<P_{T}<5.01";
0011 const char * etainterval = "|#eta|<2.5";
0012 const char * rs  = "recoTracks_rsWithMaterialTracks__RsWithMaterial";
0013 const char * ctf = "recoTracks_ctfWithMaterialTracks__FinalFits";
0014 const char * sim = "SimTracks_SimG4Object__Test";
0015 ////////////////////////////////////////
0016 vector<float> tot;
0017 vector<float> ptrmsRS;
0018 vector<float> efficRS;
0019 vector<float> ptrmsCTF;
0020 vector<float> efficCTF;
0021 ostringstream title, name, expr, cut;
0022 float binfix=0.0001;
0023 int i;
0024 //****************SIMULATION***************************************
0025 cout << "Processing SIMULATION" << endl;
0026 //compute the number of generated tracks per eta interval
0027 for (i=1;i<(bins+1);i++){
0028   name.str("");
0029   title.str("");
0030   expr.str("");
0031   cut.str("");
0032   cut << "fabs("<<sim<<".obj.momentum().pseudoRapidity())>"<< intervals[i-1] <<"&&fabs("<<sim<<".obj.momentum().pseudoRapidity())<"<<intervals[i];
0033   name << "producedtracks_vs_eta" << intervals[i];
0034   expr << ""<<sim<<".@obj.size() >> " << name.str();
0035   title << name.str();
0036   TH1F prod(name.str().c_str(),title.str().c_str(),10,0,10);
0037   tree->Draw(expr.str().c_str(),cut.str().c_str(),"goff");
0038   tot.push_back(prod.GetEntries()); 
0039 //   prod.Write();
0040 }
0041 //plot number of generated tracks versus eta
0042 name.str("etagen_distrib");
0043 title.str("");
0044 expr.str("");
0045 title <<"Eta generated "<< ptinterval;
0046 expr <<sim<<".obj.momentum().pseudoRapidity() >> etagen_distrib";
0047 TH1F etagen(name.str().c_str(),title.str().c_str(),100,-3,3);
0048 tree->Draw(expr.str().c_str(),"","goff");
0049 etagen.Write();
0050 //*****************************************************************
0051 
0052 //****************ROAD SEARCH**************************************
0053 // cout << "Processing ROAD SEARCH" << endl;
0054 // for (i=1;i<(bins+1);i++){
0055 //   //plot Pt residue distribution per eta interval
0056 //   name <<"ptresRS"<<intervals[i];
0057 //   title <<"Pt residue "<< ptinterval << " "<< intervals[i-1] << "<#eta<"<<intervals[i];
0058 //   TH1F ptres2RS(name.str().c_str(),title.str().c_str(),100,-1,1);
0059 //   expr << ""<<sim<<".obj.momentum().perp()-"<<rs<<".obj.pt() >> " << name.str();
0060 //   cut << "fabs("<<sim<<".obj.momentum().pseudoRapidity())>"<< intervals[i-1] <<"&&fabs("<<sim<<".obj.momentum().pseudoRapidity())<"<<intervals[i];
0061 //   tree->Draw(expr.str().c_str(),cut.str().c_str(),"goff");
0062 //   ptrmsRS.push_back(ptres2RS.GetRMS());
0063 //   ptres2RS.Write();
0064 
0065 //   //plot the number of missing tracks per eta interval
0066 //   //and computes the efficiency in that eta interval
0067 //   title.str("");
0068 //   name.str("");
0069 //   expr.str("");
0070 //   name <<"missingtracks_vs_eta_RS"<<intervals[i];
0071 //   title <<"(Missing Tracks)/(Total Tracks) "<< ptinterval << " "<< intervals[i-1] << "<#eta<"<<intervals[i];
0072 //   expr << "(-"<<sim<<".@obj.size()+"<<rs<<".@obj.size()) >> " << name.str();
0073 //   TH1F effRS(name.str().c_str(),title.str().c_str(),3,-1.5,1.5);
0074 //   tree->Draw(expr.str().c_str(),cut.str().c_str(),"goff");
0075 //   if (tot[i-1]!=0) effRS.Scale(1/tot[i-1]);
0076 // //   cout <<  eff.GetBinContent(eff.FindBin(0)) << endl;
0077 //   efficRS.push_back( effRS.GetBinContent( effRS.FindBin(0) ) );
0078 //   effRS.Write();
0079 // }
0080 
0081 // //plot global Pt residue distribution
0082 // name.str("ptres_distrib_RS");
0083 // title.str("");
0084 // expr.str("");
0085 // title << "Pt residue " << ptinterval << " " << etainterval;
0086 // expr <<sim<<".obj.momentum().perp()-"<<rs<<".obj.pt() >> ptres_distrib_RS";
0087 // TH1F ptresRS(name.str().c_str(),title.str().c_str(),100,-1,1);
0088 // tree->Draw(expr.str().c_str(),"","goff");
0089 // ptresRS.Write();
0090 
0091 // //plot global normalized chi2 distribution
0092 // name.str("chi2_distrib_RS");
0093 // title.str("");
0094 // expr.str("");
0095 // title <<"NChi2 distribution " << ptinterval << " " << etainterval;
0096 // expr <<rs<<".obj.chi2()/"<<rs<<".obj.ndof() >> chi2_distrib_RS";
0097 // TH1F chi2histoRS(name.str().c_str(),title.str().c_str(),100,0,10);
0098 // tree->Draw(expr.str().c_str(),"","goff");
0099 // chi2histoRS.Write();
0100 
0101 // //plot global eta residue distribution
0102 // name.str("etaresidue_distrib_RS");
0103 // title.str("");
0104 // expr.str("");
0105 // title <<"Eta residue "<< ptinterval << " " << etainterval;
0106 // expr <<sim<<".obj.momentum().pseudoRapidity()-"<<rs<<".obj.eta() >> etaresidue_distrib_RS";
0107 // TH1F etaresRS(name.str().c_str(),title.str().c_str(),100,-3,3);
0108 // tree->Draw(expr.str().c_str(),"","goff");
0109 // etaresRS.Write();
0110 
0111 // //plot global eta distribution of reconstructed tracks
0112 // name.str("etafound_distrib_RS");
0113 // title.str("");
0114 // expr.str("");
0115 // title <<"Eta found "<< ptinterval << " " << etainterval;
0116 // expr <<rs<<".obj.eta() >> etafound_distrib_RS";
0117 // TH1F etafoundRS(name.str().c_str(),title.str().c_str(),100,-3,3);
0118 // tree->Draw(expr.str().c_str(),"","goff");
0119 // etafoundRS.Write();
0120 
0121 // //plot number of hits distribution of reconstructed tracks
0122 // name.str("hits_distrib_RS");
0123 // title.str("");
0124 // expr.str("");
0125 // title <<"Hits found "<< ptinterval << " " << etainterval;
0126 // expr <<rs<<".obj.numberOfValidHits() >> hits_distrib_RS";
0127 // TH1F hitsRS(name.str().c_str(),title.str().c_str(),100,0,25);
0128 // tree->Draw(expr.str().c_str(),"","goff");
0129 // hitsRS.Write();
0130 
0131 // //plot number of hits distribution of reconstructed tracks
0132 // //when the number of reconstructed tracks is >1
0133 // name.str("hits2_distrib_RS");
0134 // title.str("");
0135 // expr.str("");
0136 // cut.str("");
0137 // title <<"Hits found (>1tracks) "<< ptinterval << " " << etainterval;
0138 // expr <<rs<<".obj.numberOfValidHits() >> hits2_distrib_RS";
0139 // cut <<rs<<".@obj.size()>1";
0140 // TH1F hits2RS(name.str().c_str(),title.str().c_str(),100,0,25);
0141 // tree->Draw(expr.str().c_str(),cut.str().c_str(),"goff");
0142 // hits2RS.Write();
0143 
0144 // //plot number of hits distribution of reconstructed tracks
0145 // //when the number of reconstructed tracks is =1
0146 // name.str("hits3_distrib_RS");
0147 // title.str("");
0148 // expr.str("");
0149 // cut.str("");
0150 // title <<"Hits found (1track) "<< ptinterval << " " << etainterval;
0151 // expr <<rs<<".obj.numberOfValidHits() >> hits3_distrib_RS";
0152 // cut <<rs<<".@obj.size()==1";
0153 // TH1F hits3RS(name.str().c_str(),title.str().c_str(),100,0,25);
0154 // tree->Draw(expr.str().c_str(),cut.str().c_str(),"goff");
0155 // hits3RS.Write();
0156 
0157 // //plot the Pt RMS per eta interval
0158 // //plot the efficiency per eta interval
0159 // name.str("PtRMS_vs_eta_RS");
0160 // title.str("");
0161 // title <<"PtRMS vs #eta "<< ptinterval;
0162 // TH1F ptrmshRS(name.str().c_str(),title.str().c_str(),bins,intervals);
0163 // name.str("eff_vs_eta_RS");
0164 // title.str("");
0165 // title <<"efficiency vs #eta "<< ptinterval;
0166 // TH1F effhRS(name.str().c_str(),title.str().c_str(),bins,intervals);
0167 // for (int i=1;i<(bins+1);i++){
0168 //   ptrmshRS.Fill(intervals[i]-binfix,ptrmsRS[i-1]);
0169 //   effhRS.Fill(intervals[i]-binfix,efficRS[i-1]);
0170 // }
0171 // effhRS.Write();
0172 // ptrmshRS.Write();
0173 //*****************************************************************
0174 
0175 
0176 //****************COMBINATORIAL TRACK FINDER***********************
0177 cout << "Processing COMBINATORIAL TRACK FINDER" << endl;
0178 for (i=1;i<(bins+1);i++){
0179   //plot Pt residue distribution per eta interval
0180   name.str("");
0181   expr.str("");
0182   title.str("");
0183   cut.str("");
0184   name <<"ptresCTF"<<intervals[i];
0185   title <<"Pt residue "<< ptinterval << " "<< intervals[i-1] << "<#eta<"<<intervals[i];
0186   TH1F ptres2CTF(name.str().c_str(),title.str().c_str(),100,-1,1);
0187   expr << ""<<sim<<".obj.momentum().perp()-"<<ctf<<".obj.pt() >> " << name.str();
0188   cut << "fabs("<<sim<<".obj.momentum().pseudoRapidity())>"<< intervals[i-1] <<"&&fabs("<<sim<<".obj.momentum().pseudoRapidity())<"<<intervals[i];
0189 tree->Draw(expr.str().c_str(),cut.str().c_str(),"goff");
0190   ptrmsCTF.push_back(ptres2CTF.GetRMS());
0191   ptres2CTF.Write();
0192 
0193   //plot the number of missing tracks per eta interval
0194   //and computes the efficiency in that eta interval
0195   name.str("");
0196   expr.str("");
0197   name <<"missingtracks_vs_eta_CTF"<<intervals[i];
0198   expr << "(-"<<sim<<".@obj.size()+"<<ctf<<".@obj.size()) >> " << name.str();
0199   TH1F effCTF(name.str().c_str(),title.str().c_str(),3,-1.5,1.5);
0200   tree->Draw(expr.str().c_str(),cut.str().c_str(),"goff");
0201   if (tot[i-1]!=0) effCTF.Scale(1/tot[i-1]);
0202   efficCTF.push_back( effCTF.GetBinContent( effCTF.FindBin(0) ) );
0203   effCTF.Write();
0204 }
0205 
0206 //plot global Pt residue distribution
0207 name.str("ptres_distrib_CTF");
0208 title.str("");
0209 expr.str("");
0210 title << "Pt residue "  << ptinterval << " " << etainterval;
0211 expr <<sim<<".obj.momentum().perp()-"<<ctf<<".obj.pt() >> ptres_distrib_CTF";
0212 TH1F ptresCTF(name.str().c_str(),title.str().c_str(),100,-1,1);
0213 tree->Draw(expr.str().c_str(),"","goff");
0214 ptresCTF.Write();
0215 
0216 
0217 //plot global normalized chi2 distribution
0218 name.str("chi2_distrib_CTF");
0219 title.str("");
0220 expr.str("");
0221 title <<"NChi2 distribution " << ptinterval << " " << etainterval;
0222 expr <<ctf<<".obj.chi2()/"<<ctf<<".obj.ndof() >> chi2_distrib_CTF";
0223 TH1F chi2histoCTF(name.str().c_str(),title.str().c_str(),100,0,10);
0224 tree->Draw(expr.str().c_str(),"","goff");
0225 chi2histoCTF.Write();
0226 
0227 //plot global eta residue distribution
0228 name.str("etaresidue_distrib_CTF");
0229 title.str("");
0230 expr.str("");
0231 title <<"Eta residue "<< ptinterval << " " << etainterval;
0232 expr <<sim<<".obj.momentum().pseudoRapidity()-"<<ctf<<".obj.eta() >> etaresidue_distrib_CTF";
0233 TH1F etaresCTF(name.str().c_str(),title.str().c_str(),100,-3,3);
0234 tree->Draw(expr.str().c_str(),"","goff");
0235 etaresCTF.Write();
0236 
0237 
0238 //plot global eta distribution of reconstructed tracks
0239 name.str("etafound_distrib_CTF");
0240 title.str("");
0241 expr.str("");
0242 title <<"Eta found "<< ptinterval << " " << etainterval;
0243 expr <<ctf<<".obj.eta() >> etafound_distrib_CTF";
0244 TH1F etafoundCTF(name.str().c_str(),title.str().c_str(),100,-3,3);
0245 tree->Draw(expr.str().c_str(),"","goff");
0246 etafoundCTF.Write();
0247 
0248 //plot number of hits distribution of reconstructed tracks
0249 name.str("hits_distrib_CTF");
0250 title.str("");
0251 expr.str("");
0252 title <<"Hits found "<< ptinterval << " " << etainterval;
0253 expr <<ctf<<".obj.numberOfValidHits() >> hits_distrib_CTF";
0254 TH1F hitsCTF(name.str().c_str(),title.str().c_str(),100,0,25);
0255 tree->Draw(expr.str().c_str(),"","goff");
0256 hitsCTF.Write();
0257 
0258 //plot number of hits distribution of reconstructed tracks
0259 //when the number of reconstructed tracks is >1
0260 name.str("hits2_distrib_CTF");
0261 title.str("");
0262 expr.str("");
0263 cut.str("");
0264 title <<"Hits found (>1tracks) "<< ptinterval << " " << etainterval;
0265 expr <<ctf<<".obj.numberOfValidHits() >> hits2_distrib_CTF";
0266 cut <<ctf<<".@obj.size()>1";
0267 TH1F hits2CTF(name.str().c_str(),title.str().c_str(),100,0,25);
0268 tree->Draw(expr.str().c_str(),cut.str().c_str(),"goff");
0269 hits2CTF.Write();
0270 
0271 //plot the Pt RMS per eta interval
0272 //plot the efficiency per eta interval
0273 name.str("PtRMS_vs_eta_CTF");
0274 title.str("");
0275 title <<"PtRMS vs #eta "<< ptinterval;
0276 TH1F ptrmshCTF(name.str().c_str(),title.str().c_str(),bins,intervals);
0277 name.str("eff_vs_eta_CTF");
0278 title.str("");
0279 title <<"efficiency vs #eta "<< ptinterval;
0280 TH1F effhCTF(name.str().c_str(),title.str().c_str(),bins,intervals);
0281 for (int i=1;i<(bins+1);i++){
0282   ptrmshCTF.Fill(intervals[i]-binfix,ptrmsCTF[i-1]);
0283   effhCTF.Fill(intervals[i]-binfix,efficCTF[i-1]);
0284 }
0285 ptrmshCTF.Write();
0286 effhCTF.Write();
0287 //*****************************************************************
0288 
0289 
0290 //hits3RS.SetLineColor(2);
0291 //hits3RS.Draw();
0292 //hits2RS.Draw("SAME");
0293 file.Close();
0294 outFile.Close();
0295 }