AliPhysics  7f2a7c4 (7f2a7c4)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliClusterContainer.cxx
Go to the documentation of this file.
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 #include <algorithm>
16 #include <iostream>
17 #include <vector>
18 #include <TClonesArray.h>
19 
20 #include "AliVEvent.h"
21 #include "AliLog.h"
22 #include "AliTLorentzVector.h"
23 
24 #include "AliClusterContainer.h"
25 
29 
34  AliEmcalContainer(),
35  fClusTimeCutLow(-10),
36  fClusTimeCutUp(10),
37  fExoticCut(kTRUE),
38  fDefaultClusterEnergy(-1),
39  fIncludePHOS(kFALSE),
40  fPhosMinNcells(0),
41  fPhosMinM02(0)
42 {
43  fBaseClassName = "AliVCluster";
44  SetClassName("AliVCluster");
45 
46  for (Int_t i = 0; i <= AliVCluster::kLastUserDefEnergy; i++) {
47  fUserDefEnergyCut[i] = 0.;
48  }
49 }
50 
56  AliEmcalContainer(name),
57  fClusTimeCutLow(-10),
58  fClusTimeCutUp(10),
59  fExoticCut(kTRUE),
60  fDefaultClusterEnergy(-1),
61  fIncludePHOS(kFALSE),
62  fPhosMinNcells(0),
63  fPhosMinM02(0)
64 {
65  fBaseClassName = "AliVCluster";
66  SetClassName("AliVCluster");
67 
68  for (Int_t i = 0; i <= AliVCluster::kLastUserDefEnergy; i++) {
69  fUserDefEnergyCut[i] = 0.;
70  }
71 }
72 
78 AliVCluster* AliClusterContainer::GetLeadingCluster(const char* opt)
79 {
80  TString option(opt);
81  option.ToLower();
82 
83  double (AliTLorentzVector::*momentum)() const = 0;
84 
85  if (option.Contains("e")) {
86  momentum = &AliTLorentzVector::E;
87  }
88  else {
89  momentum = &AliTLorentzVector::Et;
90  }
91 
93 
94  for (auto cluster : accepted_momentum()) {
95  if ((clusterMax.first.*momentum)() < (cluster.first.*momentum)()) {
96  clusterMax = cluster;
97  }
98  }
99 
100  return clusterMax.second;
101 }
102 
109 {
110  if(i<0 || i>fClArray->GetEntriesFast()) return 0;
111  AliVCluster *vp = static_cast<AliVCluster*>(fClArray->At(i));
112  return vp;
113 
114 }
115 
122 {
123  AliVCluster *vc = GetCluster(i);
124  if (!vc) return 0;
125 
126  UInt_t rejectionReason = 0;
127  if (AcceptCluster(vc, rejectionReason))
128  return vc;
129  else {
130  AliDebug(2,"Cluster not accepted.");
131  return 0;
132  }
133 }
134 
141 {
142  Int_t i = GetIndexFromLabel(lab);
143  return GetCluster(i);
144 }
145 
152 {
153  Int_t i = GetIndexFromLabel(lab);
154  return GetAcceptCluster(i);
155 }
156 
162 {
163  const Int_t n = GetNEntries();
164  AliVCluster *c = 0;
165  do {
166  fCurrentID++;
167  if (fCurrentID >= n) break;
168  c = GetAcceptCluster(fCurrentID);
169  } while (!c);
170 
171  return c;
172 }
173 
179 {
180  const Int_t n = GetNEntries();
181  AliVCluster *c = 0;
182  do {
183  fCurrentID++;
184  if (fCurrentID >= n) break;
185  c = GetCluster(fCurrentID);
186  } while (!c);
187 
188  return c;
189 }
190 
191 Bool_t AliClusterContainer::GetMomentum(TLorentzVector &mom, const AliVCluster* vc, Double_t mass) const
192 {
193  if (mass < 0) mass = 0;
194 
195  Double_t energy = 0;
196 
197  if (fDefaultClusterEnergy >= 0 && fDefaultClusterEnergy <= AliVCluster::kLastUserDefEnergy) {
198  energy = vc->GetUserDefEnergy((AliVCluster::VCluUserDefEnergy_t)fDefaultClusterEnergy);
199  }
200  else {
201  energy = vc->E();
202  }
203 
204  Double_t p = TMath::Sqrt(energy*energy - mass*mass);
205 
206  Float_t pos[3];
207  vc->GetPosition(pos);
208 
209  pos[0]-=fVertex[0];
210  pos[1]-=fVertex[1];
211  pos[2]-=fVertex[2];
212 
213  Double_t r = TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]+pos[2]*pos[2]) ;
214 
215  if (r > 1e-12) {
216  mom.SetPxPyPzE( p*pos[0]/r, p*pos[1]/r, p*pos[2]/r, energy) ;
217  }
218  else {
219  AliInfo("Null cluster radius, momentum calculation not possible");
220  return kFALSE;
221  }
222 
223  return kTRUE;
224 }
225 
226 Bool_t AliClusterContainer::GetMomentum(TLorentzVector &mom, const AliVCluster* vc) const
227 {
228  if (fMassHypothesis > 0) return GetMomentum(mom, vc, fMassHypothesis);
229 
230  if (vc) {
231  if (fDefaultClusterEnergy >= 0 && fDefaultClusterEnergy <= AliVCluster::kLastUserDefEnergy) {
232  vc->GetMomentum(mom, fVertex, (AliVCluster::VCluUserDefEnergy_t)fDefaultClusterEnergy);
233  }
234  else {
235  vc->GetMomentum(mom, fVertex);
236  }
237  }
238  else {
239  mom.SetPtEtaPhiM(0, 0, 0, 0.139);
240  return kFALSE;
241  }
242  return kTRUE;
243 }
244 
251 Bool_t AliClusterContainer::GetMomentum(TLorentzVector &mom, Int_t i) const
252 {
253  AliVCluster *vc = GetCluster(i);
254  return GetMomentum(mom, vc);
255 }
256 
263 {
264  AliVCluster *vc = GetNextCluster();
265  return GetMomentum(mom, vc);
266 }
267 
275 {
276  AliVCluster *vc = GetAcceptCluster(i);
277  return GetMomentum(mom, vc);
278 }
279 
286 {
287  AliVCluster *vc = GetNextAcceptCluster();
288  return GetMomentum(mom, vc);
289 }
290 
292 {
293  Bool_t r = ApplyClusterCuts(GetCluster(i), rejectionReason);
294  if (!r) return kFALSE;
295 
296  AliTLorentzVector mom;
297  GetMomentum(mom, i);
298 
299  return ApplyKinematicCuts(mom, rejectionReason);
300 }
301 
302 Bool_t AliClusterContainer::AcceptCluster(const AliVCluster* clus, UInt_t &rejectionReason) const
303 {
304  Bool_t r = ApplyClusterCuts(clus, rejectionReason);
305  if (!r) return kFALSE;
306 
307  AliTLorentzVector mom;
308  GetMomentum(mom, clus);
309 
310  return ApplyKinematicCuts(mom, rejectionReason);
311 }
312 
318 Bool_t AliClusterContainer::ApplyClusterCuts(const AliVCluster* clus, UInt_t &rejectionReason) const
319 {
320  if (!clus) {
321  rejectionReason |= kNullObject;
322  return kFALSE;
323  }
324 
325  Bool_t bInAcceptance = clus->IsEMCAL();
326  if (fIncludePHOS) bInAcceptance = clus->IsEMCAL() || (clus->GetType() == AliVCluster::kPHOSNeutral);
327  if (!bInAcceptance) {
328  rejectionReason |= kIsEMCalCut;
329  return kFALSE;
330  }
331 
332  if (clus->TestBits(fBitMap) != (Int_t)fBitMap) {
333  rejectionReason |= kBitMapCut;
334  return kFALSE;
335  }
336 
337  if (fMinMCLabel >= 0 && TMath::Abs(clus->GetLabel()) < fMinMCLabel) {
338  rejectionReason |= kMCLabelCut;
339  return kFALSE;
340  }
341 
342  if (fMaxMCLabel >= 0 && TMath::Abs(clus->GetLabel()) > fMaxMCLabel) {
343  rejectionReason |= kMCLabelCut;
344  return kFALSE;
345  }
346 
347  if (clus->GetTOF() > fClusTimeCutUp || clus->GetTOF() < fClusTimeCutLow) {
348  rejectionReason |= kTimeCut;
349  return kFALSE;
350  }
351 
352  if (fExoticCut && clus->GetIsExotic()) {
353  rejectionReason |= kExoticCut;
354  return kFALSE;
355  }
356 
357  for (Int_t i = 0; i <= AliVCluster::kLastUserDefEnergy; i++) {
358  if (clus->GetUserDefEnergy((VCluUserDefEnergy_t)i) < fUserDefEnergyCut[i]) {
359  rejectionReason |= kEnergyCut;
360  return kFALSE;
361  }
362  }
363 
364  if (fIncludePHOS && clus->GetType() == AliVCluster::kPHOSNeutral) {
365  if (clus->GetNCells() < fPhosMinNcells) {
366  rejectionReason |= kExoticCut;
367  return kFALSE;
368  }
369 
370  if (clus->GetM02() < fPhosMinM02) {
371  rejectionReason |= kExoticCut;
372  return kFALSE;
373  }
374  }
375 
376  return kTRUE;
377 }
378 
384 {
385  UInt_t rejectionReason = 0;
386  Int_t nClus = 0;
387  for(int iclust = 0; iclust < this->fClArray->GetEntries(); ++iclust){
388  AliVCluster *clust = this->GetCluster(iclust);
389  if(this->AcceptCluster(clust, rejectionReason)) nClus++;
390  }
391  return nClus;
392 }
393 
400 {
401  if (t >= 0 && t <= AliVCluster::kLastUserDefEnergy){
402  return fUserDefEnergyCut[t];
403  }
404  else {
405  return fMinE;
406  }
407 }
408 
415 {
416  if (t >= 0 && t <= AliVCluster::kLastUserDefEnergy){
417  fUserDefEnergyCut[t] = cut;
418  }
419  else {
420  fMinE = cut;
421  }
422 }
423 
430  return AliClusterIterableContainer(this, false);
431 }
432 
439  return AliClusterIterableContainer(this, true);
440 }
441 
448  return AliClusterIterableMomentumContainer(this, false);
449 }
450 
457  return AliClusterIterableMomentumContainer(this, true);
458 }
459 
460 const char* AliClusterContainer::GetTitle() const
461 {
462  static TString clusterString;
463 
465 
466  if (Ecut == 0) {
467  clusterString = TString::Format("%s_E0000", GetArrayName().Data());
468  }
469  else if (Ecut < 1.0) {
470  clusterString = TString::Format("%s_E0%3.0f", GetArrayName().Data(), Ecut*1000.0);
471  }
472  else {
473  clusterString = TString::Format("%s_E%4.0f", GetArrayName().Data(), Ecut*1000.0);
474  }
475 
476  return clusterString.Data();
477 }
478 
479 
480 /******************************************
481  * Unit tests *
482  ******************************************/
483 
484 int TestClusterContainerIterator(const AliClusterContainer *const cont, int iteratorType, bool verbose){
485  std::vector<AliVCluster *> reference, variation;
486  AliVCluster *test = NULL;
487  for(int iclust = 0; iclust < cont->GetNClusters(); iclust++){
488  test = cont->GetCluster(iclust);
489  if(!iteratorType){
490  UInt_t rejectionReason = 0;
491  if(!cont->AcceptCluster(test, rejectionReason)) continue;
492  }
493  reference.push_back(test);
494  }
495 
496  if(!iteratorType){
497  // test accept iterator
498  for(auto cluster : cont->accepted()){
499  variation.push_back(cluster);
500  }
501  } else {
502  // test all iterator
503  for(auto cluster : cont->all()){
504  variation.push_back(cluster);
505  }
506  }
507 
508  if(verbose){
509  // Printing cluster addresses:
510  if(reference.size() < 30){
511  std::cout << "Clusters in reference container: " << std::endl;
512  std::cout << "===========================================" << std::endl;
513  for(std::vector<AliVCluster *>::iterator refit = reference.begin(); refit != reference.end(); ++refit){
514  std::cout << "Address: " << *refit << std::endl;
515  }
516  }
517  if(variation.size() < 30){
518  std::cout << "Clusters in test container: " << std::endl;
519  std::cout << "===========================================" << std::endl;
520  for(std::vector<AliVCluster *>::iterator varit = variation.begin(); varit != variation.end(); ++varit){
521  std::cout << "Address: " << *varit << std::endl;
522  }
523  }
524  }
525 
526  int testresult = 0;
527  // compare distributions: all clusters in one vector needs to be found in the other vector and vice versa
528  bool failure = false;
529  // first test: cleck if all clusters are found by the iterator
530  for(std::vector<AliVCluster *>::iterator clit = reference.begin(); clit != reference.end(); ++clit){
531  if(std::find(variation.begin(), variation.end(), *clit) == variation.end()) {
532  if(verbose)
533  std::cout << "Could not find cluster with address " << *clit << " in test container" << std::endl;
534  failure = true;
535  break;
536  }
537  }
538  if(!failure){
539  // second test: check if there are no extra clusters found
540  for(std::vector<AliVCluster *>::iterator clit = variation.begin(); clit != variation.end(); ++clit){
541  if(std::find(reference.begin(), reference.end(), *clit) == reference.end()) {
542  if(verbose)
543  std::cout << "Could not find cluster with address " << *clit << " in reference container" << std::endl;
544  failure = true;
545  break;
546  }
547  }
548  if(failure) testresult = 2;
549  } else {
550  testresult = 1;
551  }
552 
553  if(verbose){
554  std::cout << "Unit test cluster container, iterator type " << iteratorType << std::endl;
555  std::cout << "Number of expected clusters: " << reference.size() << ", number of found clusters: " << variation.size() << std::endl;
556  std::cout << "Test result: " << testresult << std::endl;
557  std::cout << "-----------------------------------------------------------------------" << std::endl;
558  }
559  return testresult;
560 }
AliVCluster * GetAcceptClusterWithLabel(Int_t lab) const
virtual Bool_t ApplyClusterCuts(const AliVCluster *clus, UInt_t &rejectionReason) const
Bool_t GetNextAcceptMomentum(TLorentzVector &mom)
const char * GetTitle() const
double Double_t
Definition: External.C:58
EMCALIterableContainer::AliEmcalIterableContainerT< AliVCluster, EMCALIterableContainer::operator_star_object< AliVCluster > > AliClusterIterableContainer
Bool_t GetNextMomentum(TLorentzVector &mom)
const AliClusterIterableMomentumContainer accepted_momentum() const
Double_t mass
energy
Definition: HFPtSpectrum.C:44
Bool_t fIncludePHOS
whether or not to include PHOS clusters in addition to EMCal clusters
Declaration of class AliTLorentzVector.
Double_t fClusTimeCutLow
low time cut for clusters
Double_t fUserDefEnergyCut[AliVCluster::kLastUserDefEnergy+1]
cut on the energy of the cluster after higher level corrections (see AliVCluster.h) ...
TCanvas * c
Definition: TestFitELoss.C:172
void SetClusUserDefEnergyCut(Int_t t, Double_t cut)
AliVCluster * GetLeadingCluster(const char *opt="")
Double_t GetClusUserDefEnergyCut(Int_t t) const
Int_t GetDefaultClusterEnergy() const
const AliClusterIterableContainer all() const
enum AliVCluster::VCluUserDefEnergy_t VCluUserDefEnergy_t
std::pair< AliTLorentzVector, T * > momentum_object_pair
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
virtual Bool_t AcceptCluster(Int_t i, UInt_t &rejectionReason) const
AliVCluster * GetAcceptCluster(Int_t i) const
Int_t GetNClusters() const
Int_t GetNAcceptedClusters() const
const AliClusterIterableContainer accepted() const
Double_t fClusTimeCutUp
up time cut for clusters
AliVCluster * GetCluster(Int_t i) const
Bool_t fExoticCut
reject clusters marked as "exotic"
AliVCluster * GetNextCluster()
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
EMCALIterableContainer::AliEmcalIterableContainerT< AliVCluster, EMCALIterableContainer::operator_star_pair< AliVCluster > > AliClusterIterableMomentumContainer
const AliClusterIterableMomentumContainer all_momentum() const
Bool_t GetAcceptMomentum(TLorentzVector &mom, Int_t i) const
AliVCluster * GetClusterWithLabel(Int_t lab) const
Double_t fPhosMinM02
min value of M02 for phos clusters
int TestClusterContainerIterator(const AliClusterContainer *const cont, int iteratorType, bool verbose)
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
Bool_t GetMomentum(TLorentzVector &mom, const AliVCluster *vc, Double_t mass) const
bool Bool_t
Definition: External.C:53
Int_t fDefaultClusterEnergy
default cluster energy: -1 for clus->E(); otherwise clus->GetUserDefEnergy(fDefaultClusterEnergy) ...
Container structure for EMCAL clusters.
AliVCluster * GetNextAcceptCluster()
Int_t fPhosMinNcells
min number of phos cells per cluster