[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / AliJetConstituentTagCopier.cxx

// $Id: AliJetConstituentTagCopier.cxx  $
//
// Copy tags from particle level constituent to detector level
//
// Author: S. Aiola

#include "AliJetConstituentTagCopier.h"

#include <TClonesArray.h>
#include <TMath.h>
#include <TLorentzVector.h>

#include "AliNamedArrayI.h"
#include "AliVCluster.h"
#include "AliVParticle.h"
#include "AliEmcalParticle.h"
#include "AliLog.h"

ClassImp(AliJetConstituentTagCopier)

//________________________________________________________________________
AliJetConstituentTagCopier::AliJetConstituentTagCopier() : 
  AliAnalysisTaskEmcal("AliJetConstituentTagCopier", kFALSE),
  fMCParticlesName(),
  fCleanBeforeCopy(kFALSE),
  fMCParticles(0),
  fMCParticlesMap(0)
{
  // Default constructor.
}

//________________________________________________________________________
AliJetConstituentTagCopier::AliJetConstituentTagCopier(const char *name) : 
  AliAnalysisTaskEmcal(name, kFALSE),
  fMCParticlesName("MCParticles"),
  fCleanBeforeCopy(kFALSE),
  fMCParticles(0),
  fMCParticlesMap(0)
{
  // Standard constructor.
}

//________________________________________________________________________
AliJetConstituentTagCopier::~AliJetConstituentTagCopier()
{
  // Destructor
}

//________________________________________________________________________
void AliJetConstituentTagCopier::ExecOnce()
{
  // Execute once.

  if (!fMCParticles && !fMCParticlesName.IsNull()) {
    fMCParticles = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fMCParticlesName));
    if (!fMCParticles) {
      AliError(Form("%s: Could not retrieve MC particles %s!", GetName(), fMCParticlesName.Data()));
      return;
    }
    else if (!fMCParticles->GetClass()->GetBaseClass("AliVParticle")) {
      AliError(Form("%s: Collection %s does not contain AliVParticle objects!", GetName(), fMCParticlesName.Data())); 
      fMCParticles = 0;
      return;
    }
  }

  if (fMCParticles && !fMCParticlesMap) {
    fMCParticlesMap = dynamic_cast<AliNamedArrayI*>(InputEvent()->FindListObject(fMCParticlesName + "_Map"));
    // this is needed to map the MC labels with the indexes of the MC particle collection
    // if teh map is not given, the MC labels are assumed to be consistent with the indexes (which is not the case if AliEmcalMCTrackSelector is used)
    if (!fMCParticlesMap) {
      AliWarning(Form("%s: Could not retrieve map for MC particles %s! Will assume MC labels consistent with indexes...", GetName(), fMCParticlesName.Data())); 
      fMCParticlesMap = new AliNamedArrayI("tracksMap",9999);
      for (Int_t i = 0; i < 9999; i++) {
	fMCParticlesMap->AddAt(i,i);
      }
    }
  }

  AliAnalysisTaskEmcal::ExecOnce();
}

//________________________________________________________________________
Bool_t AliJetConstituentTagCopier::Run()
{
  if (fTracks) {
    if (fTracks->GetClass()->GetBaseClass("AliVParticle"))
      DoTrackLoop(fTracks);
    else if (fTracks->GetClass()->GetBaseClass("AliEmcalParticle"))
      DoEmcalParticleLoop(fTracks);
    else 
      AliError(Form("%s: Object type not recognized in collection %s. Nothing will be done.", GetName(), fTracks->GetName()));
  }

  if (fCaloClusters) {
    if (fCaloClusters->GetClass()->GetBaseClass("AliVCluster"))
      DoClusterLoop(fCaloClusters);
    else if (fCaloClusters->GetClass()->GetBaseClass("AliEmcalParticle"))
      DoEmcalParticleLoop(fCaloClusters);
    else 
      AliError(Form("%s: Object type not recognized in collection %s. Nothing will be done.", GetName(), fCaloClusters->GetName()));
  }

  return kTRUE;
}

//________________________________________________________________________
void AliJetConstituentTagCopier::DoClusterLoop(TClonesArray *array)
{
  if (fCleanBeforeCopy) {
    for (Int_t i = 0; i < array->GetEntries(); i++) {
      AliVCluster *cluster = static_cast<AliVCluster*>(array->At(i));
      if (!cluster) {
	AliError(Form("%s: Could not get cluster %d", GetName(), i));
	continue;
      }
      if (!AcceptCluster(cluster))
	continue;
      Int_t mcLabel = cluster->GetLabel();
      if (mcLabel > 0)
	cluster->SetBit(TObject::kBitMask, kFALSE);
    }
  }

  if (!fMCParticles)
    return;
  
  Double_t totalEnergy = 0;
  for (Int_t i = 0; i < array->GetEntries(); i++) {
    AliVCluster *cluster = static_cast<AliVCluster*>(array->At(i));
    if (!cluster) {
      AliError(Form("%s: Could not get cluster %d", GetName(), i));
      continue;
    }
    if (!AcceptCluster(cluster))
      continue;
    Int_t mcLabel = cluster->GetLabel();
    if (mcLabel > 0) {
      TLorentzVector vect;
      cluster->GetMomentum(vect, fVertex);
      AliDebug(2, Form("Cluster %d, pt = %f, eta = %f, phi = %f, label = %d",
		       i, cluster->E(), vect.Eta(), vect.Phi(), mcLabel));
      totalEnergy += cluster->E();
      Int_t index = -1;
      if (mcLabel < fMCParticlesMap->GetSize())
	index = fMCParticlesMap->At(mcLabel);
      if (index < 0)
	continue;
      AliVParticle *part = static_cast<AliVParticle*>(fMCParticles->At(index));
      if (!part) {
	AliError(Form("%s: Could not get MC particle %d", GetName(), index));
	continue;
      }      
      AliDebug(2, Form("Matched with particle %d, pt = %f, eta = %f, phi = %f", 
		       index, part->E(), part->Eta(), part->Phi()));
      UInt_t bits = (UInt_t)part->TestBits(TObject::kBitMask);
      cluster->SetBit(bits);
    }
  }

  AliDebug(2, Form("Total energy of MC clusters = %f", 
		   totalEnergy));
}

//________________________________________________________________________
void AliJetConstituentTagCopier::DoTrackLoop(TClonesArray *array)
{
  if (fCleanBeforeCopy) {
    for (Int_t i = 0; i < array->GetEntries(); i++) {
      AliVParticle *track = static_cast<AliVParticle*>(array->At(i));
      if (!track) {
	AliError(Form("%s: Could not get track %d", GetName(), i));
	continue;
      }
      if (!AcceptTrack(track))
      continue;
      Int_t mcLabel = TMath::Abs(track->GetLabel());
      if (mcLabel != 0) 
	track->SetBit(TObject::kBitMask, kFALSE);
    }
  }

  if (!fMCParticles)
    return;

  for (Int_t i = 0; i < array->GetEntries(); i++) {
    AliVParticle *track = static_cast<AliVParticle*>(array->At(i));
    if (!track) {
      AliError(Form("%s: Could not get track %d", GetName(), i));
      continue;
    }
    if (!AcceptTrack(track))
      continue;
    Int_t mcLabel = TMath::Abs(track->GetLabel());
    if (mcLabel != 0) {
      Int_t index = -1;
      if (mcLabel < fMCParticlesMap->GetSize())
	index = fMCParticlesMap->At(mcLabel);
      if (index < 0)
	continue;
      AliVParticle *part = static_cast<AliVParticle*>(fMCParticles->At(index));
      if (!part) {
	AliError(Form("%s: Could not get MC particle %d", GetName(), index));
	continue;
      }
      AliDebug(3, Form("Track %d, pt = %f, eta = %f, phi = %f, label = %d is matched with particle %d, pt = %f, eta = %f, phi = %f", 
		       i, track->Pt(), track->Eta(), track->Phi(), mcLabel, index, part->Pt(), part->Eta(), part->Phi()));
      UInt_t bits = (UInt_t)part->TestBits(TObject::kBitMask);
      track->SetBit(bits);
    }
  }
}

//________________________________________________________________________
void AliJetConstituentTagCopier::DoEmcalParticleLoop(TClonesArray *array)
{
  if (fCleanBeforeCopy) {
    for (Int_t i = 0; i < array->GetEntries(); i++) {
      AliEmcalParticle *emcpart = static_cast<AliEmcalParticle*>(array->At(i));
      if (!emcpart) {
	AliError(Form("%s: Could not get EmcalParticle %d", GetName(), i));
	continue;
      }
      if (!AcceptEmcalPart(emcpart))
	continue;
      AliVCluster *cluster = emcpart->GetCluster();
      AliVParticle *track = emcpart->GetTrack();
      Int_t mcLabel = 0;
      if (cluster)
	mcLabel = cluster->GetLabel();
      else if (track)
	mcLabel = TMath::Abs(track->GetLabel());
      if (mcLabel != 0) 
	emcpart->SetBit(TObject::kBitMask, kFALSE);
    }
  }

  if (!fMCParticles)
    return;

  for (Int_t i = 0; i < array->GetEntries(); i++) {
    AliEmcalParticle *emcpart = static_cast<AliEmcalParticle*>(array->At(i));
    if (!emcpart) {
      AliError(Form("%s: Could not get EmcalParticle %d", GetName(), i));
      continue;
    }
    if (!AcceptEmcalPart(emcpart))
      continue;
    AliVCluster *cluster = emcpart->GetCluster();
    AliVParticle *track = emcpart->GetTrack();
    Int_t mcLabel = 0;
    if (cluster)
      mcLabel = cluster->GetLabel();
    else if (track)
      mcLabel = TMath::Abs(track->GetLabel());
    if (mcLabel != 0) {
      Int_t index = -1;
      if (mcLabel < fMCParticlesMap->GetSize())
	index = fMCParticlesMap->At(mcLabel);
      if (index < 0)
	continue;
      AliVParticle *part = static_cast<AliVParticle*>(fMCParticles->At(index));
      if (!part) {
	AliError(Form("%s: Could not get MC particle %d", GetName(), index));
	continue;
      }
      UInt_t bits = (UInt_t)part->TestBits(TObject::kBitMask);
      emcpart->SetBit(bits);
    }
  }
}
Commit	Line	Data
8afc7c8d	1	// $Id: AliJetConstituentTagCopier.cxx $
	2	//
	3	// Copy tags from particle level constituent to detector level
	4	//
	5	// Author: S. Aiola
	6
	7	#include "AliJetConstituentTagCopier.h"
	8
	9	#include <TClonesArray.h>
5ce8ae64	10	#include <TMath.h>
4358e58a	11	#include <TLorentzVector.h>
8afc7c8d	12
5be3857d	13	#include "AliNamedArrayI.h"
8afc7c8d	14	#include "AliVCluster.h"
	15	#include "AliVParticle.h"
	16	#include "AliEmcalParticle.h"
	17	#include "AliLog.h"
	18
	19	ClassImp(AliJetConstituentTagCopier)
	20
	21	//________________________________________________________________________
	22	AliJetConstituentTagCopier::AliJetConstituentTagCopier() :
	23	AliAnalysisTaskEmcal("AliJetConstituentTagCopier", kFALSE),
	24	fMCParticlesName(),
f0c10df2	25	fCleanBeforeCopy(kFALSE),
8afc7c8d	26	fMCParticles(0),
	27	fMCParticlesMap(0)
	28	{
	29	// Default constructor.
	30	}
	31
	32	//________________________________________________________________________
	33	AliJetConstituentTagCopier::AliJetConstituentTagCopier(const char *name) :
	34	AliAnalysisTaskEmcal(name, kFALSE),
	35	fMCParticlesName("MCParticles"),
f0c10df2	36	fCleanBeforeCopy(kFALSE),
8afc7c8d	37	fMCParticles(0),
	38	fMCParticlesMap(0)
	39	{
	40	// Standard constructor.
	41	}
	42
	43	//________________________________________________________________________
	44	AliJetConstituentTagCopier::~AliJetConstituentTagCopier()
	45	{
	46	// Destructor
	47	}
	48
	49	//________________________________________________________________________
	50	void AliJetConstituentTagCopier::ExecOnce()
	51	{
	52	// Execute once.
	53
f0c10df2	54	if (!fMCParticles && !fMCParticlesName.IsNull()) {
8afc7c8d	55	fMCParticles = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fMCParticlesName));
	56	if (!fMCParticles) {
	57	AliError(Form("%s: Could not retrieve MC particles %s!", GetName(), fMCParticlesName.Data()));
	58	return;
	59	}
	60	else if (!fMCParticles->GetClass()->GetBaseClass("AliVParticle")) {
	61	AliError(Form("%s: Collection %s does not contain AliVParticle objects!", GetName(), fMCParticlesName.Data()));
	62	fMCParticles = 0;
	63	return;
	64	}
	65	}
	66
f0c10df2	67	if (fMCParticles && !fMCParticlesMap) {
5be3857d	68	fMCParticlesMap = dynamic_cast<AliNamedArrayI*>(InputEvent()->FindListObject(fMCParticlesName + "_Map"));
8afc7c8d	69	// this is needed to map the MC labels with the indexes of the MC particle collection
4358e58a	70	// if teh map is not given, the MC labels are assumed to be consistent with the indexes (which is not the case if AliEmcalMCTrackSelector is used)
8afc7c8d	71	if (!fMCParticlesMap) {
8afc7c8d	72	AliWarning(Form("%s: Could not retrieve map for MC particles %s! Will assume MC labels consistent with indexes...", GetName(), fMCParticlesName.Data()));
5be3857d	73	fMCParticlesMap = new AliNamedArrayI("tracksMap",9999);
8afc7c8d	74	for (Int_t i = 0; i < 9999; i++) {
5be3857d	75	fMCParticlesMap->AddAt(i,i);
8afc7c8d	76	}
	77	}
	78	}
	79
	80	AliAnalysisTaskEmcal::ExecOnce();
	81	}
	82
	83	//________________________________________________________________________
	84	Bool_t AliJetConstituentTagCopier::Run()
	85	{
	86	if (fTracks) {
	87	if (fTracks->GetClass()->GetBaseClass("AliVParticle"))
	88	DoTrackLoop(fTracks);
	89	else if (fTracks->GetClass()->GetBaseClass("AliEmcalParticle"))
	90	DoEmcalParticleLoop(fTracks);
	91	else
	92	AliError(Form("%s: Object type not recognized in collection %s. Nothing will be done.", GetName(), fTracks->GetName()));
	93	}
	94
	95	if (fCaloClusters) {
	96	if (fCaloClusters->GetClass()->GetBaseClass("AliVCluster"))
	97	DoClusterLoop(fCaloClusters);
	98	else if (fCaloClusters->GetClass()->GetBaseClass("AliEmcalParticle"))
	99	DoEmcalParticleLoop(fCaloClusters);
	100	else
	101	AliError(Form("%s: Object type not recognized in collection %s. Nothing will be done.", GetName(), fCaloClusters->GetName()));
	102	}
	103
	104	return kTRUE;
	105	}
	106
	107	//________________________________________________________________________
	108	void AliJetConstituentTagCopier::DoClusterLoop(TClonesArray *array)
	109	{
f0c10df2	110	if (fCleanBeforeCopy) {
	111	for (Int_t i = 0; i < array->GetEntries(); i++) {
	112	AliVCluster cluster = static_cast<AliVCluster>(array->At(i));
	113	if (!cluster) {
	114	AliError(Form("%s: Could not get cluster %d", GetName(), i));
	115	continue;
	116	}
	117	if (!AcceptCluster(cluster))
	118	continue;
	119	Int_t mcLabel = cluster->GetLabel();
	120	if (mcLabel > 0)
	121	cluster->SetBit(TObject::kBitMask, kFALSE);
	122	}
	123	}
	124
	125	if (!fMCParticles)
	126	return;
	127
4358e58a	128	Double_t totalEnergy = 0;
8afc7c8d	129	for (Int_t i = 0; i < array->GetEntries(); i++) {
	130	AliVCluster cluster = static_cast<AliVCluster>(array->At(i));
	131	if (!cluster) {
	132	AliError(Form("%s: Could not get cluster %d", GetName(), i));
	133	continue;
	134	}
	135	if (!AcceptCluster(cluster))
	136	continue;
	137	Int_t mcLabel = cluster->GetLabel();
	138	if (mcLabel > 0) {
4358e58a	139	TLorentzVector vect;
	140	cluster->GetMomentum(vect, fVertex);
	141	AliDebug(2, Form("Cluster %d, pt = %f, eta = %f, phi = %f, label = %d",
	142	i, cluster->E(), vect.Eta(), vect.Phi(), mcLabel));
	143	totalEnergy += cluster->E();
5be3857d	144	Int_t index = -1;
	145	if (mcLabel < fMCParticlesMap->GetSize())
	146	index = fMCParticlesMap->At(mcLabel);
5ce8ae64	147	if (index < 0)
5ce8ae64	148	continue;
8afc7c8d	149	AliVParticle part = static_cast<AliVParticle>(fMCParticles->At(index));
	150	if (!part) {
	151	AliError(Form("%s: Could not get MC particle %d", GetName(), index));
	152	continue;
4358e58a	153	}
	154	AliDebug(2, Form("Matched with particle %d, pt = %f, eta = %f, phi = %f",
	155	index, part->E(), part->Eta(), part->Phi()));
8afc7c8d	156	UInt_t bits = (UInt_t)part->TestBits(TObject::kBitMask);
	157	cluster->SetBit(bits);
	158	}
	159	}
4358e58a	160
	161	AliDebug(2, Form("Total energy of MC clusters = %f",
	162	totalEnergy));
8afc7c8d	163	}
	164
	165	//________________________________________________________________________
	166	void AliJetConstituentTagCopier::DoTrackLoop(TClonesArray *array)
	167	{
f0c10df2	168	if (fCleanBeforeCopy) {
	169	for (Int_t i = 0; i < array->GetEntries(); i++) {
	170	AliVParticle track = static_cast<AliVParticle>(array->At(i));
	171	if (!track) {
	172	AliError(Form("%s: Could not get track %d", GetName(), i));
	173	continue;
	174	}
	175	if (!AcceptTrack(track))
	176	continue;
	177	Int_t mcLabel = TMath::Abs(track->GetLabel());
	178	if (mcLabel != 0)
	179	track->SetBit(TObject::kBitMask, kFALSE);
	180	}
	181	}
	182
	183	if (!fMCParticles)
	184	return;
	185
8afc7c8d	186	for (Int_t i = 0; i < array->GetEntries(); i++) {
	187	AliVParticle track = static_cast<AliVParticle>(array->At(i));
	188	if (!track) {
	189	AliError(Form("%s: Could not get track %d", GetName(), i));
	190	continue;
	191	}
	192	if (!AcceptTrack(track))
	193	continue;
5ce8ae64	194	Int_t mcLabel = TMath::Abs(track->GetLabel());
8afc7c8d	195	if (mcLabel != 0) {
5be3857d	196	Int_t index = -1;
	197	if (mcLabel < fMCParticlesMap->GetSize())
	198	index = fMCParticlesMap->At(mcLabel);
5ce8ae64	199	if (index < 0)
5ce8ae64	200	continue;
8afc7c8d	201	AliVParticle part = static_cast<AliVParticle>(fMCParticles->At(index));
	202	if (!part) {
	203	AliError(Form("%s: Could not get MC particle %d", GetName(), index));
	204	continue;
	205	}
4358e58a	206	AliDebug(3, Form("Track %d, pt = %f, eta = %f, phi = %f, label = %d is matched with particle %d, pt = %f, eta = %f, phi = %f",
5ce8ae64	207	i, track->Pt(), track->Eta(), track->Phi(), mcLabel, index, part->Pt(), part->Eta(), part->Phi()));
8afc7c8d	208	UInt_t bits = (UInt_t)part->TestBits(TObject::kBitMask);
	209	track->SetBit(bits);
	210	}
	211	}
	212	}
	213
	214	//________________________________________________________________________
	215	void AliJetConstituentTagCopier::DoEmcalParticleLoop(TClonesArray *array)
	216	{
f0c10df2	217	if (fCleanBeforeCopy) {
	218	for (Int_t i = 0; i < array->GetEntries(); i++) {
	219	AliEmcalParticle emcpart = static_cast<AliEmcalParticle>(array->At(i));
	220	if (!emcpart) {
	221	AliError(Form("%s: Could not get EmcalParticle %d", GetName(), i));
	222	continue;
	223	}
	224	if (!AcceptEmcalPart(emcpart))
	225	continue;
	226	AliVCluster *cluster = emcpart->GetCluster();
	227	AliVParticle *track = emcpart->GetTrack();
	228	Int_t mcLabel = 0;
	229	if (cluster)
	230	mcLabel = cluster->GetLabel();
	231	else if (track)
	232	mcLabel = TMath::Abs(track->GetLabel());
	233	if (mcLabel != 0)
	234	emcpart->SetBit(TObject::kBitMask, kFALSE);
	235	}
	236	}
	237
	238	if (!fMCParticles)
	239	return;
	240
8afc7c8d	241	for (Int_t i = 0; i < array->GetEntries(); i++) {
	242	AliEmcalParticle emcpart = static_cast<AliEmcalParticle>(array->At(i));
	243	if (!emcpart) {
	244	AliError(Form("%s: Could not get EmcalParticle %d", GetName(), i));
	245	continue;
	246	}
	247	if (!AcceptEmcalPart(emcpart))
	248	continue;
	249	AliVCluster *cluster = emcpart->GetCluster();
	250	AliVParticle *track = emcpart->GetTrack();
	251	Int_t mcLabel = 0;
	252	if (cluster)
	253	mcLabel = cluster->GetLabel();
	254	else if (track)
5ce8ae64	255	mcLabel = TMath::Abs(track->GetLabel());
8afc7c8d	256	if (mcLabel != 0) {
5be3857d	257	Int_t index = -1;
	258	if (mcLabel < fMCParticlesMap->GetSize())
	259	index = fMCParticlesMap->At(mcLabel);
5ce8ae64	260	if (index < 0)
5ce8ae64	261	continue;
8afc7c8d	262	AliVParticle part = static_cast<AliVParticle>(fMCParticles->At(index));
	263	if (!part) {
	264	AliError(Form("%s: Could not get MC particle %d", GetName(), index));
	265	continue;
	266	}
	267	UInt_t bits = (UInt_t)part->TestBits(TObject::kBitMask);
	268	emcpart->SetBit(bits);
	269	}
	270	}
	271	}