PWGGA/CaloTrackCorrelations/AliAnaRandomTrigger.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 // Gerenate a random trigger, input for other analysis
  18 // Set flat energy distribution over acceptance of EMCAL, PHOS or CTS
  19 // Be careful, correlate only with Min Bias events this trigger
  20 //
  21 //
  22 //-- Author: Gustavo Conesa (LPSC-Grenoble)
  23 //_________________________________________________________________________
  24
  25
  26 // --- ROOT system ---
  27 #include <TH2F.h>
  28 #include <TClonesArray.h>
  29
  30 //---- AliRoot system ----
  31 #include "AliAnaRandomTrigger.h"
  32 #include "AliAODPWG4ParticleCorrelation.h"
  33 #include "AliEMCALGeometry.h"
  34
  35 ClassImp(AliAnaRandomTrigger)
  36
  37 //__________________________________________
  38 AliAnaRandomTrigger::AliAnaRandomTrigger() :
  39     AliAnaCaloTrackCorrBaseClass(),
  40     fTriggerDetector(kEMCAL),
  41     fTriggerDetectorString("EMCAL"),
  42     fRandom(0),         fNRandom(0),
  43     fMomentum(),
  44     fhE(0),             fhPt(0),
  45     fhPhi(0),           fhEta(0),
  46     fhEtaPhi(0)
  47 {
  48   //Default Ctor
  49
  50   //Initialize parameters
  51   InitParameters();
  52
  53 }
  54
  55 //_________________________________________________________________________
  56 Bool_t AliAnaRandomTrigger::ExcludeDeadBadRegions(Float_t eta, Float_t phi)
  57 {
  58   // Check if there is a dead or bad region in a detector
  59   // Now only EMCAL
  60
  61   if(fTriggerDetector!=kEMCAL) return kFALSE;
  62
  63   //-------------------------------------
  64   // Get the corresponding cell in EMCAL, check if it exists in acceptance (phi gaps, borders)
  65   //-------------------------------------
  66
  67   Int_t absId = -1;
  68   if(!GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(eta,phi, absId)) return kTRUE; // remove if out of EMCAL acceptance, phi gaps
  69
  70   Int_t icol = -1, irow = -1, iRCU = -1;
  71   Int_t sm = GetCaloUtils()->GetModuleNumberCellIndexes(absId,kEMCAL, icol, irow, iRCU);
  72
  73   //printf("eta %f, phi %f, ieta %d, iphi %d, sm %d\n",eta,phi,icol,irow,sm);
  74
  75   //-------------------------------------
  76   // Remove in case of close to border, by default always 1 but check what was set in reco utils
  77   //-------------------------------------
  78
  79   Bool_t okrow = kFALSE;
  80         Bool_t okcol = kFALSE;
  81   Int_t nborder = GetCaloUtils()->GetEMCALRecoUtils()->GetNumberOfCellsFromEMCALBorder();
  82   if (nborder<1) nborder = 1;
  83
  84   // Rows
  85   if(sm < 10)
  86   {
  87     if(irow >= nborder && irow < 24-nborder) okrow =kTRUE;
  88   }
  89   else
  90   {
  91     if((GetCaloUtils()->EMCALGeometryName()).Contains("12SM")) // 1/3 SM
  92     {
  93       if(irow >= nborder && irow < 8-nborder) okrow =kTRUE;
  94     }
  95     else // 1/2 SM
  96     {
  97       if(irow >= nborder && irow <12-nborder) okrow =kTRUE;
  98     }
  99   }
 100
 101   // Columns
 102   if(sm%2==0)
 103   {
 104     if(icol >= nborder)     okcol = kTRUE;
 105   }
 106   else
 107   {
 108     if(icol <  48-nborder)  okcol = kTRUE;
 109   }
 110
 111   //printf("okcol %d, okrow %d\n",okcol,okrow);
 112   if (!okcol || !okrow) return kTRUE;
 113
 114   //-------------------------------------
 115   // Check if the cell or those around are bad
 116   //-------------------------------------
 117
 118   if(GetCaloUtils()->GetEMCALChannelStatus(sm,icol, irow) > 0) return kTRUE ; // trigger falls into a bad channel
 119
 120   // Check if close there was a bad channel
 121 //  for(Int_t i = -1; i <= 1; i++)
 122 //  {
 123 //    for(Int_t j = -1; j <= 1; j++)
 124 //    {
 125 //      //printf("\t check icol %d, irow %d \n",icol+i, irow+j);
 126 //      if(GetCaloUtils()->GetEMCALChannelStatus(sm,icol+i, irow+j) > 0) return kTRUE ; // trigger falls into a bad channel
 127 //      //printf("\t ok\n");
 128 //    }
 129 //  }
 130
 131    //printf("\t OK\n");
 132
 133   return kFALSE;
 134
 135 }
 136
 137
 138 //__________________________________________________
 139 TObjString *  AliAnaRandomTrigger::GetAnalysisCuts()
 140 {
 141   //Save parameters used for analysis
 142   TString parList ; //this will be list of parameters used for this analysis.
 143   const Int_t buffersize = 255;
 144   char onePar[buffersize] ;
 145
 146   snprintf(onePar,buffersize,"--- AliAnaRandomTrigger ---\n") ;
 147   parList+=onePar ;
 148   snprintf(onePar,buffersize,"Detector: %s\n"    , fTriggerDetectorString.Data()) ;
 149   parList+=onePar ;
 150   snprintf(onePar,buffersize,"N per event = %d\n", fNRandom       ) ;
 151   parList+=onePar ;
 152   snprintf(onePar,buffersize,"Min E   = %3.2f - Max E   = %3.2f\n", GetMinPt(), GetMaxPt()) ;
 153   parList+=onePar ;
 154   snprintf(onePar,buffersize,"Min Eta = %3.2f - Max Eta = %3.2f\n", fEtaCut[0], fEtaCut[1]) ;
 155   parList+=onePar ;
 156   snprintf(onePar,buffersize,"Min Phi = %3.2f - Max Phi = %3.2f\n", fPhiCut[0], fPhiCut[1]) ;
 157   parList+=onePar ;
 158
 159   return new TObjString(parList) ;
 160 }
 161
 162 //_______________________________________________________
 163 TList *  AliAnaRandomTrigger::GetCreateOutputObjects()
 164 {
 165   // Create histograms to be saved in output file and
 166   // store them in fOutputContainer
 167
 168   TList * outputContainer = new TList() ;
 169   outputContainer->SetName("RandomTrigger") ;
 170
 171   Int_t nptbins  = GetHistogramRanges()->GetHistoPtBins(); Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Int_t netabins = GetHistogramRanges()->GetHistoEtaBins();
 172   Float_t ptmax  = GetHistogramRanges()->GetHistoPtMax();  Float_t phimax = GetHistogramRanges()->GetHistoPhiMax();  Float_t etamax = GetHistogramRanges()->GetHistoEtaMax();
 173   Float_t ptmin  = GetHistogramRanges()->GetHistoPtMin();  Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();  Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
 174
 175   fhE  = new TH1F ("hE","Random E distribution", nptbins,ptmin,ptmax);
 176   fhE->SetXTitle("E (GeV)");
 177   outputContainer->Add(fhE);
 178
 179   fhPt  = new TH1F ("hPt","Random p_{T} distribution", nptbins,ptmin,ptmax);
 180   fhPt->SetXTitle("p_{T} (GeV/c)");
 181   outputContainer->Add(fhPt);
 182
 183   fhPhi  = new TH2F ("hPhi","Random #phi distribution",
 184                         nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 185   fhPhi->SetYTitle("#phi (rad)");
 186   fhPhi->SetXTitle("p_{T} (GeV/c)");
 187   outputContainer->Add(fhPhi);
 188
 189   fhEta  = new TH2F ("hEta","Random #eta distribution",
 190                         nptbins,ptmin,ptmax, netabins,etamin,etamax);
 191   fhEta->SetYTitle("#eta ");
 192   fhEta->SetXTitle("p_{T} (GeV/c)");
 193   outputContainer->Add(fhEta);
 194
 195   fhEtaPhi  = new TH2F ("hEtaPhi","Random #eta vs #phi ",netabins,etamin,etamax, nphibins,phimin,phimax);
 196   fhEtaPhi->SetXTitle("#eta ");
 197   fhEtaPhi->SetYTitle("#phi (rad)");
 198   outputContainer->Add(fhEtaPhi);
 199
 200   return outputContainer;
 201
 202 }
 203
 204 //___________________________________________
 205 void AliAnaRandomTrigger::InitParameters()
 206 {
 207   //Initialize the parameters of the analysis.
 208   SetOutputAODClassName("AliAODPWG4ParticleCorrelation");
 209   SetOutputAODName("RandomTrigger");
 210
 211   AddToHistogramsName("AnaRandomTrigger_");
 212
 213   fNRandom   = 1    ;
 214   fPhiCut[0] = 0.   ;
 215   fPhiCut[1] = TMath::TwoPi() ;
 216   fEtaCut[0] =-1.   ;
 217   fEtaCut[1] = 1.   ;
 218
 219 }
 220
 221 //____________________________________________________________
 222 void AliAnaRandomTrigger::Print(const Option_t * opt) const
 223 {
 224   //Print some relevant parameters set for the analysis
 225   if(! opt)
 226     return;
 227
 228   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
 229   AliAnaCaloTrackCorrBaseClass::Print(" ");
 230
 231   printf("Detector = %s\n",  fTriggerDetectorString.Data());
 232   printf("Min E   = %3.2f - Max E   = %3.2f\n", GetMinPt(), GetMaxPt());
 233   printf("Min Eta = %3.2f - Max Eta = %3.2f\n", fEtaCut[0], fEtaCut[1]);
 234   printf("Min Phi = %3.2f - Max Phi = %3.2f\n", fPhiCut[0], fPhiCut[1]);
 235
 236 }
 237
 238 //______________________________________________
 239 void  AliAnaRandomTrigger::MakeAnalysisFillAOD()
 240 {
 241   // Do analysis and fill aods
 242   // Generate particle randomly
 243   // fNRandom particles per event
 244
 245   for(Int_t irandom = 0; irandom < fNRandom; irandom++)
 246   {
 247     // Get the random variables of the trigger
 248     Float_t pt  = fRandom.Uniform(GetMinPt(), GetMaxPt());
 249     Float_t eta = fRandom.Uniform(fEtaCut[0], fEtaCut[1]);
 250     Float_t phi = fRandom.Uniform(fPhiCut[0], fPhiCut[1]);
 251
 252     // Check if particle falls into a dead region, if inside, get new
 253     Bool_t excluded =  ExcludeDeadBadRegions(eta,phi);
 254
 255     // if excluded, generate a new trigger until accepted
 256     while (excluded)
 257     {
 258       pt  = fRandom.Uniform(GetMinPt(), GetMaxPt());
 259       eta = fRandom.Uniform(fEtaCut[0], fEtaCut[1]);
 260       phi = fRandom.Uniform(fPhiCut[0], fPhiCut[1]);
 261
 262       excluded = ExcludeDeadBadRegions(eta,phi);
 263     }
 264
 265     // Create the AOD trigger object
 266     fMomentum.SetPtEtaPhiM(pt,eta,phi,0);
 267
 268     AliAODPWG4Particle trigger = AliAODPWG4Particle(fMomentum);
 269     trigger.SetDetectorTag(fTriggerDetector);
 270
 271     if(GetDebug() > 1)
 272       printf("AliAnaRandomTrigger::MakeAnalysisFillAOD() -  iRandom %d, Trigger e %2.2f pt %2.2f, phi %2.2f, eta %2.2f \n",
 273              irandom, trigger.E(), trigger.Pt(), trigger.Phi(), trigger.Eta());
 274
 275     AddAODParticle(trigger);
 276   }
 277
 278   if(GetDebug() > 0)
 279     printf("AliAnaRandomTrigger::MakeAnalysisFillAOD() - Final aod branch entries %d\n", GetOutputAODBranch()->GetEntriesFast());
 280 }
 281
 282 //_____________________________________________________
 283 void  AliAnaRandomTrigger::MakeAnalysisFillHistograms()
 284 {
 285   // Fill histograms
 286
 287   //Loop on stored AODParticles
 288   Int_t naod = GetOutputAODBranch()->GetEntriesFast();
 289
 290   if(GetDebug() > 0)
 291     printf("AliAnaRandomTrigger::MakeAnalysisFillHistograms() - aod branch entries %d, fNRandom %d\n", naod, fNRandom);
 292
 293   for(Int_t iaod = 0; iaod < naod ; iaod++)
 294   {
 295     AliAODPWG4Particle* trigger =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
 296
 297     fhPt    ->Fill(trigger->Pt());
 298     fhE     ->Fill(trigger->E());
 299     fhPhi   ->Fill(trigger->Pt(), trigger->Phi());
 300     fhEta   ->Fill(trigger->Pt(), trigger->Eta());
 301     fhEtaPhi->Fill(trigger->Eta(),trigger->Phi());
 302
 303   }// aod branch loop
 304
 305 }
 306
 307
 308 //_________________________________________________________
 309 void AliAnaRandomTrigger::SetTriggerDetector(TString & det)
 310 {
 311   // Set the detrimeter for the analysis
 312
 313   fTriggerDetectorString = det;
 314
 315   if     (det=="EMCAL") fTriggerDetector = kEMCAL;
 316   else if(det=="PHOS" ) fTriggerDetector = kPHOS;
 317   else if(det=="CTS")   fTriggerDetector = kCTS;
 318   else if(det=="DCAL")  fTriggerDetector = kDCAL;
 319   else if(det.Contains("DCAL") && det.Contains("PHOS")) fTriggerDetector = kDCALPHOS;
 320   else AliFatal(Form("Detector < %s > not known!", det.Data()));
 321
 322 }
 323
 324 //______________________________________________________
 325 void AliAnaRandomTrigger::SetTriggerDetector(Int_t det)
 326 {
 327   // Set the detrimeter for the analysis
 328
 329   fTriggerDetector = det;
 330
 331   if     (det==kEMCAL)    fTriggerDetectorString = "EMCAL";
 332   else if(det==kPHOS )    fTriggerDetectorString = "PHOS";
 333   else if(det==kCTS)      fTriggerDetectorString = "CTS";
 334   else if(det==kDCAL)     fTriggerDetectorString = "DCAL";
 335   else if(det==kDCALPHOS) fTriggerDetectorString = "DCAL_PHOS";
 336   else AliFatal(Form("Detector < %d > not known!", det));
 337
 338 }
 339
 340