PWG4/PartCorrBase/AliFiducialCut.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 /* $Id: AliFiducialCut.cxx 21839 2007-10-29 13:49:42Z gustavo $ */
  16
  17 //_________________________________________________________________________
  18 // Class for track/cluster acceptance selection
  19 // Selection in Central barrel, EMCAL and PHOS
  20 //
  21 // Several selection regions possible for the different
  22 // detectors
  23 //
  24 //*-- Author: Gustavo Conesa (LNF-INFN)
  25 //////////////////////////////////////////////////////////////////////////////
  26
  27
  28 // --- ROOT system ---
  29 #include <TMath.h>
  30 #include <TLorentzVector.h>
  31 #include <TString.h>
  32 //#include <TArrayF.h>
  33
  34 //---- ANALYSIS system ----
  35 #include "AliFiducialCut.h"
  36
  37 ClassImp(AliFiducialCut)
  38
  39
  40 //____________________________________________________________________________
  41 AliFiducialCut::AliFiducialCut() :
  42   TObject(),
  43   fEMCALFiducialCut(0),  fPHOSFiducialCut(0),  fCTSFiducialCut(0),
  44   fCTSFidCutMinEta(0x0),fCTSFidCutMinPhi(0x0),fCTSFidCutMaxEta(0x0), fCTSFidCutMaxPhi(0x0),
  45   fEMCALFidCutMinEta(0x0),fEMCALFidCutMinPhi(0x0),fEMCALFidCutMaxEta(0x0), fEMCALFidCutMaxPhi(0x0),
  46   fPHOSFidCutMinEta(0x0),fPHOSFidCutMinPhi(0x0),fPHOSFidCutMaxEta(0x0), fPHOSFidCutMaxPhi(0x0)
  47
  48 {
  49   //Ctor
  50
  51   //Initialize parameters
  52   InitParameters();
  53 }
  54 /*
  55 //____________________________________________________________________________
  56 AliFiducialCut::AliFiducialCut(const AliFiducialCut & g) :
  57   TObject(g),
  58   fEMCALFiducialCut(g.fEMCALFiducialCut),  fPHOSFiducialCut(g.fPHOSFiducialCut),  fCTSFiducialCut(g. fCTSFiducialCut),
  59   fCTSFidCutMinEta(g.fCTSFidCutMinEta?new TArrayF(*g.fCTSFidCutMinEta):0x0),
  60   fCTSFidCutMinPhi(g.fCTSFidCutMinPhi?new TArrayF(*g.fCTSFidCutMinPhi):0x0),
  61   fCTSFidCutMaxEta(g.fCTSFidCutMaxEta?new TArrayF(*g.fCTSFidCutMaxEta):0x0),
  62   fCTSFidCutMaxPhi(g.fCTSFidCutMaxPhi?new TArrayF(*g.fCTSFidCutMaxPhi):0x0),
  63   fEMCALFidCutMinEta(g.fEMCALFidCutMinEta?new TArrayF(*g.fEMCALFidCutMinEta):0x0),
  64   fEMCALFidCutMinPhi(g.fEMCALFidCutMinPhi?new TArrayF(*g.fEMCALFidCutMinPhi):0x0),
  65   fEMCALFidCutMaxEta(g.fEMCALFidCutMaxEta?new TArrayF(*g.fEMCALFidCutMaxEta):0x0),
  66   fEMCALFidCutMaxPhi(g.fEMCALFidCutMaxPhi?new TArrayF(*g.fEMCALFidCutMaxPhi):0x0),
  67   fPHOSFidCutMinEta(g.fPHOSFidCutMinEta?new TArrayF(*g.fPHOSFidCutMinEta):0x0),
  68   fPHOSFidCutMinPhi(g.fPHOSFidCutMinPhi?new TArrayF(*g.fPHOSFidCutMinPhi):0x0),
  69   fPHOSFidCutMaxEta(g.fPHOSFidCutMaxEta?new TArrayF(*g.fPHOSFidCutMaxEta):0x0),
  70   fPHOSFidCutMaxPhi(g.fPHOSFidCutMaxPhi?new TArrayF(*g.fPHOSFidCutMaxPhi):0x0)
  71
  72 {
  73   // cpy ctor
  74
  75 }
  76
  77 //_________________________________________________________________________
  78 AliFiducialCut & AliFiducialCut::operator = (const AliFiducialCut & source)
  79 {
  80   // assignment operator
  81
  82   if(&source == this) return *this;
  83
  84   fEMCALFiducialCut = source.fEMCALFiducialCut;
  85   fPHOSFiducialCut = source.fPHOSFiducialCut;
  86   fCTSFiducialCut = source.fCTSFiducialCut;
  87
  88   fCTSFidCutMinEta = source.fCTSFidCutMinEta?new TArrayF(*source.fCTSFidCutMinEta):0x0;
  89   fCTSFidCutMinPhi = source.fCTSFidCutMinPhi?new TArrayF(*source.fCTSFidCutMinPhi):0x0;
  90   fCTSFidCutMaxEta = source.fCTSFidCutMaxEta?new TArrayF(*source.fCTSFidCutMaxEta):0x0;
  91   fCTSFidCutMaxPhi = source.fCTSFidCutMaxPhi?new TArrayF(*source.fCTSFidCutMaxPhi):0x0;
  92   fEMCALFidCutMinEta = source.fEMCALFidCutMinEta?new TArrayF(*source.fEMCALFidCutMinEta):0x0;
  93   fEMCALFidCutMinPhi = source.fEMCALFidCutMinPhi?new TArrayF(*source.fEMCALFidCutMinPhi):0x0;
  94   fEMCALFidCutMaxEta = source.fEMCALFidCutMaxEta?new TArrayF(*source.fEMCALFidCutMaxEta):0x0;
  95   fEMCALFidCutMaxPhi = source.fEMCALFidCutMaxPhi?new TArrayF(*source.fEMCALFidCutMaxPhi):0x0;
  96   fPHOSFidCutMinEta = source.fPHOSFidCutMinEta?new TArrayF(*source.fPHOSFidCutMinEta):0x0;
  97   fPHOSFidCutMinPhi = source.fPHOSFidCutMinPhi?new TArrayF(*source.fPHOSFidCutMinPhi):0x0;
  98   fPHOSFidCutMaxEta = source.fPHOSFidCutMaxEta?new TArrayF(*source.fPHOSFidCutMaxEta):0x0;
  99   fPHOSFidCutMaxPhi = source.fPHOSFidCutMaxPhi?new TArrayF(*source.fPHOSFidCutMaxPhi):0x0;
 100
 101   return *this;
 102
 103 }
 104 */
 105 //_________________________________
 106 AliFiducialCut::~AliFiducialCut() {
 107   //Dtor
 108
 109   delete fCTSFidCutMinEta ;
 110   delete fCTSFidCutMinPhi ;
 111   delete fCTSFidCutMaxEta ;
 112   delete fCTSFidCutMaxPhi ;
 113
 114   delete fEMCALFidCutMinEta ;
 115   delete fEMCALFidCutMinPhi ;
 116   delete fEMCALFidCutMaxEta ;
 117   delete fEMCALFidCutMaxPhi ;
 118
 119   delete fPHOSFidCutMinEta ;
 120   delete fPHOSFidCutMinPhi ;
 121   delete fPHOSFidCutMaxEta ;
 122   delete fPHOSFidCutMaxPhi ;
 123
 124 }
 125
 126
 127 //_______________________________________________________________
 128 Bool_t AliFiducialCut::IsInFiducialCut(const TLorentzVector momentum, const TString det) const
 129 {
 130   //Selects EMCAL or PHOS cluster or CTS track if it is inside eta-phi defined regions
 131
 132   if(det == "CTS"){
 133           if(!fCTSFiducialCut)  return kTRUE; //Fiducial cut not requested, accept all tracks
 134           else return CheckFiducialRegion(momentum, fCTSFidCutMinPhi  , fCTSFidCutMaxPhi , fCTSFidCutMinEta  , fCTSFidCutMaxEta  );
 135   }
 136   else   if(det == "EMCAL") {
 137           if(!fEMCALFiducialCut) return kTRUE; //Fiducial cut not requested, accept all clusters
 138           else return CheckFiducialRegion(momentum, fEMCALFidCutMinPhi, fEMCALFidCutMaxPhi, fEMCALFidCutMinEta, fEMCALFidCutMaxEta);
 139   }
 140   else   if(det == "PHOS")  {
 141           if(!fPHOSFiducialCut) return kTRUE; //Fiducial cut not requested, accept all clusters
 142           else return CheckFiducialRegion(momentum, fPHOSFidCutMinPhi , fPHOSFidCutMaxPhi , fPHOSFidCutMinEta , fPHOSFidCutMaxEta );
 143   }
 144   else{
 145           printf("AliFiducialCut::IsInFiducialCut() - Wrong detector name = %s\n", det.Data());
 146           return kFALSE;
 147   }
 148
 149 }
 150
 151 //_______________________________________________________________
 152 Bool_t AliFiducialCut::CheckFiducialRegion(const TLorentzVector momentum, const TArrayF* minphi, const TArrayF* maxphi, const TArrayF* mineta, const TArrayF* maxeta) const {
 153   //Given the selection regions in Eta and Phi, check if particle is in this region.
 154
 155     Double_t phi = momentum.Phi();
 156         if(phi < 0) phi+=TMath::TwoPi() ;
 157         Double_t eta =  momentum.Eta();
 158         //printf("IsInFiducialCut::Det: %s, phi = %f, eta = %f\n", det.Data(),phi*TMath::RadToDeg(), eta);
 159
 160     Int_t netaregions = maxeta->GetSize();
 161     Int_t nphiregions = maxphi->GetSize();
 162     if(netaregions !=  mineta->GetSize() || nphiregions !=  minphi->GetSize())
 163                 printf("AliFiducialCut::IsInFiducialCut() - Wrong number of fiducial cut regions: nmaxeta %d != nmineta %d; nmaxphi %d != nminphi %d\n",
 164                            netaregions, mineta->GetSize(),  nphiregions, minphi->GetSize());
 165
 166         //Eta fiducial cut
 167         Bool_t bInEtaFidCut = kFALSE;
 168         for(Int_t ieta = 0; ieta < netaregions; ieta++)
 169                 if(eta > mineta->GetAt(ieta) && eta < maxeta->GetAt(ieta)) bInEtaFidCut = kTRUE;
 170
 171         if(bInEtaFidCut){
 172                 //printf("Eta cut passed\n");
 173                 //Phi fiducial cut
 174                 Bool_t bInPhiFidCut = kFALSE;
 175                 for(Int_t iphi = 0; iphi < nphiregions; iphi++)
 176                         if(phi > minphi->GetAt(iphi) *TMath::DegToRad()&& phi < maxphi->GetAt(iphi)*TMath::DegToRad()) bInPhiFidCut = kTRUE ;
 177
 178                 if(bInPhiFidCut) {
 179                         //printf("IsInFiducialCut:: %s cluster/track accepted\n",det.Data());
 180                         return kTRUE;
 181                 }
 182                 else return kFALSE;
 183
 184         }//In eta fid cut
 185         else
 186     return kFALSE;
 187 }
 188
 189
 190 //_______________________________________________________________
 191 void AliFiducialCut::InitParameters()
 192 {
 193
 194   //Initialize the parameters of the analysis.
 195
 196   fEMCALFiducialCut = kTRUE ;
 197   fPHOSFiducialCut = kTRUE ;
 198   fCTSFiducialCut = kTRUE ;
 199
 200   fCTSFidCutMinEta = new TArrayF(1);
 201   fCTSFidCutMinEta->SetAt(-0.9,0);
 202   fCTSFidCutMaxEta = new TArrayF(1);
 203   fCTSFidCutMaxEta->SetAt(0.9,0);
 204
 205   fCTSFidCutMinPhi = new TArrayF(1);
 206   fCTSFidCutMinPhi->SetAt(0,0);
 207   fCTSFidCutMaxPhi = new TArrayF(1);
 208   fCTSFidCutMaxPhi->SetAt(360.,0);
 209
 210   fEMCALFidCutMinEta = new TArrayF(1);
 211   fEMCALFidCutMinEta->SetAt(-0.7,0);
 212   fEMCALFidCutMaxEta = new TArrayF(1);
 213   fEMCALFidCutMaxEta->SetAt(0.7,0);
 214
 215   fEMCALFidCutMinPhi = new TArrayF(1);
 216   fEMCALFidCutMinPhi->SetAt(80.,0);
 217   fEMCALFidCutMaxPhi = new TArrayF(1);
 218   fEMCALFidCutMaxPhi->SetAt(190.,0);
 219
 220   fPHOSFidCutMinEta = new TArrayF(1);
 221   fPHOSFidCutMinEta->SetAt(-0.13,0);
 222   fPHOSFidCutMaxEta = new TArrayF(1);
 223   fPHOSFidCutMaxEta->SetAt(0.13,0);
 224
 225   fPHOSFidCutMinPhi = new TArrayF(1);
 226   fPHOSFidCutMinPhi->SetAt(220.,0);
 227   fPHOSFidCutMaxPhi = new TArrayF(1);
 228   fPHOSFidCutMaxPhi->SetAt(320.,0);
 229
 230 }
 231
 232
 233 //________________________________________________________________
 234 void AliFiducialCut::Print(const Option_t * opt) const
 235 {
 236
 237   //Print some relevant parameters set for the analysis
 238   if(! opt)
 239     return;
 240
 241   printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
 242
 243   if(fCTSFiducialCut){
 244     Int_t netaregions =  fCTSFidCutMaxEta->GetSize();
 245     Int_t nphiregions =  fCTSFidCutMaxPhi->GetSize();
 246     printf(">> CTS Fiducial regions : phi %d eta %d\n", netaregions, nphiregions) ;
 247     for(Int_t ieta = 0; ieta < netaregions; ieta++)
 248       printf(" region %d : %3.2f < eta < %3.2f\n", ieta, fCTSFidCutMinEta->GetAt(ieta), fCTSFidCutMaxEta->GetAt(ieta)) ;
 249     for(Int_t iphi = 0; iphi < nphiregions; iphi++)
 250       printf(" region %d : %3.1f < phi < %3.1f\n", iphi, fCTSFidCutMinPhi->GetAt(iphi), fCTSFidCutMaxPhi->GetAt(iphi)) ;
 251   }
 252   else printf(">>No fiducial cuts in CTS\n");
 253
 254   if(fEMCALFiducialCut){
 255     Int_t netaregions =  fEMCALFidCutMaxEta->GetSize();
 256     Int_t nphiregions =  fEMCALFidCutMaxPhi->GetSize();
 257     printf(">>EMCAL Fiducial regions : phi %d eta %d\n", netaregions, nphiregions) ;
 258     for(Int_t ieta = 0; ieta < netaregions; ieta++)
 259       printf(" region %d : %3.2f < eta < %3.2f\n", ieta, fEMCALFidCutMinEta->GetAt(ieta), fEMCALFidCutMaxEta->GetAt(ieta)) ;
 260     for(Int_t iphi = 0; iphi < nphiregions; iphi++)
 261       printf(" region %d : %3.1f < phi < %3.1f\n", iphi, fEMCALFidCutMinPhi->GetAt(iphi), fEMCALFidCutMaxPhi->GetAt(iphi)) ;
 262   }
 263   else printf(">>No fiducial cuts in EMCAL\n");
 264
 265   if(fPHOSFiducialCut){
 266     Int_t netaregions =  fPHOSFidCutMaxEta->GetSize();
 267     Int_t nphiregions =  fPHOSFidCutMaxPhi->GetSize();
 268     printf(">>PHOS Fiducial regions : phi %d eta %d\n", netaregions, nphiregions) ;
 269     for(Int_t ieta = 0; ieta < netaregions; ieta++)
 270       printf(" region %d : %3.2f < eta < %3.2f\n", ieta, fPHOSFidCutMinEta->GetAt(ieta), fPHOSFidCutMaxEta->GetAt(ieta)) ;
 271     for(Int_t iphi = 0; iphi < nphiregions; iphi++)
 272       printf(" region %d : %3.1f < phi < %3.1f\n", iphi, fPHOSFidCutMinPhi->GetAt(iphi), fPHOSFidCutMaxPhi->GetAt(iphi)) ;
 273   }
 274   else printf(">>No fiducial cuts in PHOS\n");
 275   printf("    \n") ;
 276 }
 277
 278
 279 //_______________________________________________________________
 280 void AliFiducialCut::SetSimpleCTSFiducialCut(const Float_t eta, const Float_t minphi, const Float_t maxphi){
 281
 282   //Method to set simple acceptance cut to CTS
 283   fCTSFidCutMinEta->Set(1);
 284   fCTSFidCutMaxEta->Set(1);
 285   fCTSFidCutMinPhi->Set(1);
 286   fCTSFidCutMaxPhi->Set(1);
 287
 288   fCTSFidCutMinEta->SetAt(-eta,0);
 289   fCTSFidCutMaxEta->SetAt(eta,0);
 290   fCTSFidCutMinPhi->SetAt(minphi,0);
 291   fCTSFidCutMaxPhi->SetAt(maxphi,0);
 292
 293
 294 }
 295
 296
 297 //_______________________________________________________________
 298 void AliFiducialCut::SetSimpleEMCALFiducialCut(const Float_t eta, const Float_t minphi, const Float_t maxphi){
 299   //Method to set simple acceptance cut to EMCAL
 300
 301   fEMCALFidCutMinEta->Set(1);
 302   fEMCALFidCutMaxEta->Set(1);
 303   fEMCALFidCutMinPhi->Set(1);
 304   fEMCALFidCutMaxPhi->Set(1);
 305
 306   fEMCALFidCutMinEta->SetAt(-eta,0);
 307   fEMCALFidCutMaxEta->SetAt(eta,0);
 308   fEMCALFidCutMinPhi->SetAt(minphi,0);
 309   fEMCALFidCutMaxPhi->SetAt(maxphi,0);
 310
 311
 312 }
 313
 314 //_______________________________________________________________
 315 void AliFiducialCut::SetSimplePHOSFiducialCut(const Float_t eta, const Float_t minphi, const Float_t maxphi){
 316
 317   //Method to set simple acceptance cut to PHOS
 318   fPHOSFidCutMinEta->Set(1);
 319   fPHOSFidCutMaxEta->Set(1);
 320   fPHOSFidCutMinPhi->Set(1);
 321   fPHOSFidCutMaxPhi->Set(1);
 322
 323   fPHOSFidCutMinEta->SetAt(-eta,0);
 324   fPHOSFidCutMaxEta->SetAt(eta,0);
 325   fPHOSFidCutMinPhi->SetAt(minphi,0);
 326   fPHOSFidCutMaxPhi->SetAt(maxphi,0);
 327
 328
 329 }