]>
Commit | Line | Data |
---|---|---|
20963cb1 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | #include <TCanvas.h> | |
17 | #include <TChain.h> | |
18 | #include <TF1.h> | |
19 | #include <TFile.h> | |
20 | #include <TH1F.h> | |
21 | #include <TH2F.h> | |
22 | #include <TLegend.h> | |
23 | #include <TROOT.h> | |
24 | #include <TVector3.h> | |
25 | ||
6e06db1d | 26 | #include <AliESDEvent.h> |
20963cb1 | 27 | #include <AliLog.h> |
28 | #include <AliPID.h> | |
29 | ||
30 | #include <AliAnalysisTask.h> //qa() | |
31 | #include <AliAnalysisManager.h> //qa() | |
32 | #include <TBenchmark.h> //qa() | |
33 | #include <TProof.h> //qa() | |
34 | ||
35 | class AliHMPIDQaEsd : public AliAnalysisTask { | |
36 | ||
37 | public: | |
38 | AliHMPIDQaEsd() ; | |
39 | virtual ~AliHMPIDQaEsd() ; | |
40 | ||
41 | virtual void Exec(Option_t * opt = "") ; | |
42 | virtual void ConnectInputData(Option_t *); | |
43 | virtual void CreateOutputObjects(); | |
44 | virtual void Terminate(Option_t * opt = ""); | |
45 | ||
46 | private: | |
47 | TTree * fChain ; //!pointer to the analyzed TTree or TChain | |
6e06db1d | 48 | AliESDEvent * fESD ; //! Declaration of leave types |
20963cb1 | 49 | |
50 | TObjArray * fOutputContainer; //output data container | |
51 | ||
52 | ClassDef(AliHMPIDQaEsd,0); // | |
53 | }; | |
54 | ||
55 | ||
56 | ||
57 | //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |
58 | AliHMPIDQaEsd::AliHMPIDQaEsd():AliAnalysisTask("HmpidQaTask",""), fChain(0), fESD(0) | |
59 | { | |
60 | // Constructor. | |
61 | DefineInput (0,TChain::Class()); // Input slot #0 works with an Ntuple | |
62 | DefineOutput(0,TObjArray::Class()) ; // Output slot #0 writes into a TH1 container | |
63 | } | |
64 | //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |
65 | AliHMPIDQaEsd::~AliHMPIDQaEsd() | |
66 | { | |
67 | // dtor | |
68 | fOutputContainer->Clear() ; delete fOutputContainer ; | |
69 | } | |
70 | //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |
71 | void AliHMPIDQaEsd::ConnectInputData(const Option_t*) | |
72 | { | |
73 | //Virtual from AliAnalysisTask invoked by AliAnalysisTask::CheckNotify() which in turn invoked by AliAnalysisDataContainer::SetData() | |
74 | fChain = dynamic_cast<TChain *>(GetInputData(0)) ; | |
75 | if (!fChain) {AliError(Form("Input 0 for %s not found\n", GetName())); return;} | |
76 | ||
77 | // One should first check if the branch address was taken by some other task | |
78 | char ** address = (char **)GetBranchAddress(0, "ESD"); | |
79 | if (address) { | |
6e06db1d | 80 | fESD = (AliESDEvent*)(*address); |
20963cb1 | 81 | } else { |
6e06db1d | 82 | fESD = new AliESDEvent(); |
e4290ede | 83 | fESD->ReadFromTree(fChain); //clm: new ESD access works for local, need to test it for PROOF! |
84 | //SetBranchAddress(0, "esdTree", &fESD); | |
85 | //fChain->SetBranchStatus("*", 1); | |
86 | //fChain->SetBranchStatus("fTracks.*", 1); | |
20963cb1 | 87 | } |
88 | }//ConnectInputData() | |
89 | //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |
90 | void AliHMPIDQaEsd::CreateOutputObjects() | |
91 | { | |
92 | ||
93 | ||
94 | // create output container | |
95 | fOutputContainer = new TObjArray(9) ; fOutputContainer->SetName(GetName()) ; | |
96 | ||
97 | fOutputContainer->AddAt(new TH2F("CkovP" , "#theta_{c}, [rad];P, [GeV]" , 150, 0, 7 ,100, 0, 1) , 0) ; | |
98 | fOutputContainer->AddAt(new TH2F("SigP" ,"#sigma_{#theta_c} [mrad];[GeV]", 150, 0, 7 ,100, 0, 1) , 1) ; | |
99 | fOutputContainer->AddAt(new TH2F("MipXY" ,"mip position" , 260, 0,130 ,252, 0,126) , 2) ; | |
100 | fOutputContainer->AddAt(new TH2F("DifXY" ,"diff" , 200, -10, 10 ,200,-10,10) , 3) ; | |
101 | fOutputContainer->AddAt(new TH1F("PidE" ,"PID: e yellow #mu magenta" ,100,0,1) , 4) ; | |
102 | fOutputContainer->AddAt(new TH1F("PidMu","pid of #mu" ,100,0,1) , 5) ; | |
103 | fOutputContainer->AddAt(new TH1F("PidPi","PID: #pi red K green p blue",100,0,1) , 6) ; | |
104 | fOutputContainer->AddAt(new TH1F("PidK" ,"pid of K" ,100,0,1) , 7) ; | |
105 | fOutputContainer->AddAt(new TH1F("PidP" ,"pid of p" ,100,0,1) , 8) ; | |
106 | //options for drawing | |
107 | ||
108 | }//CreateOutputObjects() | |
109 | //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |
110 | void AliHMPIDQaEsd::Exec(Option_t *) | |
111 | { | |
112 | // Virtual from TTask. | |
113 | // Invoked by AliAnalysisManager::StartAnalysis()->AliAnalysisManager::ExecAnalysis()->TTask::ExecuteTask() in case of mgr->StartAnalysis("local") | |
114 | // Invoked by AliAnalysisSelector::Process()->AliAnalysisManager::ExecAnalysis()->TTask::ExecuteTask() in case of mgr->StartAnalysis("local") | |
115 | ||
116 | fChain->GetReadEntry() ; | |
117 | ||
118 | if (!fESD) { | |
119 | AliError("fESD is not connected to the input!") ; | |
120 | return ; | |
121 | } | |
122 | ||
123 | for(Int_t iTrk = 0 ; iTrk < fESD->GetNumberOfTracks() ; iTrk++){ | |
124 | AliESDtrack *pTrk = fESD->GetTrack(iTrk) ; | |
125 | ||
126 | ((TH2F*)fOutputContainer->At(0))->Fill( pTrk->GetP(), pTrk->GetHMPIDsignal() ) ; | |
127 | ((TH2F*)fOutputContainer->At(1))->Fill( pTrk->GetP(), TMath::Sqrt(pTrk->GetHMPIDchi2()) ) ; | |
128 | ||
129 | Float_t xm,ym; Int_t q,np; pTrk->GetHMPIDmip(xm,ym,q,np); //mip info | |
130 | ((TH2F*)fOutputContainer->At(2))->Fill(xm,ym); | |
131 | Float_t xRad,yRad,th,ph; pTrk->GetHMPIDtrk(xRad,yRad,th,ph); //track info at the middle of the radiator | |
132 | Float_t xPc = xRad+9.25*TMath::Tan(th)*TMath::Cos(ph); // temporar: linear extrapol (B=0!) | |
133 | Float_t yPc = yRad+9.25*TMath::Tan(th)*TMath::Sin(ph); // temporar: " | |
134 | ((TH2F*)fOutputContainer->At(3))->Fill(xm-xPc,ym-yPc); //track info | |
135 | ||
136 | Double_t pid[5] ; pTrk->GetHMPIDpid(pid) ; | |
137 | for(Int_t i = 0 ; i < 5 ; i++) ((TH1F*)fOutputContainer->At(4+i))->Fill(pid[i]) ; | |
138 | }//tracks loop | |
139 | ||
140 | PostData(0,fOutputContainer); | |
141 | }//Exec() | |
142 | //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |
143 | void AliHMPIDQaEsd::Terminate(Option_t *) | |
144 | { | |
145 | //Virual from Processing when the event loop is ended | |
146 | TObjArray *out=(TObjArray*)GetOutputData(0); | |
147 | ||
148 | TH2F *hAngP = (TH2F*)out->At(0); | |
149 | TH2F *hErrP = (TH2F*)out->At(1); | |
150 | TH2F *hMipXY = (TH2F*)out->At(2); | |
151 | TH2F *hDifXY = (TH2F*)out->At(3); | |
152 | TH1F *hProE = (TH1F*)out->At(4); | |
153 | TH1F *hProMu = (TH1F*)out->At(5); | |
154 | TH1F *hProPi = (TH1F*)out->At(6); | |
155 | TH1F *hProK = (TH1F*)out->At(7); | |
156 | TH1F *hProP = (TH1F*)out->At(8); | |
157 | ||
158 | hProE ->SetLineColor(kYellow); | |
159 | hProMu->SetLineColor(kMagenta); | |
160 | hProPi->SetLineColor(kRed); | |
161 | hProK ->SetLineColor(kGreen); | |
162 | hProP ->SetLineColor(kBlue); | |
163 | ||
164 | Float_t n = 1.292 ; //mean freon ref idx | |
165 | AliPID dummy ; //just to initialize AliPID to get the correct particle masses | |
166 | TF1* funPi = new TF1("RiPiTheo", "acos(sqrt(x*x+[0]*[0])/(x*[1]))", 1.2, 7); funPi->SetLineWidth(1); funPi->SetParameter(1,n) ; | |
167 | ||
168 | funPi->SetLineColor(kRed); funPi->SetParameter(0,AliPID::ParticleMass(AliPID::kPion)); | |
169 | TF1* funK=static_cast<TF1*>(funPi->Clone()) ; funK ->SetLineColor(kGreen) ; funK->SetParameter(0,AliPID::ParticleMass(AliPID::kKaon)) ; | |
170 | TF1* funP=static_cast<TF1*>(funPi->Clone()) ; funP ->SetLineColor(kBlue) ; funP->SetParameter(0,AliPID::ParticleMass(AliPID::kProton)) ; | |
171 | ||
172 | TCanvas * can = new TCanvas("HmpidCanvas","HMPID ESD Test"); can->SetFillColor(10) ; can->SetHighLightColor(10) ; can->Divide(3,2) ; | |
173 | ||
174 | can->cd(1);hAngP->Draw();funPi->Draw("same");funK->Draw("same");funP->Draw("same"); can->cd(2);hMipXY->Draw(); can->cd(3);hProE->Draw();hProMu->Draw("same"); | |
175 | can->cd(4);hErrP->Draw(); can->cd(5);hDifXY->Draw(); can->cd(6);hProPi->Draw();hProK->Draw("same");hProP->Draw("same"); | |
176 | }//Terminate() | |
177 | //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |
178 | void qa(Int_t mode=0) | |
179 | { | |
e4290ede | 180 | |
181 | /* | |
182 | AliAODHandler* aodHandler = new AliAODHandler(); | |
183 | mgr->SetEventHandler(aodHandler); | |
184 | */ | |
185 | ||
20963cb1 | 186 | gBenchmark->Start("HMPID QA"); |
187 | ||
188 | TChain* chain =new TChain("esdTree"); | |
e4290ede | 189 | AliAnalysisManager *mgr=new AliAnalysisManager("FunnyName"); //clm: |
190 | //AliAODHandler* aodHandler = new AliAODHandler(); | |
191 | //mgr->SetEventHandler(aodHandler); | |
192 | ||
20963cb1 | 193 | AliAnalysisTask *qa=new AliHMPIDQaEsd(); |
194 | qa->ConnectInput (0,mgr->CreateContainer("EsdChain",TChain::Class() ,AliAnalysisManager::kInputContainer)); | |
195 | qa->ConnectOutput(0,mgr->CreateContainer("HistLst",TObjArray::Class(),AliAnalysisManager::kOutputContainer)); | |
196 | ||
197 | mgr->AddTask(qa); | |
198 | if(!mgr->InitAnalysis()) return; | |
199 | mgr->PrintStatus(); | |
200 | ||
201 | switch(mode){ | |
202 | case 0: chain->Add("AliESDs.root"); | |
203 | mgr->StartAnalysis("local",chain); | |
204 | break; | |
205 | ||
206 | case 1: if(TProof::Open("proof://hmpid@lxb6046.cern.ch")==0x0) return; | |
207 | gProof->UploadPackage("ESD.par"); gProof->EnablePackage("ESD"); | |
208 | gProof->UploadPackage("ANALYSIS.par"); gProof->EnablePackage("ANALYSIS"); | |
209 | mgr->StartAnalysis("proof",chain); | |
210 | break; | |
211 | ||
212 | case 2: mgr->StartAnalysis("grid" ,chain); | |
213 | break; | |
214 | } | |
215 | gBenchmark->Show("HMPID QA"); | |
216 | } | |
217 | //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |
218 |