8 gROOT->LoadMacro(
"$VMCWORKDIR/gconfig/rootlogon.C");
10 TString inPidFile =
" evt_pid_tpc.root";
13 TFile *
inFile = TFile::Open(inSimFile,
"READ");
15 TTree *
tree=(TTree *) inFile->Get(
"pndsim") ;
16 tree->AddFriend(
"pndsim",inPidFile);
18 TClonesArray* cand_array=
new TClonesArray(
"PndPidCandidate");
19 tree->SetBranchAddress(
"PidChargedCand", &cand_array);
21 TClonesArray*
mc_array=
new TClonesArray(
"PndMCTrack");
22 tree->SetBranchAddress(
"MCTrack", &mc_array);
24 TFile *
out = TFile::Open(
"out_test_2pi_tpc.root",
"RECREATE");
25 TNtuple *nt =
new TNtuple(
"nt",
"nt",
"evt:id:mc_p:mc_theta:mc_phi:mc_pid:mc_stt:mc_tpc:mc_mvd:mc_gem:p:theta:phi:mult:stt:tpc:mvd");
26 TNtuple *ntEvt =
new TNtuple(
"ntEvt",
"ntEvt",
"evt:acc_stt:acc_tpc:eff:effct:cand:mc_fourpi:reco_fourpi");
27 Float_t mass_k = 0.13957018;
29 if (nEntries==0) nEntries = tree->GetEntriesFast();
30 for (Int_t j=0; j< nEntries; j++){
33 if ((j%100)==0) cout <<
"processing event " << j <<
"\n";
35 Float_t mc_mom = 0, mc_theta = 0, mc_phi = 0;
36 Float_t rec_mom = -1, rec_theta = -1, rec_phi = 0;
37 Int_t stt_mccount = 0, tpc_mccount = 0, reco_stt = 0, reco_tpc = 0,
reco_mvd = 0, reco_gem = 0, reco_count = 0, reco_ctcount = 0;
38 TLorentzVector mc_k[2], reco_k[2], mc_fourpi, reco_fourpi;
39 for (Int_t mc = 0; mc < mc_array->GetEntriesFast(); mc++)
47 mc_theta = mctrack->
GetMomentum().Theta()*TMath::RadToDeg();
48 mc_phi = mctrack->
GetMomentum().Phi()*TMath::RadToDeg();
51 for (Int_t pp=0; pp<cand_array->GetEntriesFast(); pp++)
55 if ( (cand_mult==0) || ((cand_mult>0) && (
fabs(rec_mom-mc_mom)>
fabs(pidCand->
GetMomentum().Mag()-mc_mom))) )
62 reco_tpc = pidCand->GetTpcHits();
70 if (cand_mult>0) reco_count++;
71 if ((cand_mult>0) && ((reco_stt>0) || (reco_tpc>0))) reco_ctcount++;
73 Float_t ntuple_nt[] = {
74 j,mc, mc_mom,mc_theta,mc_phi, mc_pid,
76 rec_mom, rec_theta*TMath::RadToDeg(), rec_phi*TMath::RadToDeg(), cand_mult, reco_stt, reco_tpc,
reco_mvd
81 mc_fourpi = mc_k[0] + mc_k[1];
82 reco_fourpi = reco_k[0] + reco_k[1];
84 Float_t ntuple_evt[] = {
85 j, stt_mccount, tpc_mccount, reco_count, reco_ctcount, cand_array->GetEntriesFast(),
86 mc_fourpi.M(), reco_fourpi.M()
88 ntEvt->Fill(ntuple_evt);
93 nt->Draw(
"mc_theta>>hMcTheta(75,0,150)");
94 nt->Draw(
"mc_theta>>hMcSttTheta(75,0,150)",
"mc_stt>0");
95 nt->Draw(
"mc_theta>>hMcSttMvdTheta(75,0,150)",
"(mc_stt>0)||(mc_mvd>0)");
96 nt->Draw(
"mc_theta>>hMcTpcTheta(75,0,150)",
"mc_tpc>0");
97 nt->Draw(
"mc_theta>>hMcTpcMvdTheta(75,0,150)",
"(mc_tpc>0)||(mc_tpc>0)");
98 nt->Draw(
"mc_theta>>hMcMvdTheta(75,0,150)",
"mc_mvd>0");
99 nt->Draw(
"mc_theta>>hMcGemTheta(75,0,150)",
"mc_gem>0");
100 nt->Draw(
"mc_theta>>hRecoTheta(75,0,150)",
"mult>0");
101 nt->Draw(
"mc_theta>>hRecoCtTheta(75,0,150)",
"mult>0&&(stt>0||tpc>0)");
102 nt->Draw(
"mc_theta>>hRecoEffTheta(75,0,150)",
"mult>0");
103 nt->Draw(
"mc_theta>>hRecoEffCtTheta(75,0,150)",
"mult>0&&(stt>0||tpc>0)");
104 hRecoEffTheta->Divide(hMcTheta);
105 hRecoEffCtTheta->Divide(hMcTheta);
106 nt->Draw(
"mc_theta>>hMcAccSttTheta(75,0,150)",
"mc_stt>0");
107 hMcAccSttTheta->Divide(hMcTheta);
108 nt->Draw(
"mc_theta>>hMcAccTpcTheta(75,0,150)",
"mc_tpc>0");
109 hMcAccTpcTheta->Divide(hMcTheta);
111 nt->Draw(
"mc_phi>>hMcPhi(50,-200,200)");
112 nt->Draw(
"mc_phi>>hMcSttPhi(50,-200,200)",
"mc_stt>0");
113 nt->Draw(
"mc_phi>>hMcSttMvdPhi(50,-200,200)",
"(mc_stt>0)||(mc_mvd>0)");
114 nt->Draw(
"mc_phi>>hMcTpcPhi(50,-200,200)",
"mc_tpc>0");
115 nt->Draw(
"mc_phi>>hMcTpcMvdPhi(50,-200,200)",
"(mc_tpc>0)||(mc_mvd>0)");
116 nt->Draw(
"mc_phi>>hMcMvdPhi(50,-200,200)",
"mc_mvd>0");
117 nt->Draw(
"mc_phi>>hMcGemPhi(50,-200,200)",
"mc_gem>0");
118 nt->Draw(
"mc_phi>>hRecoPhi(50,-200,200)",
"mult>0");
119 nt->Draw(
"mc_phi>>hRecoCtPhi(50,-200,200)",
"mult>0&&(stt>0||tpc>0)");
120 nt->Draw(
"mc_phi>>hRecoEffPhi(50,-200,200)",
"mult>0");
121 nt->Draw(
"mc_phi>>hRecoEffCtPhi(50,-200,200)",
"mult>0&&(stt>0||tpc>0)");
122 hRecoEffPhi->Divide(hMcPhi);
123 hRecoEffCtPhi->Divide(hMcPhi);
124 nt->Draw(
"mc_phi>>hMcAccSttPhi(50,-200,200)",
"mc_stt>0");
125 hMcAccSttPhi->Divide(hMcPhi);
126 nt->Draw(
"mc_phi>>hMcAccTpcPhi(50,-200,200)",
"mc_tpc>0");
127 hMcAccTpcPhi->Divide(hMcPhi);
129 nt->Draw(
"mc_p>>hMcP(100,0,5)");
130 nt->Draw(
"mc_p>>hMcSttP(100,0,5)",
"mc_stt>0");
131 nt->Draw(
"mc_p>>hMcSttMvdP(100,0,5)",
"(mc_stt>0)||(mc_mvd>0)");
132 nt->Draw(
"mc_p>>hMcTpcP(100,0,5)",
"mc_tpc>0");
133 nt->Draw(
"mc_p>>hMcTpcMvdP(100,0,5)",
"(mc_tpc>0)||(mc_mvd>0)");
134 nt->Draw(
"mc_p>>hMcMvdP(100,0,5)",
"mc_mvd>0");
135 nt->Draw(
"mc_p>>hMcGemP(100,0,5)",
"mc_gem>0");
136 nt->Draw(
"mc_p>>hRecoP(100,0,5)",
"mult>0");
137 nt->Draw(
"mc_p>>hRecoCtP(100,0,5)",
"mult>0&&(stt>0||tpc>0)");
138 nt->Draw(
"mc_p>>hRecoEffP(100,0,5)",
"mult>0");
139 nt->Draw(
"mc_p>>hRecoEffCtP(100,0,5)",
"mult>0&&(stt>0||tpc>0)");
140 hRecoEffP->Divide(hMcP);
141 hRecoEffCtP->Divide(hMcP);
142 nt->Draw(
"mc_p>>hMcAccSttP(100,0,5)",
"mc_stt>0");
143 hMcAccSttP->Divide(hMcP);
144 nt->Draw(
"mc_p>>hMcAccTpcP(100,0,5)",
"mc_tpc>0");
145 hMcAccTpcP->Divide(hMcP);
147 nt->Draw(
"(mc_p-p):mc_p>>hresp_p(100,0,5,100,-0.5,0.5)",
"mult>0");
148 nt->Draw(
"(mc_p-p):mc_theta>>hresp_theta(100,0,100,100,-0.5,0.5)",
"mult>0",
"colz");
149 nt->Draw(
"(mc_theta-theta):mc_theta>>hrestheta_theta(100,0,100,100,-2,2)",
"mult>0");
150 nt->Draw(
"(mc_theta-theta):mc_p>>hrestheta_p(100,0,5,100,-2,2)",
"mult>0");
151 nt->Draw(
"(mc_phi-phi):mc_phi>>hresphi_phi(100,-200,200,100,-2,2)",
"mult>0");
154 ntEvt->Draw(
"reco_fourpi>>hpi(100,2.00,4.00)",
"eff==2");
160 hMcTheta->Write(); hMcSttTheta->Write(); hMcSttMvdTheta->Write(); hMcMvdTheta->Write(); hMcGemTheta->Write(); hRecoTheta->Write();
161 hMcTpcTheta->Write(); hMcTpcMvdTheta->Write(); hRecoEffTheta->Write(); hMcAccSttTheta->Write(); hMcAccTpcTheta->Write();
162 hMcPhi->Write(); hMcSttPhi->Write(); hMcSttMvdPhi->Write(); hMcMvdPhi->Write(); hMcGemPhi->Write(); hRecoPhi->Write();
163 hMcTpcPhi->Write(); hMcTpcMvdPhi->Write(); hRecoEffPhi->Write(); hMcAccSttPhi->Write(); hMcAccTpcPhi->Write();
164 hMcP->Write(); hMcSttP->Write(); hMcSttMvdP->Write(); hMcMvdP->Write(); hMcGemP->Write(); hRecoP->Write();
165 hMcTpcP->Write(); hMcTpcMvdP->Write(); hRecoEffP->Write();hMcAccSttP->Write(); hMcAccTpcP->Write();
166 hRecoCtTheta->Write();
169 hRecoEffCtTheta->Write();
170 hRecoEffCtPhi->Write();
171 hRecoEffCtP->Write();
176 hresp_theta->Write();
177 hrestheta_theta->Write();
178 hrestheta_p->Write();
179 hresphi_phi->Write();
Int_t GetNPoints(DetectorId detId) const
TVector3 GetMomentum() const
TLorentzVector Get4Momentum() const
friend F32vec4 fabs(const F32vec4 &a)
Int_t GetMotherID() const
TVector3 GetMomentum() const