CORRFW/AliCFEventRecCuts.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 // Cut on the Event at reconstructed level: for the moment
  16 // the requirements on the number of charged tracks and on
  17 // the vertex position and resolution are implemented
  18 // The argument of IsSelected member function (passed object) is cast into
  19 // an AliESDEvent. In the future may be modified to use AliVEvent interface
  20 // and include more cut variables.
  21 // The class derives from AliCFCutBase
  22 // Author:S.Arcelli Silvia.Arcelli@cern.ch
  23 //
  24 //
  25 #include "TH1F.h"
  26 #include "TList.h"
  27 #include "TBits.h"
  28 #include "AliLog.h"
  29 #include "AliESDEvent.h"
  30 #include "AliESDVertex.h"
  31 #include "AliCFEventRecCuts.h"
  32 ClassImp(AliCFEventRecCuts)
  33 //____________________________________________________________________
  34 AliCFEventRecCuts::AliCFEventRecCuts() :
  35   AliCFCutBase(),
  36   fNTracksMin(-1),
  37   fNTracksMax(1000000),
  38   fRequireVtxCuts(kFALSE),
  39   fVtxXMax(1.e99),
  40   fVtxYMax(1.e99),
  41   fVtxZMax(1.e99),
  42   fVtxXMin(-1.e99),
  43   fVtxYMin(-1.e99),
  44   fVtxZMin(-1.e99),
  45   fVtxXResMax(1.e99),
  46   fVtxYResMax(1.e99),
  47   fVtxZResMax(1.e99),
  48   fVtxNCtrbMin(0),
  49   fVtxNCtrbMax((Int_t)1.e9),
  50   fVtxTPC(0),
  51   fVtxSPD(0),
  52   fBitMap(0x0)
  53 {
  54   //
  55   //ctor
  56   //
  57   fBitMap=new TBits(0);
  58   Initialise();
  59 }
  60
  61 //____________________________________________________________________
  62 AliCFEventRecCuts::AliCFEventRecCuts(Char_t* name, Char_t* title) :
  63   AliCFCutBase(name,title),
  64   fNTracksMin(-1),
  65   fNTracksMax(1000000),
  66   fRequireVtxCuts(kFALSE),
  67   fVtxXMax(1.e99),
  68   fVtxYMax(1.e99),
  69   fVtxZMax(1.e99),
  70   fVtxXMin(-1.e99),
  71   fVtxYMin(-1.e99),
  72   fVtxZMin(-1.e99),
  73   fVtxXResMax(1.e99),
  74   fVtxYResMax(1.e99),
  75   fVtxZResMax(1.e99),
  76   fVtxNCtrbMin(0),
  77   fVtxNCtrbMax((Int_t)1.e9),
  78   fVtxTPC(0),
  79   fVtxSPD(0),
  80   fBitMap(0x0)
  81  {
  82   //
  83   //ctor
  84   //
  85   fBitMap=new TBits(0);
  86   Initialise();
  87  }
  88
  89 //____________________________________________________________________
  90 AliCFEventRecCuts::AliCFEventRecCuts(const AliCFEventRecCuts& c) :
  91   AliCFCutBase(c),
  92   fNTracksMin(c.fNTracksMin),
  93   fNTracksMax(c.fNTracksMax),
  94   fRequireVtxCuts(c.fRequireVtxCuts),
  95   fVtxXMax(c.fVtxXMax),
  96   fVtxYMax(c.fVtxYMax),
  97   fVtxZMax(c.fVtxZMax),
  98   fVtxXMin(c.fVtxXMin),
  99   fVtxYMin(c.fVtxYMin),
 100   fVtxZMin(c.fVtxZMin),
 101   fVtxXResMax(c.fVtxXResMax),
 102   fVtxYResMax(c.fVtxYResMax),
 103   fVtxZResMax(c.fVtxZResMax),
 104   fVtxNCtrbMin(c.fVtxNCtrbMin),
 105   fVtxNCtrbMax(c.fVtxNCtrbMax),
 106   fVtxTPC(c.fVtxTPC),
 107   fVtxSPD(c.fVtxSPD),
 108   fBitMap(c.fBitMap)
 109 {
 110   //
 111   //copy constructor
 112   //
 113   for (Int_t i=0; i<c.kNCuts; i++){
 114     for (Int_t j=0; j<c.kNStepQA; j++){
 115       if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
 116     }
 117   }
 118
 119 }
 120
 121 //____________________________________________________________________
 122 AliCFEventRecCuts::~AliCFEventRecCuts() {
 123   //
 124   //dtor
 125   //
 126
 127   for (Int_t i=0; i<kNCuts; i++){
 128     for (Int_t j=0; j<kNStepQA; j++){
 129       if(fhQA[i][j]) delete fhQA[i][j];
 130     }
 131   }
 132
 133   if(fBitMap)delete fBitMap;
 134
 135 }
 136 //_____________________________________________________________________________
 137 void AliCFEventRecCuts::Initialise()
 138 {
 139
 140   //
 141   //initialization
 142   //
 143
 144   //
 145   // sets pointers to histos to zero
 146   //
 147
 148   for(Int_t i=0; i<kNCuts; i++){
 149     for(Int_t j =0; j<kNStepQA; j++){
 150       fhQA[i][j]=0x0;
 151     }
 152   }
 153 }
 154
 155 //____________________________________________________________________
 156 AliCFEventRecCuts& AliCFEventRecCuts::operator=(const AliCFEventRecCuts& c)
 157 {
 158   //
 159   // Assignment operator
 160   //
 161   if (this != &c) {
 162     AliCFCutBase::operator=(c) ;
 163     fNTracksMin=c.fNTracksMin;
 164     fNTracksMax=c.fNTracksMax;
 165     fRequireVtxCuts=c.fRequireVtxCuts;
 166     fVtxXMax=c.fVtxXMax;
 167     fVtxYMax=c.fVtxYMax;
 168     fVtxZMax=c.fVtxZMax;
 169     fVtxXMin=c.fVtxXMin;
 170     fVtxYMin=c.fVtxYMin;
 171     fVtxZMin=c.fVtxZMin;
 172     fVtxXResMax=c.fVtxXResMax;
 173     fVtxYResMax=c.fVtxYResMax;
 174     fVtxZResMax=c.fVtxZResMax;
 175     fVtxNCtrbMin=c.fVtxNCtrbMin;
 176     fVtxNCtrbMax=c.fVtxNCtrbMax;
 177     fVtxTPC=c.fVtxTPC;
 178     fVtxSPD=c.fVtxSPD;
 179     fBitMap=c.fBitMap;
 180   }
 181
 182   for (Int_t i=0; i<c.kNCuts; i++){
 183     for (Int_t j=0; j<c.kNStepQA; j++){
 184       if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
 185     }
 186   }
 187
 188
 189   return *this ;
 190 }
 191
 192 //____________________________________________________________________
 193 Bool_t AliCFEventRecCuts::IsSelected(TObject* obj) {
 194   //
 195   //Check if the requested cuts are passed
 196   //
 197
 198
 199   SelectionBitMap(obj);
 200
 201   if (fIsQAOn) FillHistograms(obj,0);
 202   Bool_t isSelected = kTRUE;
 203
 204   for (UInt_t icut=0; icut<fBitMap->GetNbits();icut++)
 205         if(!fBitMap->TestBitNumber(icut)) isSelected = kFALSE;
 206
 207   if (!isSelected) return kFALSE ;
 208   if (fIsQAOn) FillHistograms(obj,1);
 209   return kTRUE;
 210
 211 }
 212 //____________________________________________________________________
 213 void AliCFEventRecCuts::SelectionBitMap(TObject* obj) {
 214   //
 215   //cut on the number of charged tracks and on the event vertex.
 216   //so far specific to AliESDEvents
 217   //
 218
 219   //Check if the requested cuts are passed and return a bitmap
 220   for(Int_t j=0;j<kNCuts;j++)fBitMap->SetBitNumber(j,kFALSE);
 221   AliESDEvent* esd = dynamic_cast<AliESDEvent *>(obj);
 222   if ( !esd ) return;
 223
 224   //now start checking the cuts,
 225   //first assume the event will be accepted:
 226   for(Int_t j=0;j<kNCuts;j++)fBitMap->SetBitNumber(j,kTRUE);
 227
 228   //Number of charged tracks:
 229   Int_t nTracks = esd->GetNumberOfTracks();
 230   if(nTracks<fNTracksMin || nTracks>fNTracksMax)
 231     fBitMap->SetBitNumber(0,kFALSE);
 232
 233   if(fRequireVtxCuts){
 234     const AliESDVertex* vtxESD = 0x0;
 235     if      (fVtxTPC) vtxESD = esd->GetPrimaryVertexTPC() ;
 236     else if (fVtxSPD) vtxESD = esd->GetPrimaryVertexSPD() ;
 237     else              vtxESD = esd->GetPrimaryVertexTracks() ;
 238     if(!vtxESD){
 239       for(Int_t j=1;j<kNCuts;j++)fBitMap->SetBitNumber(j,kFALSE);
 240       AliWarning("Cannot get vertex, skipping event");
 241       return;
 242     }
 243     // Require the vertex to have been reconstructed successfully
 244     if (strcmp(vtxESD->GetName(), "default")==0){
 245       AliWarning(Form(" No reconstructed vertex found, skipping event"));
 246       for(Int_t j=1;j<kNCuts;j++)fBitMap->SetBitNumber(j,kFALSE);
 247       return;
 248     }
 249     // Pick up the position and uncertainties
 250
 251     Double_t vtxPos[3];
 252     vtxPos[0] = vtxESD->GetXv();
 253     vtxPos[1] = vtxESD->GetYv();
 254     vtxPos[2] = vtxESD->GetZv();
 255
 256     Double_t vtxRes[3];
 257     vtxRes[0] = vtxESD->GetXRes();
 258     vtxRes[1] = vtxESD->GetYRes();
 259     vtxRes[2] = vtxESD->GetZRes();
 260
 261     Int_t nCtrb = vtxESD->GetNContributors();
 262
 263     // Apply the cut
 264
 265     if (vtxPos[0]>fVtxXMax || vtxPos[0]<fVtxXMin)
 266       fBitMap->SetBitNumber(1,kFALSE);
 267     if (vtxPos[1]>fVtxYMax || vtxPos[1]<fVtxYMin)
 268       fBitMap->SetBitNumber(2,kFALSE);
 269     if (vtxPos[2]>fVtxZMax || vtxPos[2]<fVtxZMin)
 270       fBitMap->SetBitNumber(3,kFALSE);
 271     if (vtxRes[0]==0 || vtxRes[0]>fVtxXResMax)
 272       fBitMap->SetBitNumber(4,kFALSE);
 273     if (vtxRes[1]==0 || vtxRes[1]>fVtxYResMax)
 274       fBitMap->SetBitNumber(5,kFALSE);
 275     if (vtxRes[2]==0 || vtxRes[2]>fVtxZResMax)
 276       fBitMap->SetBitNumber(6,kFALSE);
 277     if (nCtrb<fVtxNCtrbMin || nCtrb>fVtxNCtrbMax)
 278       fBitMap->SetBitNumber(7,kFALSE);
 279
 280   }
 281   return;
 282 }
 283
 284 //_____________________________________________________________________________
 285 void AliCFEventRecCuts::FillHistograms(TObject* obj, Bool_t b)
 286 {
 287   //
 288   // fill the QA histograms
 289   //
 290
 291   if(!fIsQAOn) return;
 292   // cast TObject into VParticle
 293   AliESDEvent* esd = dynamic_cast<AliESDEvent *>(obj);
 294   if (!esd ) return  ;
 295
 296   // index = 0: fill histograms before cuts
 297   // index = 1: fill histograms after cuts
 298   Int_t index = -1;
 299   index = ((b) ? 1 : 0);
 300
 301
 302   //number of charged tracks:
 303   Int_t nTracks = esd->GetNumberOfTracks();
 304   fhQA[kNTracks][index]->Fill(nTracks);
 305
 306   //look at vertex parameters:
 307   const AliESDVertex* vtxESD = 0x0;
 308   if      (fVtxTPC) vtxESD = esd->GetPrimaryVertexTPC() ;
 309   else if (fVtxSPD) vtxESD = esd->GetPrimaryVertexSPD() ;
 310   else              vtxESD = esd->GetPrimaryVertexTracks() ;
 311   if(!vtxESD)return;
 312   // Require the vertex to have been reconstructed successfully
 313   if (strcmp(vtxESD->GetName(), "default")==0)return;
 314   // vertex position and uncertainties
 315   fhQA[kVtxPosX] [index]->Fill(vtxESD->GetXv());
 316   fhQA[kVtxPosY] [index]->Fill(vtxESD->GetYv());
 317   fhQA[kVtxPosZ] [index]->Fill(vtxESD->GetZv());
 318   fhQA[kVtxResX] [index]->Fill(vtxESD->GetXRes());
 319   fhQA[kVtxResY] [index]->Fill(vtxESD->GetYRes());
 320   fhQA[kVtxResZ] [index]->Fill(vtxESD->GetZRes());
 321   fhQA[kVtxNCtrb][index]->Fill(vtxESD->GetNContributors());
 322
 323 }
 324
 325 //____________________________________________________________________
 326 void AliCFEventRecCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t *bins)
 327 {
 328   //
 329   //setting x-axis bin limits of QA histogram fhQA[index]
 330   //
 331
 332   for(Int_t i=0;i<kNStepQA;i++){
 333     if(!fhQA[index][i]){AliWarning("non-existing histogram!");
 334     return;
 335     }
 336     fhQA[index][i]->GetXaxis()->Set(nbins,bins);
 337   }
 338 }
 339 //____________________________________________________________________
 340 void AliCFEventRecCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t xmin, Double_t xmax)
 341 {
 342   //
 343   //setting x-axis bins and range of QA histogram fhQA[index]
 344   //
 345
 346   for(Int_t i=0;i<kNStepQA;i++){
 347     if(!fhQA[index][i]){AliWarning("non-existing histogram!");
 348     return;
 349     }
 350     fhQA[index][i]->GetXaxis()->Set(nbins,xmin,xmax);
 351   }
 352 }
 353
 354 //_____________________________________________________________________________
 355  void AliCFEventRecCuts::DefineHistograms() {
 356   //
 357   // histograms for cut variables
 358   //
 359   Int_t color = 2;
 360
 361   if(!fIsQAOn) {
 362     AliInfo(Form("No QA histos requested, Please first set the QA flag on!"));
 363     return;
 364   }
 365
 366   // book QA histograms
 367
 368   Char_t str[5];
 369   for (Int_t i=0; i<kNStepQA; i++) {
 370     if (i==0) snprintf(str,5," ");
 371     else snprintf(str,5,"_cut");
 372
 373     fhQA[kNTracks][i]   = new  TH1F(Form("%s_NTracks%s",GetName(),str),                 "",501,-0.5,500.5);
 374     fhQA[kVtxPosX][i]   = new  TH1F(Form("%s_Vtx_Pos_X%s",GetName(),str),               "",100,-5.,5.);
 375     fhQA[kVtxPosY][i]   = new  TH1F(Form("%s_Vtx_Pos_Y%s",GetName(),str),               "",100,-5.,5.);
 376     fhQA[kVtxPosZ][i]   = new  TH1F(Form("%s_Vtx_Pos_Z%s",GetName(),str),               "",200,-50.,50.);
 377     fhQA[kVtxResX][i]   = new  TH1F(Form("%s_Vtx_Res_X%s",GetName(),str),               "",100,-1.,1.);
 378     fhQA[kVtxResY][i]   = new  TH1F(Form("%s_Vtx_Res_Y%s",GetName(),str),               "",100,-1.,1.);
 379     fhQA[kVtxResZ][i]   = new  TH1F(Form("%s_Vtx_Res_Z%s",GetName(),str),               "",100,-1.,1.);
 380     fhQA[kVtxNCtrb][i]  = new  TH1F(Form("%s_Vtx_N_Ctrb%s",GetName(),str),              "",1000,0.,1000.);
 381
 382     fhQA[kNTracks][i]   ->SetXTitle("Number of ESD tracks");
 383     fhQA[kVtxPosX][i]   ->SetXTitle("Vertex Position X (cm)");
 384     fhQA[kVtxPosY][i]   ->SetXTitle("Vertex Position Y (cm)");
 385     fhQA[kVtxPosZ][i]   ->SetXTitle("Vertex Position Z (cm)");
 386     fhQA[kVtxResX][i]   ->SetXTitle("Vertex Resolution X (cm)");
 387     fhQA[kVtxResY][i]   ->SetXTitle("Vertex Resolution Y (cm)");
 388     fhQA[kVtxResZ][i]   ->SetXTitle("Vertex Resolution Z (cm)");
 389     fhQA[kVtxNCtrb][i]  ->SetXTitle("Number of contributors");
 390   }
 391
 392   for(Int_t i=0; i<kNCuts; i++) fhQA[i][1]->SetLineColor(color);
 393
 394 }
 395
 396 //_____________________________________________________________________________
 397 void AliCFEventRecCuts::AddQAHistograms(TList *qaList) {
 398   //
 399   // saves the histograms in a TList
 400   //
 401
 402   DefineHistograms();
 403
 404   for (Int_t j=0; j<kNStepQA; j++) {
 405     for(Int_t i=0; i<kNCuts; i++)
 406         qaList->Add(fhQA[i][j]);
 407   }
 408 }