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

// $Id$
//
// Emcal jet class.
//
// Author: C.Loizides

#include "AliEmcalJet.h"

#include "AliLog.h"
#include "Riostream.h"

ClassImp(AliEmcalJet)

//__________________________________________________________________________________________________
AliEmcalJet::AliEmcalJet() : 
  AliVParticle(), 
  fPt(0), 
  fEta(0), 
  fPhi(0), 
  fM(0), 
  fNEF(0),
  fArea(0),       
  fAreaEta(0),       
  fAreaPhi(0),       
  fAreaEmc(-1), 
  fAxisInEmcal(0), 
  fFlavourTagging(0),
  fMaxCPt(0), 
  fMaxNPt(0), 
  fMCPt(0),
  fNn(0), 
  fNch(0),        
  fPtEmc(0),
  fNEmc(0),
  fClusterIDs(),
  fTrackIDs(),
  fMatched(2),
  fMatchingType(0),
  fTaggedJet(0x0),
  fTagStatus(-1),
  fPtSub(0),
  fPtVectSub(0),
  fTriggers(0)
{
  // Constructor.

  fClosestJets[0] = 0;
  fClosestJets[1] = 0; 
  fClosestJetsDist[0] = 999; 
  fClosestJetsDist[1] = 999; 
}

//__________________________________________________________________________________________________
AliEmcalJet::AliEmcalJet(Double_t px, Double_t py, Double_t pz) : 
  AliVParticle(), 
  fPt(TMath::Sqrt(px*px+py*py)), 
  fEta(TMath::ASinH(pz/fPt)),
  fPhi(0), 
  fM(0), 
  fNEF(0), 
  fArea(0), 
  fAreaEta(0),       
  fAreaPhi(0),       
  fAreaEmc(-1), 
  fAxisInEmcal(0),
  fFlavourTagging(0),
  fMaxCPt(0), 
  fMaxNPt(0), 
  fMCPt(0),
  fNn(0),
  fNch(0),
  fPtEmc(0),
  fNEmc(0),
  fClusterIDs(), 
  fTrackIDs(),
  fMatched(2),
  fMatchingType(0),
  fTaggedJet(0x0),
  fTagStatus(-1),
  fPtSub(0),
  fPtVectSub(0),
  fTriggers(0)
{    
  // Constructor.

  if (fPt != 0) {
    fPhi = TMath::ATan2(py, px);
    if (fPhi<0.) 
      fPhi += 2. * TMath::Pi();
  }

  fClosestJets[0] = 0; 
  fClosestJets[1] = 0;
  fClosestJetsDist[0] = 999; 
  fClosestJetsDist[1] = 999;
}

//_________________________________________________________________________________________________
AliEmcalJet::AliEmcalJet(Double_t pt, Double_t eta, Double_t phi, Double_t m) :
  AliVParticle(), 
  fPt(pt), 
  fEta(eta), 
  fPhi(phi), 
  fM(m), 
  fNEF(0), 
  fArea(0), 
  fAreaEta(0),       
  fAreaPhi(0),       
  fAreaEmc(-1), 
  fAxisInEmcal(0),
  fFlavourTagging(0),
  fMaxCPt(0), 
  fMaxNPt(0),
  fMCPt(0),
  fNn(0),
  fNch(0), 
  fPtEmc(0),
  fNEmc(0),
  fClusterIDs(), 
  fTrackIDs(),
  fMatched(2),
  fMatchingType(0),
  fTaggedJet(0x0),
  fTagStatus(-1),
  fPtSub(0),
  fPtVectSub(0),
  fTriggers(0)
{
  // Constructor.

  if (fPhi<0.) 
    fPhi += TMath::TwoPi();

  fClosestJets[0] = 0; 
  fClosestJets[1] = 0;
  fClosestJetsDist[0] = 999; 
  fClosestJetsDist[1] = 999;
}

//_________________________________________________________________________________________________
AliEmcalJet::AliEmcalJet(const AliEmcalJet &jet) :
  AliVParticle(jet),
  fPt(jet.fPt), 
  fEta(jet.fEta), 
  fPhi(jet.fPhi), 
  fM(jet.fM), 
  fNEF(jet.fNEF), 
  fArea(jet.fArea), 
  fAreaEta(jet.fAreaEta),       
  fAreaPhi(jet.fAreaPhi),       
  fAreaEmc(jet.fAreaEmc), 
  fAxisInEmcal(jet.fAxisInEmcal),
  fFlavourTagging(jet.fFlavourTagging),
  fMaxCPt(jet.fMaxCPt), 
  fMaxNPt(jet.fMaxNPt), 
  fMCPt(jet.fMCPt),
  fNn(jet.fNn),
  fNch(jet.fNch),
  fPtEmc(jet.fPtEmc),
  fNEmc(jet.fNEmc),
  fClusterIDs(jet.fClusterIDs), 
  fTrackIDs(jet.fTrackIDs),
  fMatched(jet.fMatched),
  fMatchingType(jet.fMatchingType),
  fTaggedJet(jet.fTaggedJet),
  fTagStatus(jet.fTagStatus),
  fPtSub(jet.fPtSub),
  fPtVectSub(jet.fPtVectSub),
  fTriggers(jet.fTriggers)
{
  // Copy constructor.

  fClosestJets[0]     = jet.fClosestJets[0]; 
  fClosestJets[1]     = jet.fClosestJets[1]; 
  fClosestJetsDist[0] = jet.fClosestJetsDist[0];  
  fClosestJetsDist[1] = jet.fClosestJetsDist[1]; 
}

//_________________________________________________________________________________________________
AliEmcalJet &AliEmcalJet::operator=(const AliEmcalJet &jet)
{
  // Assignment operator.

  if (this!=&jet) {
    AliVParticle::operator=(jet);
    fPt                 = jet.fPt;
    fEta                = jet.fEta;
    fPhi                = jet.fPhi;
    fM                  = jet.fM; 
    fNEF                = jet.fNEF;
    fArea               = jet.fArea; 
    fAreaEta            = jet.fAreaEta; 
    fAreaPhi            = jet.fAreaPhi; 
    fAreaEmc            = jet.fAreaEmc; 
    fAxisInEmcal        = jet.fAxisInEmcal; 
    fFlavourTagging     = jet.fFlavourTagging;
    fMaxCPt             = jet.fMaxCPt; 
    fMaxNPt             = jet.fMaxNPt;
    fMCPt               = jet.fMCPt;
    fNn                 = jet.fNn;
    fNch                = jet.fNch;
    fPtEmc              = jet.fPtEmc;
    fNEmc               = jet.fNEmc;
    fClusterIDs         = jet.fClusterIDs;
    fTrackIDs           = jet.fTrackIDs;
    fClosestJets[0]     = jet.fClosestJets[0]; 
    fClosestJets[1]     = jet.fClosestJets[1]; 
    fClosestJetsDist[0] = jet.fClosestJetsDist[0];  
    fClosestJetsDist[1] = jet.fClosestJetsDist[1]; 
    fMatched            = jet.fMatched;
    fTaggedJet          = jet.fTaggedJet;
    fTagStatus          = jet.fTagStatus;
    fPtSub              = jet.fPtSub;
    fPtVectSub          = jet.fPtVectSub;
    fTriggers           = jet.fTriggers;
  }

  return *this;
}

//_________________________________________________________________________________________________
Int_t AliEmcalJet::Compare(const TObject* obj) const
{
  //Return -1 if this is smaller than obj, 0 if objects are equal and 1 if this is larger than obj.

  const AliEmcalJet *jet = static_cast<const AliEmcalJet *>(obj);
  if (!obj)
    return 0;
  if (Pt()>jet->Pt())
    return -1;
  return 1;
}

//__________________________________________________________________________________________________
void AliEmcalJet::GetMom(TLorentzVector &vec) const
{
  // Return momentum as four vector.

  Double_t p = fPt *TMath::CosH(fEta);
  vec.SetPtEtaPhiE(fPt,fEta,fPhi,TMath::Sqrt(p*p+fM*fM));
}

//__________________________________________________________________________________________________
void AliEmcalJet::Print(Option_t* /*option*/) const
{
  // Print jet information.

  Printf("Jet pt=%.2f, eta=%.2f, phi=%.2f, area=%.2f, NEF=%.2f", fPt, fEta, fPhi, fArea, fNEF);
}

//__________________________________________________________________________________________________
Double_t AliEmcalJet::PtSubVect(Double_t rho) const
{
  // Return vectorial subtracted transverse momentum.

  Double_t dx = Px() - rho * fArea * TMath::Cos(fAreaPhi);
  Double_t dy = Py() - rho * fArea * TMath::Sin(fAreaPhi);
  //Double_t dz = Pz() - rho * fArea * TMath::SinH(fAreaEta);
  return TMath::Sqrt(dx*dx+dy*dy);
}

//__________________________________________________________________________________________________
void AliEmcalJet::SortConstituents()
{
  // Sort constituent by index (increasing).

  std::sort(fClusterIDs.GetArray(), fClusterIDs.GetArray() + fClusterIDs.GetSize());
  std::sort(fTrackIDs.GetArray(), fTrackIDs.GetArray() + fTrackIDs.GetSize());
}

//__________________________________________________________________________________________________
AliVParticle* AliEmcalJet::GetLeadingTrack(TClonesArray *tracks) const
{
  AliVParticle* maxTrack = 0;
  for (Int_t i = 0; i < GetNumberOfTracks(); i++) {
    AliVParticle *track = TrackAt(i, tracks);
    if (!track) {
      AliError(Form("Unable to find jet track %d in collection %s (pos in collection %d, max %d)",
		    i,tracks->GetName(),TrackAt(i),tracks->GetEntriesFast()));
      continue;
    }
    if (!maxTrack || track->Pt() > maxTrack->Pt()) 
      maxTrack = track;
  }

  return maxTrack;
}

//__________________________________________________________________________________________________
AliVCluster* AliEmcalJet::GetLeadingCluster(TClonesArray *clusters) const
{
  AliVCluster* maxCluster = 0;
  for (Int_t i = 0; i < GetNumberOfClusters(); i++) {
    AliVCluster *cluster = ClusterAt(i, clusters);
    if (!cluster) {
      AliError(Form("Unable to find jet cluster %d in collection %s (pos in collection %d, max %d)",
		    i,clusters->GetName(),ClusterAt(i),clusters->GetEntriesFast()));
      continue;
    }
    if (!maxCluster || cluster->E() > maxCluster->E()) 
      maxCluster = cluster;
  }

  return maxCluster;
}

//__________________________________________________________________________________________________
void AliEmcalJet::ResetMatching()
{
  fClosestJets[0] = 0;
  fClosestJets[1] = 0; 
  fClosestJetsDist[0] = 999; 
  fClosestJetsDist[1] = 999; 
  fMatched = 2;
}
Commit	Line	Data
7df864a3	1	// $Id$
04f6e9de	2	//
	3	// Emcal jet class.
	4	//
cd231d42	5	// Author: C.Loizides
7efbea04	6
914d486c	7	#include "AliEmcalJet.h"
7df864a3	8
59ebe296	9	#include "AliLog.h"
	10	#include "Riostream.h"
	11
914d486c	12	ClassImp(AliEmcalJet)
7efbea04	13
96919f91	14	//__________________________________________________________________________________________________
	15	AliEmcalJet::AliEmcalJet() :
	16	AliVParticle(),
	17	fPt(0),
	18	fEta(0),
	19	fPhi(0),
	20	fM(0),
	21	fNEF(0),
	22	fArea(0),
db8a0184	23	fAreaEta(0),
db8a0184	24	fAreaPhi(0),
96919f91	25	fAreaEmc(-1),
96919f91	26	fAxisInEmcal(0),
3c9775d9	27	fFlavourTagging(0),
96919f91	28	fMaxCPt(0),
	29	fMaxNPt(0),
	30	fMCPt(0),
	31	fNn(0),
	32	fNch(0),
8895e61b	33	fPtEmc(0),
8895e61b	34	fNEmc(0),
96919f91	35	fClusterIDs(),
96919f91	36	fTrackIDs(),
8895e61b	37	fMatched(2),
cfc2ac24	38	fMatchingType(0),
254370e1	39	fTaggedJet(0x0),
254370e1	40	fTagStatus(-1),
db8a0184	41	fPtSub(0),
85d48773	42	fPtVectSub(0),
85d48773	43	fTriggers(0)
96919f91	44	{
	45	// Constructor.
	46
	47	fClosestJets[0] = 0;
	48	fClosestJets[1] = 0;
	49	fClosestJetsDist[0] = 999;
	50	fClosestJetsDist[1] = 999;
96919f91	51	}
96919f91	52
7efbea04	53	//__________________________________________________________________________________________________
101cefde	54	AliEmcalJet::AliEmcalJet(Double_t px, Double_t py, Double_t pz) :
	55	AliVParticle(),
	56	fPt(TMath::Sqrt(pxpx+pypy)),
	57	fEta(TMath::ASinH(pz/fPt)),
	58	fPhi(0),
	59	fM(0),
	60	fNEF(0),
	61	fArea(0),
db8a0184	62	fAreaEta(0),
db8a0184	63	fAreaPhi(0),
101cefde	64	fAreaEmc(-1),
101cefde	65	fAxisInEmcal(0),
3c9775d9	66	fFlavourTagging(0),
101cefde	67	fMaxCPt(0),
101cefde	68	fMaxNPt(0),
96919f91	69	fMCPt(0),
	70	fNn(0),
	71	fNch(0),
8895e61b	72	fPtEmc(0),
8895e61b	73	fNEmc(0),
101cefde	74	fClusterIDs(),
96919f91	75	fTrackIDs(),
8895e61b	76	fMatched(2),
cfc2ac24	77	fMatchingType(0),
254370e1	78	fTaggedJet(0x0),
254370e1	79	fTagStatus(-1),
db8a0184	80	fPtSub(0),
85d48773	81	fPtVectSub(0),
85d48773	82	fTriggers(0)
7efbea04	83	{
	84	// Constructor.
	85
	86	if (fPt != 0) {
	87	fPhi = TMath::ATan2(py, px);
	88	if (fPhi<0.)
	89	fPhi += 2. * TMath::Pi();
	90	}
96919f91	91
	92	fClosestJets[0] = 0;
	93	fClosestJets[1] = 0;
	94	fClosestJetsDist[0] = 999;
	95	fClosestJetsDist[1] = 999;
7efbea04	96	}
7efbea04	97
9cffd2d0	98	//_________________________________________________________________________________________________
914d486c	99	AliEmcalJet::AliEmcalJet(Double_t pt, Double_t eta, Double_t phi, Double_t m) :
9cffd2d0	100	AliVParticle(),
f8087a81	101	fPt(pt),
	102	fEta(eta),
	103	fPhi(phi),
	104	fM(m),
	105	fNEF(0),
	106	fArea(0),
db8a0184	107	fAreaEta(0),
db8a0184	108	fAreaPhi(0),
101cefde	109	fAreaEmc(-1),
101cefde	110	fAxisInEmcal(0),
3c9775d9	111	fFlavourTagging(0),
f8087a81	112	fMaxCPt(0),
96919f91	113	fMaxNPt(0),
	114	fMCPt(0),
	115	fNn(0),
	116	fNch(0),
8895e61b	117	fPtEmc(0),
8895e61b	118	fNEmc(0),
f8087a81	119	fClusterIDs(),
96919f91	120	fTrackIDs(),
8895e61b	121	fMatched(2),
cfc2ac24	122	fMatchingType(0),
254370e1	123	fTaggedJet(0x0),
254370e1	124	fTagStatus(-1),
db8a0184	125	fPtSub(0),
85d48773	126	fPtVectSub(0),
85d48773	127	fTriggers(0)
9cffd2d0	128	{
	129	// Constructor.
	130
4f149872	131	if (fPhi<0.)
4f149872	132	fPhi += TMath::TwoPi();
96919f91	133
	134	fClosestJets[0] = 0;
	135	fClosestJets[1] = 0;
	136	fClosestJetsDist[0] = 999;
	137	fClosestJetsDist[1] = 999;
9cffd2d0	138	}
9cffd2d0	139
7efbea04	140	//_________________________________________________________________________________________________
914d486c	141	AliEmcalJet::AliEmcalJet(const AliEmcalJet &jet) :
7efbea04	142	AliVParticle(jet),
f8087a81	143	fPt(jet.fPt),
	144	fEta(jet.fEta),
	145	fPhi(jet.fPhi),
	146	fM(jet.fM),
	147	fNEF(jet.fNEF),
	148	fArea(jet.fArea),
db8a0184	149	fAreaEta(jet.fAreaEta),
db8a0184	150	fAreaPhi(jet.fAreaPhi),
101cefde	151	fAreaEmc(jet.fAreaEmc),
101cefde	152	fAxisInEmcal(jet.fAxisInEmcal),
3c9775d9	153	fFlavourTagging(jet.fFlavourTagging),
f8087a81	154	fMaxCPt(jet.fMaxCPt),
f8087a81	155	fMaxNPt(jet.fMaxNPt),
96919f91	156	fMCPt(jet.fMCPt),
	157	fNn(jet.fNn),
	158	fNch(jet.fNch),
8895e61b	159	fPtEmc(jet.fPtEmc),
8895e61b	160	fNEmc(jet.fNEmc),
f8087a81	161	fClusterIDs(jet.fClusterIDs),
4f149872	162	fTrackIDs(jet.fTrackIDs),
8895e61b	163	fMatched(jet.fMatched),
cfc2ac24	164	fMatchingType(jet.fMatchingType),
254370e1	165	fTaggedJet(jet.fTaggedJet),
254370e1	166	fTagStatus(jet.fTagStatus),
db8a0184	167	fPtSub(jet.fPtSub),
85d48773	168	fPtVectSub(jet.fPtVectSub),
85d48773	169	fTriggers(jet.fTriggers)
7df864a3	170	{
101cefde	171	// Copy constructor.
4f149872	172
	173	fClosestJets[0] = jet.fClosestJets[0];
	174	fClosestJets[1] = jet.fClosestJets[1];
	175	fClosestJetsDist[0] = jet.fClosestJetsDist[0];
	176	fClosestJetsDist[1] = jet.fClosestJetsDist[1];
7efbea04	177	}
	178
	179	//_________________________________________________________________________________________________
914d486c	180	AliEmcalJet &AliEmcalJet::operator=(const AliEmcalJet &jet)
7efbea04	181	{
	182	// Assignment operator.
	183
	184	if (this!=&jet) {
	185	AliVParticle::operator=(jet);
4f149872	186	fPt = jet.fPt;
	187	fEta = jet.fEta;
	188	fPhi = jet.fPhi;
	189	fM = jet.fM;
	190	fNEF = jet.fNEF;
	191	fArea = jet.fArea;
db8a0184	192	fAreaEta = jet.fAreaEta;
db8a0184	193	fAreaPhi = jet.fAreaPhi;
4f149872	194	fAreaEmc = jet.fAreaEmc;
4f149872	195	fAxisInEmcal = jet.fAxisInEmcal;
3c9775d9	196	fFlavourTagging = jet.fFlavourTagging;
4f149872	197	fMaxCPt = jet.fMaxCPt;
	198	fMaxNPt = jet.fMaxNPt;
	199	fMCPt = jet.fMCPt;
	200	fNn = jet.fNn;
	201	fNch = jet.fNch;
8895e61b	202	fPtEmc = jet.fPtEmc;
8895e61b	203	fNEmc = jet.fNEmc;
4f149872	204	fClusterIDs = jet.fClusterIDs;
	205	fTrackIDs = jet.fTrackIDs;
	206	fClosestJets[0] = jet.fClosestJets[0];
	207	fClosestJets[1] = jet.fClosestJets[1];
	208	fClosestJetsDist[0] = jet.fClosestJetsDist[0];
	209	fClosestJetsDist[1] = jet.fClosestJetsDist[1];
	210	fMatched = jet.fMatched;
254370e1	211	fTaggedJet = jet.fTaggedJet;
254370e1	212	fTagStatus = jet.fTagStatus;
8895e61b	213	fPtSub = jet.fPtSub;
f660c2d6	214	fPtVectSub = jet.fPtVectSub;
85d48773	215	fTriggers = jet.fTriggers;
7efbea04	216	}
	217
	218	return *this;
7df864a3	219	}
7efbea04	220
2103dc6a	221	//_________________________________________________________________________________________________
db8a0184	222	Int_t AliEmcalJet::Compare(const TObject* obj) const
	223	{
	224	//Return -1 if this is smaller than obj, 0 if objects are equal and 1 if this is larger than obj.
	225
	226	const AliEmcalJet jet = static_cast<const AliEmcalJet >(obj);
	227	if (!obj)
	228	return 0;
	229	if (Pt()>jet->Pt())
	230	return -1;
	231	return 1;
	232	}
	233
7efbea04	234	//__________________________________________________________________________________________________
914d486c	235	void AliEmcalJet::GetMom(TLorentzVector &vec) const
7efbea04	236	{
	237	// Return momentum as four vector.
	238
	239	Double_t p = fPt *TMath::CosH(fEta);
	240	vec.SetPtEtaPhiE(fPt,fEta,fPhi,TMath::Sqrt(pp+fMfM));
	241	}
	242
	243	//__________________________________________________________________________________________________
914d486c	244	void AliEmcalJet::Print(Option_t* /option/) const
7efbea04	245	{
	246	// Print jet information.
	247
59ebe296	248	Printf("Jet pt=%.2f, eta=%.2f, phi=%.2f, area=%.2f, NEF=%.2f", fPt, fEta, fPhi, fArea, fNEF);
7efbea04	249	}
629e03fc	250
db8a0184	251	//__________________________________________________________________________________________________
	252	Double_t AliEmcalJet::PtSubVect(Double_t rho) const
	253	{
	254	// Return vectorial subtracted transverse momentum.
	255
	256	Double_t dx = Px() - rho * fArea * TMath::Cos(fAreaPhi);
	257	Double_t dy = Py() - rho * fArea * TMath::Sin(fAreaPhi);
	258	//Double_t dz = Pz() - rho * fArea * TMath::SinH(fAreaEta);
	259	return TMath::Sqrt(dxdx+dydy);
	260	}
	261
629e03fc	262	//__________________________________________________________________________________________________
	263	void AliEmcalJet::SortConstituents()
	264	{
f8087a81	265	// Sort constituent by index (increasing).
	266
	267	std::sort(fClusterIDs.GetArray(), fClusterIDs.GetArray() + fClusterIDs.GetSize());
	268	std::sort(fTrackIDs.GetArray(), fTrackIDs.GetArray() + fTrackIDs.GetSize());
629e03fc	269	}
f660c2d6	270
	271	//__________________________________________________________________________________________________
	272	AliVParticle* AliEmcalJet::GetLeadingTrack(TClonesArray *tracks) const
	273	{
	274	AliVParticle* maxTrack = 0;
	275	for (Int_t i = 0; i < GetNumberOfTracks(); i++) {
	276	AliVParticle *track = TrackAt(i, tracks);
59ebe296	277	if (!track) {
3c9775d9	278	AliError(Form("Unable to find jet track %d in collection %s (pos in collection %d, max %d)",
3c9775d9	279	i,tracks->GetName(),TrackAt(i),tracks->GetEntriesFast()));
59ebe296	280	continue;
59ebe296	281	}
f660c2d6	282	if (!maxTrack \|\| track->Pt() > maxTrack->Pt())
	283	maxTrack = track;
	284	}
	285
	286	return maxTrack;
	287	}
	288
	289	//__________________________________________________________________________________________________
	290	AliVCluster* AliEmcalJet::GetLeadingCluster(TClonesArray *clusters) const
	291	{
	292	AliVCluster* maxCluster = 0;
	293	for (Int_t i = 0; i < GetNumberOfClusters(); i++) {
	294	AliVCluster *cluster = ClusterAt(i, clusters);
59ebe296	295	if (!cluster) {
3c9775d9	296	AliError(Form("Unable to find jet cluster %d in collection %s (pos in collection %d, max %d)",
3c9775d9	297	i,clusters->GetName(),ClusterAt(i),clusters->GetEntriesFast()));
59ebe296	298	continue;
59ebe296	299	}
f660c2d6	300	if (!maxCluster \|\| cluster->E() > maxCluster->E())
	301	maxCluster = cluster;
	302	}
	303
	304	return maxCluster;
	305	}
2103dc6a	306
	307	//__________________________________________________________________________________________________
	308	void AliEmcalJet::ResetMatching()
	309	{
	310	fClosestJets[0] = 0;
	311	fClosestJets[1] = 0;
	312	fClosestJetsDist[0] = 999;
	313	fClosestJetsDist[1] = 999;
	314	fMatched = 2;
	315	}