AliPhysics  f05a842 (f05a842)
 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 {
40  fBaseClassName = "AliVCluster";
41  SetClassName("AliVCluster");
42 
43  for (Int_t i = 0; i <= AliVCluster::kLastUserDefEnergy; i++) {
44  fUserDefEnergyCut[i] = 0.;
45  }
46 }
47 
53  AliEmcalContainer(name),
54  fClusTimeCutLow(-10),
55  fClusTimeCutUp(10),
56  fExoticCut(kTRUE),
57  fDefaultClusterEnergy(-1)
58 {
59  fBaseClassName = "AliVCluster";
60  SetClassName("AliVCluster");
61 
62  for (Int_t i = 0; i <= AliVCluster::kLastUserDefEnergy; i++) {
63  fUserDefEnergyCut[i] = 0.;
64  }
65 }
66 
72 AliVCluster* AliClusterContainer::GetLeadingCluster(const char* opt)
73 {
74  TString option(opt);
75  option.ToLower();
76 
77  double (AliTLorentzVector::*momentum)() const = 0;
78 
79  if (option.Contains("e")) {
80  momentum = &AliTLorentzVector::E;
81  }
82  else {
83  momentum = &AliTLorentzVector::Et;
84  }
85 
87 
88  for (auto cluster : accepted_momentum()) {
89  if ((clusterMax.first.*momentum)() < (cluster.first.*momentum)()) {
90  clusterMax = cluster;
91  }
92  }
93 
94  return clusterMax.second;
95 }
96 
103 {
104  if(i<0 || i>fClArray->GetEntriesFast()) return 0;
105  AliVCluster *vp = static_cast<AliVCluster*>(fClArray->At(i));
106  return vp;
107 
108 }
109 
116 {
117  AliVCluster *vc = GetCluster(i);
118  if (!vc) return 0;
119 
120  UInt_t rejectionReason = 0;
121  if (AcceptCluster(vc, rejectionReason))
122  return vc;
123  else {
124  AliDebug(2,"Cluster not accepted.");
125  return 0;
126  }
127 }
128 
135 {
136  Int_t i = GetIndexFromLabel(lab);
137  return GetCluster(i);
138 }
139 
146 {
147  Int_t i = GetIndexFromLabel(lab);
148  return GetAcceptCluster(i);
149 }
150 
156 {
157  const Int_t n = GetNEntries();
158  AliVCluster *c = 0;
159  do {
160  fCurrentID++;
161  if (fCurrentID >= n) break;
162  c = GetAcceptCluster(fCurrentID);
163  } while (!c);
164 
165  return c;
166 }
167 
173 {
174  const Int_t n = GetNEntries();
175  AliVCluster *c = 0;
176  do {
177  fCurrentID++;
178  if (fCurrentID >= n) break;
179  c = GetCluster(fCurrentID);
180  } while (!c);
181 
182  return c;
183 }
184 
185 Bool_t AliClusterContainer::GetMomentum(TLorentzVector &mom, const AliVCluster* vc, Double_t mass) const
186 {
187  if (mass < 0) mass = 0;
188 
189  Double_t energy = 0;
190 
191  if (fDefaultClusterEnergy >= 0 && fDefaultClusterEnergy <= AliVCluster::kLastUserDefEnergy) {
192  energy = vc->GetUserDefEnergy((AliVCluster::VCluUserDefEnergy_t)fDefaultClusterEnergy);
193  }
194  else {
195  energy = vc->E();
196  }
197 
198  Double_t p = TMath::Sqrt(energy*energy - mass*mass);
199 
200  Float_t pos[3];
201  vc->GetPosition(pos);
202 
203  pos[0]-=fVertex[0];
204  pos[1]-=fVertex[1];
205  pos[2]-=fVertex[2];
206 
207  Double_t r = TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]+pos[2]*pos[2]) ;
208 
209  if (r > 1e-12) {
210  mom.SetPxPyPzE( p*pos[0]/r, p*pos[1]/r, p*pos[2]/r, energy) ;
211  }
212  else {
213  AliInfo("Null cluster radius, momentum calculation not possible");
214  return kFALSE;
215  }
216 
217  return kTRUE;
218 }
219 
220 Bool_t AliClusterContainer::GetMomentum(TLorentzVector &mom, const AliVCluster* vc) const
221 {
222  if (fMassHypothesis > 0) return GetMomentum(mom, vc, fMassHypothesis);
223 
224  if (vc) {
225  if (fDefaultClusterEnergy >= 0 && fDefaultClusterEnergy <= AliVCluster::kLastUserDefEnergy) {
226  vc->GetMomentum(mom, fVertex, (AliVCluster::VCluUserDefEnergy_t)fDefaultClusterEnergy);
227  }
228  else {
229  vc->GetMomentum(mom, fVertex);
230  }
231  }
232  else {
233  mom.SetPtEtaPhiM(0, 0, 0, 0.139);
234  return kFALSE;
235  }
236  return kTRUE;
237 }
238 
245 Bool_t AliClusterContainer::GetMomentum(TLorentzVector &mom, Int_t i) const
246 {
247  AliVCluster *vc = GetCluster(i);
248  return GetMomentum(mom, vc);
249 }
250 
257 {
258  AliVCluster *vc = GetNextCluster();
259  return GetMomentum(mom, vc);
260 }
261 
269 {
270  AliVCluster *vc = GetAcceptCluster(i);
271  return GetMomentum(mom, vc);
272 }
273 
280 {
281  AliVCluster *vc = GetNextAcceptCluster();
282  return GetMomentum(mom, vc);
283 }
284 
286 {
287  Bool_t r = ApplyClusterCuts(GetCluster(i), rejectionReason);
288  if (!r) return kFALSE;
289 
290  AliTLorentzVector mom;
291  GetMomentum(mom, i);
292 
293  return ApplyKinematicCuts(mom, rejectionReason);
294 }
295 
296 Bool_t AliClusterContainer::AcceptCluster(const AliVCluster* clus, UInt_t &rejectionReason) const
297 {
298  Bool_t r = ApplyClusterCuts(clus, rejectionReason);
299  if (!r) return kFALSE;
300 
301  AliTLorentzVector mom;
302  GetMomentum(mom, clus);
303 
304  return ApplyKinematicCuts(mom, rejectionReason);
305 }
306 
312 Bool_t AliClusterContainer::ApplyClusterCuts(const AliVCluster* clus, UInt_t &rejectionReason) const
313 {
314  if (!clus) {
315  rejectionReason |= kNullObject;
316  return kFALSE;
317  }
318 
319  if (!clus->IsEMCAL()) {
320  rejectionReason |= kIsEMCalCut;
321  return kFALSE;
322  }
323 
324  if (clus->TestBits(fBitMap) != (Int_t)fBitMap) {
325  rejectionReason |= kBitMapCut;
326  return kFALSE;
327  }
328 
329  if (fMinMCLabel >= 0 && TMath::Abs(clus->GetLabel()) < fMinMCLabel) {
330  rejectionReason |= kMCLabelCut;
331  return kFALSE;
332  }
333 
334  if (fMaxMCLabel >= 0 && TMath::Abs(clus->GetLabel()) > fMaxMCLabel) {
335  rejectionReason |= kMCLabelCut;
336  return kFALSE;
337  }
338 
339  if (clus->GetTOF() > fClusTimeCutUp || clus->GetTOF() < fClusTimeCutLow) {
340  rejectionReason |= kTimeCut;
341  return kFALSE;
342  }
343 
344  if (fExoticCut && clus->GetIsExotic()) {
345  rejectionReason |= kExoticCut;
346  return kFALSE;
347  }
348 
349  for (Int_t i = 0; i <= AliVCluster::kLastUserDefEnergy; i++) {
350  if (clus->GetUserDefEnergy((VCluUserDefEnergy_t)i) < fUserDefEnergyCut[i]) {
351  rejectionReason |= kEnergyCut;
352  return kFALSE;
353  }
354  }
355 
356  return kTRUE;
357 }
358 
364 {
365  UInt_t rejectionReason = 0;
366  Int_t nClus = 0;
367  for(int iclust = 0; iclust < this->fClArray->GetEntries(); ++iclust){
368  AliVCluster *clust = this->GetCluster(iclust);
369  if(this->AcceptCluster(clust, rejectionReason)) nClus++;
370  }
371  return nClus;
372 }
373 
380 {
381  if (t >= 0 && t <= AliVCluster::kLastUserDefEnergy){
382  return fUserDefEnergyCut[t];
383  }
384  else {
385  return fMinE;
386  }
387 }
388 
395 {
396  if (t >= 0 && t <= AliVCluster::kLastUserDefEnergy){
397  fUserDefEnergyCut[t] = cut;
398  }
399  else {
400  fMinE = cut;
401  }
402 }
403 
410  return AliClusterIterableContainer(this, false);
411 }
412 
419  return AliClusterIterableContainer(this, true);
420 }
421 
428  return AliClusterIterableMomentumContainer(this, false);
429 }
430 
437  return AliClusterIterableMomentumContainer(this, true);
438 }
439 
440 const char* AliClusterContainer::GetTitle() const
441 {
442  static TString clusterString;
443 
445 
446  if (Ecut == 0) {
447  clusterString = TString::Format("%s_E0000", GetArrayName().Data());
448  }
449  else if (Ecut < 1.0) {
450  clusterString = TString::Format("%s_E0%3.0f", GetArrayName().Data(), Ecut*1000.0);
451  }
452  else {
453  clusterString = TString::Format("%s_E%4.0f", GetArrayName().Data(), Ecut*1000.0);
454  }
455 
456  return clusterString.Data();
457 }
458 
459 
460 /******************************************
461  * Unit tests *
462  ******************************************/
463 
464 int TestClusterContainerIterator(const AliClusterContainer *const cont, int iteratorType, bool verbose){
465  std::vector<AliVCluster *> reference, variation;
466  AliVCluster *test = NULL;
467  for(int iclust = 0; iclust < cont->GetNClusters(); iclust++){
468  test = cont->GetCluster(iclust);
469  if(!iteratorType){
470  UInt_t rejectionReason = 0;
471  if(!cont->AcceptCluster(test, rejectionReason)) continue;
472  }
473  reference.push_back(test);
474  }
475 
476  if(!iteratorType){
477  // test accept iterator
478  for(auto cluster : cont->accepted()){
479  variation.push_back(cluster);
480  }
481  } else {
482  // test all iterator
483  for(auto cluster : cont->all()){
484  variation.push_back(cluster);
485  }
486  }
487 
488  if(verbose){
489  // Printing cluster addresses:
490  if(reference.size() < 30){
491  std::cout << "Clusters in reference container: " << std::endl;
492  std::cout << "===========================================" << std::endl;
493  for(std::vector<AliVCluster *>::iterator refit = reference.begin(); refit != reference.end(); ++refit){
494  std::cout << "Address: " << *refit << std::endl;
495  }
496  }
497  if(variation.size() < 30){
498  std::cout << "Clusters in test container: " << std::endl;
499  std::cout << "===========================================" << std::endl;
500  for(std::vector<AliVCluster *>::iterator varit = variation.begin(); varit != variation.end(); ++varit){
501  std::cout << "Address: " << *varit << std::endl;
502  }
503  }
504  }
505 
506  int testresult = 0;
507  // compare distributions: all clusters in one vector needs to be found in the other vector and vice versa
508  bool failure = false;
509  // first test: cleck if all clusters are found by the iterator
510  for(std::vector<AliVCluster *>::iterator clit = reference.begin(); clit != reference.end(); ++clit){
511  if(std::find(variation.begin(), variation.end(), *clit) == variation.end()) {
512  if(verbose)
513  std::cout << "Could not find cluster with address " << *clit << " in test container" << std::endl;
514  failure = true;
515  break;
516  }
517  }
518  if(!failure){
519  // second test: check if there are no extra clusters found
520  for(std::vector<AliVCluster *>::iterator clit = variation.begin(); clit != variation.end(); ++clit){
521  if(std::find(reference.begin(), reference.end(), *clit) == reference.end()) {
522  if(verbose)
523  std::cout << "Could not find cluster with address " << *clit << " in reference container" << std::endl;
524  failure = true;
525  break;
526  }
527  }
528  if(failure) testresult = 2;
529  } else {
530  testresult = 1;
531  }
532 
533  if(verbose){
534  std::cout << "Unit test cluster container, iterator type " << iteratorType << std::endl;
535  std::cout << "Number of expected clusters: " << reference.size() << ", number of found clusters: " << variation.size() << std::endl;
536  std::cout << "Test result: " << testresult << std::endl;
537  std::cout << "-----------------------------------------------------------------------" << std::endl;
538  }
539  return testresult;
540 }
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
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)
energy
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
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()