X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=STEER%2FAliAODHandler.cxx;h=8f50ed62b87cde32adf91cafcdb83019f2cc6550;hb=a2fbb6ecd59f0cf8c017124df53e30d16d06d727;hp=45bf43b2c50db58f1212916f82e71c3b4b5f96e2;hpb=dce1b636c0de7c9be243ca69350eae014069766d;p=u%2Fmrichter%2FAliRoot.git diff --git a/STEER/AliAODHandler.cxx b/STEER/AliAODHandler.cxx index 45bf43b2c50..8f50ed62b87 100644 --- a/STEER/AliAODHandler.cxx +++ b/STEER/AliAODHandler.cxx @@ -1,3 +1,4 @@ + /************************************************************************** * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. * * * @@ -25,6 +26,7 @@ #include #include #include +#include #include "AliLog.h" #include "AliAODHandler.h" @@ -64,7 +66,8 @@ AliAODHandler::AliAODHandler() : fMCEventH(NULL), fTreeA(NULL), fFileA(NULL), - fFileName("") + fFileName(""), + fExtensions(NULL) { // default constructor } @@ -88,13 +91,15 @@ AliAODHandler::AliAODHandler(const char* name, const char* title): fMCEventH(NULL), fTreeA(NULL), fFileA(NULL), - fFileName("") + fFileName(""), + fExtensions(NULL) { } //______________________________________________________________________________ AliAODHandler::~AliAODHandler() { + // Destructor. delete fAODEvent; if(fFileA){ // is already handled in TerminateIO @@ -102,7 +107,7 @@ AliAODHandler::~AliAODHandler() delete fFileA; } delete fTreeA; - // destructor + if (fExtensions) delete fExtensions; } //______________________________________________________________________________ @@ -119,27 +124,28 @@ Bool_t AliAODHandler::Init(Option_t* opt) // File opening according to execution mode TString option(opt); option.ToLower(); + TDirectory *owd = gDirectory; if (option.Contains("proof")) { // proof - if (option.Contains("special")) { - // File for tree already opened on slave -> merging via files - fFileA = gFile; - CreateTree(1); - } else { - // Merging in memory - CreateTree(0); - } + // 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)); + fFileA = gFile; } else { // local and grid - TDirectory *owd = gDirectory; fFileA = new TFile(fFileName.Data(), "RECREATE"); - CreateTree(1); - owd->cd(); } + CreateTree(1); + owd->cd(); + if (fExtensions) { + TIter next(fExtensions); + AliAODExtension *ext; + while ((ext=(AliAODExtension*)next())) ext->Init(option); + } return kTRUE; } - +//______________________________________________________________________________ void AliAODHandler::StoreMCParticles(){ // @@ -184,118 +190,118 @@ 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(); - - // 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 = (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)); } - } - 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... @@ -307,13 +313,17 @@ 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 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)); } } @@ -350,6 +360,7 @@ void AliAODHandler::StoreMCParticles(){ } +//______________________________________________________________________________ Bool_t AliAODHandler::FinishEvent() { // Fill data structures @@ -357,6 +368,14 @@ Bool_t AliAODHandler::FinishEvent() fAODEvent->MakeEntriesReferencable(); StoreMCParticles(); FillTree(); + if (fExtensions) { + TIter next(fExtensions); + AliAODExtension *ext; + while ((ext=(AliAODExtension*)next())) { + ext->GetAOD()->MakeEntriesReferencable(); + ext->GetTree()->Fill(); + } + } } if (fIsStandard) fAODEvent->ResetStd(); @@ -368,20 +387,31 @@ Bool_t AliAODHandler::FinishEvent() //______________________________________________________________________________ Bool_t AliAODHandler::Terminate() { - // Terminate - AddAODtoTreeUserInfo(); - return kTRUE; + // Terminate + AddAODtoTreeUserInfo(); + if (fExtensions) { + TIter next(fExtensions); + AliAODExtension *ext; + while ((ext=(AliAODExtension*)next())) ext->GetTree()->GetUserInfo()->Add(ext->GetAOD()); + } + return kTRUE; } //______________________________________________________________________________ Bool_t AliAODHandler::TerminateIO() { - // Terminate IO - if (fFileA) { - fFileA->Close(); - delete fFileA; - } - return kTRUE; + // Terminate IO + if (fFileA) { + fFileA->Close(); + delete fFileA; + fFileA = 0; + } + if (fExtensions) { + TIter next(fExtensions); + AliAODExtension *ext; + while ((ext=(AliAODExtension*)next())) ext->TerminateIO(); + } + return kTRUE; } //______________________________________________________________________________ @@ -408,9 +438,15 @@ void AliAODHandler::AddAODtoTreeUserInfo() } //______________________________________________________________________________ -void AliAODHandler::AddBranch(const char* cname, void* addobj) +void AliAODHandler::AddBranch(const char* cname, void* addobj, const char* filename) { - // Add a new branch to the aod + // Add a new branch to the aod. Added optional filename parameter if the + // branch should be written to a separate file. + if (strlen(filename)) { + AliAODExtension *ext = AddExtension(filename); + ext->AddBranch(cname, addobj); + return; + } TDirectory *owd = gDirectory; if (fFileA) { fFileA->cd(); @@ -434,6 +470,24 @@ void AliAODHandler::AddBranch(const char* cname, void* addobj) owd->cd(); } +//______________________________________________________________________________ +AliAODExtension *AliAODHandler::AddExtension(const char *filename, const char *title) +{ +// Add an AOD extension with some branches in a different file. + TString fname(filename); + if (!fname.EndsWith(".root")) fname += ".root"; + if (!fExtensions) { + fExtensions = new TObjArray(); + fExtensions->SetOwner(); + } + AliAODExtension *ext = (AliAODExtension*)fExtensions->FindObject(fname); + if (!ext) { + ext = new AliAODExtension(fname, title); + fExtensions->Add(ext); + } + return ext; +} + //______________________________________________________________________________ void AliAODHandler::SetOutputFileName(const char* fname) { @@ -448,6 +502,7 @@ const char *AliAODHandler::GetOutputFileName() return fFileName.Data(); } +//______________________________________________________________________________ void AliAODHandler::SetMCHeaderInfo(AliAODMCHeader *mcHeader,AliGenEventHeader *genHeader){ @@ -455,17 +510,14 @@ void AliAODHandler::SetMCHeaderInfo(AliAODMCHeader *mcHeader,AliGenEventHeader // Needed since different ProcessType and ImpactParamter are not // in the base class... // We don't encode process types for event cocktails yet - // coul be done e.g. by adding offsets depnding on the generator - + // could be done e.g. by adding offsets depnding on the generator mcHeader->AddGeneratorName(genHeader->GetName()); - - - if(!genHeader)return; AliGenPythiaEventHeader *pythiaGenHeader = dynamic_cast(genHeader); if (pythiaGenHeader) { mcHeader->SetEventType(pythiaGenHeader->ProcessType()); + mcHeader->SetPtHard(pythiaGenHeader->GetPtHard()); return; } @@ -485,3 +537,92 @@ void AliAODHandler::SetMCHeaderInfo(AliAODMCHeader *mcHeader,AliGenEventHeader AliWarning(Form("MC Eventheader not known: %s",genHeader->GetName())); } + +ClassImp(AliAODExtension) + +//------------------------------------------------------------------------- +// Support class for AOD extensions. This is created by the user analysis +// that requires a separate file for some AOD branches. The name of the +// AliAODExtension object is the file name where the AOD branches will be +// stored. +//------------------------------------------------------------------------- + +//______________________________________________________________________________ +AliAODExtension::~AliAODExtension() +{ +// Destructor. + delete fAODEvent; + if(fFileE){ + // is already handled in TerminateIO + fFileE->Close(); + delete fFileE; + } + delete fTreeE; +} + +//______________________________________________________________________________ +void AliAODExtension::AddBranch(const char* cname, void* addobj) +{ + // Add a new branch to the aod + 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(); + } + char** apointer = (char**) addobj; + TObject* obj = (TObject*) *apointer; + + fAODEvent->AddObject(obj); + + const Int_t kSplitlevel = 99; // default value in TTree::Branch() + const Int_t kBufsize = 32000; // default value in TTree::Branch() + + if (!fTreeE->FindBranch(obj->GetName())) { + // Do the same as if we book via + // TTree::Branch(TCollection*) + + fTreeE->Bronch(obj->GetName(), cname, fAODEvent->GetList()->GetObjectRef(obj), + kBufsize, kSplitlevel - 1); + // fTreeA->Branch(obj->GetName(), cname, addobj); + } + owd->cd(); +} + +//______________________________________________________________________________ +Bool_t AliAODExtension::Init(Option_t *option) +{ +// Initialize IO. + if(!fAODEvent) fAODEvent = new AliAODEvent(); + TDirectory *owd = gDirectory; + TString opt(option); + opt.ToLower(); + if (opt.Contains("proof")) { + // proof + // Merging via files. Need to access analysis manager via interpreter. + gROOT->ProcessLine(Form("AliAnalysisManager::GetAnalysisManager()->OpenProofFile(\"%s\", \"RECREATE\");", fName.Data())); + fFileE = gFile; + } else { + fFileE = new TFile(GetName(), "RECREATE"); + } + fTreeE = new TTree("aodTree", "AliAOD tree"); + fTreeE->Branch(fAODEvent->GetList()); + owd->cd(); + return kTRUE; +} + +//______________________________________________________________________________ +Bool_t AliAODExtension::TerminateIO() +{ + // Terminate IO + if (fFileE) { + fFileE->Write(); + fFileE->Close(); + delete fFileE; + fFileE = 0; + } + return kTRUE; +}