X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=ANALYSIS%2FAliAnalysisTaskMCParticleFilter.cxx;h=7f6bae2200d94548dfdb49a254ff58e2f26d8fc5;hb=7fac86691d69b311589e0c594e64aa092475e568;hp=354b3a050f9dce9c6f397329d71b63f72d1904f3;hpb=3189f57dd8e1f6224ef9754556b4ec80a2bad2e7;p=u%2Fmrichter%2FAliRoot.git diff --git a/ANALYSIS/AliAnalysisTaskMCParticleFilter.cxx b/ANALYSIS/AliAnalysisTaskMCParticleFilter.cxx index 354b3a050f9..7f6bae2200d 100644 --- a/ANALYSIS/AliAnalysisTaskMCParticleFilter.cxx +++ b/ANALYSIS/AliAnalysisTaskMCParticleFilter.cxx @@ -35,6 +35,7 @@ #include "AliStack.h" #include "AliMCEvent.h" #include "AliMCEventHandler.h" +#include "AliESDInputHandler.h" #include "AliAODEvent.h" #include "AliAODHeader.h" #include "AliAODMCHeader.h" @@ -46,6 +47,7 @@ #include "AliGenHijingEventHeader.h" #include "AliGenPythiaEventHeader.h" #include "AliGenCocktailEventHeader.h" +#include "AliGenEventHeaderTunedPbPb.h" #include "AliLog.h" @@ -70,7 +72,8 @@ Bool_t AliAnalysisTaskMCParticleFilter::Notify() // Implemented Notify() to read the cross sections // from pyxsec.root // - TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree(); + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + TTree *tree = mgr->GetTree(); Double_t xsection = 0; UInt_t ntrials = 0; if(tree){ @@ -81,7 +84,11 @@ Bool_t AliAnalysisTaskMCParticleFilter::Notify() } TString fileName(curfile->GetName()); - if(fileName.Contains("AliESDs.root")){ + TString datafile = mgr->GetInputEventHandler()->GetInputFileName(); + if (fileName.Contains(datafile)) { + fileName.ReplaceAll(datafile, "pyxsec.root"); + } + else if(fileName.Contains("AliESDs.root")){ fileName.ReplaceAll("AliESDs.root", "pyxsec.root"); } else if(fileName.Contains("AliAOD.root")){ @@ -164,7 +171,7 @@ AliAnalysisTaskMCParticleFilter& AliAnalysisTaskMCParticleFilter::operator=(cons void AliAnalysisTaskMCParticleFilter::UserCreateOutputObjects() { // Create the output container - PostData(1,fHistList); + if (OutputTree()&&fTrackFilterMother) OutputTree()->GetUserInfo()->Add(fTrackFilterMother); @@ -212,17 +219,19 @@ void AliAnalysisTaskMCParticleFilter::UserCreateOutputObjects() h1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}"); fHistList->Add(h1Trials); + PostData(1,fHistList); } //____________________________________________________________________ void AliAnalysisTaskMCParticleFilter::UserExec(Option_t */*option*/) { -// Execute analysis for current event -// - -// Fill AOD tracks from Kinematic stack - + // Execute analysis for current event + // + + // Fill AOD tracks from Kinematic stack + PostData(1,fHistList); + // get AliAOD Event AliAODEvent* aod = AODEvent(); if (!aod) { @@ -269,12 +278,17 @@ void AliAnalysisTaskMCParticleFilter::UserExec(Option_t */*option*/) AliGenHijingEventHeader *hiH = 0; AliCollisionGeometry *colG = 0; AliGenDPMjetEventHeader *dpmH = 0; + AliGenEventHeaderTunedPbPb *tunedH = 0; + // it can be only one save some casts // assuming PYTHIA and HIJING are the most likely ones... if(!pyH){ hiH = dynamic_cast(mcEH); if(!hiH){ dpmH = dynamic_cast(mcEH); + if(!dpmH){ + tunedH = dynamic_cast(mcEH); + } } } @@ -290,10 +304,10 @@ void AliAnalysisTaskMCParticleFilter::UserExec(Option_t */*option*/) if (ccEH) { TList *genHeaders = ccEH->GetHeaders(); for (int imch=0; imchGetEntries(); imch++) { - if(!pyH)dynamic_cast(genHeaders->At(imch)); - if(!hiH)dynamic_cast(genHeaders->At(imch)); + if(!pyH)pyH = dynamic_cast(genHeaders->At(imch)); + if(!hiH)hiH = dynamic_cast(genHeaders->At(imch)); if(!colG)colG = dynamic_cast(genHeaders->At(imch)); - if(!dpmH)dynamic_cast(genHeaders->At(imch)); + if(!dpmH)dpmH = dynamic_cast(genHeaders->At(imch)); } } } @@ -311,13 +325,10 @@ void AliAnalysisTaskMCParticleFilter::UserExec(Option_t */*option*/) fAODMcHeader->SetReactionPlaneAngle(colG->ReactionPlaneAngle()); AliInfo(Form("Found Collision Geometry. Got Reaction Plane %lf\n", colG->ReactionPlaneAngle())); } - - - - // check varibales for charm need all daughters - static int iTaken = 0; - static int iAll = 0; - static int iCharm = 0; + else if(tunedH) { + fAODMcHeader->SetReactionPlaneAngle(tunedH->GetPsi2()); + fAODMcHeader->SetCrossSection(tunedH->GetCentrality()); + } Int_t j=0; @@ -344,7 +355,8 @@ void AliAnalysisTaskMCParticleFilter::UserExec(Option_t */*option*/) imo = mother->GetMother(); } // Select according to pseudorapidity and production point of primary ancestor - if (imo < nprim && Select(((AliMCParticle*) mcE->GetTrack(imo))->Particle(), rv, zv))write = kTRUE; + if (imo < nprim)write = kTRUE; + // if(!Select(((AliMCParticle*) mcE->GetTrack(imo))->Particle(), rv, zv))write = kFALSE; // selection on eta and phi of the mother } else if (part->GetUniqueID() == kPPair) { // Now look for pair production Int_t imo = mcpart->GetMother(); @@ -364,44 +376,49 @@ void AliAnalysisTaskMCParticleFilter::UserExec(Option_t */*option*/) write = kTRUE; } } - /* - else if (part->GetUniqueID() == 13){ - // Evaporation - // Check that we end at a primary particle - TParticle* mother = part; - Int_t imo = part->GetFirstMother(); - while((imo >= nprim) && ((mother->GetUniqueID() == 4 ) || ( mother->GetUniqueID() == 13))) { - mother = mcE->Stack()->Particle(imo); - imo = mother->GetFirstMother(); - } - // Select according to pseudorapidity and production point - if (imo < nprim && Select(mother, rv, zv)) - write = kTRUE; - } - */ + if (write) { if(mcH)mcH->SelectParticle(ip); - j++; - - // debug info to check fro charm daugthers - if((TMath::Abs(part->GetPdgCode()))==411){ - iCharm++; - Int_t d0 = mcpart->GetFirstDaughter(); - Int_t d1 = mcpart->GetLastDaughter(); - if(d0>0&&d1>0){ - for(int id = d0;id <= d1;id++){ - iAll++; - if(mcH->IsParticleSelected(id))iTaken++; + j++; + } + } + + /* + + // check varibales for charm need all daughters + static int iTaken = 0; + static int iAll = 0; + static int iCharm = 0; + + for (Int_t ip = 0; ip < np; ip++){ + AliMCParticle* mcpart = (AliMCParticle*) mcE->GetTrack(ip); + TParticle* part = mcpart->Particle(); + // debug info to check fro charm daugthers + if((TMath::Abs(part->GetPdgCode()))==411){ + iCharm++; + Int_t d0 = part->GetFirstDaughter(); + Int_t d1 = part->GetLastDaughter(); + if(d0>0&&d1>0){ + for(int id = d0;id <= d1;id++){ + iAll++; + if(mcH->IsParticleSelected(id)){ + iTaken++; + Printf("> %d/%d Taken",id,nprim); + PrintMCParticle( (AliMCParticle*)mcE->GetTrack(id),id); + } + else{ + Printf("> %d/%d NOT Taken",id,nprim); + PrintMCParticle( (AliMCParticle*)mcE->GetTrack(id),id); } } - }// if charm - } + } + }// if charm + AliInfo(Form("Taken daughters %d/%d of %d charm",iTaken,iAll,iCharm)); } + */ - AliInfo(Form("Taken daughters %d/%d of %d charm",iTaken,iAll,iCharm)); aodH->StoreMCParticles(); - PostData(1,fHistList); return; } @@ -461,3 +478,14 @@ Bool_t AliAnalysisTaskMCParticleFilter::Select(TParticle* part, Float_t rv, Floa } +void AliAnalysisTaskMCParticleFilter::PrintMCParticle(const AliMCParticle *mcp,Int_t np){ + + const TParticle *p = mcp->Particle(); + Printf("Nr.%d --- Status %d ---- Mother1 %d Mother2 %d Daughter1 %d Daughte\ +r2 %d ", + np,p->GetStatusCode(),p->GetMother(0),p->GetMother(1),p->GetDaughter \ + (0),p->GetDaughter(1)); + Printf("Eta %3.3f Phi %3.3f ",p->Eta(),p->Phi()); + p->Print(); + Printf("---------------------------------------"); +}