HMPID/AliHMPIDQaEsd.C

   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
  26 #include <AliESD.h>
  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
  48   AliESD  * fESD ;              //! Declaration of leave types
  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) {
  80     fESD = (AliESD*)(*address);
  81   } else {
  82     fESD = new AliESD();
  83     SetBranchAddress(0, "ESD", &fESD);
  84     fChain->SetBranchStatus("*", 1);
  85     fChain->SetBranchStatus("fTracks.*", 1);
  86   }
  87 }//ConnectInputData()
  88 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  89 void AliHMPIDQaEsd::CreateOutputObjects()
  90 {
  91
  92
  93   // create output container
  94   fOutputContainer = new TObjArray(9) ;  fOutputContainer->SetName(GetName()) ;
  95
  96   fOutputContainer->AddAt(new TH2F("CkovP" , "#theta_{c}, [rad];P, [GeV]"   , 150,   0,  7  ,100, 0, 1)  ,      0) ;
  97   fOutputContainer->AddAt(new TH2F("SigP"  ,"#sigma_{#theta_c} [mrad];[GeV]", 150,   0,  7  ,100, 0, 1)  ,      1) ;
  98   fOutputContainer->AddAt(new TH2F("MipXY" ,"mip position"                  , 260,   0,130  ,252, 0,126) ,      2) ;
  99   fOutputContainer->AddAt(new TH2F("DifXY" ,"diff"                          , 200, -10, 10  ,200,-10,10) ,      3) ;
 100   fOutputContainer->AddAt(new TH1F("PidE" ,"PID: e yellow #mu magenta"  ,100,0,1)                        ,      4) ;
 101   fOutputContainer->AddAt(new TH1F("PidMu","pid of #mu"                 ,100,0,1)                        ,      5) ;
 102   fOutputContainer->AddAt(new TH1F("PidPi","PID: #pi red K green p blue",100,0,1)                        ,      6) ;
 103   fOutputContainer->AddAt(new TH1F("PidK" ,"pid of K"                   ,100,0,1)                        ,      7) ;
 104   fOutputContainer->AddAt(new TH1F("PidP" ,"pid of p"                   ,100,0,1)                        ,      8) ;
 105   //options for drawing
 106
 107 }//CreateOutputObjects()
 108 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 109 void AliHMPIDQaEsd::Exec(Option_t *)
 110 {
 111 // Virtual from TTask.
 112 // Invoked by AliAnalysisManager::StartAnalysis()->AliAnalysisManager::ExecAnalysis()->TTask::ExecuteTask() in case of mgr->StartAnalysis("local")
 113 // Invoked by AliAnalysisSelector::Process()->AliAnalysisManager::ExecAnalysis()->TTask::ExecuteTask() in case of mgr->StartAnalysis("local")
 114
 115   fChain->GetReadEntry() ;
 116
 117   if (!fESD) {
 118     AliError("fESD is not connected to the input!") ;
 119     return ;
 120   }
 121
 122   for(Int_t iTrk = 0 ; iTrk < fESD->GetNumberOfTracks() ; iTrk++){
 123     AliESDtrack *pTrk = fESD->GetTrack(iTrk) ;
 124
 125     ((TH2F*)fOutputContainer->At(0))->Fill( pTrk->GetP(), pTrk->GetHMPIDsignal() ) ;
 126     ((TH2F*)fOutputContainer->At(1))->Fill( pTrk->GetP(), TMath::Sqrt(pTrk->GetHMPIDchi2()) ) ;
 127
 128      Float_t xm,ym; Int_t q,np;  pTrk->GetHMPIDmip(xm,ym,q,np);                       //mip info
 129      ((TH2F*)fOutputContainer->At(2))->Fill(xm,ym);
 130      Float_t xRad,yRad,th,ph;        pTrk->GetHMPIDtrk(xRad,yRad,th,ph);              //track info at the middle of the radiator
 131      Float_t xPc = xRad+9.25*TMath::Tan(th)*TMath::Cos(ph); // temporar: linear extrapol (B=0!)
 132      Float_t yPc = yRad+9.25*TMath::Tan(th)*TMath::Sin(ph); // temporar:          "
 133      ((TH2F*)fOutputContainer->At(3))->Fill(xm-xPc,ym-yPc); //track info
 134
 135     Double_t pid[5] ;      pTrk->GetHMPIDpid(pid) ;
 136     for(Int_t i = 0 ; i < 5 ; i++)       ((TH1F*)fOutputContainer->At(4+i))->Fill(pid[i]) ;
 137   }//tracks loop
 138
 139   PostData(0,fOutputContainer);
 140 }//Exec()
 141 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 142 void AliHMPIDQaEsd::Terminate(Option_t *)
 143 {
 144 //Virual from Processing when the event loop is ended
 145   TObjArray *out=(TObjArray*)GetOutputData(0);
 146
 147   TH2F *hAngP    = (TH2F*)out->At(0);
 148   TH2F *hErrP    = (TH2F*)out->At(1);
 149   TH2F *hMipXY   = (TH2F*)out->At(2);
 150   TH2F *hDifXY   = (TH2F*)out->At(3);
 151   TH1F *hProE    = (TH1F*)out->At(4);
 152   TH1F *hProMu   = (TH1F*)out->At(5);
 153   TH1F *hProPi   = (TH1F*)out->At(6);
 154   TH1F *hProK    = (TH1F*)out->At(7);
 155   TH1F *hProP    = (TH1F*)out->At(8);
 156
 157   hProE ->SetLineColor(kYellow);
 158   hProMu->SetLineColor(kMagenta);
 159   hProPi->SetLineColor(kRed);
 160   hProK ->SetLineColor(kGreen);
 161   hProP ->SetLineColor(kBlue);
 162
 163   Float_t n = 1.292 ; //mean freon ref idx
 164   AliPID dummy ;      //just to initialize AliPID to get the correct particle masses
 165   TF1* funPi = new TF1("RiPiTheo", "acos(sqrt(x*x+[0]*[0])/(x*[1]))", 1.2, 7); funPi->SetLineWidth(1);   funPi->SetParameter(1,n) ;
 166
 167                                                 funPi->SetLineColor(kRed);     funPi->SetParameter(0,AliPID::ParticleMass(AliPID::kPion));
 168   TF1* funK=static_cast<TF1*>(funPi->Clone()) ; funK ->SetLineColor(kGreen) ;   funK->SetParameter(0,AliPID::ParticleMass(AliPID::kKaon)) ;
 169   TF1* funP=static_cast<TF1*>(funPi->Clone()) ; funP ->SetLineColor(kBlue) ;    funP->SetParameter(0,AliPID::ParticleMass(AliPID::kProton)) ;
 170
 171   TCanvas * can = new TCanvas("HmpidCanvas","HMPID ESD Test"); can->SetFillColor(10) ;   can->SetHighLightColor(10) ;   can->Divide(3,2) ;
 172
 173   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");
 174   can->cd(4);hErrP->Draw();                                                           can->cd(5);hDifXY->Draw();  can->cd(6);hProPi->Draw();hProK->Draw("same");hProP->Draw("same");
 175 }//Terminate()
 176 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 177 void qa(Int_t mode=0)
 178 {
 179   gBenchmark->Start("HMPID QA");
 180
 181   TChain* chain =new TChain("esdTree");
 182   AliAnalysisManager *mgr=new AliAnalysisManager("FunnyName");
 183
 184   AliAnalysisTask *qa=new AliHMPIDQaEsd();
 185   qa->ConnectInput (0,mgr->CreateContainer("EsdChain",TChain::Class()   ,AliAnalysisManager::kInputContainer));
 186   qa->ConnectOutput(0,mgr->CreateContainer("HistLst",TObjArray::Class(),AliAnalysisManager::kOutputContainer));
 187
 188   mgr->AddTask(qa);
 189   if(!mgr->InitAnalysis()) return;
 190   mgr->PrintStatus();
 191
 192   switch(mode){
 193     case 0:   chain->Add("AliESDs.root");
 194               mgr->StartAnalysis("local",chain);
 195               break;
 196
 197     case 1:   if(TProof::Open("proof://hmpid@lxb6046.cern.ch")==0x0) return;
 198               gProof->UploadPackage("ESD.par"); gProof->EnablePackage("ESD");
 199               gProof->UploadPackage("ANALYSIS.par"); gProof->EnablePackage("ANALYSIS");
 200               mgr->StartAnalysis("proof",chain);
 201               break;
 202
 203     case 2:   mgr->StartAnalysis("grid" ,chain);
 204               break;
 205   }
 206   gBenchmark->Show("HMPID QA");
 207 }
 208 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 209