AliPhysics  3c9d42f (3c9d42f)
 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 
35  fClusTimeCutLow(-10),
36  fClusTimeCutUp(10),
37  fExoticCut(kTRUE),
38  fDefaultClusterEnergy(-1)
39 {
40  fClassName = "AliVCluster";
41 
42  for (Int_t i = 0; i <= AliVCluster::kLastUserDefEnergy; i++) {
43  fUserDefEnergyCut[i] = 0.;
44  }
45 }
46 
52  AliEmcalContainer(name),
53  fClusTimeCutLow(-10),
54  fClusTimeCutUp(10),
55  fExoticCut(kTRUE),
56  fDefaultClusterEnergy(-1)
57 {
58  fClassName = "AliVCluster";
59 
60  for (Int_t i = 0; i <= AliVCluster::kLastUserDefEnergy; i++) {
61  fUserDefEnergyCut[i] = 0.;
62  }
63 }
64 
70 AliVCluster* AliClusterContainer::GetLeadingCluster(const char* opt)
71 {
72  TString option(opt);
73  option.ToLower();
74 
75  Int_t tempID = fCurrentID;
77 
78  AliVCluster *clusterMax = GetNextAcceptCluster();
79  AliVCluster *cluster = 0;
80 
81  if (option.Contains("e")) {
82  while ((cluster = GetNextAcceptCluster())) {
83  if (cluster->E() > clusterMax->E()) clusterMax = cluster;
84  }
85  }
86  else {
87  Double_t et = 0;
88  Double_t etmax = 0;
89  while ((cluster = GetNextAcceptCluster())) {
90  TLorentzVector mom;
91  cluster->GetMomentum(mom,const_cast<Double_t*>(fVertex));
92  et = mom.Et();
93  if (et > etmax) {
94  clusterMax = cluster;
95  etmax = et;
96  }
97  }
98  }
99 
100  fCurrentID = tempID;
101 
102  return clusterMax;
103 }
104 
110 AliVCluster* AliClusterContainer::GetCluster(Int_t i) const
111 {
112  if(i<0 || i>fClArray->GetEntriesFast()) return 0;
113  AliVCluster *vp = static_cast<AliVCluster*>(fClArray->At(i));
114  return vp;
115 
116 }
117 
123 AliVCluster* AliClusterContainer::GetAcceptCluster(Int_t i) const
124 {
125  AliVCluster *vc = GetCluster(i);
126  if (!vc) return 0;
127 
128  UInt_t rejectionReason = 0;
129  if (AcceptCluster(vc, rejectionReason))
130  return vc;
131  else {
132  AliDebug(2,"Cluster not accepted.");
133  return 0;
134  }
135 }
136 
142 AliVCluster* AliClusterContainer::GetClusterWithLabel(Int_t lab) const
143 {
144  Int_t i = GetIndexFromLabel(lab);
145  return GetCluster(i);
146 }
147 
154 {
155  Int_t i = GetIndexFromLabel(lab);
156  return GetAcceptCluster(i);
157 }
158 
164 {
165  const Int_t n = GetNEntries();
166  AliVCluster *c = 0;
167  do {
168  fCurrentID++;
169  if (fCurrentID >= n) break;
171  } while (!c);
172 
173  return c;
174 }
175 
181 {
182  const Int_t n = GetNEntries();
183  AliVCluster *c = 0;
184  do {
185  fCurrentID++;
186  if (fCurrentID >= n) break;
187  c = GetCluster(fCurrentID);
188  } while (!c);
189 
190  return c;
191 }
192 
193 Bool_t AliClusterContainer::GetMomentum(TLorentzVector &mom, const AliVCluster* vc, Double_t mass) const
194 {
195  if (mass < 0) mass = 0;
196 
197  Double_t energy = 0;
198 
199  if (fDefaultClusterEnergy >= 0 && fDefaultClusterEnergy <= AliVCluster::kLastUserDefEnergy) {
200  energy = vc->GetUserDefEnergy((AliVCluster::VCluUserDefEnergy_t)fDefaultClusterEnergy);
201  }
202  else {
203  energy = vc->E();
204  }
205 
206  Double_t p = TMath::Sqrt(energy*energy - mass*mass);
207 
208  Float_t pos[3];
209  vc->GetPosition(pos);
210 
211  pos[0]-=fVertex[0];
212  pos[1]-=fVertex[1];
213  pos[2]-=fVertex[2];
214 
215  Double_t r = TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]+pos[2]*pos[2]) ;
216 
217  if (r > 1e-12) {
218  mom.SetPxPyPzE( p*pos[0]/r, p*pos[1]/r, p*pos[2]/r, energy) ;
219  }
220  else {
221  AliInfo("Null cluster radius, momentum calculation not possible");
222  return kFALSE;
223  }
224 
225  return kTRUE;
226 }
227 
228 Bool_t AliClusterContainer::GetMomentum(TLorentzVector &mom, const AliVCluster* vc) const
229 {
230  if (fMassHypothesis > 0) return GetMomentum(mom, vc, fMassHypothesis);
231 
232  if (vc) {
233  if (fDefaultClusterEnergy >= 0 && fDefaultClusterEnergy <= AliVCluster::kLastUserDefEnergy) {
234  vc->GetMomentum(mom, fVertex, (AliVCluster::VCluUserDefEnergy_t)fDefaultClusterEnergy);
235  }
236  else {
237  vc->GetMomentum(mom, fVertex);
238  }
239  }
240  else {
241  mom.SetPtEtaPhiM(0, 0, 0, 0.139);
242  return kFALSE;
243  }
244  return kTRUE;
245 }
246 
253 Bool_t AliClusterContainer::GetMomentum(TLorentzVector &mom, Int_t i) const
254 {
255  AliVCluster *vc = GetCluster(i);
256  return GetMomentum(mom, vc);
257 }
258 
264 Bool_t AliClusterContainer::GetNextMomentum(TLorentzVector &mom)
265 {
266  AliVCluster *vc = GetNextCluster();
267  return GetMomentum(mom, vc);
268 }
269 
276 Bool_t AliClusterContainer::GetAcceptMomentum(TLorentzVector &mom, Int_t i) const
277 {
278  AliVCluster *vc = GetAcceptCluster(i);
279  return GetMomentum(mom, vc);
280 }
281 
288 {
289  AliVCluster *vc = GetNextAcceptCluster();
290  return GetMomentum(mom, vc);
291 }
292 
293 Bool_t AliClusterContainer::AcceptCluster(Int_t i, UInt_t &rejectionReason) const
294 {
295  Bool_t r = ApplyClusterCuts(GetCluster(i), rejectionReason);
296  if (!r) return kFALSE;
297 
298  AliTLorentzVector mom;
299  GetMomentum(mom, i);
300 
301  return ApplyKinematicCuts(mom, rejectionReason);
302 }
303 
304 Bool_t AliClusterContainer::AcceptCluster(const AliVCluster* clus, UInt_t &rejectionReason) const
305 {
306  Bool_t r = ApplyClusterCuts(clus, rejectionReason);
307  if (!r) return kFALSE;
308 
309  AliTLorentzVector mom;
310  GetMomentum(mom, clus);
311 
312  return ApplyKinematicCuts(mom, rejectionReason);
313 }
314 
320 Bool_t AliClusterContainer::ApplyClusterCuts(const AliVCluster* clus, UInt_t &rejectionReason) const
321 {
322  if (!clus) {
323  rejectionReason |= kNullObject;
324  return kFALSE;
325  }
326 
327  if (!clus->IsEMCAL()) {
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  return kTRUE;
365 }
366 
372 {
373  UInt_t rejectionReason = 0;
374  Int_t nClus = 0;
375  for(int iclust = 0; iclust < this->fClArray->GetEntries(); ++iclust){
376  AliVCluster *clust = this->GetCluster(iclust);
377  if(this->AcceptCluster(clust, rejectionReason)) nClus++;
378  }
379  return nClus;
380 }
381 
383 {
384  if (t >= 0 && t <= AliVCluster::kLastUserDefEnergy){
385  return fUserDefEnergyCut[t];
386  }
387  else {
388  return fMinE;
389  }
390 }
391 
393 {
394  if (t >= 0 && t <= AliVCluster::kLastUserDefEnergy){
395  fUserDefEnergyCut[t] = cut;
396  }
397  else {
398  fMinE = cut;
399  }
400 }
401 
402 
407 void AliClusterContainer::SetClassName(const char *clname)
408 {
409  TClass cls(clname);
410  if (cls.InheritsFrom("AliVCluster")) fClassName = clname;
411  else AliError(Form("Unable to set class name %s for a AliClusterContainer, it must inherits from AliVCluster!",clname));
412 }
413 
414 const char* AliClusterContainer::GetTitle() const
415 {
416  static TString clusterString;
417 
419 
420  if (Ecut == 0) {
421  clusterString = TString::Format("%s_E0000", GetArrayName().Data());
422  }
423  else if (Ecut < 1.0) {
424  clusterString = TString::Format("%s_E0%3.0f", GetArrayName().Data(), Ecut*1000.0);
425  }
426  else {
427  clusterString = TString::Format("%s_E%4.0f", GetArrayName().Data(), Ecut*1000.0);
428  }
429 
430  return clusterString.Data();
431 }
432 
439  return all_iterator(this, 0, true);
440 }
441 
448  return all_iterator(this, this->GetNClusters(), true);
449 }
450 
457  return all_iterator(this, this->GetNClusters()-1, false);
458 }
459 
466  return all_iterator(this, -1, false);
467 }
468 
475  return accept_iterator(this, 0, true);
476 }
477 
484  return accept_iterator(this, this->GetNAcceptedClusters(), true);
485 }
486 
493  return accept_iterator(this, this->GetNAcceptedClusters()-1, false);
494 }
495 
502  return accept_iterator(this, -1, false);
503 }
504 
505 
511 
528  const AliClusterContainer *cont,
529  int start,
530  bool forward):
531  fkContainer(cont),
532  fAcceptIndices(),
533  fCurrentPos(start),
534  fForward(forward)
535 {
537  UInt_t rejection = 0;
538  Int_t mappos = 0;
539  for(int i = 0; i < fkContainer->GetNEntries(); i++){
540  if(fkContainer->AcceptCluster(i, rejection)){
541  fAcceptIndices[mappos++] = i;
542  }
543  }
544 }
545 
552  fkContainer(other.fkContainer),
553  fAcceptIndices(other.fAcceptIndices),
554  fCurrentPos(other.fCurrentPos),
555  fForward(other.fForward)
556 {
557 }
558 
566  if(this != &other){
567  fkContainer = other.fkContainer;
568  fAcceptIndices = other.fAcceptIndices;
569  fCurrentPos = other.fCurrentPos;
570  fForward = other.fForward;
571  }
572  return *this;
573 }
574 
582  return fCurrentPos != other.fCurrentPos;
583 }
584 
591  if(fForward){
592  fCurrentPos++;
593  } else {
594  fCurrentPos--;
595  }
596  return *this;
597 }
598 
607  this->operator++();
608  return result;
609 }
610 
617  if(fForward){
618  fCurrentPos--;
619  } else {
620  fCurrentPos++;
621  }
622  return *this;
623 }
624 
633  this->operator++();
634  return result;
635 }
636 
643  if(fCurrentPos < 0 || fCurrentPos >= fAcceptIndices.GetSize()){
644  return NULL;
645  }
646  return fkContainer->GetCluster(fAcceptIndices[fCurrentPos]);
647 }
648 
654 
669  const AliClusterContainer *cont,
670  int startpos,
671  bool forward):
672  fkContainer(cont),
673  fCurrentPos(startpos),
674  fForward(forward)
675 {
676 }
677 
684  fkContainer(other.fkContainer),
685  fCurrentPos(other.fCurrentPos),
686  fForward(other.fForward)
687 {
688 }
689 
697  if(&other != this){
698  fkContainer = other.fkContainer;
699  fCurrentPos = other.fCurrentPos;
700  fForward = other.fForward;
701  }
702  return *this;
703 }
704 
712  return other.fCurrentPos != fCurrentPos;
713 }
714 
721  if(fForward){
722  fCurrentPos++;
723  } else {
724  fCurrentPos--;
725  }
726  return *this;
727 }
728 
736  all_iterator tmp(*this);
737  operator++();
738  return tmp;
739 }
740 
747  if(fForward){
748  fCurrentPos--;
749  } else {
750  fCurrentPos++;
751  }
752  return *this;
753 }
754 
762  all_iterator tmp(*this);
763  operator--();
764  return tmp;
765 }
766 
773  if(fCurrentPos >= 0 && fCurrentPos < fkContainer->GetNClusters())
774  return fkContainer->GetCluster(fCurrentPos);
775  return NULL;
776 }
777 
778 /******************************************
779  * Unit tests *
780  ******************************************/
781 
782 int TestClusterContainerIterator(const AliClusterContainer *const cont, int iteratorType, bool verbose){
783  std::vector<AliVCluster *> reference, variation;
784  AliVCluster *test = NULL;
785  for(int iclust = 0; iclust < cont->GetNClusters(); iclust++){
786  test = cont->GetCluster(iclust);
787  if(!iteratorType){
788  UInt_t rejectionReason = 0;
789  if(!cont->AcceptCluster(test, rejectionReason)) continue;
790  }
791  reference.push_back(test);
792  }
793 
794  if(!iteratorType){
795  // test accept iterator
796  for(AliClusterContainer::accept_iterator iter = cont->accept_begin(); iter != cont->accept_end(); ++iter){
797  variation.push_back(*iter);
798  }
799  } else {
800  // test all iterator
801  for(AliClusterContainer::all_iterator iter = cont->begin(); iter != cont->end(); ++iter){
802  variation.push_back(*iter);
803  }
804  }
805 
806  if(verbose){
807  // Printing cluster addresses:
808  if(reference.size() < 30){
809  std::cout << "Clusters in reference container: " << std::endl;
810  std::cout << "===========================================" << std::endl;
811  for(std::vector<AliVCluster *>::iterator refit = reference.begin(); refit != reference.end(); ++refit){
812  std::cout << "Address: " << *refit << std::endl;
813  }
814  }
815  if(variation.size() < 30){
816  std::cout << "Clusters in test container: " << std::endl;
817  std::cout << "===========================================" << std::endl;
818  for(std::vector<AliVCluster *>::iterator varit = variation.begin(); varit != variation.end(); ++varit){
819  std::cout << "Address: " << *varit << std::endl;
820  }
821  }
822  }
823 
824  int testresult = 0;
825  // compare distributions: all clusters in one vector needs to be found in the other vector and vice versa
826  bool failure = false;
827  // first test: cleck if all clusters are found by the iterator
828  for(std::vector<AliVCluster *>::iterator clit = reference.begin(); clit != reference.end(); ++clit){
829  if(std::find(variation.begin(), variation.end(), *clit) == variation.end()) {
830  if(verbose)
831  std::cout << "Could not find cluster with address " << *clit << " in test container" << std::endl;
832  failure = true;
833  break;
834  }
835  }
836  if(!failure){
837  // second test: check if there are no extra clusters found
838  for(std::vector<AliVCluster *>::iterator clit = variation.begin(); clit != variation.end(); ++clit){
839  if(std::find(reference.begin(), reference.end(), *clit) == reference.end()) {
840  if(verbose)
841  std::cout << "Could not find cluster with address " << *clit << " in reference container" << std::endl;
842  failure = true;
843  break;
844  }
845  }
846  if(failure) testresult = 2;
847  } else {
848  testresult = 1;
849  }
850 
851  if(verbose){
852  std::cout << "Unit test cluster container, iterator type " << iteratorType << std::endl;
853  std::cout << "Number of expected clusters: " << reference.size() << ", number of found clusters: " << variation.size() << std::endl;
854  std::cout << "Test result: " << testresult << std::endl;
855  std::cout << "-----------------------------------------------------------------------" << std::endl;
856  }
857  return testresult;
858 }
AliVCluster * GetAcceptClusterWithLabel(Int_t lab) const
bool fForward
Direction, expressed in forward direction.
bool fForward
Direction, expressed in forward direction.
virtual Bool_t ApplyClusterCuts(const AliVCluster *clus, UInt_t &rejectionReason) const
Bool_t GetNextAcceptMomentum(TLorentzVector &mom)
virtual Bool_t ApplyKinematicCuts(const AliTLorentzVector &mom, UInt_t &rejectionReason) const
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
stl iterator over accepted clusters
const char * GetTitle() const
all_iterator begin() const
accept_iterator accept_rbegin() const
Double_t fMinE
Max. cut on particle energy.
Bool_t GetNextMomentum(TLorentzVector &mom)
Int_t GetIndexFromLabel(Int_t lab) const
Double_t mass
bool operator!=(const accept_iterator &other) const
accept_iterator accept_end() const
UInt_t fBitMap
bitmap mask
Double_t fClusTimeCutLow
low time cut for clusters
TString fClassName
name of the class in the TClonesArray
Double_t fUserDefEnergyCut[AliVCluster::kLastUserDefEnergy+1]
cut on the energy of the cluster after higher level corrections (see AliVCluster.h) ...
const AliClusterContainer * fkContainer
Container iterated over.
accept_iterator accept_begin() const
all_iterator end() const
int fCurrentPos
Current position inside the container.
void SetClusUserDefEnergyCut(Int_t t, Double_t cut)
AliVCluster * GetLeadingCluster(const char *opt="")
Cell time cut not passed.
Double_t GetClusUserDefEnergyCut(Int_t t) const
Int_t GetDefaultClusterEnergy() const
accept_iterator & operator=(const accept_iterator &other)
enum AliVCluster::VCluUserDefEnergy_t VCluUserDefEnergy_t
Cluster is exotic cluster.
Energy below threshold.
stl iterator over all clusters
all_iterator & operator=(const all_iterator &other)
Base class for container structures within the EMCAL framework.
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
int fCurrentPos
Current position inside the container.
all_iterator rend() 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 operator!=(const all_iterator &other) const
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)
void SetClassName(const char *clname)
energy
Int_t fMinMCLabel
minimum MC label
accept_iterator accept_rend() const
const TString & GetArrayName() const
Double_t fMassHypothesis
if < 0 it will use a PID mass when available
Cluster not in the EMCAL.
Double_t fVertex[3]
! event vertex array
Bool_t GetAcceptMomentum(TLorentzVector &mom, Int_t i) const
AliVCluster * GetClusterWithLabel(Int_t lab) const
int TestClusterContainerIterator(const AliClusterContainer *const cont, int iteratorType, bool verbose)
TArrayI fAcceptIndices
Indices of accepted clusters.
Bool_t GetMomentum(TLorentzVector &mom, const AliVCluster *vc, Double_t mass) const
Int_t fMaxMCLabel
maximum MC label
TClonesArray * fClArray
! Pointer to array in input event
Int_t GetNEntries() const
Int_t fDefaultClusterEnergy
default cluster energy: -1 for clus->E(); otherwise clus->GetUserDefEnergy(fDefaultClusterEnergy) ...
const AliClusterContainer * fkContainer
Container iterated over.
all_iterator rbegin() const
Int_t fCurrentID
! current ID for automatic loops
Container structure for EMCAL clusters.
AliVCluster * GetNextAcceptCluster()
void ResetCurrentID(Int_t i=-1)