X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=STEER%2FAliAODHandler.cxx;h=16955922af49d823d31f89b0cb757f691e2af426;hb=e21df713c517d602f8df755137c5a742f33ec698;hp=a97464c2b9c7c1c2b0988c35c51e2ae0dbaae30f;hpb=a36d34b58a270a35103f38d302e4ab4bbe46b8bc;p=u%2Fmrichter%2FAliRoot.git diff --git a/STEER/AliAODHandler.cxx b/STEER/AliAODHandler.cxx index a97464c2b9c..16955922af4 100644 --- a/STEER/AliAODHandler.cxx +++ b/STEER/AliAODHandler.cxx @@ -1,3 +1,4 @@ + /************************************************************************** * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. * * * @@ -23,7 +24,6 @@ #include #include -#include #include #include #include @@ -61,13 +61,15 @@ AliAODHandler::AliAODHandler() : fNeedsJetsBranchReplication(kFALSE), fNeedsFMDClustersBranchReplication(kFALSE), fNeedsCaloClustersBranchReplication(kFALSE), + fNeedsMCParticlesBranchReplication(kFALSE), fAODIsReplicated(kFALSE), fAODEvent(NULL), fMCEventH(NULL), fTreeA(NULL), fFileA(NULL), fFileName(""), - fExtensions(NULL) + fExtensions(NULL), + fFilters(NULL) { // default constructor } @@ -86,28 +88,33 @@ AliAODHandler::AliAODHandler(const char* name, const char* title): fNeedsJetsBranchReplication(kFALSE), fNeedsFMDClustersBranchReplication(kFALSE), fNeedsCaloClustersBranchReplication(kFALSE), + fNeedsMCParticlesBranchReplication(kFALSE), fAODIsReplicated(kFALSE), fAODEvent(NULL), fMCEventH(NULL), fTreeA(NULL), fFileA(NULL), fFileName(""), - fExtensions(NULL) + fExtensions(NULL), + fFilters(NULL) { +// Normal constructor. } //______________________________________________________________________________ AliAODHandler::~AliAODHandler() { // Destructor. - delete fAODEvent; + if (fAODEvent) delete fAODEvent; if(fFileA){ // is already handled in TerminateIO fFileA->Close(); delete fFileA; + fTreeA = 0; } - delete fTreeA; - if (fExtensions) delete fExtensions; + if (fTreeA) delete fTreeA; + if (fExtensions) {fExtensions->Delete(); delete fExtensions;} + if (fFilters) {fFilters->Delete(); delete fFilters;} } //______________________________________________________________________________ @@ -128,8 +135,8 @@ Bool_t AliAODHandler::Init(Option_t* opt) if (option.Contains("proof")) { // proof // Merging via files. Need to access analysis manager via interpreter. - gROOT->ProcessLine(Form("AliAnalysisManager::GetAnalysisManager()->OpenProofFile(\"%s\", \"RECREATE\");", fFileName.Data())); - gROOT->ProcessLine(Form("AliAnalysisManager::GetAnalysisManager()->GetCommonOutputContainer()->SetFile((TFile*)0x%lx);", gFile)); + gROOT->ProcessLine(Form("AliAnalysisDataContainer *c_common_out = AliAnalysisManager::GetAnalysisManager()->GetCommonOutputContainer();")); + gROOT->ProcessLine(Form("AliAnalysisManager::GetAnalysisManager()->OpenProofFile(c_common_out, \"RECREATE\");")); fFileA = gFile; } else { // local and grid @@ -142,6 +149,14 @@ Bool_t AliAODHandler::Init(Option_t* opt) AliAODExtension *ext; while ((ext=(AliAODExtension*)next())) ext->Init(option); } + if (fFilters) { + TIter nextf(fFilters); + AliAODExtension *filteredAOD; + while ((filteredAOD=(AliAODExtension*)nextf())) { + filteredAOD->SetEvent(fAODEvent); + filteredAOD->Init(option); + } + } return kTRUE; } @@ -190,118 +205,119 @@ void AliAODHandler::StoreMCParticles(){ // AliHeader* header = fMCEventH->MCEvent()->Header(); - if (!header)return; - - // get the MC vertex - AliGenEventHeader* genHeader = header->GenEventHeader(); - TArrayF vtxMC(3); - genHeader->PrimaryVertex(vtxMC); - mcHeader->SetVertex(vtxMC[0],vtxMC[1],vtxMC[2]); - - // we search the MCEventHeaders first - // Two cases, cocktail or not... - AliGenCocktailEventHeader* genCocktailHeader = dynamic_cast(genHeader); - if(genCocktailHeader){ - // we have a coktail header - mcHeader->AddGeneratorName(genHeader->GetName()); - // Loop from the back so that the first one sets the process type - TList* headerList = genCocktailHeader->GetHeaders(); - for(int i = headerList->GetEntries()-1;i>=0;--i){ - AliGenEventHeader *headerEntry = dynamic_cast(headerList->At(i)); - SetMCHeaderInfo(mcHeader,headerEntry); - } - } - else{ - // No Cocktail just take the first one - SetMCHeaderInfo(mcHeader,genHeader); + // get the MC vertex + AliGenEventHeader* genHeader = 0; + if (header) genHeader = header->GenEventHeader(); + if (genHeader) { + TArrayF vtxMC(3); + genHeader->PrimaryVertex(vtxMC); + mcHeader->SetVertex(vtxMC[0],vtxMC[1],vtxMC[2]); + + // we search the MCEventHeaders first + // Two cases, cocktail or not... + AliGenCocktailEventHeader* genCocktailHeader = dynamic_cast(genHeader); + if(genCocktailHeader){ + // we have a coktail header + mcHeader->AddGeneratorName(genHeader->GetName()); + // Loop from the back so that the first one sets the process type + TList* headerList = genCocktailHeader->GetHeaders(); + for(int i = headerList->GetEntries()-1;i>=0;--i){ + AliGenEventHeader *headerEntry = dynamic_cast(headerList->At(i)); + SetMCHeaderInfo(mcHeader,headerEntry); + } + } + else{ + // No Cocktail just take the first one + SetMCHeaderInfo(mcHeader,genHeader); + } } + // 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(iIsPhysicalPrimary(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(); + for(int i = 0; i < np; ++i){ + if(fMCEventH->IsParticleSelected(i)){ + Int_t flag = 0; + AliMCParticle* mcpart = (AliMCParticle*) mcEvent->GetTrack(i); + if(iIsPhysicalPrimary(i))flag |= AliAODMCParticle::kPhysicalPrim; + + if(fMCEventH->GetNewLabel(i)!=j){ + AliError(Form("MISMATCH New label %d j: %d",fMCEventH->GetNewLabel(i),j)); + } - // 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); + AliAODMCParticle mcpart_tmp(mcpart,i,flag); + + mcpart_tmp.SetStatus(mcpart->Particle()->GetStatusCode()); + // + 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)); } - } - 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"); + 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... @@ -313,13 +329,21 @@ void AliAODHandler::StoreMCParticles(){ 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 sign = 1; + Int_t label = track->GetLabel(); + if(label<0){ // preserve the sign for later usage + label *= -1; + sign = -1; } - if(fMCEventH->GetNewLabel(track->GetLabel())==0){ - AliWarning(Form("New label not found for %d",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(label) == 0){ + AliWarning(Form("New label not found for %5d (%5d)",track->GetLabel(), label)); } - track->SetLabel(fMCEventH->GetNewLabel(track->GetLabel())); + track->SetLabel(sign*fMCEventH->GetNewLabel(label)); } } @@ -362,18 +386,22 @@ Bool_t AliAODHandler::FinishEvent() // Fill data structures if(fFillAOD){ fAODEvent->MakeEntriesReferencable(); - StoreMCParticles(); + // StoreMCParticles(); FillTree(); if (fExtensions) { TIter next(fExtensions); AliAODExtension *ext; - while ((ext=(AliAODExtension*)next())) { - ext->GetAOD()->MakeEntriesReferencable(); - ext->GetTree()->Fill(); + while ((ext=(AliAODExtension*)next())) ext->FinishEvent(); + } + if (fFilters) { + TIter nextf(fFilters); + AliAODExtension *ext; + while ((ext=(AliAODExtension*)nextf())) { +// ext->SetEvent(fAODEvent); + ext->FinishEvent(); } } } - if (fIsStandard) fAODEvent->ResetStd(); // Reset AOD replication flag fAODIsReplicated = kFALSE; @@ -401,12 +429,19 @@ Bool_t AliAODHandler::TerminateIO() fFileA->Close(); delete fFileA; fFileA = 0; + // When closing the file, the tree is also deleted. + fTreeA = 0; } if (fExtensions) { TIter next(fExtensions); AliAODExtension *ext; while ((ext=(AliAODExtension*)next())) ext->TerminateIO(); } + if (fFilters) { + TIter nextf(fFilters); + AliAODExtension *ext; + while ((ext=(AliAODExtension*)nextf())) ext->TerminateIO(); + } return kTRUE; } @@ -416,7 +451,6 @@ void AliAODHandler::CreateTree(Int_t flag) // Creates the AOD Tree fTreeA = new TTree("aodTree", "AliAOD tree"); fTreeA->Branch(fAODEvent->GetList()); - TRef junk = (TObject*)fTreeA->BranchRef(); if (flag == 0) fTreeA->SetDirectory(0); } @@ -432,6 +466,8 @@ void AliAODHandler::AddAODtoTreeUserInfo() { // Add aod event to tree user info fTreeA->GetUserInfo()->Add(fAODEvent); + // Now the tree owns our fAODEvent... + fAODEvent = 0; } //______________________________________________________________________________ @@ -484,6 +520,38 @@ AliAODExtension *AliAODHandler::AddExtension(const char *filename, const char *t } return ext; } + +//______________________________________________________________________________ +AliAODExtension *AliAODHandler::GetExtension(const char *filename) const +{ +// Getter for AOD extensions via file name. + if (!fExtensions) return NULL; + return (AliAODExtension*)fExtensions->FindObject(filename); +} + +//______________________________________________________________________________ +AliAODExtension *AliAODHandler::AddFilteredAOD(const char *filename, const char *filtername) +{ +// Add an AOD extension that can write only AOD events that pass a user filter. + if (!fFilters) { + fFilters = new TObjArray(); + fFilters->SetOwner(); + } + AliAODExtension *filter = (AliAODExtension*)fFilters->FindObject(filename); + if (!filter) { + filter = new AliAODExtension(filename, filtername, kTRUE); + fFilters->Add(filter); + } + return filter; +} + +//______________________________________________________________________________ +AliAODExtension *AliAODHandler::GetFilteredAOD(const char *filename) const +{ +// Getter for AOD filters via file name. + if (!fFilters) return NULL; + return (AliAODExtension*)fFilters->FindObject(filename); +} //______________________________________________________________________________ void AliAODHandler::SetOutputFileName(const char* fname) @@ -514,6 +582,7 @@ void AliAODHandler::SetMCHeaderInfo(AliAODMCHeader *mcHeader,AliGenEventHeader AliGenPythiaEventHeader *pythiaGenHeader = dynamic_cast(genHeader); if (pythiaGenHeader) { mcHeader->SetEventType(pythiaGenHeader->ProcessType()); + mcHeader->SetPtHard(pythiaGenHeader->GetPtHard()); return; } @@ -543,24 +612,50 @@ ClassImp(AliAODExtension) // stored. //------------------------------------------------------------------------- +//______________________________________________________________________________ +AliAODExtension::AliAODExtension(const char* name, const char* title, Bool_t isfilter) + :TNamed(name,title), + fAODEvent(0), + fTreeE(0), + fFileE(0), + fNtotal(0), + fNpassed(0), + fSelected(kFALSE) +{ +// Constructor. + if (isfilter) { + TObject::SetBit(kFilteredAOD); + printf("####### Added AOD filter %s\n", name); + } else printf("####### Added AOD extension %s\n", name); +} + //______________________________________________________________________________ AliAODExtension::~AliAODExtension() { // Destructor. - delete fAODEvent; if(fFileE){ // is already handled in TerminateIO fFileE->Close(); delete fFileE; + fTreeE = 0; + fAODEvent = 0; } - delete fTreeE; + if (fTreeE) delete fTreeE; } //______________________________________________________________________________ void AliAODExtension::AddBranch(const char* cname, void* addobj) { // Add a new branch to the aod - if (!fAODEvent) Init(""); + if (IsFilteredAOD()) { + Error("AddBranch", "Not allowed to add branched to filtered AOD's."); + return; + } + if (!fAODEvent) { + char type[20]; + gROOT->ProcessLine(Form("TString s_tmp; AliAnalysisManager::GetAnalysisManager()->GetAnalysisTypeString(s_tmp); sprintf((char*)0x%lx, \"%%s\", s_tmp.Data());", type)); + Init(type); + } TDirectory *owd = gDirectory; if (fFileE) { fFileE->cd(); @@ -584,6 +679,24 @@ void AliAODExtension::AddBranch(const char* cname, void* addobj) owd->cd(); } +//______________________________________________________________________________ +Bool_t AliAODExtension::FinishEvent() +{ +// Fill current event. + fNtotal++; + if (!IsFilteredAOD()) { + fAODEvent->MakeEntriesReferencable(); + fTreeE->Fill(); + return kTRUE; + } + // Filtered AOD. Fill only if event is selected. + if (!fSelected) return kTRUE; + fNpassed++; + fTreeE->Fill(); + fSelected = kFALSE; // so that next event will not be selected unless demanded + return kTRUE; +} + //______________________________________________________________________________ Bool_t AliAODExtension::Init(Option_t *option) { @@ -602,20 +715,40 @@ Bool_t AliAODExtension::Init(Option_t *option) } fTreeE = new TTree("aodTree", "AliAOD tree"); fTreeE->Branch(fAODEvent->GetList()); - TRef junk = (TObject*)fTreeE->BranchRef(); owd->cd(); return kTRUE; } +//______________________________________________________________________________ +void AliAODExtension::SetEvent(AliAODEvent *event) +{ +// Connects to an external event + if (!IsFilteredAOD()) { + Error("SetEvent", "Not allowed to set external event for filtered AOD's"); + return; + } + // Use the copy constructor or assignment operator to synchronize with external event. +// AliAODEvent &other = *event; +// if (!fAODEvent) fAODEvent = new AliAODEvent(other); +// else if (fSelected) *fAODEvent = other; + fAODEvent = event; +} + //______________________________________________________________________________ Bool_t AliAODExtension::TerminateIO() { // Terminate IO + if (TObject::TestBit(kFilteredAOD)) + printf("AOD Filter %s: events processed: %d passed: %d\n", GetName(), fNtotal, fNpassed); + else + printf("AOD extension %s: events processed: %d\n", GetName(), fNtotal); if (fFileE) { fFileE->Write(); fFileE->Close(); delete fFileE; fFileE = 0; + fTreeE = 0; + fAODEvent = 0; } return kTRUE; }