[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / AliRsnLoopEffPair.cxx

//
// Class AliRsnLoopEffPair
//
// Inherits from basic AliRsnLoopEff for efficiency,
// and computed efficiencies for single-tracks
//
// author: Alberto Pulvirenti (alberto.pulvirenti@ct.infn.it)
//

#include <Riostream.h>

#include "AliStack.h"
#include "AliGenEventHeader.h"
#include "AliAODMCHeader.h"

#include "AliRsnEvent.h"
#include "AliRsnCutManager.h"
#include "AliRsnDaughterDef.h"
#include "AliRsnPairDef.h"
#include "AliRsnLoopEffPair.h"

ClassImp(AliRsnLoopEffPair)

//_____________________________________________________________________________
AliRsnLoopEffPair::AliRsnLoopEffPair(const char *name, AliRsnPairDef *def) :
   AliRsnLoopEff(name, 1),
   fDef(def),
   fMother()
{
//
// Default constructor.
// Do not repeat 'DefineOutput' since it is done in base class and we don't add new ones.
//
}

//_____________________________________________________________________________
AliRsnLoopEffPair::AliRsnLoopEffPair(const AliRsnLoopEffPair& copy) :
   AliRsnLoopEff(copy),
   fDef(copy.fDef),
   fMother(copy.fMother)
{
//
// Copy constrtuctor.
//
}

//_____________________________________________________________________________
AliRsnLoopEffPair& AliRsnLoopEffPair::operator=(const AliRsnLoopEffPair& copy)
{
//
// Assignment operator.
// Owned data members are meaningless for this operator.
//

   AliRsnLoopEff::operator=(copy);
   fDef = copy.fDef;

   return (*this);
}

//__________________________________________________________________________________________________
Bool_t AliRsnLoopEffPair::AssignMotherAndDaughters(AliRsnEvent *rsnEvent, Int_t ipart)
{
//
// Calls the appropriate assignment method
//

   // setup pointers
   fMother.SetDaughter(0, &fDaughter[0]);
   fMother.SetDaughter(1, &fDaughter[1]);
   fMother.SetRefEvent(rsnEvent);
   fDaughter[0].SetOwnerEvent(rsnEvent);
   fDaughter[1].SetOwnerEvent(rsnEvent);

   if (rsnEvent->IsESD())
      return AssignMotherAndDaughtersESD(rsnEvent, ipart);
   else if (rsnEvent->IsAOD())
      return AssignMotherAndDaughtersAOD(rsnEvent, ipart);
   else {
      AliError("Unrecognized input event");
      return kFALSE;
   }
}

//__________________________________________________________________________________________________
Bool_t AliRsnLoopEffPair::AssignMotherAndDaughtersESD(AliRsnEvent *rsnEvent, Int_t ipart)
{
//
// Gets a particle in the MC event and try to assign it to the mother.
// If it has two daughters, retrieve them and assign also them.
// NOTE: assignment is done only for MC, since reconstructed match is assigned in the same way
//       for ESD and AOD, if available
// ---
// Implementation for ESD inputs
//

   AliMCEvent    *mc      = rsnEvent->GetRefMCESD();
   AliStack      *stack   = mc->Stack();
   AliMCParticle *mother  = (AliMCParticle*)mc->GetTrack(ipart);
   TParticle     *motherP = mother->Particle();
   Int_t          ntracks = stack->GetNtrack();
   
   // check PDG code and exit if it is wrong
   if (TMath::Abs(motherP->GetPdgCode()) != fDef->GetMotherPDG()) return kFALSE;
   
   // check number of daughters and exit if it is not 2
   if (motherP->GetNDaughters() < 2) return kFALSE;
   
   // check distance from primary vertex
   TLorentzVector vprod;
   motherP->ProductionVertex(vprod);
   if (!CheckDistanceFromPV(vprod.X(), vprod.Y(), vprod.Z())) {
      AliDebugClass(1, "Distant production vertex");
      return kFALSE;
   }
   
   // get the daughters and check their PDG code and charge:
   // if they match one of the pair daughter definitions, 
   // assign them as MC reference of the 'fDaughter' objects
   fDaughter[0].Reset();
   fDaughter[1].Reset();
   Int_t index[2] = {motherP->GetFirstDaughter(), motherP->GetLastDaughter()};
   Int_t i, pdg;
   Short_t charge;
   AliMCParticle *daughter = 0x0;
   for (i = 0; i < 2; i++) {
      // check index for stack
      if (index[i] < 0 || index[i] > ntracks) {
         AliError(Form("Index %d overflow: value = %d, max = %d", i, index[i], ntracks));
         return kFALSE;
      }
      // get daughter and its PDG and charge
      daughter = (AliMCParticle*)mc->GetTrack(index[i]);
      pdg      = TMath::Abs(daughter->Particle()->GetPdgCode());
      charge   = (Short_t)(daughter->Particle()->GetPDG()->Charge() / 3);
      // check if it matches one definition
      if (fDef->GetDef1().MatchesPDG(pdg) && fDef->GetDef1().MatchesChargeS(charge)) {
         fDaughter[0].SetGood();
         fDaughter[0].SetRefMC(daughter);
         fDaughter[0].SetLabel(index[i]);
      } else if (fDef->GetDef2().MatchesPDG(pdg) && fDef->GetDef2().MatchesChargeS(charge)) {
         fDaughter[1].SetGood();
         fDaughter[1].SetRefMC(daughter);
         fDaughter[1].SetLabel(index[i]);
      }
   }
   
   // return success if both daughters were assigned
   if (fDaughter[0].IsOK() && fDaughter[1].IsOK()) {
      return kTRUE;
   } else {
      fDaughter[0].Reset();
      fDaughter[1].Reset();
      return kFALSE;
   }
}

//__________________________________________________________________________________________________
Bool_t AliRsnLoopEffPair::AssignMotherAndDaughtersAOD(AliRsnEvent *rsnEvent, Int_t ipart)
{
   AliAODEvent      *aod     = rsnEvent->GetRefAOD();
   TClonesArray     *listAOD = (TClonesArray*)(aod->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
   AliAODMCParticle *mother  = (AliAODMCParticle*)listAOD->At(ipart);
   Int_t             ntracks = listAOD->GetEntries();
   
   // check PDG code and exit if it is wrong
   if (TMath::Abs(mother->GetPdgCode()) != fDef->GetMotherPDG()) return kFALSE;
   
   // check number of daughters and exit if it is not 2
   if (mother->GetNDaughters() < 2) return kFALSE;
   
   // check distance from primary vertex
   Double_t vprod[3] = {(Double_t)mother->Xv(), (Double_t)mother->Yv(), (Double_t)mother->Zv()};
   Double_t dv = DistanceFromPV(vprod[0], vprod[1], vprod[2]);
   if (dv > fMaxDistPV) {
      AliDebugClass(1, "Distant production vertex");
      return kFALSE;
   }
   
   // get the daughters and check their PDG code and charge:
   // if they match one of the pair daughter definitions, 
   // assign them as MC reference of the 'fDaughter' objects
   fDaughter[0].Reset();
   fDaughter[1].Reset();
   Int_t index[2] = {(Int_t)mother->GetDaughter(0), (Int_t)mother->GetDaughter(1)};
   Int_t i, pdg;
   Short_t charge;
   AliAODMCParticle *daughter = 0x0;
   for (i = 0; i < 2; i++) {
      // check index for stack
      if (index[i] < 0 || index[i] > ntracks) {
         AliError(Form("Index %d overflow: value = %d, max = %d", i, index[i], ntracks));
         return kFALSE;
      }
      // get daughter and its PDG and charge
      daughter = (AliAODMCParticle*)listAOD->At(index[i]);
      pdg      = TMath::Abs(daughter->GetPdgCode());
      charge   = (Short_t)(daughter->Charge() / 3);
      // check if it matches one definition
      if (fDef->GetDef1().MatchesPDG(pdg) && fDef->GetDef1().MatchesChargeS(charge)) {
         fDaughter[0].SetGood();
         fDaughter[0].SetRefMC(daughter);
         fDaughter[0].SetLabel(index[i]);
      } else if (fDef->GetDef2().MatchesPDG(pdg) && fDef->GetDef2().MatchesChargeS(charge)) {
         fDaughter[1].SetGood();
         fDaughter[1].SetRefMC(daughter);
         fDaughter[1].SetLabel(index[i]);
      }
   }
   
   // return success if both daughters were assigned
   if (fDaughter[0].IsOK() && fDaughter[1].IsOK()) {
      return kTRUE;
   } else {
      fDaughter[0].Reset();
      fDaughter[1].Reset();
      return kFALSE;
   }
}

//_____________________________________________________________________________
Int_t AliRsnLoopEffPair::DoLoop(AliRsnEvent *rsn, AliRsnDaughterSelector*, AliRsnEvent*, AliRsnDaughterSelector*)
{
//
// Loop on event and fill containers
//

   // check event cuts
   if (!OkEvent(rsn)) return 0;
   
   // retrieve output
   fOutput = (AliRsnListOutput*)fOutputs[0];
   
   // check presence of MC reference
   if (!rsn->GetRefMC()) {
      AliError("Need a MC to compute efficiency");
      return 0;
   }
   
   // check presence of event
   if (!rsn->GetRef()) {
      AliError("Need an event to compute efficiency");
      return 0;
   }
   
   // check event type:
   // must be ESD or AOD, and then use a bool to know in the rest
   if (!rsn->IsESD() && !rsn->IsAOD()) {
      AliError("Need to process ESD or AOD input");
      return 0;
   }
   
   // retrieve the MC primary vertex position
   // and do some additional coherence checks
   Int_t i, npart = 0;
   if (rsn->IsESD()) {
      TArrayF mcVertex(3);
      AliGenEventHeader *genEH = rsn->GetRefMCESD()->GenEventHeader();
      genEH->PrimaryVertex(mcVertex);
      for (i = 0; i < 3; i++) fVertex[i] = (Double_t)mcVertex[i];
      npart = rsn->GetRefMCESD()->GetNumberOfTracks();
   } else {
      for (i = 0; i < 3; i++) fVertex[i] = 0.0;
      AliAODEvent *aod = rsn->GetRefMCAOD();
      TClonesArray *listAOD = (TClonesArray*)(aod->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
      if (listAOD) npart = listAOD->GetEntries();
      AliAODMCHeader *mcH = static_cast<AliAODMCHeader*>(aod->FindListObject(AliAODMCHeader::StdBranchName()));
      if (mcH) mcH->GetVertex(fVertex);
   }
   
   // check number of particles
   if (!npart) {
      AliInfo("Empty event");
      return 0;
   }
   
   // utility variables
   Int_t ipart, istep, count = 0, nsteps = fSteps.GetEntries();
   Int_t ntracks = rsn->GetAbsoluteSum();
   AliRsnDaughter check;
   
   // loop over particles
   for (ipart = 0; ipart < npart; ipart++) {
      // check i-th particle
      if (!AssignMotherAndDaughters(rsn, ipart)) continue;
      // if assignment was successful, for first step, use MC info
      fDaughter[0].SetRef(fDaughter[0].GetRefMC());
      fDaughter[1].SetRef(fDaughter[1].GetRefMC());
      fMother.SetRefEvent(rsn);
      fMother.ComputeSum(fDef->GetDef1().GetMass(), fDef->GetDef2().GetMass(), fDef->GetMotherMass());
      fOutput->Fill(&fMother, 0);
      count++;
      // for each further step, try to find two tracks which pass the related cuts
      for (istep = 0; istep < nsteps; istep++) {
         AliRsnCutManager *cuts = (AliRsnCutManager*)fSteps[istep];
         fDaughter[0].SetBad();
         fDaughter[1].SetBad();
         for (i = 0; i < ntracks; i++) {
            rsn->SetDaughter(check, i);
            if (!cuts->PassCommonDaughterCuts(&check)) continue;
            if (TMath::Abs(check.GetLabel()) == fDaughter[0].GetLabel()) {
               if (!cuts->PassDaughter1Cuts(&check)) continue;
               fDaughter[0].SetRef(check.GetRef());
               fDaughter[0].SetGood();
            } else if (TMath::Abs(check.GetLabel()) == fDaughter[1].GetLabel()) {
               if (!cuts->PassDaughter2Cuts(&check)) continue;
               fDaughter[1].SetRef(check.GetRef());
               fDaughter[1].SetGood();
            }
            if (fDaughter[0].IsOK() && fDaughter[1].IsOK()) {
               fMother.ComputeSum(fDef->GetDef1().GetMass(), fDef->GetDef2().GetMass(), fDef->GetMotherMass());
               if (cuts->PassMotherCuts(&fMother)) {
                  fOutput->Fill(&fMother, 1 + istep);
                  break;
               }
            }
         }
      }
   }
   
   return count;
}
Commit	Line	Data
c865cb1d	1	//
	2	// Class AliRsnLoopEffPair
	3	//
	4	// Inherits from basic AliRsnLoopEff for efficiency,
	5	// and computed efficiencies for single-tracks
	6	//
	7	// author: Alberto Pulvirenti (alberto.pulvirenti@ct.infn.it)
	8	//
	9
f34f960b	10	#include <Riostream.h>
f34f960b	11
c865cb1d	12	#include "AliStack.h"
f34f960b	13	#include "AliGenEventHeader.h"
f34f960b	14	#include "AliAODMCHeader.h"
c865cb1d	15
e187bd70	16	#include "AliRsnEvent.h"
c865cb1d	17	#include "AliRsnCutManager.h"
	18	#include "AliRsnDaughterDef.h"
	19	#include "AliRsnPairDef.h"
	20	#include "AliRsnLoopEffPair.h"
	21
	22	ClassImp(AliRsnLoopEffPair)
	23
	24	//_____________________________________________________________________________
	25	AliRsnLoopEffPair::AliRsnLoopEffPair(const char name, AliRsnPairDef def) :
f34f960b	26	AliRsnLoopEff(name, 1),
c865cb1d	27	fDef(def),
	28	fMother()
	29	{
	30	//
	31	// Default constructor.
	32	// Do not repeat 'DefineOutput' since it is done in base class and we don't add new ones.
	33	//
	34	}
	35
	36	//_____________________________________________________________________________
	37	AliRsnLoopEffPair::AliRsnLoopEffPair(const AliRsnLoopEffPair& copy) :
	38	AliRsnLoopEff(copy),
	39	fDef(copy.fDef),
	40	fMother(copy.fMother)
	41	{
	42	//
	43	// Copy constrtuctor.
	44	//
	45	}
	46
	47	//_____________________________________________________________________________
	48	AliRsnLoopEffPair& AliRsnLoopEffPair::operator=(const AliRsnLoopEffPair& copy)
	49	{
	50	//
	51	// Assignment operator.
	52	// Owned data members are meaningless for this operator.
	53	//
	54
	55	AliRsnLoopEff::operator=(copy);
	56	fDef = copy.fDef;
	57
	58	return (*this);
	59	}
	60
f34f960b	61	//__________________________________________________________________________________________________
	62	Bool_t AliRsnLoopEffPair::AssignMotherAndDaughters(AliRsnEvent *rsnEvent, Int_t ipart)
	63	{
	64	//
	65	// Calls the appropriate assignment method
	66	//
	67
	68	// setup pointers
	69	fMother.SetDaughter(0, &fDaughter[0]);
	70	fMother.SetDaughter(1, &fDaughter[1]);
	71	fMother.SetRefEvent(rsnEvent);
	72	fDaughter[0].SetOwnerEvent(rsnEvent);
	73	fDaughter[1].SetOwnerEvent(rsnEvent);
	74
	75	if (rsnEvent->IsESD())
	76	return AssignMotherAndDaughtersESD(rsnEvent, ipart);
	77	else if (rsnEvent->IsAOD())
	78	return AssignMotherAndDaughtersAOD(rsnEvent, ipart);
	79	else {
	80	AliError("Unrecognized input event");
	81	return kFALSE;
	82	}
	83	}
	84
	85	//__________________________________________________________________________________________________
	86	Bool_t AliRsnLoopEffPair::AssignMotherAndDaughtersESD(AliRsnEvent *rsnEvent, Int_t ipart)
	87	{
	88	//
	89	// Gets a particle in the MC event and try to assign it to the mother.
	90	// If it has two daughters, retrieve them and assign also them.
	91	// NOTE: assignment is done only for MC, since reconstructed match is assigned in the same way
	92	// for ESD and AOD, if available
	93	// ---
	94	// Implementation for ESD inputs
	95	//
	96
	97	AliMCEvent *mc = rsnEvent->GetRefMCESD();
	98	AliStack *stack = mc->Stack();
	99	AliMCParticle mother = (AliMCParticle)mc->GetTrack(ipart);
	100	TParticle *motherP = mother->Particle();
	101	Int_t ntracks = stack->GetNtrack();
	102
	103	// check PDG code and exit if it is wrong
	104	if (TMath::Abs(motherP->GetPdgCode()) != fDef->GetMotherPDG()) return kFALSE;
	105
	106	// check number of daughters and exit if it is not 2
	107	if (motherP->GetNDaughters() < 2) return kFALSE;
	108
	109	// check distance from primary vertex
	110	TLorentzVector vprod;
	111	motherP->ProductionVertex(vprod);
	112	if (!CheckDistanceFromPV(vprod.X(), vprod.Y(), vprod.Z())) {
	113	AliDebugClass(1, "Distant production vertex");
	114	return kFALSE;
	115	}
	116
	117	// get the daughters and check their PDG code and charge:
	118	// if they match one of the pair daughter definitions,
	119	// assign them as MC reference of the 'fDaughter' objects
	120	fDaughter[0].Reset();
	121	fDaughter[1].Reset();
	122	Int_t index[2] = {motherP->GetFirstDaughter(), motherP->GetLastDaughter()};
	123	Int_t i, pdg;
	124	Short_t charge;
125	AliMCParticle *daughter = 0x0;
126	for (i = 0; i < 2; i++) {
127	// check index for stack
128	if (index[i] < 0 \|\| index[i] > ntracks) {
129	AliError(Form("Index %d overflow: value = %d, max = %d", i, index[i], ntracks));
130	return kFALSE;
131	}
132	// get daughter and its PDG and charge
133	daughter = (AliMCParticle*)mc->GetTrack(index[i]);
134	pdg = TMath::Abs(daughter->Particle()->GetPdgCode());
135	charge = (Short_t)(daughter->Particle()->GetPDG()->Charge() / 3);
136	// check if it matches one definition
137	if (fDef->GetDef1().MatchesPDG(pdg) && fDef->GetDef1().MatchesChargeS(charge)) {
138	fDaughter[0].SetGood();
139	fDaughter[0].SetRefMC(daughter);
140	fDaughter[0].SetLabel(index[i]);
141	} else if (fDef->GetDef2().MatchesPDG(pdg) && fDef->GetDef2().MatchesChargeS(charge)) {
142	fDaughter[1].SetGood();
143	fDaughter[1].SetRefMC(daughter);
144	fDaughter[1].SetLabel(index[i]);
145	}
146	}
147
148	// return success if both daughters were assigned
149	if (fDaughter[0].IsOK() && fDaughter[1].IsOK()) {
150	return kTRUE;
151	} else {
152	fDaughter[0].Reset();
153	fDaughter[1].Reset();
154	return kFALSE;
155	}
156	}
157
158	//__________________________________________________________________________________________________
159	Bool_t AliRsnLoopEffPair::AssignMotherAndDaughtersAOD(AliRsnEvent *rsnEvent, Int_t ipart)
160	{
161	AliAODEvent *aod = rsnEvent->GetRefAOD();
162	TClonesArray listAOD = (TClonesArray)(aod->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
163	AliAODMCParticle mother = (AliAODMCParticle)listAOD->At(ipart);
164	Int_t ntracks = listAOD->GetEntries();
165
166	// check PDG code and exit if it is wrong
167	if (TMath::Abs(mother->GetPdgCode()) != fDef->GetMotherPDG()) return kFALSE;
168
169	// check number of daughters and exit if it is not 2
170	if (mother->GetNDaughters() < 2) return kFALSE;
171
172	// check distance from primary vertex
173	Double_t vprod[3] = {(Double_t)mother->Xv(), (Double_t)mother->Yv(), (Double_t)mother->Zv()};
174	Double_t dv = DistanceFromPV(vprod[0], vprod[1], vprod[2]);
175	if (dv > fMaxDistPV) {
176	AliDebugClass(1, "Distant production vertex");
177	return kFALSE;
178	}
179
180	// get the daughters and check their PDG code and charge:
181	// if they match one of the pair daughter definitions,
182	// assign them as MC reference of the 'fDaughter' objects
183	fDaughter[0].Reset();
184	fDaughter[1].Reset();
185	Int_t index[2] = {(Int_t)mother->GetDaughter(0), (Int_t)mother->GetDaughter(1)};
186	Int_t i, pdg;
187	Short_t charge;
188	AliAODMCParticle *daughter = 0x0;
189	for (i = 0; i < 2; i++) {
190	// check index for stack
191	if (index[i] < 0 \|\| index[i] > ntracks) {
192	AliError(Form("Index %d overflow: value = %d, max = %d", i, index[i], ntracks));
193	return kFALSE;
194	}
195	// get daughter and its PDG and charge
196	daughter = (AliAODMCParticle*)listAOD->At(index[i]);
197	pdg = TMath::Abs(daughter->GetPdgCode());
198	charge = (Short_t)(daughter->Charge() / 3);
199	// check if it matches one definition
200	if (fDef->GetDef1().MatchesPDG(pdg) && fDef->GetDef1().MatchesChargeS(charge)) {
201	fDaughter[0].SetGood();
202	fDaughter[0].SetRefMC(daughter);
203	fDaughter[0].SetLabel(index[i]);
204	} else if (fDef->GetDef2().MatchesPDG(pdg) && fDef->GetDef2().MatchesChargeS(charge)) {
205	fDaughter[1].SetGood();
206	fDaughter[1].SetRefMC(daughter);
207	fDaughter[1].SetLabel(index[i]);
208	}
209	}
210
211	// return success if both daughters were assigned
212	if (fDaughter[0].IsOK() && fDaughter[1].IsOK()) {
213	return kTRUE;
214	} else {
215	fDaughter[0].Reset();
216	fDaughter[1].Reset();
217	return kFALSE;
218	}
219	}
220
c865cb1d	221	//_____________________________________________________________________________
	222	Int_t AliRsnLoopEffPair::DoLoop(AliRsnEvent rsn, AliRsnDaughterSelector, AliRsnEvent, AliRsnDaughterSelector)
	223	{
	224	//
	225	// Loop on event and fill containers
	226	//
	227
	228	// check event cuts
	229	if (!OkEvent(rsn)) return 0;
	230
f34f960b	231	// retrieve output
	232	fOutput = (AliRsnListOutput*)fOutputs[0];
	233
c865cb1d	234	// check presence of MC reference
	235	if (!rsn->GetRefMC()) {
	236	AliError("Need a MC to compute efficiency");
	237	return 0;
	238	}
	239
f34f960b	240	// check presence of event
	241	if (!rsn->GetRef()) {
	242	AliError("Need an event to compute efficiency");
	243	return 0;
	244	}
	245
c865cb1d	246	// check event type:
	247	// must be ESD or AOD, and then use a bool to know in the rest
	248	if (!rsn->IsESD() && !rsn->IsAOD()) {
	249	AliError("Need to process ESD or AOD input");
	250	return 0;
	251	}
c865cb1d	252
f34f960b	253	// retrieve the MC primary vertex position
	254	// and do some additional coherence checks
	255	Int_t i, npart = 0;
e187bd70	256	if (rsn->IsESD()) {
f34f960b	257	TArrayF mcVertex(3);
	258	AliGenEventHeader *genEH = rsn->GetRefMCESD()->GenEventHeader();
	259	genEH->PrimaryVertex(mcVertex);
	260	for (i = 0; i < 3; i++) fVertex[i] = (Double_t)mcVertex[i];
	261	npart = rsn->GetRefMCESD()->GetNumberOfTracks();
c865cb1d	262	} else {
f34f960b	263	for (i = 0; i < 3; i++) fVertex[i] = 0.0;
	264	AliAODEvent *aod = rsn->GetRefMCAOD();
	265	TClonesArray listAOD = (TClonesArray)(aod->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
	266	if (listAOD) npart = listAOD->GetEntries();
	267	AliAODMCHeader mcH = static_cast<AliAODMCHeader>(aod->FindListObject(AliAODMCHeader::StdBranchName()));
	268	if (mcH) mcH->GetVertex(fVertex);
c865cb1d	269	}
	270
	271	// check number of particles
	272	if (!npart) {
	273	AliInfo("Empty event");
	274	return 0;
	275	}
c865cb1d	276
c865cb1d	277	// utility variables
f34f960b	278	Int_t ipart, istep, count = 0, nsteps = fSteps.GetEntries();
	279	Int_t ntracks = rsn->GetAbsoluteSum();
	280	AliRsnDaughter check;
c865cb1d	281
	282	// loop over particles
	283	for (ipart = 0; ipart < npart; ipart++) {
f34f960b	284	// check i-th particle
	285	if (!AssignMotherAndDaughters(rsn, ipart)) continue;
	286	// if assignment was successful, for first step, use MC info
	287	fDaughter[0].SetRef(fDaughter[0].GetRefMC());
	288	fDaughter[1].SetRef(fDaughter[1].GetRefMC());
	289	fMother.SetRefEvent(rsn);
	290	fMother.ComputeSum(fDef->GetDef1().GetMass(), fDef->GetDef2().GetMass(), fDef->GetMotherMass());
	291	fOutput->Fill(&fMother, 0);
c865cb1d	292	count++;
f34f960b	293	// for each further step, try to find two tracks which pass the related cuts
	294	for (istep = 0; istep < nsteps; istep++) {
	295	AliRsnCutManager cuts = (AliRsnCutManager)fSteps[istep];
	296	fDaughter[0].SetBad();
	297	fDaughter[1].SetBad();
	298	for (i = 0; i < ntracks; i++) {
	299	rsn->SetDaughter(check, i);
	300	if (!cuts->PassCommonDaughterCuts(&check)) continue;
	301	if (TMath::Abs(check.GetLabel()) == fDaughter[0].GetLabel()) {
	302	if (!cuts->PassDaughter1Cuts(&check)) continue;
	303	fDaughter[0].SetRef(check.GetRef());
	304	fDaughter[0].SetGood();
	305	} else if (TMath::Abs(check.GetLabel()) == fDaughter[1].GetLabel()) {
	306	if (!cuts->PassDaughter2Cuts(&check)) continue;
	307	fDaughter[1].SetRef(check.GetRef());
	308	fDaughter[1].SetGood();
	309	}
	310	if (fDaughter[0].IsOK() && fDaughter[1].IsOK()) {
	311	fMother.ComputeSum(fDef->GetDef1().GetMass(), fDef->GetDef2().GetMass(), fDef->GetMotherMass());
	312	if (cuts->PassMotherCuts(&fMother)) {
	313	fOutput->Fill(&fMother, 1 + istep);
	314	break;
	315	}
	316	}
c865cb1d	317	}
	318	}
	319	}
	320
	321	return count;
	322	}