Main use case: analysis of merged events.
return;
}
- AliStack* stack = mcE->Stack();
Int_t np = mcE->GetNumberOfTracks();
- Int_t nprim = stack->GetNprimary();
+ Int_t nprim = mcE->GetNumberOfPrimaries();
// TODO ADD MC VERTEX
// We take all real primaries
} else if (part->GetUniqueID() == 4) {
// Particles from decay
// Check that the decay chain ends at a primary particle
- TParticle* mother = part;
- Int_t imo = part->GetFirstMother();
+ AliMCParticle* mother = mcpart;
+ Int_t imo = mcpart->GetMother();
while((imo >= nprim) && (mother->GetUniqueID() == 4)) {
- mother = mcE->Stack()->Particle(imo);
- imo = mother->GetFirstMother();
+ mother = mcE->GetTrack(imo);
+ imo = mother->GetMother();
}
// Select according to pseudorapidity and production point of primary ancestor
- if (imo < nprim && Select(mcE->Stack()->Particle(imo), rv, zv))write = kTRUE;
+ if (imo < nprim && Select(mcE->GetTrack(imo)->Particle(), rv, zv))write = kTRUE;
} else if (part->GetUniqueID() == 5) {
// Now look for pair production
- Int_t imo = part->GetFirstMother();
+ Int_t imo = mcpart->GetMother();
if (imo < nprim) {
// Select, if the gamma is a primary
write = kTRUE;
} else {
// Check if the gamma comes from the decay chain of a primary particle
- TParticle* mother = mcE->Stack()->Particle(imo);
- imo = mother->GetFirstMother();
+ AliMCParticle* mother = mcE->GetTrack(imo);
+ imo = mother->GetMother();
while((imo >= nprim) && (mother->GetUniqueID() == 4)) {
- mother = mcE->Stack()->Particle(imo);
- imo = mother->GetFirstMother();
+ mother = mcE->GetTrack(imo);
+ imo = mother->GetMother();
}
// Select according to pseudorapidity and production point
- if (imo < nprim && Select(mother, rv, zv))
+ if (imo < nprim && Select(mother->Particle(), rv, zv))
write = kTRUE;
}
}
// cases
iCharm++;
Printf("Decay Mother %s",part->GetPDG()->GetName());
- Int_t d0 = part->GetFirstDaughter();
- Int_t d1 = part->GetLastDaughter();
+ Int_t d0 = mcpart->GetFirstDaughter();
+ Int_t d1 = mcpart->GetLastDaughter();
if(d0>0&&d1>0){
for(int id = d0;id <= d1;id++){
- TParticle* daughter = mcE->Stack()->Particle(id);
- Printf("Decay Daughter %s",daughter->GetPDG()->GetName());
+ AliMCParticle* daughter = mcE->GetTrack(id);
+ Printf("Decay Daughter %s",daughter->Particle()->GetPDG()->GetName());
iAll++;
if(mcH->IsParticleSelected(id))iTaken++;
}
// Store the AliAODParticlesMC
-
- Int_t np = pStack->GetNtrack();
- Int_t nprim = pStack->GetNprimary();
+ AliMCEvent* mcEvent = fMCEventH->MCEvent();
+
+ Int_t np = mcEvent->GetNumberOfTracks();
+ Int_t nprim = mcEvent->GetNumberOfPrimaries();
Int_t j = 0;
TClonesArray& l = *mcarray;
- for(int i = 0;i < np;++i){
- if(fMCEventH->IsParticleSelected(i)){
-
- Int_t flag = 0;
- TParticle *part = pStack->Particle(i);
- if(i<nprim)flag |= AliAODMCParticle::kPrimary;
- if(pStack->IsPhysicalPrimary(i))flag |= AliAODMCParticle::kPhysicalPrim;
-
- if(fMCEventH->GetNewLabel(i)!=j){
- AliError(Form("MISMATCH New label %d j: %d",fMCEventH->GetNewLabel(i),j));
- }
- AliAODMCParticle mcpart_tmp(part,i,flag);
-
- //
- Int_t d0 = mcpart_tmp.GetDaughter(0);
- Int_t d1 = mcpart_tmp.GetDaughter(1);
- Int_t m = mcpart_tmp.GetMother();
-
- // other than for the track labels, negative values mean
- // no daughter/mother so preserve it
-
- if(d0<0 && d1<0){
- // no first daughter -> no second daughter
- // nothing to be done
- // second condition not needed just for sanity check at the end
- mcpart_tmp.SetDaughter(0,d0);
- mcpart_tmp.SetDaughter(1,d1);
- }
- else if(d1 < 0 && d0 >= 0){
- // Only one daughter
- // second condition not needed just for sanity check at the end
- if(fMCEventH->IsParticleSelected(d0)){
- mcpart_tmp.SetDaughter(0,fMCEventH->GetNewLabel(d0));
- }
- else{
- mcpart_tmp.SetDaughter(0,-1);
- }
- mcpart_tmp.SetDaughter(1,d1);
- }
- else if (d0 > 0 && d1 > 0 ){
- // we have two or more daughters loop on the stack to see if they are
- // selected
- Int_t d0_tmp = -1;
- Int_t d1_tmp = -1;
- for(int id = d0; id<=d1;++id){
- if(fMCEventH->IsParticleSelected(id)){
- if(d0_tmp==-1){
- // first time
- d0_tmp = fMCEventH->GetNewLabel(id);
- d1_tmp = d0_tmp; // this is to have the same schema as on the stack i.e. with one daugther d0 and d1 are the same
- }
- else d1_tmp = fMCEventH->GetNewLabel(id);
+ for(int i = 0; i < np; ++i){
+ if(fMCEventH->IsParticleSelected(i)){
+ Int_t flag = 0;
+ AliMCParticle* mcpart = mcEvent->GetTrack(i);
+ TParticle *part = mcpart->Particle();
+ if(i<nprim)flag |= AliAODMCParticle::kPrimary;
+
+ if(mcEvent->IsPhysicalPrimary(i))flag |= AliAODMCParticle::kPhysicalPrim;
+
+ if(fMCEventH->GetNewLabel(i)!=j){
+ AliError(Form("MISMATCH New label %d j: %d",fMCEventH->GetNewLabel(i),j));
}
- }
- mcpart_tmp.SetDaughter(0,d0_tmp);
- mcpart_tmp.SetDaughter(1,d1_tmp);
- }
- else{
- AliError(Form("Unxpected indices %d %d",d0,d1));
- }
- if(m<0){
- mcpart_tmp.SetMother(m);
- }
- else{
- if(fMCEventH->IsParticleSelected(m))mcpart_tmp.SetMother(fMCEventH->GetNewLabel(m));
- else AliError("PROBLEM Mother not selected");
+ AliAODMCParticle mcpart_tmp(mcpart,i,flag);
+
+ //
+ Int_t d0 = mcpart_tmp.GetDaughter(0);
+ Int_t d1 = mcpart_tmp.GetDaughter(1);
+ Int_t m = mcpart_tmp.GetMother();
+
+ // other than for the track labels, negative values mean
+ // no daughter/mother so preserve it
+
+ if(d0<0 && d1<0){
+ // no first daughter -> no second daughter
+ // nothing to be done
+ // second condition not needed just for sanity check at the end
+ mcpart_tmp.SetDaughter(0,d0);
+ mcpart_tmp.SetDaughter(1,d1);
+ } else if(d1 < 0 && d0 >= 0) {
+ // Only one daughter
+ // second condition not needed just for sanity check at the end
+ if(fMCEventH->IsParticleSelected(d0)){
+ mcpart_tmp.SetDaughter(0,fMCEventH->GetNewLabel(d0));
+ } else {
+ mcpart_tmp.SetDaughter(0,-1);
+ }
+ mcpart_tmp.SetDaughter(1,d1);
+ }
+ else if (d0 > 0 && d1 > 0 ){
+ // we have two or more daughters loop on the stack to see if they are
+ // selected
+ Int_t d0_tmp = -1;
+ Int_t d1_tmp = -1;
+ for(int id = d0; id<=d1;++id){
+ if(fMCEventH->IsParticleSelected(id)){
+ if(d0_tmp==-1){
+ // first time
+ d0_tmp = fMCEventH->GetNewLabel(id);
+ d1_tmp = d0_tmp; // this is to have the same schema as on the stack i.e. with one daugther d0 and d1 are the same
+ }
+ else d1_tmp = fMCEventH->GetNewLabel(id);
+ }
+ }
+ mcpart_tmp.SetDaughter(0,d0_tmp);
+ mcpart_tmp.SetDaughter(1,d1_tmp);
+ } else {
+ AliError(Form("Unxpected indices %d %d",d0,d1));
+ }
+
+ if(m<0){
+ mcpart_tmp.SetMother(m);
+ } else {
+ if(fMCEventH->IsParticleSelected(m))mcpart_tmp.SetMother(fMCEventH->GetNewLabel(m));
+ else AliError(Form("PROBLEM Mother not selected %d \n", m));
+ }
+
+ new (l[j++]) AliAODMCParticle(mcpart_tmp);
+
}
-
- new (l[j++]) AliAODMCParticle(mcpart_tmp);
-
- }
}
- AliInfo(Form("AliAODHandler::StoreMCParticles: Selected %d (Primaries %d / total %d) after validation",
+ AliInfo(Form("AliAODHandler::StoreMCParticles: Selected %d (Primaries %d / total %d) after validation \n",
j,nprim,np));
// Set the labels in the AOD output...
for(int it = 0; it < fAODEvent->GetNTracks();++it){
AliAODTrack *track = fAODEvent->GetTrack(it);
- if(TMath::Abs(track->GetLabel())>np||track->GetLabel()==0){
- AliWarning(Form("Wrong ESD track label %d",track->GetLabel()));
+ Int_t label = TMath::Abs(track->GetLabel());
+ if (label >= AliMCEvent::BgLabelOffset()) label = mcEvent->BgLabelToIndex(label);
+
+
+ if(label > np || track->GetLabel() == 0){
+ AliWarning(Form("Wrong ESD track label %5d (%5d)",track->GetLabel(), label));
}
- if(fMCEventH->GetNewLabel(track->GetLabel())==0){
- AliWarning(Form("New label not found for %d",track->GetLabel()));
+ if(fMCEventH->GetNewLabel(label) == 0){
+ AliWarning(Form("New label not found for %5d (%5d)",track->GetLabel(), label));
}
- track->SetLabel(fMCEventH->GetNewLabel(track->GetLabel()));
+ track->SetLabel(fMCEventH->GetNewLabel(label));
}
}
}
-AliAODMCParticle::AliAODMCParticle(TParticle* part, Int_t label,Int_t flag):
+AliAODMCParticle::AliAODMCParticle(AliMCParticle* mcpart, Int_t label,Int_t flag):
AliVParticle(),
- fPdgCode(part->GetPdgCode()),
+ fPdgCode(mcpart->Particle()->GetPdgCode()),
fFlag(flag),
fLabel(label),
- fMother(part->GetMother(0)),
- fPx(part->Px()),
- fPy(part->Py()),
- fPz(part->Pz()),
- fE(part->Energy()),
- fVx(part->Vx()),
- fVy(part->Vy()),
- fVz(part->Vz())
+ fMother(mcpart->GetMother()),
+ fPx(mcpart->Particle()->Px()),
+ fPy(mcpart->Particle()->Py()),
+ fPz(mcpart->Particle()->Pz()),
+ fE(mcpart->Particle()->Energy()),
+ fVx(mcpart->Particle()->Vx()),
+ fVy(mcpart->Particle()->Vy()),
+ fVz(mcpart->Particle()->Vz())
{
- fDaughter[0] = part->GetDaughter(0);
- fDaughter[1] = part->GetDaughter(1);
+ fDaughter[0] = mcpart->GetFirstDaughter();
+ fDaughter[1] = mcpart->GetLastDaughter();
// Set unique id
- TObject::SetUniqueID(part->GetUniqueID());
+ TObject::SetUniqueID(mcpart->Particle()->GetUniqueID());
}
#include "AliTrackReference.h"
#include "AliVParticle.h"
+#include "AliMCParticle.h"
class AliAODEvent;
class TParticle;
class AliAODMCParticle: public AliVParticle {
public:
AliAODMCParticle();
- AliAODMCParticle(TParticle* part, Int_t label=0,Int_t flag = 0);
+ AliAODMCParticle(AliMCParticle* part, Int_t label=0,Int_t flag = 0);
virtual ~AliAODMCParticle(){};
AliAODMCParticle(const AliAODMCParticle& mcPart);
AliAODMCParticle& operator=(const AliAODMCParticle& mcPart);
#include <TParticle.h>
#include <TClonesArray.h>
#include <TRefArray.h>
+#include <TList.h>
#include "AliLog.h"
#include "AliMCEvent.h"
#include "AliGenEventHeader.h"
+Int_t AliMCEvent::fgkBgLabelOffset(10000000);
+
+
AliMCEvent::AliMCEvent():
AliVEvent(),
fStack(0),
fTmpTreeTR(0),
fTmpFileTR(0),
fNprimaries(-1),
- fNparticles(-1)
+ fNparticles(-1),
+ fSubsidiaryEvents(0),
+ fPrimaryOffset(0),
+ fSecondaryOffset(0)
{
// Default constructor
}
AliMCEvent::AliMCEvent(const AliMCEvent& mcEvnt) :
AliVEvent(mcEvnt),
- fStack(0),
- fMCParticles(0),
- fMCParticleMap(0),
- fHeader(new AliHeader()),
- fTRBuffer(0),
- fTrackReferences(0),
- fTreeTR(0),
- fTmpTreeTR(0),
- fTmpFileTR(0),
- fNprimaries(-1),
- fNparticles(-1)
+ fStack(mcEvnt.fStack),
+ fMCParticles(mcEvnt.fMCParticles),
+ fMCParticleMap(mcEvnt.fMCParticleMap),
+ fHeader(mcEvnt.fHeader),
+ fTRBuffer(mcEvnt.fTRBuffer),
+ fTrackReferences(mcEvnt.fTrackReferences),
+ fTreeTR(mcEvnt.fTreeTR),
+ fTmpTreeTR(mcEvnt.fTmpTreeTR),
+ fTmpFileTR(mcEvnt.fTmpFileTR),
+ fNprimaries(mcEvnt.fNprimaries),
+ fNparticles(mcEvnt.fNparticles),
+ fSubsidiaryEvents(0),
+ fPrimaryOffset(0),
+ fSecondaryOffset(0)
{
// Copy constructor
}
fNparticles = -1;
fNprimaries = -1;
fStack = 0;
-
+// fSubsidiaryEvents->Clear();
+ fSubsidiaryEvents = 0;
}
AliMCParticle* AliMCEvent::GetTrack(Int_t i) const
{
// Get MC Particle i
+ //
+ //
+ // Check first if this explicitely accesses the subsidiary event
+ if (i >= BgLabelOffset()) {
+ if (fSubsidiaryEvents) {
+ AliMCEvent* bgEvent = (AliMCEvent*) (fSubsidiaryEvents->At(1));
+ return (bgEvent->GetTrack(i - BgLabelOffset()));
+ } else {
+ return 0;
+ }
+ }
+
+ //
AliMCParticle *mcParticle = 0;
TParticle *particle = 0;
TClonesArray *trefs = 0;
Int_t ntref = 0;
TRefArray *rarray = 0;
+
+
+
// Out of range check
if (i < 0 || i >= fNparticles) {
AliWarning(Form("AliMCEvent::GetEntry: Index out of range"));
mcParticle = 0;
return (mcParticle);
}
+
+ if (fSubsidiaryEvents) {
+ AliMCEvent* mc;
+ Int_t idx = FindIndexAndEvent(i, mc);
+ return (mc->GetTrack(idx));
+ }
//
- // First check of the MC Particle has been already cached
+ // First check If the MC Particle has been already cached
if(!fMCParticleMap->At(i)) {
// Get particle from the stack
particle = fStack->Particle(i);
new ((*fMCParticles)[nentries]) AliMCParticle(particle, rarray, i);
mcParticle = dynamic_cast<AliMCParticle*>((*fMCParticles)[nentries]);
fMCParticleMap->AddAt(mcParticle, i);
+
+ TParticle* part = mcParticle->Particle();
+ Int_t imo = part->GetFirstMother();
+ Int_t id1 = part->GetFirstDaughter();
+ Int_t id2 = part->GetLastDaughter();
+ if (fPrimaryOffset > 0 || fSecondaryOffset > 0) {
+ // Remapping of the mother and daughter indices
+ if (imo < fNprimaries) {
+ mcParticle->SetMother(imo + fPrimaryOffset);
+ } else {
+ mcParticle->SetMother(imo + fSecondaryOffset - fNprimaries);
+ }
+
+ if (id1 < fNprimaries) {
+ mcParticle->SetFirstDaughter(id1 + fPrimaryOffset);
+ mcParticle->SetLastDaughter (id2 + fPrimaryOffset);
+ } else {
+ mcParticle->SetFirstDaughter(id1 + fSecondaryOffset - fNprimaries);
+ mcParticle->SetLastDaughter (id2 + fSecondaryOffset - fNprimaries);
+ }
+
+
+ if (i > fNprimaries) {
+ mcParticle->SetLabel(i + fPrimaryOffset);
+ } else {
+ mcParticle->SetLabel(i + fSecondaryOffset - fNprimaries);
+ }
+
+
+ } else {
+ mcParticle->SetFirstDaughter(id1);
+ mcParticle->SetLastDaughter (id2);
+ mcParticle->SetMother (imo);
+ }
+
} else {
mcParticle = dynamic_cast<AliMCParticle*>(fMCParticleMap->At(i));
}
+
+
return mcParticle;
}
- AliGenEventHeader* AliMCEvent::GenEventHeader() {return (fHeader->GenEventHeader());}
+ AliGenEventHeader* AliMCEvent::GenEventHeader() {return (fHeader->GenEventHeader());}
+
+
+void AliMCEvent::AddSubsidiaryEvent(AliMCEvent* event)
+{
+ // Add a subsidiary event to the list; for example merged background event.
+ if (!fSubsidiaryEvents) {
+ TList* events = new TList();
+ events->Add(new AliMCEvent(*this));
+ fSubsidiaryEvents = events;
+ }
+
+ fSubsidiaryEvents->Add(event);
+}
+
+Int_t AliMCEvent::FindIndexAndEvent(Int_t oldidx, AliMCEvent*& event) const
+{
+ // Find the index and event in case of composed events like signal + background
+ TIter next(fSubsidiaryEvents);
+ next.Reset();
+ if (oldidx < fNprimaries) {
+ while((event = (AliMCEvent*)next())) {
+ if (oldidx < (event->GetPrimaryOffset() + event->GetNumberOfPrimaries())) break;
+ }
+ return (oldidx - event->GetPrimaryOffset());
+ } else {
+ while((event = (AliMCEvent*)next())) {
+ if (oldidx < (event->GetSecondaryOffset() + (event->GetNumberOfTracks() - event->GetNumberOfPrimaries()))) break;
+ }
+ return (oldidx - event->GetSecondaryOffset() + event->GetNumberOfPrimaries());
+ }
+}
+
+Int_t AliMCEvent::BgLabelToIndex(Int_t label)
+{
+ // Convert a background label to an absolute index
+ if (fSubsidiaryEvents) {
+ AliMCEvent* bgEvent = (AliMCEvent*) (fSubsidiaryEvents->At(1));
+ label -= BgLabelOffset();
+ if (label < bgEvent->GetNumberOfPrimaries()) {
+ label += bgEvent->GetPrimaryOffset();
+ } else {
+ label += (bgEvent->GetSecondaryOffset() - fNprimaries);
+ }
+ }
+ return (label);
+}
+
+
+Bool_t AliMCEvent::IsPhysicalPrimary(Int_t i)
+{
+//
+// Delegate to subevent if necesarry
+
+
+ if (!fSubsidiaryEvents) {
+ return fStack->IsPhysicalPrimary(i);
+ } else {
+ AliMCEvent* evt = 0;
+ Int_t idx = FindIndexAndEvent(i, evt);
+ return (evt->IsPhysicalPrimary(idx));
+ }
+}
+
+void AliMCEvent::InitEvent()
+{
+ if (fSubsidiaryEvents) {
+ TIter next(fSubsidiaryEvents);
+ AliMCEvent* evt;
+ fNprimaries = 0;
+ fNparticles = 0;
+
+ while((evt = (AliMCEvent*)next())) {
+ fNprimaries += evt->GetNumberOfPrimaries();
+ fNparticles += evt->GetNumberOfTracks();
+ }
+
+ Int_t ioffp = 0;
+ Int_t ioffs = fNprimaries;
+ next.Reset();
+
+ while((evt = (AliMCEvent*)next())) {
+ evt->SetPrimaryOffset(ioffp);
+ evt->SetSecondaryOffset(ioffs);
+ ioffp += evt->GetNumberOfPrimaries();
+ ioffs += (evt->GetNumberOfTracks() - evt->GetNumberOfPrimaries());
+ }
+ }
+}
+
ClassImp(AliMCEvent)
class AliGenEventHeader;
class TClonesArray;
+class TList;
class AliMCEvent : public AliVEvent {
virtual void ConnectTreeK (TTree* tree);
virtual void ConnectTreeTR(TTree* tree);
virtual void Clean();
+ virtual void InitEvent();
virtual void FinishEvent();
virtual Int_t GetParticleAndTR(Int_t i, TParticle*& particle, TClonesArray*& trefs);
virtual void DrawCheck(Int_t i, Int_t search);
+ virtual void AddSubsidiaryEvent(AliMCEvent* event);
+ virtual Int_t GetNumberOfPrimaries() {return fNprimaries;}
+ virtual Int_t GetPrimaryOffset() const {return fPrimaryOffset;}
+ virtual Int_t GetSecondaryOffset() const {return fSecondaryOffset;}
+ virtual void SetPrimaryOffset(Int_t ioff) {fPrimaryOffset = ioff;}
+ virtual void SetSecondaryOffset(Int_t ioff) {fSecondaryOffset = ioff;}
+ virtual Bool_t IsPhysicalPrimary(Int_t i);
+ virtual Int_t BgLabelToIndex(Int_t label);
+ static Int_t BgLabelOffset() {return fgkBgLabelOffset;}
+
private:
virtual void ReorderAndExpandTreeTR();
-
-private:
- AliStack *fStack; // Current pointer to stack
- TClonesArray *fMCParticles; // Pointer to list of particles
- TRefArray *fMCParticleMap; // Map of MC Particles
- AliHeader *fHeader; // Current pointer to header
- TClonesArray *fTRBuffer; // Track reference buffer
- TClonesArray *fTrackReferences; // Array of track references
- TTree *fTreeTR; // Pointer to Track Reference Tree
- TTree *fTmpTreeTR; // Temporary tree TR to read old format
- TFile *fTmpFileTR; // Temporary file with TreeTR to read old format
- Int_t fNprimaries; // Number of primaries
- Int_t fNparticles; // Number of particles
+ virtual Int_t FindIndexAndEvent(Int_t oldidx, AliMCEvent*& event) const;
+private:
+ AliStack *fStack; // Current pointer to stack
+ TClonesArray *fMCParticles; // Pointer to list of particles
+ TRefArray *fMCParticleMap; // Map of MC Particles
+ AliHeader *fHeader; // Current pointer to header
+ TClonesArray *fTRBuffer; // Track reference buffer
+ TClonesArray *fTrackReferences; // Array of track references
+ TTree *fTreeTR; // Pointer to Track Reference Tree
+ TTree *fTmpTreeTR; // Temporary tree TR to read old format
+ TFile *fTmpFileTR; // Temporary file with TreeTR to read old format
+ Int_t fNprimaries; // Number of primaries
+ Int_t fNparticles; // Number of particles
+ TList *fSubsidiaryEvents; // List of possible subsidiary events (for example merged underlying event)
+ Int_t fPrimaryOffset; // Offset for primaries
+ Int_t fSecondaryOffset; // Offset for secondaries
+ static Int_t fgkBgLabelOffset; // Standard branch name
ClassDef(AliMCEvent, 1) // AliVEvent realisation for MC data
};
-/************************************************************************* * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved *
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved *
* *
* Author: The ALICE Off-line Project. *
* Contributors are mentioned in the code where appropriate. *
#include "AliMCEventHandler.h"
#include "AliMCEvent.h"
+#include "AliMCParticle.h"
#include "AliPDG.h"
#include "AliTrackReference.h"
#include "AliHeader.h"
#include <TTree.h>
#include <TFile.h>
+#include <TList.h>
#include <TParticle.h>
#include <TString.h>
#include <TClonesArray.h>
fFileNumber(0),
fEventsPerFile(0),
fReadTR(kTRUE),
- fInitOk(kFALSE)
+ fInitOk(kFALSE),
+ fSubsidiaryHandlers(0),
+ fEventsInContainer(0)
{
//
// Default constructor
fFileNumber(0),
fEventsPerFile(0),
fReadTR(kTRUE),
- fInitOk(kFALSE)
+ fInitOk(kFALSE),
+ fSubsidiaryHandlers(0),
+ fEventsInContainer(0)
{
//
// Constructor
fFileNumber = 0;
AliInfo(Form("Number of events in this directory %5d \n", fNEvent));
fInitOk = kTRUE;
+
+
+ if (fSubsidiaryHandlers) {
+ TIter next(fSubsidiaryHandlers);
+ AliMCEventHandler *handler;
+ while((handler = (AliMCEventHandler*)next())) {
+ handler->Init(opt);
+ handler->SetNumberOfEventsInContainer(fNEvent);
+ }
+ }
+
return kTRUE;
}
// Connect TR to MCEvent
fMCEvent->ConnectTreeTR(fTreeTR);
}
+
//
return kTRUE;
}
}
fInitOk = kTRUE;
+
return kTRUE;
}
fParticleSelected.Delete();
fLabelMap.Delete();
// Read the next event
+
+ if (fEventsInContainer != 0) {
+ entry = (Long64_t) entry * Float_t(fNEvent) / Float_t (fEventsInContainer);
+ }
+
+
if (entry == -1) {
fEvent++;
entry = fEvent;
}
if (entry >= fNEvent) {
- AliWarning(Form("AliMCEventHandler: Event number out of range %5d %5d\n", entry,fNEvent));
+ AliWarning(Form("AliMCEventHandler: Event number out of range %5d %5d\n", entry, fNEvent));
return kFALSE;
}
- return GetEvent(entry);
+
+ Bool_t result = GetEvent(entry);
+
+ if (fSubsidiaryHandlers) {
+ TIter next(fSubsidiaryHandlers);
+ AliMCEventHandler *handler;
+ while((handler = (AliMCEventHandler*)next())) {
+ handler->BeginEvent(entry);
+ }
+ next.Reset();
+ while((handler = (AliMCEventHandler*)next())) {
+ fMCEvent->AddSubsidiaryEvent(handler->MCEvent());
+ }
+ fMCEvent->InitEvent();
+ }
+ return result;
+
}
void AliMCEventHandler::SelectParticle(Int_t i){
// taking the absolute values here, need to take care
// of negative daughter and mother
// IDs when setting!
- if(!IsParticleSelected(TMath::Abs(i)))fParticleSelected.Add(TMath::Abs(i),1);
+ if (TMath::Abs(i) >= AliMCEvent::BgLabelOffset()) i = fMCEvent->BgLabelToIndex(TMath::Abs(i));
+ if(!IsParticleSelected(TMath::Abs(i)))fParticleSelected.Add(TMath::Abs(i),1);
}
Bool_t AliMCEventHandler::IsParticleSelected(Int_t i) {
}
VerifySelectedParticles();
- AliStack *pStack = fMCEvent->Stack();
Int_t iNew = 0;
- for(int i = 0;i < pStack->GetNtrack();++i){
+ for(int i = 0;i < fMCEvent->GetNumberOfTracks();++i){
if(IsParticleSelected(i)){
fLabelMap.Add(i,iNew);
iNew++;
}
Int_t AliMCEventHandler::GetNewLabel(Int_t i) {
- // Gets the labe from the new created Map
+ // Gets the label from the new created Map
// Call CreatLabelMap before
// otherwise only 0 returned
return fLabelMap.GetValue(TMath::Abs(i));
// Private, should be only called by CreateLabelMap
if(!fMCEvent){
- fParticleSelected.Delete();
- return;
+ fParticleSelected.Delete();
+ return;
}
- AliStack *pStack = fMCEvent->Stack();
- Int_t nprim = pStack->GetNprimary();
-
- for(int i = 0;i < pStack->GetNtrack();++i){
- if(i<nprim){
- SelectParticle(i);// take all primaries
- continue;
- }
- if(!IsParticleSelected(i))continue;
- TParticle *part = pStack->Particle(i);
- Int_t imo = part->GetFirstMother();
- while((imo >= nprim)&&!IsParticleSelected(imo)){
- // Mother not yet selected
- SelectParticle(imo);
- TParticle *mother = pStack->Particle(imo);
- imo = mother->GetFirstMother();
- }
- // after last step we may have a unselected primary
+ Int_t nprim = fMCEvent->GetNumberOfPrimaries();
+
+ for(int i = 0;i < fMCEvent->GetNumberOfTracks(); ++i){
+ if(i < nprim){
+ SelectParticle(i);// take all primaries
+ continue;
+ }
+
+ if(!IsParticleSelected(i))continue;
+
+ AliMCParticle* mcpart = fMCEvent->GetTrack(i);
+ Int_t imo = mcpart->GetMother();
+ while((imo >= nprim)&&!IsParticleSelected(imo)){
+ // Mother not yet selected
+ SelectParticle(imo);
+ AliMCParticle* mcpart = fMCEvent->GetTrack(imo);
+ imo = mcpart->GetMother();
+ }
+ // after last step we may have an unselected primary
// mother
if(imo>=0){
if(!IsParticleSelected(imo))
else if (fileName.BeginsWith("root:")) {
fileName.Append("?ZIP=");
}
+
*fPathName = fileName;
- AliInfo(Form("Notify() Path: %s\n", fPathName->Data()));
+ AliInfo(Form("Path: -%s-\n", fPathName->Data()));
ResetIO();
InitIO("");
+// Handle subsidiary handlers
+ if (fSubsidiaryHandlers) {
+ TIter next(fSubsidiaryHandlers);
+ AliMCEventHandler *handler;
+ while((handler = (AliMCEventHandler*) next())) {
+ TString* spath = handler->GetInputPath();
+ if (spath->Contains("merged")) {
+ if (! fPathName->IsNull()) {
+ handler->Notify(Form("%s/../.", fPathName->Data()));
+ } else {
+ handler->Notify("../");
+ }
+ }
+ }
+ }
+
return kTRUE;
}
if (fFileTR) {delete fFileTR; fFileTR = 0;}
fExtension="";
fInitOk = kFALSE;
+
+ if (fSubsidiaryHandlers) {
+ TIter next(fSubsidiaryHandlers);
+ AliMCEventHandler *handler;
+ while((handler = (AliMCEventHandler*)next())) {
+ handler->ResetIO();
+ }
+ }
+
}
delete fDirTR; fDirTR = 0;
delete fDirK; fDirK = 0;
if (fInitOk) fMCEvent->FinishEvent();
+
+ if (fSubsidiaryHandlers) {
+ TIter next(fSubsidiaryHandlers);
+ AliMCEventHandler *handler;
+ while((handler = (AliMCEventHandler*)next())) {
+ handler->FinishEvent();
+ }
+ }
+
return kTRUE;
}
delete fPathName;
fPathName = new TString(fname);
}
+
+void AliMCEventHandler::AddSubsidiaryHandler(AliMCEventHandler* handler)
+{
+ // Add a subsidiary handler. For example for background events
+
+ if (!fSubsidiaryHandlers) fSubsidiaryHandlers = new TList();
+ fSubsidiaryHandlers->Add(handler);
+}
class TFile;
class TTree;
+class TList;
+
class TParticle;
class TString;
class TClonesArray;
virtual void ResetIO();
virtual Bool_t GetEvent(Int_t iev);
virtual void SetReadTR(Bool_t flag) { fReadTR = flag; }
+ virtual void AddSubsidiaryHandler(AliMCEventHandler* handler);
+ virtual void SetNumberOfEventsInContainer(Int_t nev) {fEventsInContainer = nev;}
//
AliMCEvent* MCEvent() const {return fMCEvent;}
TTree* TreeTR() const {return fTreeTR;}
AliMCEventHandler(const AliMCEventHandler& handler);
AliMCEventHandler& operator=(const AliMCEventHandler& handler);
private:
- AliMCEvent *fMCEvent; //! MC Event
- TFile *fFileE; //! File with TreeE
- TFile *fFileK; //! File with TreeK
- TFile *fFileTR; //! File with TreeTR
- TTree *fTreeE; //! TreeE (Event Headers)
- TTree *fTreeK; //! TreeK (kinematics tree)
- TTree *fTreeTR; //! TreeTR (track references tree)
- TDirectoryFile *fDirK; //! Directory for Kine Tree
- TDirectoryFile *fDirTR; //! Directory for TR Tree
- TExMap fParticleSelected; //! List of selected MC particles for t
- TExMap fLabelMap; //! Stores the Map of MC (ESDLabel,AODlabel)
- Int_t fNEvent; //! Number of events
- Int_t fEvent; //! Current event
- TString *fPathName; //! Input file path
- const Char_t *fExtension; //! File name extension
- Int_t fFileNumber; //! Input file number
- Int_t fEventsPerFile; //! Number of events per file
- Bool_t fReadTR; // determines if TR shall be read
- Bool_t fInitOk; // Initialization ok
+ AliMCEvent *fMCEvent; //! MC Event
+ TFile *fFileE; //! File with TreeE
+ TFile *fFileK; //! File with TreeK
+ TFile *fFileTR; //! File with TreeTR
+ TTree *fTreeE; //! TreeE (Event Headers)
+ TTree *fTreeK; //! TreeK (kinematics tree)
+ TTree *fTreeTR; //! TreeTR (track references tree)
+ TDirectoryFile *fDirK; //! Directory for Kine Tree
+ TDirectoryFile *fDirTR; //! Directory for TR Tree
+ TExMap fParticleSelected; //! List of selected MC particles for t
+ TExMap fLabelMap; //! Stores the Map of MC (ESDLabel,AODlabel)
+ Int_t fNEvent; //! Number of events
+ Int_t fEvent; //! Current event
+ TString *fPathName; //! Input file path
+ const Char_t *fExtension; //! File name extension
+ Int_t fFileNumber; //! Input file number
+ Int_t fEventsPerFile; //! Number of events per file
+ Bool_t fReadTR; // determines if TR shall be read
+ Bool_t fInitOk; // Initialization ok
+ TList *fSubsidiaryHandlers; //! List of subsidiary MC handlers (for example for Background)
+ Int_t fEventsInContainer; //! Number of events in container class
ClassDef(AliMCEventHandler,1) //MC Truth EventHandler class
};
#endif
fParticle(0),
fTrackReferences(0),
fNTrackRef(0),
- fLabel(-1)
+ fLabel(-1),
+ fMother(-1),
+ fFirstDaughter(-1),
+ fLastDaughter(-1)
{
// Constructor
}
fParticle(part),
fTrackReferences(rarray),
fNTrackRef(0),
- fLabel(index)
+ fLabel(index),
+ fMother(-1),
+ fFirstDaughter(-1),
+ fLastDaughter(-1)
{
// Constructor
if (rarray != 0) {
fParticle(0),
fTrackReferences(0),
fNTrackRef(0),
- fLabel(-1)
+ fLabel(-1),
+ fMother(-1),
+ fFirstDaughter(-1),
+ fLastDaughter(-1)
{
// Copy constructor
}
// "Trackable" criteria
Float_t GetTPCTrackLength(Float_t bz, Float_t ptmin, Int_t &counter, Float_t deadWidth);
-
+ // Navigation
+ Int_t GetMother() const {return fMother;}
+ Int_t GetFirstDaughter() const {return fFirstDaughter;}
+ Int_t GetLastDaughter() const {return fLastDaughter;}
+ void SetMother(Int_t idx) {fMother = idx;}
+ void SetFirstDaughter(Int_t idx) {fFirstDaughter = idx;}
+ void SetLastDaughter(Int_t idx) {fLastDaughter = idx;}
+ void SetLabel(Int_t label) {fLabel = label;}
+
private:
TParticle *fParticle; // The wrapped TParticle
TRefArray *fTrackReferences; // Reference array to track references
Int_t fNTrackRef; // Number of track references
Int_t fLabel; // fParticle Label in the Stack
+ Int_t fMother; // Mother particles
+ Int_t fFirstDaughter; // First daughter
+ Int_t fLastDaughter; // LastDaughter
ClassDef(AliMCParticle,0) // AliVParticle realisation for MCParticles
};
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff