From 9b0cb3c31a546ef910fd5387cb313d8149d21fdc Mon Sep 17 00:00:00 2001 From: mfloris Date: Thu, 3 Jun 2010 13:29:56 +0000 Subject: [PATCH] First version of a class to compute the number of collisions corresponding to a data sample (AliCollisionNormalization) and of a corresponding analysis task, demonstrating usage (AliCollisionNormalizationTask). Macros to run the task are also provided. Minor modification fo AliPhysiscsSelection: the enum of stat table is now public. --- ANALYSIS/ANALYSISaliceLinkDef.h | 3 + ANALYSIS/AliCollisionNormalization.cxx | 626 +++++++++++++++++++++ ANALYSIS/AliCollisionNormalization.h | 105 ++++ ANALYSIS/AliCollisionNormalizationTask.cxx | 220 ++++++++ ANALYSIS/AliCollisionNormalizationTask.h | 47 ++ ANALYSIS/AliPhysicsSelection.h | 2 +- ANALYSIS/libANALYSISalice.pkg | 3 +- ANALYSIS/macros/runProofNormalization.C | 66 +++ ANALYSIS/macros/runTaskNormalization.C | 80 +++ 9 files changed, 1150 insertions(+), 2 deletions(-) create mode 100644 ANALYSIS/AliCollisionNormalization.cxx create mode 100644 ANALYSIS/AliCollisionNormalization.h create mode 100644 ANALYSIS/AliCollisionNormalizationTask.cxx create mode 100644 ANALYSIS/AliCollisionNormalizationTask.h create mode 100644 ANALYSIS/macros/runProofNormalization.C create mode 100644 ANALYSIS/macros/runTaskNormalization.C diff --git a/ANALYSIS/ANALYSISaliceLinkDef.h b/ANALYSIS/ANALYSISaliceLinkDef.h index 098dd1eb52d..1878e3deb6e 100644 --- a/ANALYSIS/ANALYSISaliceLinkDef.h +++ b/ANALYSIS/ANALYSISaliceLinkDef.h @@ -31,6 +31,9 @@ #pragma link C++ class AliAnalysisFilter+; #pragma link C++ class AliAnalysisCuts+; +#pragma link C++ class AliCollisionNormalization+; +#pragma link C++ class AliCollisionNormalizationTask+; + #ifdef WITHXML #pragma link C++ class AliTagAnalysis+; #pragma link C++ class AliXMLCollection+; diff --git a/ANALYSIS/AliCollisionNormalization.cxx b/ANALYSIS/AliCollisionNormalization.cxx new file mode 100644 index 00000000000..69928f49e05 --- /dev/null +++ b/ANALYSIS/AliCollisionNormalization.cxx @@ -0,0 +1,626 @@ +#include "AliCollisionNormalization.h" +#include "AliPhysicsSelection.h" +#include "AliLog.h" +#include "TFile.h" +#include "TCanvas.h" +#include "AliGenPythiaEventHeader.h" +#include "AliGenDPMjetEventHeader.h" +#include "AliGenEventHeader.h" +#include "AliMCEvent.h" + +ClassImp(AliCollisionNormalization) + +const char * AliCollisionNormalization::fProcLabel[] = {"SD","DD","ND", "Unknown"}; + +AliCollisionNormalization::AliCollisionNormalization() : + TObject(), + fNbinsVz(0), + fMinVz(0), + fMaxVz(0), + fZRange(9.99), + fIsMC(0), + fReferenceXS(0), + fVerbose(0), + fHistVzData (0), + fHistProcTypes (0), + fHistStatBin0 (0) +{ + + // ctor + + fNbinsVz = 30; + fMinVz = -15; + fMaxVz = +15; + + for(Int_t iproc = 0; iproc < kNProcs; iproc++){ + fHistVzMCGen[iproc] = 0; + fHistVzMCRec[iproc] = 0; + fHistVzMCTrg[iproc] = 0; + } + + + BookAllHistos(); +} + +AliCollisionNormalization::AliCollisionNormalization(Int_t nbinz, Float_t minz, Float_t maxz): + TObject(), + fNbinsVz(0), + fMinVz(0), + fMaxVz(0), + fZRange(9.99), + fIsMC(0), + fReferenceXS(0), + fVerbose(0), + fHistVzData (0), + fHistProcTypes (0), + fHistStatBin0 (0) +{ + + // ctor, allows setting binning + fNbinsVz = nbinz; + fMinVz = minz; + fMaxVz = maxz; + + for(Int_t iproc = 0; iproc < kNProcs; iproc++){ + fHistVzMCGen[iproc] = 0; + fHistVzMCRec[iproc] = 0; + fHistVzMCTrg[iproc] = 0; + } + + BookAllHistos(); +} + +AliCollisionNormalization::AliCollisionNormalization(const char * dataFile, const char * dataListName, + const char * mcFile, const char * mcListName, + const char * eventStatFile) : + TObject(), + fNbinsVz(0), + fMinVz(0), + fMaxVz(0), + fZRange(9.99), + fIsMC(0), + fReferenceXS(0), + fVerbose(0), + fHistVzData (0), + fHistProcTypes (0), + fHistStatBin0 (0) +{ + + // ctor, loads histograms from file + for(Int_t iproc = 0; iproc < kNProcs; iproc++){ + fHistVzMCGen[iproc] = 0; + fHistVzMCRec[iproc] = 0; + fHistVzMCTrg[iproc] = 0; + } + + + TFile * fdata = new TFile (dataFile); + TFile * fmc = new TFile (mcFile ); + + TList * ldata = (TList*) fdata->Get(dataListName); + TList * lmc = (TList*) fmc ->Get(mcListName ); + + AliCollisionNormalization * cndata = (AliCollisionNormalization*) ldata->FindObject("AliCollisionNormalization"); + AliCollisionNormalization * cnmc = (AliCollisionNormalization*) lmc ->FindObject("AliCollisionNormalization"); + + + // Assign or book all histos + for(Int_t iproc = 0; iproc < kNProcs; iproc++){ + fHistVzMCGen[iproc]= cnmc->GetVzMCGen(iproc); + fHistVzMCRec[iproc]= cnmc->GetVzMCRec(iproc); + fHistVzMCTrg[iproc]= cnmc->GetVzMCTrg(iproc); + } + fHistVzData = cndata->GetVzData(); + fHistProcTypes = cnmc->GetHistProcTypes(); + + TFile * fstat = new TFile(eventStatFile); + fHistStatBin0 = (TH1F*) fstat->Get("fHistStatistics_Bin0"); + +} + + +AliCollisionNormalization::~AliCollisionNormalization(){ + + // dtor + for(Int_t iproc = 0; iproc < kNProcs; iproc++){ + if(fHistVzMCGen[iproc]) { delete fHistVzMCGen[iproc] ; fHistVzMCGen[iproc] =0;} + if(fHistVzMCRec[iproc]) { delete fHistVzMCRec[iproc] ; fHistVzMCRec[iproc] =0;} + if(fHistVzMCTrg[iproc]) { delete fHistVzMCTrg[iproc] ; fHistVzMCTrg[iproc] =0;} + } + + if(fHistVzData ) { delete fHistVzData ; fHistVzData =0;} + if(fHistStatBin0 ) { delete fHistStatBin0 ; fHistStatBin0 =0;} + if(fHistProcTypes ) { delete fHistProcTypes ; fHistProcTypes =0;} + +} + +void AliCollisionNormalization::BookAllHistos(){ + // books all histos + // Book histos of vz distributions vs multiplicity + // if vzOnly == kTRUE, it returns a 1d histo with vz dist only + + // Do not attach histos to the current directory + Bool_t oldStatus = TH1::AddDirectoryStatus(); + TH1::AddDirectory(kFALSE); + + for(Int_t iproc = 0; iproc < kNProcs; iproc++){ + fHistVzMCGen [iproc] = (TH2F*) BookVzHisto(TString("fHistVzMCGen")+ fProcLabel[iproc] ,"Vz distribution of generated events vs rec multiplicity "); + fHistVzMCRec [iproc] = (TH2F*) BookVzHisto(TString("fHistVzMCRec")+ fProcLabel[iproc] ,"Vz distribution of reconstructed events vs rec multiplicity"); + fHistVzMCTrg [iproc] = (TH2F*) BookVzHisto(TString("fHistVzMCTrg")+ fProcLabel[iproc] ,"Vz distribution of triggered events vs rec multiplicity "); + } + fHistVzData = (TH2F*) BookVzHisto("fHistVzData" ,"Vz distribution of triggered events vs rec multiplicity "); + fHistProcTypes = new TH1F ("fHistProcTypes", "Number of events in the different process classes", kNProcs, -0.5 , kNProcs-0.5); + + fHistProcTypes->GetXaxis()->SetBinLabel(kProcSD+1,"SD"); + fHistProcTypes->GetXaxis()->SetBinLabel(kProcND+1,"ND"); + fHistProcTypes->GetXaxis()->SetBinLabel(kProcDD+1,"DD"); + fHistProcTypes->GetXaxis()->SetBinLabel(kProcUnknown+1,"Unknown"); + + TH1::AddDirectory(oldStatus); + +} + +TH1 * AliCollisionNormalization::BookVzHisto(const char * name , const char * title, Bool_t vzOnly) { + + // Book histos of vz distributions vs multiplicity + // if vzOnly == kTRUE, it returns a 1d histo with vz dist only + + + // Do not attach histos to the current directory + Bool_t oldStatus = TH1::AddDirectoryStatus(); + TH1::AddDirectory(kFALSE); + + TH1 * h; + Double_t binLimitsVtx[] = {-30,-25,-20,-15,-10,-7,-5.5,-4,-3,-2,-1,0,1,2,3,4,5.5,7,10,15,20,25,30}; + if (vzOnly) { + h = new TH1F(name,title,22,binLimitsVtx); + } else { + h = new TH2F(name,title,22,binLimitsVtx,100,Double_t(-0.5),Double_t(99.5)); + } + + h->SetXTitle("V_{z} (cm)"); + h->SetYTitle("n_{trk}"); + h->Sumw2(); + + TH1::AddDirectory(oldStatus); + + return h; + +} + +TH2F * AliCollisionNormalization::GetVzMCGen (Int_t procType) { + + // returns MC gen histo. If proc type < 0 sums over all proc types, reweighting the XS + + if(procType>=0) return fHistVzMCGen[procType] ; + + TH2F * sum = (TH2F*) fHistVzMCGen[kProcSD]->Clone(); + sum->Reset(); + + for(Int_t iproc = 0; iproc < kProcUnknown; iproc++){ + sum->Add(fHistVzMCGen[iproc],GetProcessWeight(iproc)); + } + + return sum; +} +TH2F * AliCollisionNormalization::GetVzMCRec (Int_t procType) { + + // returns MC rec histo. If proc type < 0 sums over all proc types, reweighting the XS + + if(procType>=0) return fHistVzMCRec[procType] ; + + TH2F * sum = (TH2F*) fHistVzMCRec[kProcSD]->Clone(); + sum->Reset(); + + for(Int_t iproc = 0; iproc < kProcUnknown; iproc++){ + sum->Add(fHistVzMCRec[iproc],GetProcessWeight(iproc)); + } + + return sum; + +} + + +TH2F * AliCollisionNormalization::GetVzMCTrg (Int_t procType) { + + // returns MC trg histo. If proc type < 0 sums over all proc types, reweighting the XS + + if(procType>=0) return fHistVzMCTrg[procType] ; + + TH2F * sum = (TH2F*) fHistVzMCTrg[kProcSD]->Clone(); + sum->Reset(); + + for(Int_t iproc = 0; iproc < kProcUnknown; iproc++){ + sum->Add(fHistVzMCTrg[iproc],GetProcessWeight(iproc)); + } + + return sum; + +} + +Double_t AliCollisionNormalization::ComputeNint() { + + // Compute the number of collisions applying all corrections + // TODO: check error propagation + + TH1 * hVzData = fHistVzData->ProjectionX("_px",2,-1,"E"); // Skip zero bin + Int_t allEventsWithVertex = (Int_t) fHistVzData->Integral(0, fHistVzData->GetNbinsX()+1, + 2, fHistVzData->GetNbinsY()+1); // include under/overflow!, skip 0 bin + + // Assign histos reweighted with measured XS + TH2F * histVzMCGen = GetVzMCGen(-1); + TH2F * histVzMCTrg = GetVzMCTrg(-1); + TH2F * histVzMCRec = GetVzMCRec(-1); + + // Get or compute BG. This assumes the CINT1B suite + Double_t triggeredEventsWith0MultWithBG = fHistVzData->Integral(0, fHistVzData->GetNbinsX()+1,1, 1); + Double_t bg = 0; // This will include beam gas + accidentals + if (fHistStatBin0->GetNbinsY() > 4) { // FIXME: we need a better criterion to decide... + AliInfo("Using BG computed by Physics Selection"); + bg = fHistStatBin0->GetBinContent(fHistStatBin0->GetNbinsX(),AliPhysicsSelection::kStatRowBG); + bg += fHistStatBin0->GetBinContent(fHistStatBin0->GetNbinsX(),AliPhysicsSelection::kStatRowAcc); + } else { + AliInfo("Computing BG using CINT1A/B/C/E, ignoring intensities"); + Int_t icol = fHistStatBin0->GetNbinsX(); + Int_t cint1B = (Int_t) fHistStatBin0->GetBinContent(icol,1); + Int_t cint1A = (Int_t) fHistStatBin0->GetBinContent(icol,2); + Int_t cint1C = (Int_t) fHistStatBin0->GetBinContent(icol,3); + Int_t cint1E = (Int_t) fHistStatBin0->GetBinContent(icol,4); + bg = cint1A + cint1C-2*cint1E ; + if (cint1B != triggeredEventsWith0MultWithBG) { + AliWarning(Form("Events in bin0 from physics selection and local counter not consistent: %d - %d", cint1B, triggeredEventsWith0MultWithBG)); + } + } + + Double_t triggeredEventsWith0Mult = triggeredEventsWith0MultWithBG - bg; + if(fVerbose > 0) AliInfo(Form("Measured events with vertex: %d",allEventsWithVertex)); + Double_t bin0 = fHistVzData->Integral(0, fHistVzData->GetNbinsX()+1, + 1, 1); + if(fVerbose > 0) AliInfo(Form("Zero Bin, Meas: %2.2f, BG: %2.2f, Meas - BG: %2.2f", bin0, bg, triggeredEventsWith0Mult)); + + // This pointers are here so that I use the same names used in Jan Fiete's code (allows for faster/easier comparison) + + TH2* eTrig = histVzMCTrg; + TH2* eTrigVtx = histVzMCRec; + TH1* eTrigVtx_projx = eTrigVtx->ProjectionX("eTrigVtx_projx", 2, eTrigVtx->GetNbinsX()+1); + + // compute trigger and vertex efficiency + TH2 * effVtxTrig = (TH2*) histVzMCRec->Clone("effVtxTrig"); + effVtxTrig->Reset(); + effVtxTrig->Divide(histVzMCRec,histVzMCGen,1,1,"B"); + // Apply correction to data to get corrected events + TH2 * correctedEvents = (TH2*) fHistVzData->Clone("correctedEvents"); + correctedEvents->Divide(effVtxTrig); + + // TH1 * correctedEvents = fHistVzData->ProjectionX("eTrigVtx_projx", 2, eTrigVtx->GetNbinsX()+1); + + // loop over vertex bins + for (Int_t i = 1; i <= fHistVzData->GetNbinsX(); i++) + { + Double_t alpha = (Double_t) hVzData->GetBinContent(i) / allEventsWithVertex; + Double_t events = alpha * triggeredEventsWith0Mult; + + if (histVzMCRec->GetBinContent(i,1) == 0) + continue; + + Double_t fZ = eTrigVtx_projx->Integral(0, eTrigVtx_projx->GetNbinsX()+1) / eTrigVtx_projx->GetBinContent(i) * + eTrig->GetBinContent(i, 1) / eTrig->Integral(0, eTrig->GetNbinsX()+1, 1, 1); + + events *= fZ; + + // multiply with trigger correction + Double_t correction = histVzMCGen->GetBinContent(i,1)/histVzMCTrg->GetBinContent(i,1); + events *= correction; + + if (fVerbose > 1) Printf(" Bin %d, alpha is %.2f%%, fZ is %.3f, number of events with 0 mult.: %.2f (MC comparison: %.2f)", i, alpha * 100., fZ, events, + histVzMCRec->GetBinContent(i,1)); + correctedEvents->SetBinContent(i, 1, events); + } + + + + // Integrate correctedEvents over full z range + Double_t allEvents = correctedEvents->Integral(0, correctedEvents->GetNbinsX()+1,0, correctedEvents->GetNbinsY()+1); + // Integrate correctedEvents over needed z range + Double_t allEventsZrange = correctedEvents->Integral(correctedEvents->GetXaxis()->FindBin(-fZRange), correctedEvents->GetXaxis()->FindBin(fZRange), + 0, correctedEvents->GetNbinsY()+1); + + if(fVerbose > 1) AliInfo(Form("Results in |Vz| < %3.3f",fZRange)); + if(fVerbose > 1) AliInfo(Form(" Events in Bin0: %2.2f, With > 1 track: %2.2f, All corrected: %2.2f", + correctedEvents->Integral(correctedEvents->GetXaxis()->FindBin(-fZRange), correctedEvents->GetXaxis()->FindBin(fZRange),1,1), + correctedEvents->Integral(correctedEvents->GetXaxis()->FindBin(-fZRange), correctedEvents->GetXaxis()->FindBin(fZRange), + 2,correctedEvents->GetNbinsX()+1), + allEventsZrange)); + if(fVerbose > 1) { + Int_t nbin = histVzMCRec->GetNbinsX(); + AliInfo(Form("Efficiency in the zero bin: %3.3f", histVzMCRec->Integral(0,nbin+1,1,1)/histVzMCGen->Integral(0,nbin+1,1,1) )); + } + + AliInfo(Form("Number of collisions in full phase space: %2.2f", allEvents)); + + return allEvents; +} + +Long64_t AliCollisionNormalization::Merge(TCollection* list) +{ + // Merge a list of AliCollisionNormalization objects with this + // (needed for PROOF). + // Returns the number of merged objects (including this). + + if (!list) + return 0; + + if (list->IsEmpty()) + return 1; + + TIterator* iter = list->MakeIterator(); + TObject* obj; + + // collections of all histograms + const Int_t nHists = kNProcs*3+4; + TList collections[nHists]; + + Int_t count = 0; + while ((obj = iter->Next())) { + + AliCollisionNormalization* entry = dynamic_cast (obj); + if (entry == 0) + continue; + Int_t ihist = -1; + + for(Int_t iproc = 0; iproc < kNProcs; iproc++){ + if (entry->fHistVzMCGen[iproc] ) collections[++ihist].Add(entry->fHistVzMCGen[iproc] ); + if (entry->fHistVzMCRec[iproc] ) collections[++ihist].Add(entry->fHistVzMCRec[iproc] ); + if (entry->fHistVzMCTrg[iproc] ) collections[++ihist].Add(entry->fHistVzMCTrg[iproc] ); + } + if (entry->fHistVzData ) collections[++ihist].Add(entry->fHistVzData ); + if (entry->fHistProcTypes ) collections[++ihist].Add(entry->fHistProcTypes ); + if (entry->fHistStatBin0 ) collections[++ihist].Add(entry->fHistStatBin0 ); + + count++; + } + + Int_t ihist = -1; + for(Int_t iproc = 0; iproc < kNProcs; iproc++){ + if (fHistVzMCGen[iproc] ) fHistVzMCGen[iproc] ->Merge(&collections[++ihist]); + if (fHistVzMCRec[iproc] ) fHistVzMCRec[iproc] ->Merge(&collections[++ihist]); + if (fHistVzMCTrg[iproc] ) fHistVzMCTrg[iproc] ->Merge(&collections[++ihist]); + + } + if (fHistVzData ) fHistVzData ->Merge(&collections[++ihist]); + if (fHistProcTypes ) fHistProcTypes ->Merge(&collections[++ihist]); + if (fHistStatBin0 ) fHistStatBin0 ->Merge(&collections[++ihist]); + + + delete iter; + + return count+1; +} + +void AliCollisionNormalization::FillVzMCGen(Float_t vz, Int_t ntrk, AliMCEvent * mcEvt) { + + // Fill MC gen histo and the process types statistics + + Int_t evtype = GetProcessType(mcEvt); + // When I fill gen histos, I also fill statistics of process types (used to detemine ratios afterwards). + fHistProcTypes->Fill(evtype); + fHistVzMCGen[evtype]->Fill(vz,ntrk); +} + +void AliCollisionNormalization::FillVzMCRec(Float_t vz, Int_t ntrk, AliMCEvent * mcEvt){ + + // Fill MC rec histo + Int_t evtype = GetProcessType(mcEvt); + fHistVzMCRec[evtype]->Fill(vz,ntrk); + +} +void AliCollisionNormalization::FillVzMCTrg(Float_t vz, Int_t ntrk, AliMCEvent * mcEvt) { + + // Fill MC trg histo + Int_t evtype = GetProcessType(mcEvt); + fHistVzMCTrg[evtype]->Fill(vz,ntrk); +} + + +Int_t AliCollisionNormalization::GetProcessType(AliMCEvent * mcEvt) { + + // Determine if the event was generated with pythia or phojet and return the process type + + // Check if mcEvt is fine + if (!mcEvt) { + AliFatal("NULL mc event"); + } + + // Determine if it was a pythia or phojet header, and return the correct process type + AliGenPythiaEventHeader * headPy = 0; + AliGenDPMjetEventHeader * headPho = 0; + AliGenEventHeader * htmp = mcEvt->GenEventHeader(); + if(!htmp) { + AliFatal("Cannot Get MC Header!!"); + } + if( TString(htmp->IsA()->GetName()) == "AliGenPythiaEventHeader") { + headPy = (AliGenPythiaEventHeader*) htmp; + } else if (TString(htmp->IsA()->GetName()) == "AliGenDPMjetEventHeader") { + headPho = (AliGenDPMjetEventHeader*) htmp; + } else { + AliError("Unknown header"); + } + + // Determine process type + if(headPy) { + if(headPy->ProcessType() == 92 || headPy->ProcessType() == 93) { + // single difractive + return kProcSD; + } else if (headPy->ProcessType() == 94) { + // double diffractive + return kProcDD; + } + else if(headPy->ProcessType() != 92 && headPy->ProcessType() != 93 && headPy->ProcessType() != 94) { + // non difractive + return kProcND; + } + } else if (headPho) { + if(headPho->ProcessType() == 5 || headPho->ProcessType() == 6 ) { + // single difractive + return kProcSD; + } else if (headPho->ProcessType() == 7) { + // double diffractive + return kProcDD; + } else if(headPho->ProcessType() != 5 && headPho->ProcessType() != 6 && headPho->ProcessType() != 7 ) { + // non difractive + return kProcND; + } + } + + + // no process type found? + AliError(Form("Unknown header: %s", htmp->IsA()->GetName())); + return kProcUnknown; +} + +Double_t AliCollisionNormalization::GetProcessWeight(Int_t proc) { + + // Return a weight to be used when combining the MC histos to + // compute efficiency, defined as the ratio XS measured / XS in + // generator + + Float_t ref_SD, ref_DD, ref_ND, error_SD, error_DD, error_ND; + GetRelativeFractions(fReferenceXS,ref_SD, ref_DD, ref_ND, error_SD, error_DD, error_ND); + + static Double_t total = fHistProcTypes->Integral(); + if (fHistProcTypes->GetBinContent(fHistProcTypes->FindBin(kProcUnknown)) > 0) { + AliError("There were unknown processes!!!"); + } + static Double_t SD = fHistProcTypes->GetBinContent(fHistProcTypes->FindBin(kProcSD))/total; + static Double_t DD = fHistProcTypes->GetBinContent(fHistProcTypes->FindBin(kProcDD))/total; + static Double_t ND = fHistProcTypes->GetBinContent(fHistProcTypes->FindBin(kProcND))/total; + + if (fVerbose > 2) { + AliInfo(Form("Total MC evts: %f",total)); + AliInfo(Form(" Frac SD %4.4f", SD)); + AliInfo(Form(" Frac DD %4.4f", DD)); + AliInfo(Form(" Frac ND %4.4f", ND)); + } + + switch(proc) { + case kProcSD: + return ref_SD/SD; + break; + case kProcDD: + return ref_DD/DD; + break; + case kProcND: + return ref_ND/ND; + break; + default: + AliError("Unknown process"); + } + + return 0; + +} + +void AliCollisionNormalization::GetRelativeFractions(Int_t origin, Float_t& ref_SD, Float_t& ref_DD, Float_t& ref_ND, Float_t& error_SD, Float_t& error_DD, Float_t& error_ND) +{ + // Returns fraction of XS (SD, ND, DD) and corresponding error + // Stolen from Jan Fiete's drawSystematics macro + + // origin: + // -1 = Pythia (test) + // 0 = UA5 + // 1 = Data 1.8 TeV + // 2 = Tel-Aviv + // 3 = Durham + // + + switch (origin) + { + case -10: // Pythia default at 7 GeV, 50% error + AliInfo("PYTHIA x-sections"); + ref_SD = 0.192637; error_SD = ref_SD * 0.5; + ref_DD = 0.129877; error_DD = ref_DD * 0.5; + ref_ND = 0.677486; error_ND = 0; + break; + + case -1: // Pythia default at 900 GeV, as test + AliInfo("PYTHIA x-sections"); + ref_SD = 0.223788; + ref_DD = 0.123315; + ref_ND = 0.652897; + break; + + case 0: // UA5 + AliInfo("UA5 x-sections a la first paper"); + ref_SD = 0.153; error_SD = 0.05; + ref_DD = 0.080; error_DD = 0.05; + ref_ND = 0.767; error_ND = 0; + break; + + case 10: // UA5 + AliInfo("UA5 x-sections hadron level definition for Pythia"); + // Fractions in Pythia with UA5 cuts selection for SD + // ND: 0.688662 + // SD: 0.188588 --> this should be 15.3 + // DD: 0.122750 + ref_SD = 0.224 * 0.153 / 0.189; error_SD = 0.023 * 0.224 / 0.189; + ref_DD = 0.095; error_DD = 0.06; + ref_ND = 1.0 - ref_SD - ref_DD; error_ND = 0; + break; + + case 11: // UA5 + AliInfo("UA5 x-sections hadron level definition for Phojet"); + // Fractions in Phojet with UA5 cuts selection for SD + // ND: 0.783573 + // SD: 0.151601 --> this should be 15.3 + // DD: 0.064827 + ref_SD = 0.191 * 0.153 / 0.152; error_SD = 0.023 * 0.191 / 0.152; + ref_DD = 0.095; error_DD = 0.06; + ref_ND = 1.0 - ref_SD - ref_DD; error_ND = 0; + break; + + case 20: // E710, 1.8 TeV + AliInfo("E710 x-sections hadron level definition for Pythia"); + // ND: 0.705709 + // SD: 0.166590 --> this should be 15.9 + // DD: 0.127701 + ref_SD = 0.217 * 0.159 / 0.167; error_SD = 0.024 * 0.217 / 0.167; + ref_DD = 0.075 * 1.43; error_DD = 0.02 * 1.43; + ref_ND = 1.0 - ref_SD - ref_DD; error_ND = 0; + break; + + case 21: // E710, 1.8 TeV + AliInfo("E710 x-sections hadron level definition for Phojet"); + // ND: 0.817462 + // SD: 0.125506 --> this should be 15.9 + // DD: 0.057032 + ref_SD = 0.161 * 0.159 / 0.126; error_SD = 0.024 * 0.161 / 0.126; + ref_DD = 0.075 * 1.43; error_DD = 0.02 * 1.43; + ref_ND = 1.0 - ref_SD - ref_DD; error_ND = 0; + break; + + case 1: // data 1.8 TeV + AliInfo("??? x-sections"); + ref_SD = 0.152; + ref_DD = 0.092; + ref_ND = 1 - ref_SD - ref_DD; + break; + + case 2: // tel-aviv model + AliInfo("Tel-aviv model x-sections"); + ref_SD = 0.171; + ref_DD = 0.094; + ref_ND = 1 - ref_SD - ref_DD; + break; + + case 3: // durham model + AliInfo("Durham model x-sections"); + ref_SD = 0.190; + ref_DD = 0.125; + ref_ND = 1 - ref_SD - ref_DD; + break; + + default: + AliFatal(Form("Unknown origin %d", origin)); + } +} + diff --git a/ANALYSIS/AliCollisionNormalization.h b/ANALYSIS/AliCollisionNormalization.h new file mode 100644 index 00000000000..a53b643d9de --- /dev/null +++ b/ANALYSIS/AliCollisionNormalization.h @@ -0,0 +1,105 @@ + +#ifndef ALICOLLISIONNORMALIZATION_H +#define ALICOLLISIONNORMALIZATION_H + +/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * + * See cxx source for full Copyright notice */ + +//------------------------------------------------------------------------- +// Implementation of Class AliCollisionNormalization +// +// This class is used to store the vertex ditributions in the data +// and in Monte Carlo, needed to compute the real number of +// collisions a given sample is corresponding to. +// The strategy matches what described in CERN-THESIS-2009-033 p 119. +// +// Author: Michele Floris, CERN +//------------------------------------------------------------------------- + +#include "TH2F.h" + +class TH1F; +class TH1I; +class AliMCEvent; + +class AliCollisionNormalization : public TObject + +{ + + +public: + enum { kNevBin0, kNevCollisions, kNevNbin }; + typedef enum { kProcSD, kProcDD, kProcND, kProcUnknown, kNProcs } ProcType_t; + AliCollisionNormalization(); + AliCollisionNormalization(Int_t nbinz, Float_t minz, Float_t maxz); + AliCollisionNormalization(const char * dataFile, const char * dataListName, + const char * mcFile, const char * mcListName, + const char * eventStatFile); + + ~AliCollisionNormalization(); + + void SetMC(Bool_t flag = kTRUE) { fIsMC = flag;} + + void BookAllHistos(); + TH1 * BookVzHisto(const char * name , const char * title, Bool_t vzOnly=kFALSE); + + void FillVzMCGen(Float_t vz, Int_t ntrk, AliMCEvent * mcEvt); + void FillVzMCRec(Float_t vz, Int_t ntrk, AliMCEvent * mcEvt); + void FillVzMCTrg(Float_t vz, Int_t ntrk, AliMCEvent * mcEvt); + void FillVzData (Float_t vz, Int_t ntrk) {fHistVzData ->Fill(vz,ntrk);} + + TH2F * GetVzMCGen (Int_t procType) ; + TH2F * GetVzMCRec (Int_t procType) ; + TH2F * GetVzMCTrg (Int_t procType) ; + TH2F * GetVzData () { return fHistVzData ; } + TH1F * GetStatBin0 () { return fHistStatBin0 ; } + TH1F * GetHistProcTypes () { return fHistProcTypes ; } + + + Int_t GetProcessType(AliMCEvent * mcEvt) ; + Double_t GetProcessWeight(Int_t proctype); + + void SetReferencsXS(Int_t ref) { fReferenceXS = ref;} + + Double_t ComputeNint(); + void SetZRange(Float_t zrange) { fZRange = zrange ;} + + void SetReferenceXS(Int_t ref) { fReferenceXS = ref ;} + void GetRelativeFractions(Int_t origin, Float_t& ref_SD, Float_t& ref_DD, Float_t& ref_ND, Float_t& error_SD, Float_t& error_DD, Float_t& error_ND); + + void SetVerbose(Int_t lev) { fVerbose = lev ;} + + Long64_t Merge(TCollection* list); + +protected: + + Int_t fNbinsVz; // number of z bins in the vz histo + Float_t fMinVz ; // lowest Z + Float_t fMaxVz ; // highest Z + + Float_t fZRange; // max |Z| vertex to be considered + + Bool_t fIsMC; // True if processing MC + + Int_t fReferenceXS; // index of reference cross section to be used to rescale process types in the calculation of the efficiency + + Int_t fVerbose; // Determines the ammount of printout + + TH2F * fHistVzMCGen[kNProcs] ; // Vz distribution of generated events vs rec multiplicity + TH2F * fHistVzMCRec[kNProcs] ; // Vz distribution of reconstructed events vs rec multiplicity + TH2F * fHistVzMCTrg[kNProcs] ; // Vz distribution of triggered events vs rec multiplicity + TH2F * fHistVzData ; // Vz distribution of triggered events vs rec multiplicity + TH1F * fHistProcTypes ; // Number of evts for different Process types + + TH1F * fHistStatBin0 ; // event stat histogram, created by physiscs selection; used in ComputeNint; + + static const char * fProcLabel[] ; // labels of the different process types + + ClassDef(AliCollisionNormalization, 1); + +private: + AliCollisionNormalization(const AliCollisionNormalization&); + AliCollisionNormalization& operator=(const AliCollisionNormalization&); +}; + +#endif diff --git a/ANALYSIS/AliCollisionNormalizationTask.cxx b/ANALYSIS/AliCollisionNormalizationTask.cxx new file mode 100644 index 00000000000..8344f9d70f2 --- /dev/null +++ b/ANALYSIS/AliCollisionNormalizationTask.cxx @@ -0,0 +1,220 @@ +/* $Id$ */ + + +// Simple task to test the collision normalization. Can be called for +// both MC and data and shows how to fill the collision normalization +// class +// +// Author: Michele Floris +// CERN + + +#include "AliCollisionNormalizationTask.h" + +#include +#include +#include + +#include +#include +#include + +#include "AliCollisionNormalization.h" +#include "AliAnalysisManager.h" +#include "AliInputEventHandler.h" + +//#include "AliBackgroundSelection.h" +#include "AliMultiplicity.h" +#include "AliMCEvent.h" + +ClassImp(AliCollisionNormalizationTask) + +AliCollisionNormalizationTask::AliCollisionNormalizationTask() : + AliAnalysisTaskSE("AliCollisionNormalizationTask"), + fOutput(0), + fIsMC(0), + fCollisionNormalization(0) +{ + // + // Default event handler + // + + // Define input and output slots here + DefineOutput(1, TList::Class()); + +} + +AliCollisionNormalizationTask::AliCollisionNormalizationTask(const char* name) : + AliAnalysisTaskSE(name), + fOutput(0), + fIsMC(0), + fCollisionNormalization(new AliCollisionNormalization()) +{ + // + // Constructor. Initialization of pointers + // + + // Define input and output slots here + DefineOutput(1, TList::Class()); + + // AliLog::SetClassDebugLevel("AliCollisionNormalizationTask", AliLog::kWarning); +} + +AliCollisionNormalizationTask::~AliCollisionNormalizationTask() +{ + // + // Destructor + // + + // histograms are in the output list and deleted when the output + // list is deleted by the TSelector dtor + + if (fOutput) { + delete fOutput; + fOutput = 0; + } +} + +void AliCollisionNormalizationTask::UserCreateOutputObjects() +{ + // create result objects and add to output list + + Printf("AliCollisionNormalizationTask::CreateOutputObjects"); + + fOutput = new TList; + fOutput->SetOwner(); + + if (!fCollisionNormalization) + fCollisionNormalization = new AliCollisionNormalization; + + fOutput->Add(fCollisionNormalization); +// fOutput->Add(fCollisionNormalization->GetVzCorrZeroBin ()); +// fOutput->Add(fCollisionNormalization->GetVzMCGen ()); +// fOutput->Add(fCollisionNormalization->GetVzMCRec ()); +// fOutput->Add(fCollisionNormalization->GetVzMCTrg ()); +// fOutput->Add(fCollisionNormalization->GetVzData ()); +// fOutput->Add(fCollisionNormalization->GetNEvents ()); +// fOutput->Add(fCollisionNormalization->GetStatBin0 ()); + +} + +void AliCollisionNormalizationTask::UserExec(Option_t*) +{ + // process the event + + + PostData(1, fOutput); + + // Get the ESD + AliESDEvent * aESD = dynamic_cast(fInputEvent); + if (strcmp(aESD->ClassName(),"AliESDEvent")) { + AliFatal("Not processing ESDs"); + } + // Get MC event, if needed + AliMCEvent* mcEvent = fIsMC ? MCEvent() : 0; + if (!mcEvent && fIsMC){ + AliFatal("Running on MC but no MC handler available"); + } + + // Physics selection. At least in the case of MC we cannot use + // yourTask->SelectCollisionCandidates();, because we also need to + // fill the "generated" histogram + // NB never call IsEventSelected more than once per event + // (statistics histogram would be altered) + + Bool_t isSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected(); + + // Get the Multiplicity cut + const AliMultiplicity* mult = aESD->GetMultiplicity(); + if (!mult){ + AliError("Can't get mult object"); + return; + } + + // Check if the event is "reconstructed" according to some quality + // cuts. Tipically, we mean "it has a good-eonugh vertex" + + Int_t ntracklet = mult->GetNumberOfTracklets(); + const AliESDVertex * vtxESD = aESD->GetPrimaryVertexSPD(); + if(vtxESD) { + // If there is a vertex from vertexer z with delta phi > 0.02 we + // don't consider it rec (we keep the event in bin0). If quality + // is good eneough we check the number of tracklets + // A similar selection is applied in IsEventInBinZero + // If the vertex is too far away, we don't consider it reconstructed either. + if ((vtxESD->IsFromVertexerZ() && vtxESD->GetDispersion() > 0.02) || TMath::Abs(vtxESD->GetZ()) > 15) { + vtxESD = 0; + ntracklet = 0; // Don't trust reconstructed tracklets if you don't trust clusters + } + } + + if (ntracklet > 0 && !vtxESD) { + AliError("No vertex but reconstructed tracklets?"); + } + + // assign vz. For MC we use generated vz + Float_t vz = 0; + if (!fIsMC) vz = vtxESD ? vtxESD->GetZ() : 0; // FIXME : is zv used anywhere in Gen? + else vz = mcEvent->GetPrimaryVertex()->GetZ(); + + if (fIsMC) { + // Monte Carlo: we fill 3 histos + fCollisionNormalization->FillVzMCGen(vz, ntracklet, mcEvent); + // If triggered == passing the physics selection + if (isSelected) fCollisionNormalization->FillVzMCTrg(vz, ntracklet, mcEvent); + // If reconstructer == good enough vertex + if (vtxESD) fCollisionNormalization->FillVzMCRec(vz, ntracklet, mcEvent); + } else { + if (isSelected) { + // Passing the trigger + fCollisionNormalization->FillVzData(vz,ntracklet); + } + } + +} + +void AliCollisionNormalizationTask::Terminate(Option_t *) +{ + // The Terminate() function is the last function to be called during + // a query. It always runs on the client, it can be used to present + // the results graphically or save the results to file. + + fOutput = dynamic_cast (GetOutputData(1)); + if (!fOutput) + Printf("ERROR: fOutput not available"); + + +} + +Bool_t AliCollisionNormalizationTask::IsEventInBinZero() { + + // Returns true if an event is to be assigned to the zero bin + // + // You should have your own version of this method in the class: in + // general, the definition of "reconstructed" event is subjective. + + Bool_t isZeroBin = kTRUE; + const AliESDEvent* esd= dynamic_cast(fInputEvent); + const AliMultiplicity* mult = esd->GetMultiplicity(); + if (!mult){ + Printf("AliAnalysisTaskBGvsTime::IsBinZero: Can't get mult object"); + return kFALSE; + } + Int_t ntracklet = mult->GetNumberOfTracklets(); + const AliESDVertex * vtxESD = esd->GetPrimaryVertexSPD(); + if(vtxESD) { + // If there is a vertex from vertexer z with delta phi > 0.02 we + // don't consider it rec (we keep the event in bin0). If quality + // is good eneough we check the number of tracklets + // if the vertex is more than 15 cm away, this is autamatically bin0 + if( TMath::Abs(vtxESD->GetZ()) <= 15 ) { + if (vtxESD->IsFromVertexerZ()) { + if (vtxESD->GetDispersion()<=0.02 ) { + if(ntracklet>0) isZeroBin = kFALSE; + } + } else if(ntracklet>0) isZeroBin = kFALSE; // if the event is not from Vz we chek the n of tracklets + } + } + return isZeroBin; + +} diff --git a/ANALYSIS/AliCollisionNormalizationTask.h b/ANALYSIS/AliCollisionNormalizationTask.h new file mode 100644 index 00000000000..c88b1ce3e77 --- /dev/null +++ b/ANALYSIS/AliCollisionNormalizationTask.h @@ -0,0 +1,47 @@ +/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * + * See cxx source for full Copyright notice */ + + +/* $Id$ */ + +#ifndef ALICOLLISIONNORMALIZATIONTASK_H +#define ALICOLLISIONNORMALIZATIONTASK_H + +#include "AliAnalysisTaskSE.h" + +class AliCollisionNormalization; + +class AliCollisionNormalizationTask : public AliAnalysisTaskSE { + public: + AliCollisionNormalizationTask(); + AliCollisionNormalizationTask(const char* name); + + virtual ~AliCollisionNormalizationTask(); + + virtual void UserCreateOutputObjects(); + virtual void UserExec(Option_t*); + virtual void Terminate(Option_t*); + + Bool_t IsEventInBinZero(); // returns true if the event has to be put in the bin0. + void SetMC(Bool_t flag = kTRUE) { fIsMC = flag; } + + // void SetOption(const char* opt) { fOption = opt; } + + // void SetCollisionNormalization(AliCollisionNormalization* physicsSelection) { fCollisionNormalization = physicsSelection; } + AliCollisionNormalization* GetCollisionNormalization() const { return fCollisionNormalization; } + + protected: + TList* fOutput; //! list send on output slot 1 + // TString fOption; // option string + Bool_t fIsMC; + + AliCollisionNormalization* fCollisionNormalization; // collision normalization class + + private: + AliCollisionNormalizationTask(const AliCollisionNormalizationTask&); + AliCollisionNormalizationTask& operator=(const AliCollisionNormalizationTask&); + + ClassDef(AliCollisionNormalizationTask, 1); +}; + +#endif diff --git a/ANALYSIS/AliPhysicsSelection.h b/ANALYSIS/AliPhysicsSelection.h index 02fe48582f0..6fbc6067819 100644 --- a/ANALYSIS/AliPhysicsSelection.h +++ b/ANALYSIS/AliPhysicsSelection.h @@ -31,6 +31,7 @@ class AliAnalysisTaskSE; class AliPhysicsSelection : public AliAnalysisCuts { +public: enum {kStatTriggerClass=1,kStatHWTrig,kStatV0ABG,kStatV0CBG,kStatMB1,kStatMB1Prime,kStatFMD,kStatFO1,kStatFO2,kStatV0A,kStatV0C,kStatSSD1,kStatFO1AndV0,kStatV0,kStatAny2Hits,kStatOffline,kStatBG,kStatAccepted}; @@ -43,7 +44,6 @@ class AliPhysicsSelection : public AliAnalysisCuts enum {kStatIdxAll=0,kStatIdxBin0=1}; -public: typedef Bool_t (*Bin0Callback_t)(const AliESDEvent *); AliPhysicsSelection(); diff --git a/ANALYSIS/libANALYSISalice.pkg b/ANALYSIS/libANALYSISalice.pkg index 2813d02efdb..e7af22f63e0 100644 --- a/ANALYSIS/libANALYSISalice.pkg +++ b/ANALYSIS/libANALYSISalice.pkg @@ -9,7 +9,8 @@ SRCS = AliAnalysisTaskSE.cxx AliAnalysisTaskME.cxx \ AliTriggerAnalysis.cxx \ AliPhysicsSelection.cxx AliBackgroundSelection.cxx \ AliPhysicsSelectionTask.cxx \ - AliAnalysisFilter.cxx AliAnalysisCuts.cxx + AliAnalysisFilter.cxx AliAnalysisCuts.cxx \ + AliCollisionNormalization.cxx AliCollisionNormalizationTask.cxx ifeq (yes,$(CHECKALIEN)) diff --git a/ANALYSIS/macros/runProofNormalization.C b/ANALYSIS/macros/runProofNormalization.C new file mode 100644 index 00000000000..c754bc616c8 --- /dev/null +++ b/ANALYSIS/macros/runProofNormalization.C @@ -0,0 +1,66 @@ + +void runProofNormalization(const char * dataset = "LHC09b12_7TeV_0.5T", TString dataSetPath ="/PWG0/jgrosseo/",const char * filename = "LHC09b12_7TeV_0.5T_norm.root", Bool_t isMC = 1,Int_t nev =123456789) { + + gEnv->SetValue("XSec.GSI.DelegProxy","2"); + TProof::Open("alicecaf"); + + // gSystem->AddIncludePath("-I${ALICE_ROOT}/include/ -I${ALICE_ROOT}/PWG0/ -I${ALICE_ROOT}/PWG0/dNdEta/"); + gSystem->AddIncludePath("-I${ALICE_ROOT}/include/"); + gProof->UploadPackage("$ALICE_ROOT/STEERBase"); + gProof->EnablePackage("$ALICE_ROOT/STEERBase"); + gProof->UploadPackage("$ALICE_ROOT/ESD"); + gProof->EnablePackage("$ALICE_ROOT/ESD"); + gProof->UploadPackage("$ALICE_ROOT/AOD"); + gProof->EnablePackage("$ALICE_ROOT/AOD"); + gProof->UploadPackage("$ALICE_ROOT/ANALYSIS"); + gProof->EnablePackage("$ALICE_ROOT/ANALYSIS"); + gProof->UploadPackage("$ALICE_ROOT/ANALYSISalice"); + gProof->EnablePackage("$ALICE_ROOT/ANALYSISalice"); + gProof->UploadPackage("$ALICE_ROOT/CORRFW"); + gProof->EnablePackage("$ALICE_ROOT/CORRFW"); + + + // Make the analysis manager + AliAnalysisManager *mgr = new AliAnalysisManager("TestManager"); + // mgr->SetDebugLevel(3); + // Add ESD handler + AliESDInputHandler* esdH = new AliESDInputHandler; + + mgr->SetInputEventHandler(esdH); + + if(isMC) { + AliMCEventHandler *mc = new AliMCEventHandler(); + mc->SetReadTR(kFALSE); + mgr->SetMCtruthEventHandler(mc); + } + // assign simple task +// gProof->Load("AliCollisionNormalization.cxx++g"); +// gProof->Load("AliCollisionNormalizationTask.cxx++g"); + //____________________________________________// + // assign simple task + AliCollisionNormalizationTask * task = new AliCollisionNormalizationTask("TaskNormalization"); + // task->SetMC(); + task->SetMC(isMC); + mgr->AddTask(task); + + + gROOT->LoadMacro("$(ALICE_ROOT)/ANALYSIS/macros/AddTaskPhysicsSelection.C"); + AliPhysicsSelectionTask* physSelTask = AddTaskPhysicsSelection(isMC,1,!isMC); // Use Physics Selection. Enable computation of BG if is not MC + // task->SelectCollisionCandidates(); /// This should be disabled, at least for MC: we need all the events + physSelTask->GetPhysicsSelection()->SetBin0Callback("TaskNormalization"); + + AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer(); + mgr->ConnectInput(task,0,cinput1); + + + + // Attach output + cOutput = mgr->CreateContainer("Norm", TList::Class(), AliAnalysisManager::kOutputContainer,filename); + mgr->ConnectOutput(task, 1, cOutput); + + if (!mgr->InitAnalysis()) return; + + mgr->PrintStatus(); + mgr->StartAnalysis("proof",dataSetPath+dataset+"#esdTree",nev); + +} diff --git a/ANALYSIS/macros/runTaskNormalization.C b/ANALYSIS/macros/runTaskNormalization.C new file mode 100644 index 00000000000..179e3b28dd3 --- /dev/null +++ b/ANALYSIS/macros/runTaskNormalization.C @@ -0,0 +1,80 @@ + +void runTaskNormalization(const char * incollection,const char * filename = "LHC09b12_7TeV_0.5T_norm.root", Bool_t isMC = 1,Int_t nev =123456789) { + + + // Load libraries + gSystem->Load("libANALYSIS") ; + gSystem->Load("libANALYSISalice") ; + gSystem->Load("libCORRFW") ; + gSystem->Load("libITSbase") ; + gSystem->Load("libPWG0base") ; + + + // chain + TChain* analysisChain = 0; + analysisChain = new TChain("esdTree"); + if (TString(incollection).Contains(".root")){ + analysisChain->Add(incollection); + } + else if (TString(incollection).Contains("xml")){ + TGrid::Connect("alien://"); + TAlienCollection * coll = TAlienCollection::Open (incollection); + while(coll->Next()){ + analysisChain->Add(TString("alien://")+coll->GetLFN()); + } + } else { + ifstream file_collect(incollection); + TString line; + while (line.ReadLine(file_collect) ) { + analysisChain->Add(line.Data()); + } + } + analysisChain->GetListOfFiles()->Print(); + + + // Make the analysis manager + AliAnalysisManager *mgr = new AliAnalysisManager("TestManager"); + // mgr->SetDebugLevel(3); + // Add ESD handler + AliESDInputHandler* esdH = new AliESDInputHandler; + + mgr->SetInputEventHandler(esdH); + + if(isMC) { + AliMCEventHandler *mc = new AliMCEventHandler(); + mc->SetReadTR(kFALSE); + mgr->SetMCtruthEventHandler(mc); + } + // assign simple task +// gROOT->LoadMacro("AliCollisionNormalization.cxx++g"); +// gROOT->LoadMacro("AliCollisionNormalizationTask.cxx++g"); + //____________________________________________// + // Physics selection + gROOT->LoadMacro("$(ALICE_ROOT)/ANALYSIS/macros/AddTaskPhysicsSelection.C"); + AliPhysicsSelectionTask* physSelTask = AddTaskPhysicsSelection(isMC,1,!isMC); // Use Physics Selection. Enable computation of BG if is not MC + // task->SelectCollisionCandidates(); /// This should be disabled, at least for MC: we need all the events + physSelTask->GetPhysicsSelection()->SetBin0Callback("TaskNormalization"); + + // assign simple task + AliCollisionNormalizationTask * task = new AliCollisionNormalizationTask("TaskNormalization"); + // task->SetMC(); + task->SetMC(isMC); + mgr->AddTask(task); + + + + AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer(); + mgr->ConnectInput(task,0,cinput1); + + + + // Attach output + cOutput = mgr->CreateContainer("Norm", TList::Class(), AliAnalysisManager::kOutputContainer,filename); + mgr->ConnectOutput(task, 1, cOutput); + + if (!mgr->InitAnalysis()) return; + + mgr->PrintStatus(); + mgr->StartAnalysis("local",analysisChain,nev); + +} -- 2.43.0