From: dainese Date: Wed, 24 Jun 2009 15:49:30 +0000 (+0000) Subject: New task for D0 mass distributions (Chiara B) X-Git-Url: http://git.uio.no/git/?p=u%2Fmrichter%2FAliRoot.git;a=commitdiff_plain;h=490611769248c13d91bd2fb4de6a215067354faa New task for D0 mass distributions (Chiara B) --- diff --git a/PWG3/CMake_libPWG3vertexingHF.txt b/PWG3/CMake_libPWG3vertexingHF.txt index 53de11d74c3..0e9fd8847fd 100644 --- a/PWG3/CMake_libPWG3vertexingHF.txt +++ b/PWG3/CMake_libPWG3vertexingHF.txt @@ -10,6 +10,7 @@ set(SRCS vertexingHF/AliAnalysisTaskSESelectHF.cxx vertexingHF/AliAnalysisTaskSECompareHF.cxx vertexingHF/AliAnalysisTaskSEDplus.cxx + vertexingHF/AliAnalysisTaskSED0Mass.cxx vertexingHF/AliAnalysisTaskCharmFraction.cxx vertexingHF/AliAnalysisTaskSEBkgLikeSignJPSI.cxx vertexingHF/AliAnalysisTaskSEBkgLikeSignD0.cxx diff --git a/PWG3/PWG3vertexingHFLinkDef.h b/PWG3/PWG3vertexingHFLinkDef.h index b8e86ccac41..8429f7cbce3 100644 --- a/PWG3/PWG3vertexingHFLinkDef.h +++ b/PWG3/PWG3vertexingHFLinkDef.h @@ -15,6 +15,7 @@ #pragma link C++ class AliAnalysisTaskSESelectHF+; #pragma link C++ class AliAnalysisTaskSECompareHF+; #pragma link C++ class AliAnalysisTaskSEDplus+; +#pragma link C++ class AliAnalysisTaskSED0Mass+; #pragma link C++ class AliAnalysisTaskCharmFraction+; #pragma link C++ class AliCFHeavyFlavourTask+; #pragma link C++ class AliCFHeavyFlavourTaskMultiVar+; diff --git a/PWG3/libPWG3vertexingHF.pkg b/PWG3/libPWG3vertexingHF.pkg index 14a37f51a2c..52df6731334 100644 --- a/PWG3/libPWG3vertexingHF.pkg +++ b/PWG3/libPWG3vertexingHF.pkg @@ -9,6 +9,7 @@ SRCS:= vertexingHF/AliAODRecoDecayHF.cxx \ vertexingHF/AliAnalysisTaskSESelectHF.cxx \ vertexingHF/AliAnalysisTaskSECompareHF.cxx \ vertexingHF/AliAnalysisTaskSEDplus.cxx \ + vertexingHF/AliAnalysisTaskSED0Mass.cxx \ vertexingHF/AliAnalysisTaskCharmFraction.cxx \ vertexingHF/AliCFHeavyFlavourTask.cxx \ vertexingHF/AliCFHeavyFlavourTaskMultiVar.cxx \ diff --git a/PWG3/vertexingHF/AddTaskD0Mass.C b/PWG3/vertexingHF/AddTaskD0Mass.C new file mode 100644 index 00000000000..66d77716d19 --- /dev/null +++ b/PWG3/vertexingHF/AddTaskD0Mass.C @@ -0,0 +1,37 @@ +AliAnalysisTaskSED0Mass *AddTaskD0Mass() +{ + // + // Test macro for the AliAnalysisTaskSE for D0 candidates + // invariant mass histogram and association with MC truth + // (using MC info in AOD) + // C.Bianchin chiara.bianchin@pd.infn.it + // + + + // Get the pointer to the existing analysis manager via the static access method. + //============================================================================== + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + if (!mgr) { + ::Error("AddTaskD0Mass", "No analysis manager to connect to."); + return NULL; + } + + + // Aanalysis task + AliAnalysisTaskSED0Mass *massD0Task = new AliAnalysisTaskSED0Mass("D0MassAnalysis"); + massD0Task->SetDebugLevel(2); + mgr->AddTask(massD0Task); + + // + // Create containers for input/output + AliAnalysisDataContainer *cinputmassD0 = mgr->CreateContainer("cinputmassD0",TChain::Class(), + AliAnalysisManager::kInputContainer); + AliAnalysisDataContainer *coutputmassD0 = mgr->CreateContainer("coutputmassD0",TList::Class(), + AliAnalysisManager::kOutputContainer, + "D0InvMass.root"); + mgr->ConnectInput(massD0Task,0,mgr->GetCommonInputContainer()); + + mgr->ConnectOutput(massD0Task,1,coutputmassD0); + + return massD0Task; +} diff --git a/PWG3/vertexingHF/AliAnalysisTaskSED0Mass.cxx b/PWG3/vertexingHF/AliAnalysisTaskSED0Mass.cxx new file mode 100644 index 00000000000..ede921ecef9 --- /dev/null +++ b/PWG3/vertexingHF/AliAnalysisTaskSED0Mass.cxx @@ -0,0 +1,328 @@ +/************************************************************************** + * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. * + * * + * Author: The ALICE Off-line Project. * + * Contributors are mentioned in the code where appropriate. * + * * + * Permission to use, copy, modify and distribute this software and its * + * documentation strictly for non-commercial purposes is hereby granted * + * without fee, provided that the above copyright notice appears in all * + * copies and that both the copyright notice and this permission notice * + * appear in the supporting documentation. The authors make no claims * + * about the suitability of this software for any purpose. It is * + * provided "as is" without express or implied warranty. * + **************************************************************************/ + +///////////////////////////////////////////////////////////// +// +// AliAnalysisTaskSE for D0 candidates invariant mass histogram +// and comparison with the MC truth. +// +// Authors: A.Dainese, andrea.dainese@lnl.infn.it +// and Chiara Bianchin, chiara.bianchin@pd.infn.it +///////////////////////////////////////////////////////////// + +#include +#include +#include +#include +#include + + +#include "AliAODEvent.h" +#include "AliAODVertex.h" +#include "AliAODTrack.h" +#include "AliAODMCHeader.h" +#include "AliAODMCParticle.h" +#include "AliAODRecoDecayHF2Prong.h" +#include "AliAODRecoCascadeHF.h" +#include "AliAnalysisVertexingHF.h" +#include "AliAnalysisTaskSE.h" +#include "AliAnalysisTaskSED0Mass.h" + +ClassImp(AliAnalysisTaskSED0Mass) + + +//________________________________________________________________________ +AliAnalysisTaskSED0Mass::AliAnalysisTaskSED0Mass(): +AliAnalysisTaskSE(), +fOutput(0), +//fNtupleD0Cmp(0), +fhistMass(0), +fhistSgn(0), +fhistBkg(0), +fVHF(0) +{ + // Default constructor +} + +//________________________________________________________________________ +AliAnalysisTaskSED0Mass::AliAnalysisTaskSED0Mass(const char *name): +AliAnalysisTaskSE(name), +fOutput(0), +//fNtupleD0Cmp(0), +fhistMass(0), +fhistSgn(0), +fhistBkg(0), +fVHF(0) +{ + // Default constructor + + // Output slot #1 writes into a TList container + DefineOutput(1,TList::Class()); //My private output +} + +//________________________________________________________________________ +AliAnalysisTaskSED0Mass::~AliAnalysisTaskSED0Mass() +{ + // Destructor + + if (fhistMass){ + delete fhistMass; + fhistMass=0; + } + + if (fhistSgn){ + delete fhistSgn; + fhistSgn=0; + } + + if (fhistBkg){ + delete fhistBkg; + fhistBkg=0; + } + + if (fOutput) { + delete fOutput; + fOutput = 0; + } + if (fVHF) { + delete fVHF; + fVHF = 0; + } + +} + +//________________________________________________________________________ +void AliAnalysisTaskSED0Mass::Init() +{ + // Initialization + + if(fDebug > 1) printf("AnalysisTaskSED0Mass::Init() \n"); + + gROOT->LoadMacro("$ALICE_ROOT/PWG3/vertexingHF/ConfigVertexingHF.C"); + + fVHF = (AliAnalysisVertexingHF*)gROOT->ProcessLine("ConfigVertexingHF()"); + // set dedidcated cuts + fVHF->SetD0toKpiCuts(0.7,0.03,0.8,0.06,0.06,0.05,0.05,-0.0002,0.6); //2.p-p vertex reconstructed + fVHF->PrintStatus(); + + return; +} + +//________________________________________________________________________ +void AliAnalysisTaskSED0Mass::UserCreateOutputObjects() +{ + + // Create the output container + // + if(fDebug > 1) printf("AnalysisTaskSED0Mass::UserCreateOutputObjects() \n"); + + // Several histograms are more conveniently managed in a TList + fOutput = new TList(); + fOutput->SetOwner(); + + const Int_t nhist=3; + fhistMass=new TClonesArray("TH1F",nhist); + fhistMass->SetName("fhistMass"); + fhistSgn =new TClonesArray("TH1F",nhist); + fhistSgn->SetName("fhistSgn"); + fhistBkg =new TClonesArray("TH1F",nhist); + fhistBkg->SetName("fhistBkg"); + + TString nameMass, nameSgn, nameBkg; + + for(Int_t i=0;iUncheckedAt(i)->GetName(),fhistSgn->UncheckedAt(i)->GetName(),fhistBkg->UncheckedAt(i)->GetName()); + } + + fOutput->Add(fhistMass); + fOutput->Add(fhistSgn); + fOutput->Add(fhistBkg); + + return; +} + +//________________________________________________________________________ +void AliAnalysisTaskSED0Mass::UserExec(Option_t */*option*/) +{ + // Execute analysis for current event: + // heavy flavor candidates association to MC truth + + AliAODEvent *aod = dynamic_cast (InputEvent()); + + // load D0->Kpi candidates + TClonesArray *inputArrayD0toKpi = + (TClonesArray*)aod->GetList()->FindObject("D0toKpi"); + if(!inputArrayD0toKpi) { + printf("AliAnalysisTaskSECompareHFpt::UserExec: D0toKpi branch not found!\n"); + return; + } + + // AOD primary vertex + AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex(); + //vtx1->Print(); + + // load MC particles + TClonesArray *mcArray = + (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName()); + if(!mcArray) { + printf("AliAnalysisTaskSED0Mass::UserExec: MC particles branch not found!\n"); + return; + } + + // load MC header + AliAODMCHeader *mcHeader = + (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName()); + if(!mcHeader) { + printf("AliAnalysisTaskSED0Mass::UserExec: MC header branch not found!\n"); + return; + } + + + // loop over D0->Kpi candidates + Int_t nInD0toKpi = inputArrayD0toKpi->GetEntriesFast(); + printf("Number of D0->Kpi: %d\n",nInD0toKpi); + + for (Int_t iD0toKpi = 0; iD0toKpi < nInD0toKpi; iD0toKpi++) { + + AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)inputArrayD0toKpi->UncheckedAt(iD0toKpi); + Bool_t unsetvtx=kFALSE; + if(!d->GetOwnPrimaryVtx()) { + d->SetOwnPrimaryVtx(vtx1); // needed to compute all variables + unsetvtx=kTRUE; + } + + // cout<<"Cuts applied: "; + // Double_t *cutsapplied=(Double_t*)fVHF->GetD0toKpiCuts(); + // for (Int_t i=0;i<9;i++){ + // printf(cutsapplied[i]\t); + // } + // cout<Pt(); + Int_t ptbin=0; + + //cout<<"P_t = "<0. && pt<=1.) { + ptbin=1; + //printf("I'm in the bin %d\n",ptbin); + } + else { + if(pt>1. && pt<=3.) { + ptbin=2; + //printf("I'm in the bin %d\n",ptbin); + } + else { + ptbin=3; + //printf("I'm in the bin %d\n",ptbin); + }//if(pt>3) + } + + FillHists(ptbin,d,mcArray); + + if(unsetvtx) d->UnsetOwnPrimaryVtx(); + + } + + + + // Post the data + PostData(1,fOutput); + + + return; +} +//____________________________________________________________________________* +void AliAnalysisTaskSED0Mass::FillHists(Int_t ptbin, AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC){ + // + // function used in UserExec: + // + Int_t okD0=0,okD0bar=0; + + if(part->SelectD0(fVHF->GetD0toKpiCuts(),okD0,okD0bar)) {//selected + Double_t invmassD0 = part->InvMassD0(), invmassD0bar = part->InvMassD0bar(); + //printf("SELECTED\n"); + Int_t labD0 = part->MatchToMC(421,arrMC); //return MC particle label if the array corresponds to a D0, -1 if not (cf. AliAODRecoDecay.cxx) + //printf("labD0 %d",labD0); + + if(labD0>=0) { + + AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0); + + Int_t pdgD0 = partD0->GetPdgCode(); + + //printf(" pdgD0 %d\n",pdgD0); + + if (pdgD0==421){ //D0 + //cout<<"Fill S with D0"<UncheckedAt(ptbin-1)))->Fill(invmassD0); + } + else {//D0bar + //printf("Fill S with D0bar"); + ((TH1F*)(fhistSgn->UncheckedAt(ptbin-1)))->Fill(invmassD0bar); + } + + } + else {//background + //printf("Fill background"); + ((TH1F*)(fhistBkg->UncheckedAt(ptbin-1)))->Fill(invmassD0); + ((TH1F*)(fhistBkg->UncheckedAt(ptbin-1)))->Fill(invmassD0bar); + } + + //no MC info, just cut selection + if (okD0==1) { + //printf("Fill mass with D0"); + ((TH1F*)(fhistMass->UncheckedAt(ptbin-1)))->Fill(invmassD0); + + } + if (okD0bar==1) { + ((TH1F*)fhistMass->UncheckedAt(ptbin-1))->Fill(invmassD0bar); + //printf("Fill mass with D0bar"); + } + + } //else cout<<"NOT SELECTED"< 1) printf("AnalysisTaskSED0Mass: Terminate() \n"); + + fOutput = dynamic_cast (GetOutputData(1)); + if (!fOutput) { + printf("ERROR: fOutput not available\n"); + return; + } + + + fhistMass = dynamic_cast(fOutput->FindObject("fhistMass")); + fhistSgn = dynamic_cast(fOutput->FindObject("fhistSgn")); + fhistBkg = dynamic_cast(fOutput->FindObject("fhistBkg")); + + return; +} + diff --git a/PWG3/vertexingHF/AliAnalysisTaskSED0Mass.h b/PWG3/vertexingHF/AliAnalysisTaskSED0Mass.h new file mode 100644 index 00000000000..631d5c7bd7a --- /dev/null +++ b/PWG3/vertexingHF/AliAnalysisTaskSED0Mass.h @@ -0,0 +1,57 @@ +#ifndef ALIANALYSISTASKSED0MASS_H +#define ALIANALYSISTASKSED0MASS_H + +/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. * + * See cxx source for full Copyright notice */ + +//************************************************************************* +// Class AliAnalysisTaskSED0Mass +// AliAnalysisTaskSE for D0 candidates invariant mass histogram +// and comparison to MC truth (kinematics stored in the AOD) +// Authors: A.Dainese, andrea.dainese@ln.infn.it +// and C.Bianchin, chiara.bianchin@pd.infn.it +//************************************************************************* + +#include +#include +#include +#include + +#include "AliAnalysisTaskSE.h" +#include "AliAnalysisVertexingHF.h" + +class AliAnalysisTaskSED0Mass : public AliAnalysisTaskSE +{ + public: + + AliAnalysisTaskSED0Mass(); + AliAnalysisTaskSED0Mass(const char *name); + virtual ~AliAnalysisTaskSED0Mass(); + + + // Implementation of interface methods + virtual void UserCreateOutputObjects(); + virtual void Init(); + virtual void LocalInit() {Init();} + virtual void UserExec(Option_t *option); + virtual void Terminate(Option_t *option); + + + private: + + AliAnalysisTaskSED0Mass(const AliAnalysisTaskSED0Mass &source); + AliAnalysisTaskSED0Mass& operator=(const AliAnalysisTaskSED0Mass& source); + void FillHists(Int_t ptbin, AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC); + TList *fOutput; //! list send on output slot 0 + //TNtuple *fNtupleD0Cmp; // output ntuple + TClonesArray *fhistMass; // output total invariant mass histogram - no MC truth + TClonesArray *fhistSgn; // output signal invariant mass histogram - use cuts + TClonesArray *fhistBkg; // output background invariant mass histogram - use cuts + + AliAnalysisVertexingHF *fVHF; // Vertexer heavy flavour (used to pass the cuts) + + ClassDef(AliAnalysisTaskSED0Mass,1); // AliAnalysisTaskSE for the MC association of heavy-flavour decay candidates +}; + +#endif +