PWGPP/HMPID/AliHMPIDTaskQA.cxx

   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 //==============================================================================
  17 // AliHMPIDTaskQA - Class representing a quality check tool of HMPID
  18 // A set of histograms is created.
  19 //==============================================================================
  20
  21
  22 #ifndef AliHMPIDTASKQA_CXX
  23 #define AliHMPIDTASKQA_CXX
  24
  25
  26 #include "TH1.h"
  27 #include "TH2.h"
  28 #include "TFile.h"
  29 #include "TCanvas.h"
  30 #include "TGraphErrors.h"
  31 #include "AliPID.h"
  32 #include "AliVEvent.h"
  33 #include "AliVParticle.h"
  34 #include "AliVTrack.h"
  35 #include "AliESDtrackCuts.h"
  36 #include "AliAnalysisFilter.h"
  37 #include "AliAnalysisManager.h"
  38 #include "AliESDInputHandler.h"
  39 #include "AliMCEventHandler.h"
  40 #include "AliMCEvent.h"
  41 #include "AliESDtrack.h"
  42 #include "AliPID.h"
  43 #include "AliLog.h"
  44 #include "AliHMPIDTaskQA.h"
  45
  46 ClassImp(AliHMPIDTaskQA)
  47
  48 //__________________________________________________________________________
  49 AliHMPIDTaskQA::AliHMPIDTaskQA() :
  50   fESD(0x0),fMC(0x0),fUseMC(kTRUE),
  51   fHmpHistList(0x0),
  52   fHmpNevents(0x0),
  53   fZvertex(0x0),
  54   fTrackCuts(0x0),
  55   fTrackFilter(0x0),
  56   fTree(0x0)
  57 {
  58   //
  59   //Default ctor
  60   //
  61   for (Int_t i=0; i<23; i++) fVar[i]=0;
  62 }
  63
  64 //___________________________________________________________________________
  65 AliHMPIDTaskQA::AliHMPIDTaskQA(const Char_t* name) :
  66   AliAnalysisTaskSE(name),
  67   fESD(0x0), fMC(0x0), fUseMC(kTRUE),
  68   fHmpHistList(0x0),
  69   fHmpNevents(0x0),
  70   fZvertex(0x0),
  71   fTrackCuts(0x0),
  72   fTrackFilter(0x0),
  73   fTree(0x0)
  74 {
  75   //
  76   // Constructor. Initialization of Inputs and Outputs
  77   //
  78   for (Int_t i=0; i<23; i++) fVar[i]=0;
  79
  80   DefineOutput(1,TList::Class());
  81   DefineOutput(2,TTree::Class());
  82 }
  83
  84 //___________________________________________________________________________
  85 AliHMPIDTaskQA& AliHMPIDTaskQA::operator=(const AliHMPIDTaskQA& c)
  86 {
  87   //
  88   // Assignment operator
  89   //
  90   if (this!=&c) {
  91     AliAnalysisTaskSE::operator=(c);
  92     fESD             = c.fESD;
  93     fMC              = c.fMC;
  94     fUseMC           = c.fUseMC;
  95     fHmpHistList     = c.fHmpHistList;
  96     fHmpNevents      = c.fHmpNevents;
  97     fZvertex         = c.fZvertex;
  98     fTrackCuts       = c.fTrackCuts;
  99     fTrackFilter     = c.fTrackFilter;
 100     fTree            = c.fTree;
 101     for (Int_t i=0; i<23; i++) fVar[i]=c.fVar[i];
 102   }
 103   return *this;
 104 }
 105
 106 //___________________________________________________________________________
 107 AliHMPIDTaskQA::AliHMPIDTaskQA(const AliHMPIDTaskQA& c) :
 108   AliAnalysisTaskSE(c),
 109   fESD(c.fESD),fMC(c.fMC),fUseMC(c.fUseMC),
 110   fHmpHistList(c.fHmpHistList),
 111   fHmpNevents(c.fHmpNevents),
 112   fZvertex(c.fZvertex),
 113   fTrackCuts(c.fTrackCuts),
 114   fTrackFilter(c.fTrackFilter),
 115   fTree(c.fTree)
 116 {
 117   //
 118   // Copy Constructor
 119   //
 120   for (Int_t i=0; i<23; i++) fVar[i]=c.fVar[i];
 121 }
 122
 123 //___________________________________________________________________________
 124 AliHMPIDTaskQA::~AliHMPIDTaskQA() {
 125   //
 126   //destructor
 127   //
 128   Info("~AliHMPIDTaskQA","Calling Destructor");
 129   if (fHmpHistList && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fHmpHistList;
 130 }
 131
 132 //___________________________________________________________________________
 133 void AliHMPIDTaskQA::ConnectInputData(Option_t *)
 134 {
 135   // Connect ESD here
 136
 137   AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
 138   if (!esdH) {
 139     AliDebug(2,Form("ERROR: Could not get ESDInputHandler"));
 140   } else
 141     fESD = esdH->GetEvent();
 142
 143   if (fUseMC){
 144     // Connect MC
 145     AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
 146     if (!mcH) {
 147       AliDebug(2,Form("ERROR: Could not get MCEventHandler"));
 148       fUseMC = kFALSE;
 149     } else
 150       fMC = mcH->MCEvent();
 151       if (!fMC) AliDebug(2,Form("ERROR: Could not get MCEvent"));
 152   }
 153
 154   fTrackCuts = new AliESDtrackCuts("fTrackCuts", "Standard");
 155   fTrackCuts->SetAcceptKinkDaughters(kFALSE);
 156   fTrackCuts->SetMinNClustersTPC(70);  // changed default value 80 -> 70 (mtangaro)
 157   fTrackCuts->SetMaxChi2PerClusterTPC(4);
 158   fTrackCuts->SetMaxDCAToVertexXY(3);
 159   fTrackCuts->SetMaxDCAToVertexZ(3);
 160   fTrackCuts->SetRequireTPCRefit(kTRUE);
 161   fTrackCuts->SetRequireITSRefit(kTRUE);
 162   fTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
 163
 164   fTrackFilter = new AliAnalysisFilter("trackFilter");
 165   fTrackFilter->AddCuts(fTrackCuts);
 166
 167 }
 168 //***************************************************************************************************************************************************************************
 169 void AliHMPIDTaskQA::UserExec(Option_t *)
 170 {
 171
 172   fHmpNevents->Fill(0);
 173
 174   const AliESDVertex *vertex = fESD->GetPrimaryVertexTracks();
 175
 176   if(!vertex || !vertex->GetStatus() || vertex->GetNContributors()<1) {
 177
 178     // SPD vertex
 179     vertex = fESD->GetPrimaryVertexSPD();
 180     if(!vertex) return;
 181     if(!vertex->GetStatus()) return;
 182     if(vertex->GetNContributors()<1) return; // no good vertex, skip event
 183
 184   }
 185
 186   fHmpNevents->Fill(1);
 187
 188   Double_t vtxPos[3] = {999, 999, 999};
 189   if(vertex) vertex->GetXYZ(vtxPos);
 190   fZvertex->Fill(vtxPos[2]);
 191
 192   // vertex cut
 193   if (TMath::Abs(vtxPos[2]) > 10.) return;
 194
 195   AliESDtrack *track = 0;
 196
 197   Double_t ktol = 0.001;
 198
 199   //
 200   // Main loop function, executed on Event basis
 201   //
 202   for (Int_t iTrack = 0; iTrack<fESD->GetNumberOfTracks(); iTrack++) {
 203
 204     track = fESD->GetTrack(iTrack);
 205     if(!track) continue;
 206
 207     // pT min
 208     if ( track->Pt() < 0.2 ) continue;
 209
 210     // TPC clust
 211     if( track->GetNcls(1) < 70 ) continue;
 212
 213     // DCA cut
 214     Float_t b[2];
 215     Float_t bCov[3];
 216     track->GetImpactParameters(b,bCov);
 217     if( TMath::Sqrt(b[0]*b[0] + b[1]*b[1]) > 2 ) continue;
 218
 219     // ITS-TPC refit
 220     if( !track->IsOn(AliESDtrack::kITSrefit) ) continue;
 221     if( !track->IsOn(AliESDtrack::kTPCrefit) ) continue;
 222
 223     //HMPID cuts
 224
 225     if(Equal(track->GetHMPIDsignal(),-20.,ktol)) continue;
 226     if(track->GetHMPIDcluIdx() < 0) continue;
 227
 228     Int_t q, nph;
 229     Float_t x, y;
 230     Float_t xpc, ypc, th, ph;
 231     track->GetHMPIDmip(x,y,q,nph);
 232     track->GetHMPIDtrk(xpc,ypc,th,ph);
 233
 234     if(Equal(x,0.,ktol) && Equal(y,0.,ktol) && Equal(xpc,0.,ktol) && Equal(ypc,0.,ktol)) continue;
 235
 236     Double_t pHmp[3] = {0}, pHmp3 = 0;
 237     if (track->GetOuterHmpPxPyPz(pHmp)) pHmp3 = TMath::Sqrt(pHmp[0]*pHmp[0]+pHmp[1]*pHmp[1]+pHmp[2]*pHmp[2]);
 238
 239     fVar[0] = track->GetHMPIDcluIdx()/1000000;
 240     fVar[1] = pHmp3;
 241     fVar[2] = (Float_t)track->P();
 242     fVar[3] = xpc;
 243     fVar[4] = ypc;
 244     fVar[5] = x;
 245     fVar[6] = y;
 246     fVar[7] = (Float_t)track->GetHMPIDsignal();
 247     fVar[8] = q;
 248     fVar[9] = th;
 249     fVar[10] = ph;
 250     fVar[11] = (Float_t)track->GetSign();
 251     fVar[12] = (Float_t)nph;
 252     fVar[13] = track->GetHMPIDcluIdx()%1000000/1000; // cluster size
 253     fVar[14] = (Float_t)track->Eta();
 254     fVar[15] = (Float_t)track->Phi();
 255     fVar[16] = (Float_t)pHmp[0];
 256     fVar[17] = (Float_t)pHmp[1];
 257     fVar[18] = (Float_t)pHmp[2];
 258     fVar[19] = (Float_t)track->Px();
 259     fVar[20] = (Float_t)track->Py();
 260     fVar[21] = (Float_t)track->Pz();
 261     fVar[22] = (Float_t)track->GetHMPIDchi2();
 262
 263     fTree->Fill();
 264
 265   }//track loop
 266
 267   /* PostData(0) is taken care of by AliAnalysisTaskSE */
 268   PostData(1,fHmpHistList);
 269   PostData(2,fTree);
 270 }
 271 //___________________________________________________________________________
 272 void AliHMPIDTaskQA::Terminate(Option_t*)
 273 {
 274   // The Terminate() function is the last function to be called during
 275   // a query. It always runs on the client, it can be used to present
 276   // the results graphically or save the results to file.
 277
 278   Info("Terminate"," ");
 279
 280   if (!fUseMC) return;
 281
 282   fHmpHistList = dynamic_cast<TList*> (GetOutputData(1));
 283
 284   if (!fHmpHistList) {
 285     AliError("Histogram List is not available");
 286     return;
 287   }
 288
 289
 290   AliAnalysisTaskSE::Terminate();
 291
 292 }
 293 //___________________________________________________________________________
 294 void AliHMPIDTaskQA::UserCreateOutputObjects() {
 295   //
 296   //HERE ONE CAN CREATE OUTPUT OBJECTS
 297   //TO BE SET BEFORE THE EXECUTION OF THE TASK
 298   //
 299
 300   //slot #1
 301 //   OpenFile(1);
 302    fHmpHistList = new TList();
 303    fHmpHistList->SetOwner();
 304
 305    fHmpNevents = new TH1F("fHmpNevents","Number of events",2,0,2);
 306    fHmpHistList->Add(fHmpNevents);
 307
 308    fZvertex = new TH1F("fZvertex","Z primary vertex distribution",4000,-20,20);
 309    fHmpHistList->Add(fZvertex);
 310
 311 //   OpenFile(2);
 312    fTree = new TTree("Tree","Tree with data");
 313    fTree->Branch("Chamber",&fVar[0]);
 314    fTree->Branch("pHmp3",&fVar[1]);
 315    fTree->Branch("P",&fVar[2]);
 316    fTree->Branch("Xpc",&fVar[3]);
 317    fTree->Branch("Ypc",&fVar[4]);
 318    fTree->Branch("X",&fVar[5]);
 319    fTree->Branch("Y",&fVar[6]);
 320    fTree->Branch("HMPIDsignal",&fVar[7]);
 321    fTree->Branch("Charge",&fVar[8]);
 322    fTree->Branch("Theta",&fVar[9]);
 323    fTree->Branch("Phi",&fVar[10]);
 324    fTree->Branch("Sign",&fVar[11]);
 325    fTree->Branch("NumPhotons",&fVar[12]);
 326    fTree->Branch("ClustSize",&fVar[13]);
 327    fTree->Branch("Eta",&fVar[14]);
 328    fTree->Branch("PhiTrack",&fVar[15]);
 329    fTree->Branch("pHmpX",&fVar[16]);
 330    fTree->Branch("pHmpY",&fVar[17]);
 331    fTree->Branch("pHmpZ",&fVar[18]);
 332    fTree->Branch("Px",&fVar[19]);
 333    fTree->Branch("Py",&fVar[20]);
 334    fTree->Branch("Pz",&fVar[21]);
 335    fTree->Branch("HmpSigma",&fVar[22]);
 336
 337    PostData(1,fHmpHistList);
 338    PostData(2,fTree);
 339 }
 340
 341 //____________________________________________________________________________________________________________________________________
 342 Bool_t AliHMPIDTaskQA::Equal(Double_t x, Double_t y, Double_t tolerance)
 343 {
 344  return abs(x - y) <= tolerance ;
 345 }
 346
 347 #endif