#include "AliCFManager.h"
#include "AliCFContainer.h"
#include "AliLog.h"
+#include "AliInputEventHandler.h"
#include "AliAnalysisManager.h"
#include "AliAODHandler.h"
#include "AliAODEvent.h"
#include "AliRDHFCutsDstoKKpi.h"
#include "AliRDHFCutsLctopKpi.h"
#include "AliRDHFCutsD0toKpipipi.h"
+#include "AliRDHFCutsLctoV0.h"
#include "AliCFVertexingHF2Prong.h"
#include "AliCFVertexingHF3Prong.h"
#include "AliCFVertexingHFCascade.h"
+#include "AliCFVertexingHFLctoV0bachelor.h"
#include "AliCFVertexingHF.h"
#include "AliVertexingHFUtils.h"
#include "AliAnalysisDataSlot.h"
#include "AliAnalysisDataContainer.h"
+#include "AliPIDResponse.h"
//__________________________________________________________________________
AliCFTaskVertexingHF::AliCFTaskVertexingHF() :
- AliAnalysisTaskSE(),
- fCFManager(0x0),
- fHistEventsProcessed(0x0),
- fCorrelation(0x0),
- fCountMC(0),
- fCountAcc(0),
- fCountVertex(0),
- fCountRefit(0),
- fCountReco(0),
- fCountRecoAcc(0),
- fCountRecoITSClusters(0),
- fCountRecoPPR(0),
- fCountRecoPID(0),
- fEvents(0),
- fDecayChannel(0),
- fFillFromGenerated(kFALSE),
- fOriginDselection(0),
- fAcceptanceUnf(kTRUE),
- fCuts(0),
- fUseWeight(kFALSE),
- fWeight(1.),
- fUseFlatPtWeight(kFALSE),
- fUseZWeight(kFALSE),
- fUseNchWeight(kFALSE),
- fNvar(0),
- fPartName(""),
- fDauNames(""),
- fSign(2),
- fCentralitySelection(kTRUE),
- fFakeSelection(0),
- fRejectIfNoQuark(kTRUE),
- fUseMCVertex(kFALSE),
- fDsOption(1),
- fGenDsOption(3),
- fConfiguration(kCheetah), // by default, setting the fast configuration
- fFuncWeight(0x0),
- fHistoMeasNch(0x0),
- fHistoMCNch(0x0),
- fResonantDecay(0)
+ AliAnalysisTaskSE(),
+ fCFManager(0x0),
+ fHistEventsProcessed(0x0),
+ fCorrelation(0x0),
+ fCountMC(0),
+ fCountAcc(0),
+ fCountVertex(0),
+ fCountRefit(0),
+ fCountReco(0),
+ fCountRecoAcc(0),
+ fCountRecoITSClusters(0),
+ fCountRecoPPR(0),
+ fCountRecoPID(0),
+ fEvents(0),
+ fDecayChannel(0),
+ fFillFromGenerated(kFALSE),
+ fOriginDselection(0),
+ fAcceptanceUnf(kTRUE),
+ fCuts(0),
+ fUseWeight(kFALSE),
+ fWeight(1.),
+ fUseFlatPtWeight(kFALSE),
+ fUseZWeight(kFALSE),
+ fUseNchWeight(kFALSE),
+ fNvar(0),
+ fPartName(""),
+ fDauNames(""),
+ fSign(2),
+ fCentralitySelection(kTRUE),
+ fFakeSelection(0),
+ fRejectIfNoQuark(kTRUE),
+ fUseMCVertex(kFALSE),
+ fDsOption(1),
+ fGenDsOption(3),
+ fConfiguration(kCheetah), // by default, setting the fast configuration
+ fFuncWeight(0x0),
+ fHistoMeasNch(0x0),
+ fHistoMCNch(0x0),
+ fResonantDecay(0),
+ fLctoV0bachelorOption(1),
+ fGenLctoV0bachelorOption(0)
{
- //
- //Default ctor
- //
+ //
+ //Default ctor
+ //
}
//___________________________________________________________________________
AliCFTaskVertexingHF::AliCFTaskVertexingHF(const Char_t* name, AliRDHFCuts* cuts, TF1* func) :
- AliAnalysisTaskSE(name),
- fCFManager(0x0),
- fHistEventsProcessed(0x0),
- fCorrelation(0x0),
- fCountMC(0),
- fCountAcc(0),
- fCountVertex(0),
- fCountRefit(0),
- fCountReco(0),
- fCountRecoAcc(0),
- fCountRecoITSClusters(0),
- fCountRecoPPR(0),
- fCountRecoPID(0),
- fEvents(0),
- fDecayChannel(0),
- fFillFromGenerated(kFALSE),
- fOriginDselection(0),
- fAcceptanceUnf(kTRUE),
- fCuts(cuts),
- fUseWeight(kFALSE),
- fWeight(1.),
- fUseFlatPtWeight(kFALSE),
- fUseZWeight(kFALSE),
- fUseNchWeight(kFALSE),
- fNvar(0),
- fPartName(""),
- fDauNames(""),
- fSign(2),
- fCentralitySelection(kTRUE),
- fFakeSelection(0),
- fRejectIfNoQuark(kTRUE),
- fUseMCVertex(kFALSE),
- fDsOption(1),
- fGenDsOption(3),
- fConfiguration(kCheetah), // by default, setting the fast configuration
- fFuncWeight(func),
- fHistoMeasNch(0x0),
- fHistoMCNch(0x0),
- fResonantDecay(0)
+ AliAnalysisTaskSE(name),
+ fCFManager(0x0),
+ fHistEventsProcessed(0x0),
+ fCorrelation(0x0),
+ fCountMC(0),
+ fCountAcc(0),
+ fCountVertex(0),
+ fCountRefit(0),
+ fCountReco(0),
+ fCountRecoAcc(0),
+ fCountRecoITSClusters(0),
+ fCountRecoPPR(0),
+ fCountRecoPID(0),
+ fEvents(0),
+ fDecayChannel(0),
+ fFillFromGenerated(kFALSE),
+ fOriginDselection(0),
+ fAcceptanceUnf(kTRUE),
+ fCuts(cuts),
+ fUseWeight(kFALSE),
+ fWeight(1.),
+ fUseFlatPtWeight(kFALSE),
+ fUseZWeight(kFALSE),
+ fUseNchWeight(kFALSE),
+ fNvar(0),
+ fPartName(""),
+ fDauNames(""),
+ fSign(2),
+ fCentralitySelection(kTRUE),
+ fFakeSelection(0),
+ fRejectIfNoQuark(kTRUE),
+ fUseMCVertex(kFALSE),
+ fDsOption(1),
+ fGenDsOption(3),
+ fConfiguration(kCheetah), // by default, setting the fast configuration
+ fFuncWeight(func),
+ fHistoMeasNch(0x0),
+ fHistoMCNch(0x0),
+ fResonantDecay(0),
+ fLctoV0bachelorOption(1),
+ fGenLctoV0bachelorOption(0)
{
- //
- // Constructor. Initialization of Inputs and Outputs
- //
- /*
- DefineInput(0) and DefineOutput(0)
- are taken care of by AliAnalysisTaskSE constructor
- */
- DefineOutput(1,TH1I::Class());
- DefineOutput(2,AliCFContainer::Class());
- DefineOutput(3,THnSparseD::Class());
- DefineOutput(4,AliRDHFCuts::Class());
+ //
+ // Constructor. Initialization of Inputs and Outputs
+ //
+ /*
+ DefineInput(0) and DefineOutput(0)
+ are taken care of by AliAnalysisTaskSE constructor
+ */
+ DefineOutput(1,TH1I::Class());
+ DefineOutput(2,AliCFContainer::Class());
+ DefineOutput(3,THnSparseD::Class());
+ DefineOutput(4,AliRDHFCuts::Class());
- fCuts->PrintAll();
+ fCuts->PrintAll();
}
//___________________________________________________________________________
AliCFTaskVertexingHF& AliCFTaskVertexingHF::operator=(const AliCFTaskVertexingHF& c)
{
- //
- // Assignment operator
- //
- if (this!=&c) {
- AliAnalysisTaskSE::operator=(c) ;
- fCFManager = c.fCFManager;
- fHistEventsProcessed = c.fHistEventsProcessed;
- fCuts = c.fCuts;
- fFuncWeight = c.fFuncWeight;
- fHistoMeasNch = c.fHistoMeasNch;
- fHistoMCNch = c.fHistoMCNch;
- }
- return *this;
+ //
+ // Assignment operator
+ //
+ if (this!=&c) {
+ AliAnalysisTaskSE::operator=(c) ;
+ fCFManager = c.fCFManager;
+ fHistEventsProcessed = c.fHistEventsProcessed;
+ fCuts = c.fCuts;
+ fFuncWeight = c.fFuncWeight;
+ fHistoMeasNch = c.fHistoMeasNch;
+ fHistoMCNch = c.fHistoMCNch;
+ }
+ return *this;
}
//___________________________________________________________________________
AliCFTaskVertexingHF::AliCFTaskVertexingHF(const AliCFTaskVertexingHF& c) :
- AliAnalysisTaskSE(c),
- fCFManager(c.fCFManager),
- fHistEventsProcessed(c.fHistEventsProcessed),
- fCorrelation(c.fCorrelation),
- fCountMC(c.fCountMC),
- fCountAcc(c.fCountAcc),
- fCountVertex(c.fCountVertex),
- fCountRefit(c.fCountRefit),
- fCountReco(c.fCountReco),
- fCountRecoAcc(c.fCountRecoAcc),
- fCountRecoITSClusters(c.fCountRecoITSClusters),
- fCountRecoPPR(c.fCountRecoPPR),
- fCountRecoPID(c.fCountRecoPID),
- fEvents(c.fEvents),
- fDecayChannel(c.fDecayChannel),
- fFillFromGenerated(c.fFillFromGenerated),
- fOriginDselection(c.fOriginDselection),
- fAcceptanceUnf(c.fAcceptanceUnf),
- fCuts(c.fCuts),
- fUseWeight(c.fUseWeight),
- fWeight(c.fWeight),
- fUseFlatPtWeight(c.fUseFlatPtWeight),
- fUseZWeight(c.fUseZWeight),
- fUseNchWeight(c.fUseNchWeight),
- fNvar(c.fNvar),
- fPartName(c.fPartName),
- fDauNames(c.fDauNames),
- fSign(c.fSign),
- fCentralitySelection(c.fCentralitySelection),
- fFakeSelection(c.fFakeSelection),
- fRejectIfNoQuark(c.fRejectIfNoQuark),
- fUseMCVertex(c.fUseMCVertex),
- fDsOption(c.fDsOption),
- fGenDsOption(c.fGenDsOption),
- fConfiguration(c.fConfiguration),
- fFuncWeight(c.fFuncWeight),
- fHistoMeasNch(c.fHistoMeasNch),
- fHistoMCNch(c.fHistoMCNch),
- fResonantDecay(c.fResonantDecay)
+ AliAnalysisTaskSE(c),
+ fCFManager(c.fCFManager),
+ fHistEventsProcessed(c.fHistEventsProcessed),
+ fCorrelation(c.fCorrelation),
+ fCountMC(c.fCountMC),
+ fCountAcc(c.fCountAcc),
+ fCountVertex(c.fCountVertex),
+ fCountRefit(c.fCountRefit),
+ fCountReco(c.fCountReco),
+ fCountRecoAcc(c.fCountRecoAcc),
+ fCountRecoITSClusters(c.fCountRecoITSClusters),
+ fCountRecoPPR(c.fCountRecoPPR),
+ fCountRecoPID(c.fCountRecoPID),
+ fEvents(c.fEvents),
+ fDecayChannel(c.fDecayChannel),
+ fFillFromGenerated(c.fFillFromGenerated),
+ fOriginDselection(c.fOriginDselection),
+ fAcceptanceUnf(c.fAcceptanceUnf),
+ fCuts(c.fCuts),
+ fUseWeight(c.fUseWeight),
+ fWeight(c.fWeight),
+ fUseFlatPtWeight(c.fUseFlatPtWeight),
+ fUseZWeight(c.fUseZWeight),
+ fUseNchWeight(c.fUseNchWeight),
+ fNvar(c.fNvar),
+ fPartName(c.fPartName),
+ fDauNames(c.fDauNames),
+ fSign(c.fSign),
+ fCentralitySelection(c.fCentralitySelection),
+ fFakeSelection(c.fFakeSelection),
+ fRejectIfNoQuark(c.fRejectIfNoQuark),
+ fUseMCVertex(c.fUseMCVertex),
+ fDsOption(c.fDsOption),
+ fGenDsOption(c.fGenDsOption),
+ fConfiguration(c.fConfiguration),
+ fFuncWeight(c.fFuncWeight),
+ fHistoMeasNch(c.fHistoMeasNch),
+ fHistoMCNch(c.fHistoMCNch),
+ fResonantDecay(c.fResonantDecay),
+ fLctoV0bachelorOption(c.fLctoV0bachelorOption),
+ fGenLctoV0bachelorOption(c.fGenLctoV0bachelorOption)
{
- //
- // Copy Constructor
- //
+ //
+ // Copy Constructor
+ //
}
//___________________________________________________________________________
AliCFTaskVertexingHF::~AliCFTaskVertexingHF()
{
- //
- //destructor
- //
- if (fCFManager) delete fCFManager ;
- if (fHistEventsProcessed) delete fHistEventsProcessed ;
- if (fCorrelation) delete fCorrelation ;
- if (fCuts) delete fCuts;
- if (fFuncWeight) delete fFuncWeight;
- if (fHistoMeasNch) delete fHistoMeasNch;
- if (fHistoMCNch) delete fHistoMCNch;
+ //
+ //destructor
+ //
+ if (fCFManager) delete fCFManager ;
+ if (fHistEventsProcessed) delete fHistEventsProcessed ;
+ if (fCorrelation) delete fCorrelation ;
+ if (fCuts) delete fCuts;
+ if (fFuncWeight) delete fFuncWeight;
+ if (fHistoMeasNch) delete fHistoMeasNch;
+ if (fHistoMCNch) delete fHistoMCNch;
}
//_________________________________________________________________________-
void AliCFTaskVertexingHF::Init()
{
- //
- // Initialization
- //
+ //
+ // Initialization
+ //
- if (fDebug>1) printf("AliCFTaskVertexingHF::Init()");
- if(fUseWeight && fUseZWeight) { AliFatal("Can not use at the same time pt and z-vtx weights, please choose"); return; }
- if(fUseWeight && fUseNchWeight) { AliFatal("Can not use at the same time pt and Nch weights, please choose"); return; }
- if(fUseNchWeight && !fHistoMCNch) { AliFatal("Need to pass the MC Nch distribution to use Nch weights"); return; }
- if(fUseNchWeight) CreateMeasuredNchHisto();
-
- AliRDHFCuts *copyfCuts = 0x0;
- if (!fCuts){
- AliFatal("No cuts defined - Exiting...");
- return;
- }
+ if (fDebug>1) printf("AliCFTaskVertexingHF::Init()");
+ if(fUseWeight && fUseZWeight) { AliFatal("Can not use at the same time pt and z-vtx weights, please choose"); return; }
+ if(fUseWeight && fUseNchWeight) { AliFatal("Can not use at the same time pt and Nch weights, please choose"); return; }
+ if(fUseNchWeight && !fHistoMCNch) { AliFatal("Need to pass the MC Nch distribution to use Nch weights"); return; }
+ if(fUseNchWeight) CreateMeasuredNchHisto();
- switch (fDecayChannel){
- case 2:{
- copyfCuts = new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fCuts)));
- switch (fConfiguration) {
- case kSnail: // slow configuration: all variables in
- fNvar = 16;
- break;
- case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
- fNvar = 8;
- break;
- }
- fPartName="D0";
- fDauNames="K+pi";
- break;
- }
- case 21:{
- copyfCuts = new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fCuts)));
- switch (fConfiguration) {
- case kSnail: // slow configuration: all variables in
- fNvar = 16;
- break;
- case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
- fNvar = 8;
- break;
- }
- fPartName="Dstar";
- fDauNames="K+pi+pi";
- break;
- }
- case 31:{
- copyfCuts = new AliRDHFCutsDplustoKpipi(*(static_cast<AliRDHFCutsDplustoKpipi*>(fCuts)));
- switch (fConfiguration) {
- case kSnail: // slow configuration: all variables in
- fNvar = 14;
- break;
- case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
- fNvar = 8;
- break;
- }
- fPartName="Dplus";
- fDauNames="K+pi+pi";
- break;
- }
- case 32:{
- copyfCuts = new AliRDHFCutsLctopKpi(*(static_cast<AliRDHFCutsLctopKpi*>(fCuts)));
- switch (fConfiguration) {
- case kSnail: // slow configuration: all variables in
- fNvar = 18;
- break;
- case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
- fNvar = 8;
- break;
- }
- fPartName="Lambdac";
- fDauNames="p+K+pi";
- break;
- }
- case 33:{
- copyfCuts = new AliRDHFCutsDstoKKpi(*(static_cast<AliRDHFCutsDstoKKpi*>(fCuts)));
- switch (fConfiguration) {
- case kSnail: // slow configuration: all variables in
- fNvar = 14;
- break;
- case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
- fNvar = 8;
- break;
- }
- fPartName="Ds";
- fDauNames="K+K+pi";
- break;
- }
- case 4:{
- copyfCuts = new AliRDHFCutsD0toKpipipi(*(static_cast<AliRDHFCutsD0toKpipipi*>(fCuts)));
- switch (fConfiguration) {
- case kSnail: // slow configuration: all variables in
- fNvar = 16;
- break;
- case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
- fNvar = 8;
- break;
- }
- fPartName="D0";
- fDauNames="K+pi+pi+pi";
- break;
- }
- default:
- AliFatal("The decay channel MUST be defined according to AliCFVertexing::DecayChannel - Exiting...");
- break;
- }
+ AliRDHFCuts *copyfCuts = 0x0;
+ if (!fCuts){
+ AliFatal("No cuts defined - Exiting...");
+ return;
+ }
+
+ switch (fDecayChannel){
+ case 2:{
+ copyfCuts = new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fCuts)));
+ switch (fConfiguration) {
+ case kSnail: // slow configuration: all variables in
+ fNvar = 16;
+ break;
+ case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
+ fNvar = 8;
+ break;
+ }
+ fPartName="D0";
+ fDauNames="K+pi";
+ break;
+ }
+ case 21:{
+ copyfCuts = new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fCuts)));
+ switch (fConfiguration) {
+ case kSnail: // slow configuration: all variables in
+ fNvar = 16;
+ break;
+ case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
+ fNvar = 8;
+ break;
+ }
+ fPartName="Dstar";
+ fDauNames="K+pi+pi";
+ break;
+ }
+ case 22:{
+ copyfCuts = new AliRDHFCutsLctoV0(*(static_cast<AliRDHFCutsLctoV0*>(fCuts)));
+ switch (fConfiguration) {
+ case kSnail: // slow configuration: all variables in
+ fNvar = 16;
+ break;
+ case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
+ fNvar = 8;
+ break;
+ }
+ fPartName="Lambdac";
+ fDauNames="V0+bachelor";
+ break;
+ }
+ case 31:{
+ copyfCuts = new AliRDHFCutsDplustoKpipi(*(static_cast<AliRDHFCutsDplustoKpipi*>(fCuts)));
+ switch (fConfiguration) {
+ case kSnail: // slow configuration: all variables in
+ fNvar = 14;
+ break;
+ case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
+ fNvar = 8;
+ break;
+ }
+ fPartName="Dplus";
+ fDauNames="K+pi+pi";
+ break;
+ }
+ case 32:{
+ copyfCuts = new AliRDHFCutsLctopKpi(*(static_cast<AliRDHFCutsLctopKpi*>(fCuts)));
+ switch (fConfiguration) {
+ case kSnail: // slow configuration: all variables in
+ fNvar = 18;
+ break;
+ case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
+ fNvar = 8;
+ break;
+ }
+ fPartName="Lambdac";
+ fDauNames="p+K+pi";
+ break;
+ }
+ case 33:{
+ copyfCuts = new AliRDHFCutsDstoKKpi(*(static_cast<AliRDHFCutsDstoKKpi*>(fCuts)));
+ switch (fConfiguration) {
+ case kSnail: // slow configuration: all variables in
+ fNvar = 14;
+ break;
+ case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
+ fNvar = 8;
+ break;
+ }
+ fPartName="Ds";
+ fDauNames="K+K+pi";
+ break;
+ }
+ case 4:{
+ copyfCuts = new AliRDHFCutsD0toKpipipi(*(static_cast<AliRDHFCutsD0toKpipipi*>(fCuts)));
+ switch (fConfiguration) {
+ case kSnail: // slow configuration: all variables in
+ fNvar = 16;
+ break;
+ case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
+ fNvar = 8;
+ break;
+ }
+ fPartName="D0";
+ fDauNames="K+pi+pi+pi";
+ break;
+ }
+ default:
+ AliFatal("The decay channel MUST be defined according to AliCFVertexing::DecayChannel - Exiting...");
+ break;
+ }
- const char* nameoutput=GetOutputSlot(4)->GetContainer()->GetName();
- if (copyfCuts){
- copyfCuts->SetName(nameoutput);
+ const char* nameoutput=GetOutputSlot(4)->GetContainer()->GetName();
+ if (copyfCuts){
+ copyfCuts->SetName(nameoutput);
- //Post the data
- PostData(4, copyfCuts);
- }
- else{
- AliFatal("Failing initializing AliRDHFCuts object - Exiting...");
- }
+ //Post the data
+ PostData(4, copyfCuts);
+ }
+ else{
+ AliFatal("Failing initializing AliRDHFCuts object - Exiting...");
+ }
- return;
+ return;
}
//_________________________________________________
void AliCFTaskVertexingHF::UserExec(Option_t *)
{
- //
- // Main loop function
- //
+ //
+ // Main loop function
+ //
- PostData(1,fHistEventsProcessed) ;
- PostData(2,fCFManager->GetParticleContainer()) ;
- PostData(3,fCorrelation) ;
+ PostData(1,fHistEventsProcessed) ;
+ PostData(2,fCFManager->GetParticleContainer()) ;
+ PostData(3,fCorrelation) ;
- if (fFillFromGenerated){
- AliWarning("Flag to fill container with generated value ON ---> dca, d0pi, d0K, d0xd0, cosPointingAngle will be set as dummy!");
- }
+ if (fFillFromGenerated){
+ AliWarning("Flag to fill container with generated value ON ---> dca, d0pi, d0K, d0xd0, cosPointingAngle will be set as dummy!");
+ }
- if (!fInputEvent) {
- Error("UserExec","NO EVENT FOUND!");
- return;
- }
+ if (!fInputEvent) {
+ Error("UserExec","NO EVENT FOUND!");
+ return;
+ }
- AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
+ AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
- TClonesArray *arrayBranch=0;
+ TClonesArray *arrayBranch=0;
- if(!aodEvent && AODEvent() && IsStandardAOD()) {
- // In case there is an AOD handler writing a standard AOD, use the AOD
- // event in memory rather than the input (ESD) event.
- aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());
- // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
- // have to taken from the AOD event hold by the AliAODExtension
- AliAODHandler* aodHandler = (AliAODHandler*)
- ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
- if(aodHandler->GetExtensions()) {
- AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
- AliAODEvent *aodFromExt = ext->GetAOD();
+ if(!aodEvent && AODEvent() && IsStandardAOD()) {
+ // In case there is an AOD handler writing a standard AOD, use the AOD
+ // event in memory rather than the input (ESD) event.
+ aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());
+ // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
+ // have to taken from the AOD event hold by the AliAODExtension
+ AliAODHandler* aodHandler = (AliAODHandler*)
+ ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
+ if(aodHandler->GetExtensions()) {
+ AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
+ AliAODEvent *aodFromExt = ext->GetAOD();
- switch (fDecayChannel){
- case 2:{
- arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
- break;
- }
- case 21:{
- arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
- break;
- }
- case 31:
- case 32:
- case 33:{
- arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
- break;
- }
- case 4:{
- arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm4Prong");
- break;
- }
- default:
- break;
- }
- }
- }
- else {
- switch (fDecayChannel){
- case 2:{
- arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("D0toKpi");
- break;
- }
- case 21:{
- arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Dstar");
- break;
- }
- case 31:
- case 32:
- case 33:{
- arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Charm3Prong");
- break;
- }
- case 4:{
- arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Charm4Prong");
- break;
- }
- default:
- break;
- }
- }
+ switch (fDecayChannel){
+ case 2:{
+ arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
+ break;
+ }
+ case 21:{
+ arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
+ break;
+ }
+ case 22:{
+ arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("CascadesHF");
+ break;
+ }
+ case 31:
+ case 32:
+ case 33:{
+ arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
+ break;
+ }
+ case 4:{
+ arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm4Prong");
+ break;
+ }
+ default:
+ break;
+ }
+ }
+ }
+ else {
+ switch (fDecayChannel){
+ case 2:{
+ arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("D0toKpi");
+ break;
+ }
+ case 21:{
+ arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Dstar");
+ break;
+ }
+ case 22:{
+ arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("CascadesHF");
+ break;
+ }
+ case 31:
+ case 32:
+ case 33:{
+ arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Charm3Prong");
+ break;
+ }
+ case 4:{
+ arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Charm4Prong");
+ break;
+ }
+ default:
+ break;
+ }
+ }
- AliAODVertex *aodVtx = (AliAODVertex*)aodEvent->GetPrimaryVertex();
- if (!aodVtx) return;
+ AliAODVertex *aodVtx = (AliAODVertex*)aodEvent->GetPrimaryVertex();
+ if (!aodVtx) return;
- if (!arrayBranch) {
- AliError("Could not find array of HF vertices");
- return;
- }
+ if (!arrayBranch) {
+ AliError("Could not find array of HF vertices");
+ return;
+ }
- fEvents++;
+ fEvents++;
- fCFManager->SetRecEventInfo(aodEvent);
- fCFManager->SetMCEventInfo(aodEvent);
+ fCFManager->SetRecEventInfo(aodEvent);
+ fCFManager->SetMCEventInfo(aodEvent);
- //******** DEFINE number of variables of the container***** for now set at 13, in the future in the config macro.
+ //******** DEFINE number of variables of the container***** for now set at 13, in the future in the config macro.
- //loop on the MC event
+ //loop on the MC event
- TClonesArray* mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));
- if (!mcArray) {
- AliError("Could not find Monte-Carlo in AOD");
- return;
- }
- Int_t icountMC = 0;
- Int_t icountAcc = 0;
- Int_t icountReco = 0;
- Int_t icountVertex = 0;
- Int_t icountRefit = 0;
- Int_t icountRecoAcc = 0;
- Int_t icountRecoITSClusters = 0;
- Int_t icountRecoPPR = 0;
- Int_t icountRecoPID = 0;
- Int_t cquarks = 0;
+ TClonesArray* mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));
+ if (!mcArray) {
+ AliError("Could not find Monte-Carlo in AOD");
+ return;
+ }
+ Int_t icountMC = 0;
+ Int_t icountAcc = 0;
+ Int_t icountReco = 0;
+ Int_t icountVertex = 0;
+ Int_t icountRefit = 0;
+ Int_t icountRecoAcc = 0;
+ Int_t icountRecoITSClusters = 0;
+ Int_t icountRecoPPR = 0;
+ Int_t icountRecoPID = 0;
+ Int_t cquarks = 0;
- AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
- if (!mcHeader) {
- AliError("Could not find MC Header in AOD");
- return;
- }
+ AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
+ if (!mcHeader) {
+ AliError("Could not find MC Header in AOD");
+ return;
+ }
- Double_t* containerInput = new Double_t[fNvar];
- Double_t* containerInputMC = new Double_t[fNvar];
+ Double_t* containerInput = new Double_t[fNvar];
+ Double_t* containerInputMC = new Double_t[fNvar];
- AliCFVertexingHF* cfVtxHF=0x0;
- switch (fDecayChannel){
- case 2:{
- cfVtxHF = new AliCFVertexingHF2Prong(mcArray, fOriginDselection);
- break;
- }
- case 21:{
- cfVtxHF = new AliCFVertexingHFCascade(mcArray, fOriginDselection);
- break;
- }
- case 31:
-// case 32:
- case 33:{
- cfVtxHF = new AliCFVertexingHF3Prong(mcArray, fOriginDselection, fDecayChannel);
- if(fDecayChannel==33){
- ((AliCFVertexingHF3Prong*)cfVtxHF)->SetGeneratedDsOption(fGenDsOption);
- }
- break;
- }
- case 32:{
- cfVtxHF = new AliCFVertexingHF3Prong(mcArray, fOriginDselection, fDecayChannel,fResonantDecay);
- }
- case 4:{
- //cfVtxHF = new AliCFVertexingHF4Prong(mcArray, originDselection); // not there yet
- break;
- }
- default:
- break;
- }
- if (!cfVtxHF){
- AliError("No AliCFVertexingHF initialized");
- delete[] containerInput;
- delete[] containerInputMC;
- return;
- }
+ AliCFVertexingHF* cfVtxHF=0x0;
+ switch (fDecayChannel){
+ case 2:{
+ cfVtxHF = new AliCFVertexingHF2Prong(mcArray, fOriginDselection);
+ break;
+ }
+ case 21:{
+ cfVtxHF = new AliCFVertexingHFCascade(mcArray, fOriginDselection);
+ break;
+ }
+ case 22:{
+ cfVtxHF = new AliCFVertexingHFLctoV0bachelor(mcArray, fOriginDselection,fGenLctoV0bachelorOption); // Lc -> K0S+proton
+ break;
+ }
+ case 31:
+ // case 32:
+ case 33:{
+ cfVtxHF = new AliCFVertexingHF3Prong(mcArray, fOriginDselection, fDecayChannel);
+ if(fDecayChannel==33){
+ ((AliCFVertexingHF3Prong*)cfVtxHF)->SetGeneratedDsOption(fGenDsOption);
+ }
+ break;
+ }
+ case 32:{
+ cfVtxHF = new AliCFVertexingHF3Prong(mcArray, fOriginDselection, fDecayChannel,fResonantDecay);
+ }
+ case 4:{
+ //cfVtxHF = new AliCFVertexingHF4Prong(mcArray, originDselection); // not there yet
+ break;
+ }
+ default:
+ break;
+ }
+ if (!cfVtxHF){
+ AliError("No AliCFVertexingHF initialized");
+ delete[] containerInput;
+ delete[] containerInputMC;
+ return;
+ }
- Double_t zPrimVertex = aodVtx ->GetZ();
- Double_t zMCVertex = mcHeader->GetVtxZ();
- Int_t runnumber = aodEvent->GetRunNumber();
- fWeight=1.;
- if(fUseZWeight) fWeight *= GetZWeight(zMCVertex,runnumber);
- if(fUseNchWeight){
- Int_t nChargedMCPhysicalPrimary=AliVertexingHFUtils::GetGeneratedPhysicalPrimariesInEtaRange(mcArray,-1.0,1.0);
- fWeight *= GetNchWeight(nChargedMCPhysicalPrimary);
- AliDebug(2,Form("Using Nch weights, Mult=%d Weight=%f\n",nChargedMCPhysicalPrimary,fWeight));
- }
+ Double_t zPrimVertex = aodVtx ->GetZ();
+ Double_t zMCVertex = mcHeader->GetVtxZ();
+ Int_t runnumber = aodEvent->GetRunNumber();
+ fWeight=1.;
+ if(fUseZWeight) fWeight *= GetZWeight(zMCVertex,runnumber);
+ if(fUseNchWeight){
+ Int_t nChargedMCPhysicalPrimary=AliVertexingHFUtils::GetGeneratedPhysicalPrimariesInEtaRange(mcArray,-1.0,1.0);
+ fWeight *= GetNchWeight(nChargedMCPhysicalPrimary);
+ AliDebug(2,Form("Using Nch weights, Mult=%d Weight=%f\n",nChargedMCPhysicalPrimary,fWeight));
+ }
- if (TMath::Abs(zMCVertex) > fCuts->GetMaxVtxZ()){
- AliDebug(3,Form("z coordinate of MC vertex = %f, it was required to be within [-%f, +%f], skipping event", zMCVertex, fCuts->GetMaxVtxZ(), fCuts->GetMaxVtxZ()));
- delete[] containerInput;
- delete[] containerInputMC;
- return;
- }
+ if (TMath::Abs(zMCVertex) > fCuts->GetMaxVtxZ()){
+ AliDebug(3,Form("z coordinate of MC vertex = %f, it was required to be within [-%f, +%f], skipping event", zMCVertex, fCuts->GetMaxVtxZ(), fCuts->GetMaxVtxZ()));
+ delete[] containerInput;
+ delete[] containerInputMC;
+ return;
+ }
- AliESDtrackCuts** trackCuts = new AliESDtrackCuts*[cfVtxHF->GetNProngs()];
- if (fDecayChannel == 21){
- // for the D*, setting the third element of the array of the track cuts to those for the soft pion
- for (Int_t iProng = 0; iProng<cfVtxHF->GetNProngs()-1; iProng++){
- trackCuts[iProng]=fCuts->GetTrackCuts();
- }
- trackCuts[2] = fCuts->GetTrackCutsSoftPi();
- }
- else {
- for (Int_t iProng = 0; iProng<cfVtxHF->GetNProngs(); iProng++){
- trackCuts[iProng]=fCuts->GetTrackCuts();
- }
- }
+ AliESDtrackCuts** trackCuts = new AliESDtrackCuts*[cfVtxHF->GetNProngs()];
+ if (fDecayChannel == 21){
+ // for the D*, setting the third element of the array of the track cuts to those for the soft pion
+ for (Int_t iProng = 0; iProng<cfVtxHF->GetNProngs()-1; iProng++){
+ trackCuts[iProng]=fCuts->GetTrackCuts();
+ }
+ trackCuts[2] = fCuts->GetTrackCutsSoftPi();
+ }
+ else if (fDecayChannel == 22) {
+ // for the Lc->V0+bachelor, setting the second and third elements of the array of the track cuts to those for the V0 daughters
+ trackCuts[0]=fCuts->GetTrackCuts();
+ trackCuts[1]=fCuts->GetTrackCutsV0daughters();
+ trackCuts[2]=fCuts->GetTrackCutsV0daughters();
+ }
+ else {
+ for (Int_t iProng = 0; iProng<cfVtxHF->GetNProngs(); iProng++){
+ trackCuts[iProng]=fCuts->GetTrackCuts();
+ }
+ }
- //General settings: vertex, feed down and fill reco container with generated values.
- cfVtxHF->SetRecoPrimVertex(zPrimVertex);
- cfVtxHF->SetMCPrimaryVertex(zMCVertex);
- cfVtxHF->SetFillFromGenerated(fFillFromGenerated);
- cfVtxHF->SetNVar(fNvar);
- cfVtxHF->SetFakeSelection(fFakeSelection);
- cfVtxHF->SetRejectCandidateIfNotFromQuark(fRejectIfNoQuark);
- cfVtxHF->SetConfiguration(fConfiguration);
-
- // switch-off the trigger class selection (doesn't work for MC)
- fCuts->SetTriggerClass("");
-
- // MC vertex, to be used, in case, for pp
- if (fUseMCVertex) fCuts->SetUseMCVertex();
-
- if (fCentralitySelection){ // keep only the requested centrality
- if(fCuts->IsEventSelectedInCentrality(aodEvent)!=0) {
- delete[] containerInput;
- delete[] containerInputMC;
- delete [] trackCuts;
- return;
- }
- } else { // keep all centralities
- fCuts->SetMinCentrality(0.);
- fCuts->SetMaxCentrality(100.);
- }
+ //General settings: vertex, feed down and fill reco container with generated values.
+ cfVtxHF->SetRecoPrimVertex(zPrimVertex);
+ cfVtxHF->SetMCPrimaryVertex(zMCVertex);
+ cfVtxHF->SetFillFromGenerated(fFillFromGenerated);
+ cfVtxHF->SetNVar(fNvar);
+ cfVtxHF->SetFakeSelection(fFakeSelection);
+ cfVtxHF->SetRejectCandidateIfNotFromQuark(fRejectIfNoQuark);
+ cfVtxHF->SetConfiguration(fConfiguration);
+
+ // switch-off the trigger class selection (doesn't work for MC)
+ fCuts->SetTriggerClass("");
+
+ // MC vertex, to be used, in case, for pp
+ if (fUseMCVertex) fCuts->SetUseMCVertex();
+
+ if (fCentralitySelection){ // keep only the requested centrality
+ if(fCuts->IsEventSelectedInCentrality(aodEvent)!=0) {
+ delete[] containerInput;
+ delete[] containerInputMC;
+ delete [] trackCuts;
+ return;
+ }
+ } else { // keep all centralities
+ fCuts->SetMinCentrality(0.);
+ fCuts->SetMaxCentrality(100.);
+ }
- Float_t centValue = fCuts->GetCentrality(aodEvent);
- cfVtxHF->SetCentralityValue(centValue);
+ Float_t centValue = fCuts->GetCentrality(aodEvent);
+ cfVtxHF->SetCentralityValue(centValue);
- // number of tracklets - multiplicity estimator
- Double_t multiplicity = (Double_t)(aodEvent->GetTracklets()->GetNumberOfTracklets()); // casted to double because the CF is filled with doubles
- cfVtxHF->SetMultiplicity(multiplicity);
+ // number of tracklets - multiplicity estimator
+ Double_t multiplicity = (Double_t)(aodEvent->GetTracklets()->GetNumberOfTracklets()); // casted to double because the CF is filled with doubles
+ cfVtxHF->SetMultiplicity(multiplicity);
- for (Int_t iPart=0; iPart<mcArray->GetEntriesFast(); iPart++) {
- AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(mcArray->At(iPart));
- if (!mcPart){
- AliError("Failed casting particle from MC array!, Skipping particle");
- continue;
- }
- // check the MC-level cuts, must be the desidered particle
- if (!fCFManager->CheckParticleCuts(0, mcPart)) {
- continue; // 0 stands for MC level
- }
- cfVtxHF->SetMCCandidateParam(iPart);
+ for (Int_t iPart=0; iPart<mcArray->GetEntriesFast(); iPart++) {
+ AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(mcArray->At(iPart));
+ if (!mcPart){
+ AliError("Failed casting particle from MC array!, Skipping particle");
+ continue;
+ }
+ // check the MC-level cuts, must be the desidered particle
+ if (!fCFManager->CheckParticleCuts(0, mcPart)) {
+ AliDebug(2,"Check the MC-level cuts - not desidered particle");
+ continue; // 0 stands for MC level
+ }
+ cfVtxHF->SetMCCandidateParam(iPart);
- //counting c quarks
- cquarks += cfVtxHF->MCcquarkCounting(mcPart);
+ //counting c quarks
+ cquarks += cfVtxHF->MCcquarkCounting(mcPart);
- if (!(cfVtxHF->SetLabelArray())){
- AliDebug(2,Form("Impossible to set the label array (decaychannel = %d)",fDecayChannel));
- continue;
- }
+ if (!(cfVtxHF->SetLabelArray())){
+ AliDebug(2,Form("Impossible to set the label array (decaychannel = %d)",fDecayChannel));
+ continue;
+ }
- //check the candiate family at MC level
- if (!(cfVtxHF->CheckMCPartFamily(mcPart, mcArray))) {
- AliDebug(2,Form("Check on the family wrong!!! (decaychannel = %d)",fDecayChannel));
- continue;
- }
- else{
- AliDebug(2,Form("Check on the family OK!!! (decaychannel = %d)",fDecayChannel));
- }
+ //check the candiate family at MC level
+ if (!(cfVtxHF->CheckMCPartFamily(mcPart, mcArray))) {
+ AliDebug(2,Form("Check on the family wrong!!! (decaychannel = %d)",fDecayChannel));
+ continue;
+ }
+ else{
+ AliInfo(Form("Check on the family OK!!! (decaychannel = %d)",fDecayChannel));
+ }
- //Fill the MC container
- Bool_t mcContainerFilled = cfVtxHF -> FillMCContainer(containerInputMC);
- if (mcContainerFilled) {
- if (fUseWeight){
- if (fFuncWeight){ // using user-defined function
- AliDebug(2,"Using function");
- fWeight = fFuncWeight->Eval(containerInputMC[0]);
- }
- else{ // using FONLL
- AliDebug(2,"Using FONLL");
- fWeight = GetWeight(containerInputMC[0]);
- }
- AliDebug(2,Form("pt = %f, weight = %f",containerInputMC[0], fWeight));
- }
- if (!fCuts->IsInFiducialAcceptance(containerInputMC[0],containerInputMC[1])) continue;
- //MC Limited Acceptance
- if (TMath::Abs(containerInputMC[1]) < 0.5) {
- fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepGeneratedLimAcc, fWeight);
- AliDebug(3,"MC Lim Acc container filled\n");
- }
+ //Fill the MC container
+ Bool_t mcContainerFilled = cfVtxHF -> FillMCContainer(containerInputMC);
+ AliDebug(2,Form("mcContainerFilled = %d)",mcContainerFilled));
+ if (mcContainerFilled) {
+ if (fUseWeight){
+ if (fFuncWeight){ // using user-defined function
+ AliDebug(2,"Using function");
+ fWeight = fFuncWeight->Eval(containerInputMC[0]);
+ }
+ else{ // using FONLL
+ AliDebug(2,"Using FONLL");
+ fWeight = GetWeight(containerInputMC[0]);
+ }
+ AliDebug(2,Form("pt = %f, weight = %f",containerInputMC[0], fWeight));
+ }
+ if (!fCuts->IsInFiducialAcceptance(containerInputMC[0],containerInputMC[1])) continue;
+ //MC Limited Acceptance
+ if (TMath::Abs(containerInputMC[1]) < 0.5) {
+ fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepGeneratedLimAcc, fWeight);
+ AliDebug(3,"MC Lim Acc container filled\n");
+ }
- //MC
- fCFManager->GetParticleContainer()->Fill(containerInputMC, kStepGenerated, fWeight);
- icountMC++;
- AliDebug(3,"MC cointainer filled \n");
+ //MC
+ fCFManager->GetParticleContainer()->Fill(containerInputMC, kStepGenerated, fWeight);
+ icountMC++;
+ AliDebug(3,"MC cointainer filled \n");
- // MC in acceptance
- // check the MC-Acceptance level cuts
- // since standard CF functions are not applicable, using Kine Cuts on daughters
- Bool_t mcAccepStep = cfVtxHF-> MCAcceptanceStep();
- if (mcAccepStep){
- fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepAcceptance, fWeight);
- AliDebug(3,"MC acceptance cut passed\n");
- icountAcc++;
+ // MC in acceptance
+ // check the MC-Acceptance level cuts
+ // since standard CF functions are not applicable, using Kine Cuts on daughters
+ Bool_t mcAccepStep = cfVtxHF-> MCAcceptanceStep();
+ if (mcAccepStep){
+ fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepAcceptance, fWeight);
+ AliDebug(3,"MC acceptance cut passed\n");
+ icountAcc++;
- //MC Vertex step
- if (fCuts->IsEventSelected(aodEvent)){
- // filling the container if the vertex is ok
- fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepVertex, fWeight) ;
- AliDebug(3,"Vertex cut passed and container filled\n");
- icountVertex++;
+ //MC Vertex step
+ if (fCuts->IsEventSelected(aodEvent)){
+ // filling the container if the vertex is ok
+ fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepVertex, fWeight) ;
+ AliDebug(3,"Vertex cut passed and container filled\n");
+ icountVertex++;
- //mc Refit requirement
- Bool_t mcRefitStep = cfVtxHF->MCRefitStep(aodEvent, &trackCuts[0]);
- if (mcRefitStep){
- fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepRefit, fWeight);
- AliDebug(3,"MC Refit cut passed and container filled\n");
- icountRefit++;
- }
- else{
- AliDebug(3,"MC Refit cut not passed\n");
- continue;
- }
- }
- else{
- AliDebug (3, "MC vertex step not passed\n");
- continue;
- }
- }
- else{
- AliDebug (3, "MC in acceptance step not passed\n");
- continue;
- }
- }
- else {
- AliDebug (3, "MC container not filled\n");
- }
+ //mc Refit requirement
+ Bool_t mcRefitStep = cfVtxHF->MCRefitStep(aodEvent, &trackCuts[0]);
+ if (mcRefitStep){
+ fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepRefit, fWeight);
+ AliDebug(3,"MC Refit cut passed and container filled\n");
+ icountRefit++;
+ }
+ else{
+ AliDebug(3,"MC Refit cut not passed\n");
+ continue;
+ }
}
+ else{
+ AliDebug (3, "MC vertex step not passed\n");
+ continue;
+ }
+ }
+ else{
+ AliDebug (3, "MC in acceptance step not passed\n");
+ continue;
+ }
+ }
+ else {
+ AliDebug (3, "MC container not filled\n");
+ }
+ }
- if (cquarks<2) AliDebug(2,Form("Event with %d c-quarks", cquarks));
- AliDebug(2,Form("Found %i MC particles that are %s!!",icountMC,fPartName.Data()));
- AliDebug(2,Form("Found %i MC particles that are %s and satisfy Acc cuts!!",icountAcc,fPartName.Data()));
- AliDebug(2,Form("Found %i MC particles that are %s and satisfy Vertex cuts!!",icountVertex,fPartName.Data()));
- AliDebug(2,Form("Found %i MC particles that are %s and satisfy Refit cuts!!",icountRefit,fPartName.Data()));
-
- // Now go to rec level
- fCountMC += icountMC;
- fCountAcc += icountAcc;
- fCountVertex+= icountVertex;
- fCountRefit+= icountRefit;
-
- AliDebug(2,Form("Found %d vertices for decay channel %d",arrayBranch->GetEntriesFast(),fDecayChannel));
+ if (cquarks<2) AliDebug(2,Form("Event with %d c-quarks", cquarks));
+ AliDebug(2,Form("Found %i MC particles that are %s!!",icountMC,fPartName.Data()));
+ AliDebug(2,Form("Found %i MC particles that are %s and satisfy Acc cuts!!",icountAcc,fPartName.Data()));
+ AliDebug(2,Form("Found %i MC particles that are %s and satisfy Vertex cuts!!",icountVertex,fPartName.Data()));
+ AliDebug(2,Form("Found %i MC particles that are %s and satisfy Refit cuts!!",icountRefit,fPartName.Data()));
+
+ // Now go to rec level
+ fCountMC += icountMC;
+ fCountAcc += icountAcc;
+ fCountVertex+= icountVertex;
+ fCountRefit+= icountRefit;
+
+ AliDebug(2,Form("Found %d vertices for decay channel %d",arrayBranch->GetEntriesFast(),fDecayChannel));
- for(Int_t iCandid = 0; iCandid<arrayBranch->GetEntriesFast();iCandid++){
- AliAODRecoDecayHF* charmCandidate=0x0;
- switch (fDecayChannel){
- case 2:{
- charmCandidate = (AliAODRecoDecayHF2Prong*)arrayBranch->At(iCandid);
- break;
- }
- case 21:{
- charmCandidate = (AliAODRecoCascadeHF*)arrayBranch->At(iCandid);
- break;
- }
- case 31:
- case 32:
- case 33:{
- charmCandidate = (AliAODRecoDecayHF3Prong*)arrayBranch->At(iCandid);
- break;
- }
- case 4:{
- charmCandidate = (AliAODRecoDecayHF4Prong*)arrayBranch->At(iCandid);
- break;
- }
- default:
- break;
- }
+ for(Int_t iCandid = 0; iCandid<arrayBranch->GetEntriesFast();iCandid++){
+ AliAODRecoDecayHF* charmCandidate=0x0;
+ switch (fDecayChannel){
+ case 2:{
+ charmCandidate = (AliAODRecoDecayHF2Prong*)arrayBranch->At(iCandid);
+ break;
+ }
+ case 21:
+ case 22:{
+ charmCandidate = (AliAODRecoCascadeHF*)arrayBranch->At(iCandid);
+ break;
+ }
+ case 31:
+ case 32:
+ case 33:{
+ charmCandidate = (AliAODRecoDecayHF3Prong*)arrayBranch->At(iCandid);
+ break;
+ }
+ case 4:{
+ charmCandidate = (AliAODRecoDecayHF4Prong*)arrayBranch->At(iCandid);
+ break;
+ }
+ default:
+ break;
+ }
- Bool_t unsetvtx=kFALSE;
- if(!charmCandidate->GetOwnPrimaryVtx()) {
- charmCandidate->SetOwnPrimaryVtx(aodVtx); // needed to compute all variables
- unsetvtx=kTRUE;
- }
+ Bool_t unsetvtx=kFALSE;
+ if(!charmCandidate->GetOwnPrimaryVtx()) {
+ charmCandidate->SetOwnPrimaryVtx(aodVtx); // needed to compute all variables
+ unsetvtx=kTRUE;
+ }
- Bool_t signAssociation = cfVtxHF->SetRecoCandidateParam((AliAODRecoDecayHF*)charmCandidate);
- if (!signAssociation){
- charmCandidate = 0x0;
- continue;
- }
+ Bool_t signAssociation = cfVtxHF->SetRecoCandidateParam((AliAODRecoDecayHF*)charmCandidate);
+ if (!signAssociation){
+ charmCandidate = 0x0;
+ continue;
+ }
- Int_t isPartOrAntipart = cfVtxHF->CheckReflexion(fSign);
- if (isPartOrAntipart == 0){
- AliDebug(2, Form("The candidate pdg code doesn't match the requirement set in the task (fSign = %d)",fSign));
- continue;
- }
+ Int_t isPartOrAntipart = cfVtxHF->CheckReflexion(fSign);
+ if (isPartOrAntipart == 0){
+ AliDebug(2, Form("The candidate pdg code doesn't match the requirement set in the task (fSign = %d)",fSign));
+ continue;
+ }
+ AliDebug(3,Form("iCandid=%d - signAssociation=%d, isPartOrAntipart=%d",iCandid, signAssociation, isPartOrAntipart));
- Bool_t recoContFilled = cfVtxHF->FillRecoContainer(containerInput);
- if (recoContFilled){
+ Bool_t recoContFilled = cfVtxHF->FillRecoContainer(containerInput);
+ if (recoContFilled){
- // weight according to pt
- if (fUseWeight){
- if (fFuncWeight){ // using user-defined function
- AliDebug(2, "Using function");
- fWeight = fFuncWeight->Eval(containerInput[0]);
- }
- else{ // using FONLL
- AliDebug(2, "Using FONLL");
- fWeight = GetWeight(containerInput[0]);
- }
- AliDebug(2, Form("pt = %f, weight = %f",containerInput[0], fWeight));
- }
+ // weight according to pt
+ if (fUseWeight){
+ if (fFuncWeight){ // using user-defined function
+ AliDebug(2, "Using function");
+ fWeight = fFuncWeight->Eval(containerInput[0]);
+ }
+ else{ // using FONLL
+ AliDebug(2, "Using FONLL");
+ fWeight = GetWeight(containerInput[0]);
+ }
+ AliDebug(2, Form("pt = %f, weight = %f",containerInput[0], fWeight));
+ }
- if (!fCuts->IsInFiducialAcceptance(containerInput[0],containerInput[1])) continue;
+ if (!fCuts->IsInFiducialAcceptance(containerInput[0],containerInput[1])) continue;
- //Reco Step
- Bool_t recoStep = cfVtxHF->RecoStep();
- Bool_t vtxCheck = fCuts->IsEventSelected(aodEvent);
+ //Reco Step
+ Bool_t recoStep = cfVtxHF->RecoStep();
+ Bool_t vtxCheck = fCuts->IsEventSelected(aodEvent);
- if (recoStep && recoContFilled && vtxCheck){
- fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructed, fWeight) ;
- icountReco++;
- AliDebug(3,"Reco step passed and container filled\n");
+
+ if (recoStep && recoContFilled && vtxCheck){
+ fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructed, fWeight) ;
+ icountReco++;
+ AliDebug(3,"Reco step passed and container filled\n");
- //Reco in the acceptance -- take care of UNFOLDING!!!!
- Bool_t recoAcceptanceStep = cfVtxHF->RecoAcceptStep(&trackCuts[0]);
- if (recoAcceptanceStep) {
- fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoAcceptance, fWeight) ;
- icountRecoAcc++;
- AliDebug(3,"Reco acceptance cut passed and container filled\n");
+ //Reco in the acceptance -- take care of UNFOLDING!!!!
+ Bool_t recoAcceptanceStep = cfVtxHF->RecoAcceptStep(&trackCuts[0]);
+ if (recoAcceptanceStep) {
+ fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoAcceptance, fWeight) ;
+ icountRecoAcc++;
+ AliDebug(3,"Reco acceptance cut passed and container filled\n");
- if(fAcceptanceUnf){
- Double_t fill[4]; //fill response matrix
- Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fill);
- if (bUnfolding) fCorrelation->Fill(fill);
- }
+ if(fAcceptanceUnf){
+ Double_t fill[4]; //fill response matrix
+ Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fill);
+ if (bUnfolding) fCorrelation->Fill(fill);
+ }
- //Number of ITS cluster requirements
- Int_t recoITSnCluster = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kTracks);
- if (recoITSnCluster){
- fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoITSClusters, fWeight) ;
- icountRecoITSClusters++;
- AliDebug(3,"Reco n ITS cluster cut passed and container filled\n");
+ //Number of ITS cluster requirements
+ Int_t recoITSnCluster = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kTracks);
+ if (recoITSnCluster){
+ fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoITSClusters, fWeight) ;
+ icountRecoITSClusters++;
+ AliDebug(3,"Reco n ITS cluster cut passed and container filled\n");
- Bool_t iscutsusingpid = fCuts->GetIsUsePID();
- Int_t recoAnalysisCuts = -1, recoPidSelection = -1;
- fCuts->SetUsePID(kFALSE);
- recoAnalysisCuts = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kCandidate, aodEvent);
-
- if (fDecayChannel==33){ // Ds case, where more possibilities are considered
- Bool_t keepDs=ProcessDs(recoAnalysisCuts);
- if(keepDs) recoAnalysisCuts=3;
- }
+ Bool_t iscutsusingpid = fCuts->GetIsUsePID();
+ Int_t recoAnalysisCuts = -1, recoPidSelection = -1;
+ fCuts->SetUsePID(kFALSE);
+ recoAnalysisCuts = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kCandidate, aodEvent);
+
+ if (fDecayChannel==33){ // Ds case, where more possibilities are considered
+ Bool_t keepDs=ProcessDs(recoAnalysisCuts);
+ if(keepDs) recoAnalysisCuts=3;
+ }
+ else if (fDecayChannel==22){ // Lc->V0+bachelor case, where more possibilities are considered
+ Bool_t keepLctoV0bachelor=ProcessLctoV0Bachelor(recoAnalysisCuts);
+ if (keepLctoV0bachelor) recoAnalysisCuts=3;
+ }
+
- fCuts->SetUsePID(iscutsusingpid); //restoring usage of the PID from the cuts object
- Bool_t tempAn=(recoAnalysisCuts == 3 || recoAnalysisCuts == isPartOrAntipart);
- if (fDecayChannel == 32) tempAn=(recoAnalysisCuts >0 || recoAnalysisCuts == isPartOrAntipart);
+ fCuts->SetUsePID(iscutsusingpid); //restoring usage of the PID from the cuts object
+ Bool_t tempAn=(recoAnalysisCuts == 3 || recoAnalysisCuts == isPartOrAntipart);
+ if (fDecayChannel == 32) tempAn=(recoAnalysisCuts >0 || recoAnalysisCuts == isPartOrAntipart);
- if (tempAn){
- fCFManager->GetParticleContainer()->Fill(containerInput, kStepRecoPPR, fWeight);
- icountRecoPPR++;
- AliDebug(3,"Reco Analysis cuts passed and container filled \n");
- //pid selection
- //recoPidSelection = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kPID);
- //if((fCuts->CombineSelectionLevels(3,recoAnalysisCuts,recoPidSelection)==isPartOrAntipart)||(fCuts->CombineSelectionLevels(3,recoAnalysisCuts,recoPidSelection)==3)){
- recoPidSelection = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kCandidate, aodEvent);
-
- if (fDecayChannel==33){ // Ds case, where more possibilities are considered
- Bool_t keepDs=ProcessDs(recoPidSelection);
- if(keepDs) recoPidSelection=3;
- }
- Bool_t tempPid=(recoPidSelection == 3 || recoPidSelection == isPartOrAntipart);
- if (fDecayChannel == 32) tempPid=(recoPidSelection >0 || recoPidSelection == isPartOrAntipart);
-
- if (tempPid){
- fCFManager->GetParticleContainer()->Fill(containerInput, kStepRecoPID, fWeight);
- icountRecoPID++;
- AliDebug(3,"Reco PID cuts passed and container filled \n");
- if(!fAcceptanceUnf){
- Double_t fill[4]; //fill response matrix
- Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fill);
- if (bUnfolding) fCorrelation->Fill(fill);
- }
- }
- else {
- AliDebug(3, "Analysis Cuts step not passed \n");
- continue;
- }
- }
- else {
- AliDebug(3, "PID selection not passed \n");
- continue;
- }
- }
- else{
- AliDebug(3, "Number of ITS cluster step not passed\n");
- continue;
- }
- }
- else{
- AliDebug(3, "Reco acceptance step not passed\n");
- continue;
- }
- }
- else {
- AliDebug(3, "Reco step not passed\n");
- continue;
- }
+ if (tempAn){
+ fCFManager->GetParticleContainer()->Fill(containerInput, kStepRecoPPR, fWeight);
+ icountRecoPPR++;
+ AliDebug(3,"Reco Analysis cuts passed and container filled \n");
+ //pid selection
+ //recoPidSelection = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kPID);
+ //if((fCuts->CombineSelectionLevels(3,recoAnalysisCuts,recoPidSelection)==isPartOrAntipart)||(fCuts->CombineSelectionLevels(3,recoAnalysisCuts,recoPidSelection)==3)){
+ recoPidSelection = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kCandidate, aodEvent);
+
+ if (fDecayChannel==33){ // Ds case, where more possibilities are considered
+ Bool_t keepDs=ProcessDs(recoPidSelection);
+ if(keepDs) recoPidSelection=3;
+ } else if (fDecayChannel==22){ // Lc->V0+bachelor case, where more possibilities are considered
+ Bool_t keepLctoV0bachelor=ProcessLctoV0Bachelor(recoAnalysisCuts);
+ if (keepLctoV0bachelor) recoAnalysisCuts=3;
+ }
+
+ Bool_t tempPid=(recoPidSelection == 3 || recoPidSelection == isPartOrAntipart);
+ if (fDecayChannel == 32) tempPid=(recoPidSelection >0 || recoPidSelection == isPartOrAntipart);
+
+ if (tempPid){
+ fCFManager->GetParticleContainer()->Fill(containerInput, kStepRecoPID, fWeight);
+ icountRecoPID++;
+ AliDebug(3,"Reco PID cuts passed and container filled \n");
+ if(!fAcceptanceUnf){
+ Double_t fill[4]; //fill response matrix
+ Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fill);
+ if (bUnfolding) fCorrelation->Fill(fill);
}
+ }
+ else {
+ AliDebug(3, "Analysis Cuts step not passed \n");
+ continue;
+ }
+ }
+ else {
+ AliDebug(3, "PID selection not passed \n");
+ continue;
+ }
+ }
+ else{
+ AliDebug(3, "Number of ITS cluster step not passed\n");
+ continue;
+ }
+ }
+ else{
+ AliDebug(3, "Reco acceptance step not passed\n");
+ continue;
+ }
+ }
+ else {
+ AliDebug(3, "Reco step not passed\n");
+ continue;
+ }
+ }
- if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
- } // end loop on candidate
+ if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
+ } // end loop on candidate
- fCountReco+= icountReco;
- fCountRecoAcc+= icountRecoAcc;
- fCountRecoITSClusters+= icountRecoITSClusters;
- fCountRecoPPR+= icountRecoPPR;
- fCountRecoPID+= icountRecoPID;
+ fCountReco+= icountReco;
+ fCountRecoAcc+= icountRecoAcc;
+ fCountRecoITSClusters+= icountRecoITSClusters;
+ fCountRecoPPR+= icountRecoPPR;
+ fCountRecoPID+= icountRecoPID;
- fHistEventsProcessed->Fill(0);
-
- delete[] containerInput;
- delete[] containerInputMC;
- delete cfVtxHF;
- if (trackCuts){
- // for (Int_t i=0; i<cfVtxHF->GetNProngs(); i++){
- // delete [] trackCuts[i];
- // }
- delete [] trackCuts;
- }
+ fHistEventsProcessed->Fill(0);
+
+ delete[] containerInput;
+ delete[] containerInputMC;
+ delete cfVtxHF;
+ if (trackCuts){
+ // for (Int_t i=0; i<cfVtxHF->GetNProngs(); i++){
+ // delete [] trackCuts[i];
+ // }
+ delete [] trackCuts;
+ }
}
//___________________________________________________________________________
void AliCFTaskVertexingHF::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.
+ // 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.
- AliAnalysisTaskSE::Terminate();
-
- AliInfo(Form("Found %i MC particles that are %s in MC, in %d events",fCountMC,fPartName.Data(),fEvents));
- AliInfo(Form("Found %i MC particles that are %s in MC and satisfy Acc cuts, in %d events",fCountAcc,fPartName.Data(),fEvents));
- AliInfo(Form("Found %i MC particles that are %s in MC and satisfy Acc cuts, and satisfy Vertex requirement in %d events",fCountVertex,fPartName.Data(),fEvents));
- AliInfo(Form("Found %i MC particles that are %s in MC and satisfy Acc cuts, and satisfy ITS+TPC refit requirementin %d events",fCountRefit,fPartName.Data(),fEvents));
- AliInfo(Form("Found %i reco %s that are decaying in %s, in %d events",fCountReco,fPartName.Data(),fDauNames.Data(),fEvents));
- AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and are in the requested acceptance, in %d events",fCountRecoAcc,fPartName.Data(),fDauNames.Data(),fEvents));
- AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and have at least 5 clusters in ITS, in %d events",fCountRecoITSClusters,fPartName.Data(),fDauNames.Data(),fEvents));
- AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and satisfy PPR cuts, in %d events",fCountRecoPPR,fPartName.Data(),fDauNames.Data(),fEvents));
- AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and satisfy PPR+PID cuts, in %d events",fCountRecoPID,fPartName.Data(),fDauNames.Data(),fEvents));
+ AliAnalysisTaskSE::Terminate();
+
+ AliInfo(Form("Found %i MC particles that are %s in MC, in %d events",fCountMC,fPartName.Data(),fEvents));
+ AliInfo(Form("Found %i MC particles that are %s in MC and satisfy Acc cuts, in %d events",fCountAcc,fPartName.Data(),fEvents));
+ AliInfo(Form("Found %i MC particles that are %s in MC and satisfy Acc cuts, and satisfy Vertex requirement in %d events",fCountVertex,fPartName.Data(),fEvents));
+ AliInfo(Form("Found %i MC particles that are %s in MC and satisfy Acc cuts, and satisfy ITS+TPC refit requirementin %d events",fCountRefit,fPartName.Data(),fEvents));
+ AliInfo(Form("Found %i reco %s that are decaying in %s, in %d events",fCountReco,fPartName.Data(),fDauNames.Data(),fEvents));
+ AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and are in the requested acceptance, in %d events",fCountRecoAcc,fPartName.Data(),fDauNames.Data(),fEvents));
+ AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and have at least 5 clusters in ITS, in %d events",fCountRecoITSClusters,fPartName.Data(),fDauNames.Data(),fEvents));
+ AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and satisfy PPR cuts, in %d events",fCountRecoPPR,fPartName.Data(),fDauNames.Data(),fEvents));
+ AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and satisfy PPR+PID cuts, in %d events",fCountRecoPID,fPartName.Data(),fDauNames.Data(),fEvents));
- // draw some example plots....
- AliCFContainer *cont= dynamic_cast<AliCFContainer*> (GetOutputData(2));
- if(!cont) {
- printf("CONTAINER NOT FOUND\n");
- return;
- }
- // projecting the containers to obtain histograms
- // first argument = variable, second argument = step
-
- TH1D** h = new TH1D*[3];
- if (fConfiguration == kSnail){
- //h = new TH1D[3][12];
- for (Int_t ih = 0; ih<3; ih++){
- h[ih] = new TH1D[12];
- }
- for(Int_t iC=1;iC<12; iC++){
- // MC-level
- h[0][iC] = *(cont->ShowProjection(iC,0));
- // MC-Acceptance level
- h[1][iC] = *(cont->ShowProjection(iC,1));
- // Reco-level
- h[2][iC] = *(cont->ShowProjection(iC,4));
- }
- }
- else{
- //h = new TH1D[3][12];
- for (Int_t ih = 0; ih<3; ih++){
- h[ih] = new TH1D[8];
- }
- for(Int_t iC=0;iC<8; iC++){
- // MC-level
- h[0][iC] = *(cont->ShowProjection(iC,0));
- // MC-Acceptance level
- h[1][iC] = *(cont->ShowProjection(iC,1));
- // Reco-level
- h[2][iC] = *(cont->ShowProjection(iC,4));
- }
- }
- TString* titles;
- Int_t nvarToPlot = 0;
- if (fConfiguration == kSnail){
- nvarToPlot = 12;
- titles = new TString[nvarToPlot];
- if(fDecayChannel==31){
- titles[0]="pT_Dplus (GeV/c)";
- titles[1]="rapidity";
- titles[2]="phi (rad)";
- titles[3]="cT (#mum)";
- titles[4]="cosPointingAngle";
- titles[5]="pT_1 (GeV/c)";
- titles[6]="pT_2 (GeV/c)";
- titles[7]="pT_3 (GeV/c)";
- titles[8]="d0_1 (#mum)";
- titles[9]="d0_2 (#mum)";
- titles[10]="d0_3 (#mum)";
- titles[11]="zVertex (cm)";
- }else{
- titles[0]="pT_D0 (GeV/c)";
- titles[1]="rapidity";
- titles[2]="cosThetaStar";
- titles[3]="pT_pi (GeV/c)";
- titles[4]="pT_K (Gev/c)";
- titles[5]="cT (#mum)";
- titles[6]="dca (#mum)";
- titles[7]="d0_pi (#mum)";
- titles[8]="d0_K (#mum)";
- titles[9]="d0xd0 (#mum^2)";
- titles[10]="cosPointingAngle";
- titles[11]="phi (rad)";
-
- }
- }
- else{
- nvarToPlot = 8;
- titles = new TString[nvarToPlot];
- titles[0]="pT_candidate (GeV/c)";
- titles[1]="rapidity";
- titles[2]="cT (#mum)";
- titles[3]="phi";
- titles[4]="z_{vtx}";
- titles[5]="centrality";
- titles[6]="fake";
- titles[7]="multiplicity";
- }
+ // draw some example plots....
+ AliCFContainer *cont= dynamic_cast<AliCFContainer*> (GetOutputData(2));
+ if(!cont) {
+ printf("CONTAINER NOT FOUND\n");
+ return;
+ }
+ // projecting the containers to obtain histograms
+ // first argument = variable, second argument = step
- Int_t markers[12]={20,24,21,25,27,28,
- 20,24,21,25,27,28};
- Int_t colors[3]={2,8,4};
- for(Int_t iC=0;iC<nvarToPlot; iC++){
- for(Int_t iStep=0;iStep<3;iStep++){
- h[iStep][iC].SetTitle(titles[iC].Data());
- h[iStep][iC].GetXaxis()->SetTitle(titles[iC].Data());
- Double_t maxh=h[iStep][iC].GetMaximum();
- h[iStep][iC].GetYaxis()->SetRangeUser(0,maxh*1.2);
- h[iStep][iC].SetMarkerStyle(markers[iC]);
- h[iStep][iC].SetMarkerColor(colors[iStep]);
- }
- }
+ TH1D** h = new TH1D*[3];
+ Int_t nvarToPlot = 0;
+ if (fConfiguration == kSnail){
+ //h = new TH1D[3][12];
+ for (Int_t ih = 0; ih<3; ih++){
+ if(fDecayChannel==22){
+ nvarToPlot = 16;
+ } else {
+ nvarToPlot = 12;
+ }
+ h[ih] = new TH1D[nvarToPlot];
+ }
+ for(Int_t iC=1;iC<nvarToPlot; iC++){
+ // MC-level
+ h[0][iC] = *(cont->ShowProjection(iC,0));
+ // MC-Acceptance level
+ h[1][iC] = *(cont->ShowProjection(iC,1));
+ // Reco-level
+ h[2][iC] = *(cont->ShowProjection(iC,4));
+ }
+ }
+ else {
+ //h = new TH1D[3][12];
+ nvarToPlot = 8;
+ for (Int_t ih = 0; ih<3; ih++){
+ h[ih] = new TH1D[nvarToPlot];
+ }
+ for(Int_t iC=0;iC<nvarToPlot; iC++){
+ // MC-level
+ h[0][iC] = *(cont->ShowProjection(iC,0));
+ // MC-Acceptance level
+ h[1][iC] = *(cont->ShowProjection(iC,1));
+ // Reco-level
+ h[2][iC] = *(cont->ShowProjection(iC,4));
+ }
+ }
+ TString* titles;
+ //Int_t nvarToPlot = 0;
+ if (fConfiguration == kSnail){
+ if(fDecayChannel==31){
+ //nvarToPlot = 12;
+ titles = new TString[nvarToPlot];
+ titles[0]="pT_Dplus (GeV/c)";
+ titles[1]="rapidity";
+ titles[2]="phi (rad)";
+ titles[3]="cT (#mum)";
+ titles[4]="cosPointingAngle";
+ titles[5]="pT_1 (GeV/c)";
+ titles[6]="pT_2 (GeV/c)";
+ titles[7]="pT_3 (GeV/c)";
+ titles[8]="d0_1 (#mum)";
+ titles[9]="d0_2 (#mum)";
+ titles[10]="d0_3 (#mum)";
+ titles[11]="zVertex (cm)";
+ } else if (fDecayChannel==22) {
+ //nvarToPlot = 16;
+ titles = new TString[nvarToPlot];
+ titles[0]="pT_Lc (GeV/c)";
+ titles[1]="rapidity";
+ titles[2]="phi (rad)";
+ titles[3]="cosPAV0";
+ titles[4]="onTheFlyStatusV0";
+ titles[5]="centrality";
+ titles[6]="fake";
+ titles[7]="multiplicity";
+ titles[8]="pT_bachelor (GeV/c)";
+ titles[9]="pT_V0pos (GeV/c)";
+ titles[10]="pT_V0neg (GeV/c)";
+ titles[11]="invMassV0 (GeV/c2)";
+ titles[12]="dcaV0 (nSigma)";
+ titles[13]="c#tauV0 (#mum)";
+ titles[14]="c#tau (#mum)";
+ titles[15]="cosPA";
+ } else {
+ //nvarToPlot = 12;
+ titles = new TString[nvarToPlot];
+ titles[0]="pT_D0 (GeV/c)";
+ titles[1]="rapidity";
+ titles[2]="cosThetaStar";
+ titles[3]="pT_pi (GeV/c)";
+ titles[4]="pT_K (Gev/c)";
+ titles[5]="cT (#mum)";
+ titles[6]="dca (#mum)";
+ titles[7]="d0_pi (#mum)";
+ titles[8]="d0_K (#mum)";
+ titles[9]="d0xd0 (#mum^2)";
+ titles[10]="cosPointingAngle";
+ titles[11]="phi (rad)";
+ }
+ }
+ else {
+ //nvarToPlot = 8;
+ titles = new TString[nvarToPlot];
+ if (fDecayChannel==22) {
+ titles[0]="pT_candidate (GeV/c)";
+ titles[1]="rapidity";
+ titles[2]="phi (rad)";
+ titles[3]="cosPAV0";
+ titles[4]="onTheFlyStatusV0";
+ titles[5]="centrality";
+ titles[6]="fake";
+ titles[7]="multiplicity";
+ } else {
+ titles[0]="pT_candidate (GeV/c)";
+ titles[1]="rapidity";
+ titles[2]="cT (#mum)";
+ titles[3]="phi";
+ titles[4]="z_{vtx}";
+ titles[5]="centrality";
+ titles[6]="fake";
+ titles[7]="multiplicity";
+ }
+ }
+
+ Int_t markers[16]={20,24,21,25,27,28,
+ 20,24,21,25,27,28,
+ 20,24,21,25};
+ Int_t colors[3]={2,8,4};
+ for(Int_t iC=0;iC<nvarToPlot; iC++){
+ for(Int_t iStep=0;iStep<3;iStep++){
+ h[iStep][iC].SetTitle(titles[iC].Data());
+ h[iStep][iC].GetXaxis()->SetTitle(titles[iC].Data());
+ Double_t maxh=h[iStep][iC].GetMaximum();
+ h[iStep][iC].GetYaxis()->SetRangeUser(0,maxh*1.2);
+ h[iStep][iC].SetMarkerStyle(markers[iC]);
+ h[iStep][iC].SetMarkerColor(colors[iStep]);
+ }
+ }
- gStyle->SetCanvasColor(0);
- gStyle->SetFrameFillColor(0);
- gStyle->SetTitleFillColor(0);
- gStyle->SetStatColor(0);
+ gStyle->SetCanvasColor(0);
+ gStyle->SetFrameFillColor(0);
+ gStyle->SetTitleFillColor(0);
+ gStyle->SetStatColor(0);
- // drawing in 2 separate canvas for a matter of clearity
- TCanvas * c1 =new TCanvas(Form("c1New_%d",fDecayChannel),"Vars 0, 1, 2, 3",1100,1200);
- c1->Divide(3,4);
- Int_t iPad=1;
- for(Int_t iVar=0; iVar<4; iVar++){
- c1->cd(iPad++);
- h[0][iVar].DrawCopy("p");
- c1->cd(iPad++);
- h[1][iVar].DrawCopy("p");
- c1->cd(iPad++);
- h[2][iVar].DrawCopy("p");
- }
+ // drawing in 2 separate canvas for a matter of clearity
+ TCanvas * c1 =new TCanvas(Form("c1New_%d",fDecayChannel),"Vars 0, 1, 2, 3",1100,1200);
+ c1->Divide(3,4);
+ Int_t iPad=1;
+ for(Int_t iVar=0; iVar<4; iVar++){
+ c1->cd(iPad++);
+ h[0][iVar].DrawCopy("p");
+ c1->cd(iPad++);
+ h[1][iVar].DrawCopy("p");
+ c1->cd(iPad++);
+ h[2][iVar].DrawCopy("p");
+ }
- TCanvas * c2 =new TCanvas(Form("c2New_%d",fDecayChannel),"Vars 4, 5, 6, 7",1100,1200);
- c2->Divide(3,4);
- iPad=1;
- for(Int_t iVar=4; iVar<8; iVar++){
- c2->cd(iPad++);
- h[0][iVar].DrawCopy("p");
- c2->cd(iPad++);
- h[1][iVar].DrawCopy("p");
- c2->cd(iPad++);
- h[2][iVar].DrawCopy("p");
- }
+ TCanvas * c2 =new TCanvas(Form("c2New_%d",fDecayChannel),"Vars 4, 5, 6, 7",1100,1200);
+ c2->Divide(3,4);
+ iPad=1;
+ for(Int_t iVar=4; iVar<8; iVar++){
+ c2->cd(iPad++);
+ h[0][iVar].DrawCopy("p");
+ c2->cd(iPad++);
+ h[1][iVar].DrawCopy("p");
+ c2->cd(iPad++);
+ h[2][iVar].DrawCopy("p");
+ }
- if (fConfiguration == kSnail){
- TCanvas * c3 =new TCanvas(Form("c3New_%d",fDecayChannel),"Vars 8, 9, 10, 11",1100,1200);
- c3->Divide(3,4);
- iPad=1;
- for(Int_t iVar=8; iVar<12; iVar++){
- c3->cd(iPad++);
- h[0][iVar].DrawCopy("p");
- c3->cd(iPad++);
- h[1][iVar].DrawCopy("p");
- c3->cd(iPad++);
- h[2][iVar].DrawCopy("p");
- }
- }
+ if (fConfiguration == kSnail){
+ TCanvas * c3 =new TCanvas(Form("c3New_%d",fDecayChannel),"Vars 8, 9, 10, 11",1100,1200);
+ c3->Divide(3,4);
+ iPad=1;
+ for(Int_t iVar=8; iVar<12; iVar++){
+ c3->cd(iPad++);
+ h[0][iVar].DrawCopy("p");
+ c3->cd(iPad++);
+ h[1][iVar].DrawCopy("p");
+ c3->cd(iPad++);
+ h[2][iVar].DrawCopy("p");
+ }
+ if (fDecayChannel==22) {
+ TCanvas * c4 =new TCanvas(Form("c4New_%d",fDecayChannel),"Vars 12, 13, 14, 15",1100,1200);
+ c4->Divide(3,4);
+ iPad=1;
+ for(Int_t iVar=12; iVar<16; iVar++){
+ c4->cd(iPad++);
+ h[0][iVar].DrawCopy("p");
+ c4->cd(iPad++);
+ h[1][iVar].DrawCopy("p");
+ c4->cd(iPad++);
+ h[2][iVar].DrawCopy("p");
+ }
+ }
+ }
- /*
- THnSparseD* hcorr = dynamic_cast<THnSparseD*> (GetOutputData(3));
+ /*
+ THnSparseD* hcorr = dynamic_cast<THnSparseD*> (GetOutputData(3));
- TH2D* corr1 =hcorr->Projection(0,2);
- TH2D* corr2 = hcorr->Projection(1,3);
+ TH2D* corr1 = hcorr->Projection(0,2);
+ TH2D* corr2 = hcorr->Projection(1,3);
- TCanvas * c7 =new TCanvas(Form("c7New_%d",fDecayChannel),"",800,400);
- c7->Divide(2,1);
- c7->cd(1);
- corr1->DrawCopy("text");
- c7->cd(2);
- corr2->DrawCopy("text");
- */
- TFile* file_projection = new TFile("CFtaskHFprojectionNew.root","RECREATE");
+ TCanvas * c7 =new TCanvas(Form("c7New_%d",fDecayChannel),"",800,400);
+ c7->Divide(2,1);
+ c7->cd(1);
+ corr1->DrawCopy("text");
+ c7->cd(2);
+ corr2->DrawCopy("text");
+ */
+ TFile* file_projection = new TFile("CFtaskHFprojectionNew.root","RECREATE");
- // corr1->Write();
- // corr2->Write();
+ // corr1->Write();
+ // corr2->Write();
- for(Int_t iC=0;iC<nvarToPlot; iC++){
- for(Int_t iStep=0;iStep<3;iStep++){
- h[iStep][iC].Write(Form("Step%d_%s",iStep,titles[iC].Data()));
- }
- }
- file_projection->Close();
- for (Int_t ih = 0; ih<3; ih++) delete [] h[ih];
- delete [] h;
+ for(Int_t iC=0;iC<nvarToPlot; iC++){
+ for(Int_t iStep=0;iStep<3;iStep++){
+ h[iStep][iC].Write(Form("Step%d_%s",iStep,titles[iC].Data()));
+ }
+ }
+ file_projection->Close();
+ for (Int_t ih = 0; ih<3; ih++) delete [] h[ih];
+ delete [] h;
}
//___________________________________________________________________________
void AliCFTaskVertexingHF::UserCreateOutputObjects()
{
- //HERE ONE CAN CREATE OUTPUT OBJECTS, IN PARTICULAR IF THE OBJECT PARAMETERS DON'T NEED
- //TO BE SET BEFORE THE EXECUTION OF THE TASK
- //
- Info("UserCreateOutputObjects","CreateOutputObjects of task %s\n", GetName());
+ //HERE ONE CAN CREATE OUTPUT OBJECTS, IN PARTICULAR IF THE OBJECT PARAMETERS DON'T NEED
+ //TO BE SET BEFORE THE EXECUTION OF THE TASK
+ //
+ Info("UserCreateOutputObjects","CreateOutputObjects of task %s\n", GetName());
- //slot #1
- OpenFile(1);
- const char* nameoutput=GetOutputSlot(1)->GetContainer()->GetName();
- fHistEventsProcessed = new TH1I(nameoutput,"",1,0,1) ;
+ AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
+ AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
+ AliPIDResponse *localPIDResponse = (AliPIDResponse*)inputHandler->GetPIDResponse();
+
+ if (fCuts->GetIsUsePID() && fDecayChannel==22) {
+
+ fCuts->GetPidHF()->SetPidResponse(localPIDResponse);
+ (dynamic_cast<AliRDHFCutsLctoV0*>(fCuts))->GetPidV0pos()->SetPidResponse(localPIDResponse);
+ (dynamic_cast<AliRDHFCutsLctoV0*>(fCuts))->GetPidV0neg()->SetPidResponse(localPIDResponse);
+ fCuts->GetPidHF()->SetOldPid(kFALSE);
+ (dynamic_cast<AliRDHFCutsLctoV0*>(fCuts))->GetPidV0pos()->SetOldPid(kFALSE);
+ (dynamic_cast<AliRDHFCutsLctoV0*>(fCuts))->GetPidV0neg()->SetOldPid(kFALSE);
+ }
- PostData(1,fHistEventsProcessed) ;
- PostData(2,fCFManager->GetParticleContainer()) ;
- PostData(3,fCorrelation) ;
+ //slot #1
+ OpenFile(1);
+ const char* nameoutput=GetOutputSlot(1)->GetContainer()->GetName();
+ fHistEventsProcessed = new TH1I(nameoutput,"",1,0,1) ;
+
+ PostData(1,fHistEventsProcessed) ;
+ PostData(2,fCFManager->GetParticleContainer()) ;
+ PostData(3,fCorrelation) ;
}
//_________________________________________________________________________
Double_t AliCFTaskVertexingHF::GetWeight(Float_t pt)
{
- //
- // calculating the weight to fill the container
- //
+ //
+ // calculating the weight to fill the container
+ //
- // FNOLL central:
- // p0 = 1.63297e-01 --> 0.322643
- // p1 = 2.96275e+00
- // p2 = 2.30301e+00
- // p3 = 2.50000e+00
-
- // PYTHIA
- // p0 = 1.85906e-01 --> 0.36609
- // p1 = 1.94635e+00
- // p2 = 1.40463e+00
- // p3 = 2.50000e+00
-
- Double_t func1[4] = {0.322643,2.96275,2.30301,2.5};
- Double_t func2[4] = {0.36609,1.94635,1.40463,2.5};
-
- Double_t dndpt_func1 = dNdptFit(pt,func1);
- if(fUseFlatPtWeight) dndpt_func1 = 1./30.;
- Double_t dndpt_func2 = dNdptFit(pt,func2);
- AliDebug(2,Form("pt = %f, FONLL = %f, Pythia = %f, ratio = %f",pt,dndpt_func1,dndpt_func2,dndpt_func1/dndpt_func2));
- return dndpt_func1/dndpt_func2;
+ // FNOLL central:
+ // p0 = 1.63297e-01 --> 0.322643
+ // p1 = 2.96275e+00
+ // p2 = 2.30301e+00
+ // p3 = 2.50000e+00
+
+ // PYTHIA
+ // p0 = 1.85906e-01 --> 0.36609
+ // p1 = 1.94635e+00
+ // p2 = 1.40463e+00
+ // p3 = 2.50000e+00
+
+ Double_t func1[4] = {0.322643,2.96275,2.30301,2.5};
+ Double_t func2[4] = {0.36609,1.94635,1.40463,2.5};
+
+ Double_t dndpt_func1 = dNdptFit(pt,func1);
+ if(fUseFlatPtWeight) dndpt_func1 = 1./30.;
+ Double_t dndpt_func2 = dNdptFit(pt,func2);
+ AliDebug(2,Form("pt = %f, FONLL = %f, Pythia = %f, ratio = %f",pt,dndpt_func1,dndpt_func2,dndpt_func1/dndpt_func2));
+ return dndpt_func1/dndpt_func2;
}
//__________________________________________________________________________________________________
Double_t AliCFTaskVertexingHF::dNdptFit(Float_t pt, Double_t* par)
{
- //
- // calculating dNdpt
- //
+ //
+ // calculating dNdpt
+ //
- Double_t denom = TMath::Power((pt/par[1]), par[3] );
- Double_t dNdpt = par[0]*pt/TMath::Power(1.+denom, par[2]);
+ Double_t denom = TMath::Power((pt/par[1]), par[3] );
+ Double_t dNdpt = par[0]*pt/TMath::Power(1.+denom, par[2]);
- return dNdpt;
+ return dNdpt;
}
//__________________________________________________________________________________________________
}
}
+//__________________________________________________________________________________________________
Bool_t AliCFTaskVertexingHF::ProcessDs(Int_t recoAnalysisCuts) const{
// processes output of Ds is selected
Bool_t keep=kFALSE;
}
return keep;
}
+//__________________________________________________________________________________________________
+Bool_t AliCFTaskVertexingHF::ProcessLctoV0Bachelor(Int_t recoAnalysisCuts) const{
+ // processes output of Lc->V0+bnachelor is selected
+ Bool_t keep=kFALSE;
+ if(recoAnalysisCuts > 0){
+
+ Int_t isK0Sp = recoAnalysisCuts&1;
+ Int_t isLambdaBarpi = recoAnalysisCuts&2;
+ Int_t isLambdapi = recoAnalysisCuts&4;
+
+ if(fLctoV0bachelorOption==1){
+ if(isK0Sp) keep=kTRUE;
+ }
+ else if(fLctoV0bachelorOption==2){
+ if(isLambdaBarpi) keep=kTRUE;
+ }
+ else if(fLctoV0bachelorOption==4){
+ if(isLambdapi) keep=kTRUE;
+ }
+ else if(fLctoV0bachelorOption==7) {
+ if (isK0Sp || isLambdaBarpi || isLambdapi) keep=kTRUE;
+ }
+ }
+ return keep;
+}
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2011, 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. *
+ **************************************************************************/
+
+/* $Id$ */
+
+//-----------------------------------------------------------------------
+// Class for HF corrections as a function of many variables and steps
+// For Lc->V0+bachelor
+//
+// - Based on AliCFVertexingHFCascade -
+//
+// Contact : A.De Caro - decaro@sa.infn.it
+// Centro 'E.Fermi' - Rome (Italy)
+// INFN and University of Salerno (Italy)
+//
+//-----------------------------------------------------------------------
+
+#include "AliAODRecoDecayHF2Prong.h"
+#include "AliAODMCParticle.h"
+#include "AliAODEvent.h"
+#include "TClonesArray.h"
+#include "AliCFVertexingHF.h"
+#include "AliESDtrack.h"
+#include "TDatabasePDG.h"
+#include "AliAODRecoCascadeHF.h"
+#include "AliAODv0.h"
+#include "AliCFVertexingHFLctoV0bachelor.h"
+#include "AliCFContainer.h"
+#include "AliCFTaskVertexingHF.h"
+
+#include <Riostream.h>
+
+using std::cout;
+using std::endl;
+
+ClassImp(AliCFVertexingHFLctoV0bachelor)
+
+//_________________________________________
+AliCFVertexingHFLctoV0bachelor::AliCFVertexingHFLctoV0bachelor():
+fGenLcOption(0)
+{
+ // standard constructor
+}
+
+//_____________________________________
+AliCFVertexingHFLctoV0bachelor::AliCFVertexingHFLctoV0bachelor(TClonesArray *mcArray, UShort_t originDselection, Int_t lcDecay):
+AliCFVertexingHF(mcArray, originDselection),
+ fGenLcOption(lcDecay)
+{
+ // standard constructor
+
+ SetNProngs(3);
+ fPtAccCut=new Float_t[fProngs];
+ fEtaAccCut=new Float_t[fProngs];
+ for(Int_t iP=0; iP<fProngs; iP++){
+ fPtAccCut[iP]=0.1;
+ fEtaAccCut[iP]=0.9;
+ }
+
+}
+
+
+//_____________________________________
+AliCFVertexingHFLctoV0bachelor& AliCFVertexingHFLctoV0bachelor::operator=(const AliCFVertexingHFLctoV0bachelor& c)
+{
+ // operator =
+
+ if (this != &c) {
+ AliCFVertexingHF::operator=(c);
+ }
+
+ return *this;
+
+}
+
+//__________________________________________
+Bool_t AliCFVertexingHFLctoV0bachelor::SetRecoCandidateParam(AliAODRecoDecayHF *recoCand)
+{
+ // set the AliAODRecoDecay candidate
+
+ Bool_t bSignAssoc = kFALSE;
+
+ fRecoCandidate = recoCand;
+ if (!fRecoCandidate) {
+ AliError("fRecoCandidate not found, problem in assignement\n");
+ return bSignAssoc;
+ }
+
+ if (fRecoCandidate->GetPrimaryVtx()) AliDebug(3,"fReco Candidate has a pointer to PrimVtx\n");
+
+ AliAODRecoCascadeHF* lcV0bachelor = (AliAODRecoCascadeHF*)fRecoCandidate;
+
+ if ( !(lcV0bachelor->Getv0()) ) return bSignAssoc;
+
+ Int_t nentries = fmcArray->GetEntriesFast();
+ AliDebug(3,Form("nentries = %d\n", nentries));
+
+ Int_t pdgCand = 4122;
+ Int_t mcLabel = -1;
+ Int_t mcLabelK0S = -1;
+ Int_t mcLabelLambda = -1;
+
+ // Lc->K0S+p and cc
+ Int_t pdgDgLctoV0bachelor[2]={310,2212};
+ Int_t pdgDgV0toDaughters[2]={211,211};
+ mcLabelK0S = lcV0bachelor->MatchToMC(pdgCand,pdgDgLctoV0bachelor[0],pdgDgLctoV0bachelor,pdgDgV0toDaughters,fmcArray,kTRUE);
+ // Lc->Lambda+pi and cc
+ pdgDgLctoV0bachelor[0]=3122, pdgDgLctoV0bachelor[1]=211;
+ pdgDgV0toDaughters[0]=2212, pdgDgV0toDaughters[1]=211;
+ mcLabelLambda = lcV0bachelor->MatchToMC(pdgCand,pdgDgLctoV0bachelor[0],pdgDgLctoV0bachelor,pdgDgV0toDaughters,fmcArray,kTRUE);
+
+ if (mcLabelK0S!=-1 && mcLabelLambda!=-1)
+ AliInfo("Strange: current Lc->V0+bachelor candidate has two MC different labels!");
+
+ if (fGenLcOption==kCountAllLctoV0) {
+ if (mcLabelK0S!=-1) mcLabel=mcLabelK0S;
+ if (mcLabelLambda!=-1) mcLabel=mcLabelLambda;
+ }
+ else if (fGenLcOption==kCountK0Sp) {
+ if (mcLabelK0S!=-1) mcLabel=mcLabelK0S;
+ if (mcLabelLambda!=-1) {
+ mcLabel=-1;
+ fFake = 0; // fake candidate
+ if (fFakeSelection==1) return bSignAssoc;
+ }
+ }
+ else if (fGenLcOption==kCountLambdapi) {
+ if (mcLabelLambda!=-1) mcLabel=mcLabelLambda;
+ if (mcLabelK0S!=-1) {
+ mcLabel=-1;
+ fFake = 0; // fake candidate
+ if (fFakeSelection==1) return bSignAssoc;
+ }
+ }
+
+ if (mcLabel==-1) return bSignAssoc;
+
+ if (fRecoCandidate->NumberOfFakeDaughters()>0){
+ fFake = 0; // fake candidate
+ if (fFakeSelection==1) return bSignAssoc;
+ }
+ if (fRecoCandidate->NumberOfFakeDaughters()==0){
+ fFake = 2; // non-fake candidate
+ if (fFakeSelection==2) return bSignAssoc;
+ }
+
+ SetMCLabel(mcLabel);
+ fmcPartCandidate = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fmcLabel));
+
+ if (!fmcPartCandidate){
+ AliDebug(3,"No part candidate");
+ return bSignAssoc;
+ }
+
+ bSignAssoc = kTRUE;
+ return bSignAssoc;
+}
+
+//______________________________________________
+Bool_t AliCFVertexingHFLctoV0bachelor::GetGeneratedValuesFromMCParticle(Double_t* vectorMC)
+{
+ //
+ // collecting all the necessary info (pt, y, invMassV0, cosPAwrtPVxV0, onTheFlyStatusV0) from MC particle
+ // (additional infos: pTbachelor, pTV0pos, pTV0neg, phi, dcaV0, cTV0, cT, cosPA)
+ //
+
+ Bool_t bGenValues = kFALSE;
+
+ if (fmcPartCandidate->GetNDaughters()!=2) return bGenValues;
+
+ Int_t daughter0lc = fmcPartCandidate->GetDaughter(0);
+ Int_t daughter1lc = fmcPartCandidate->GetDaughter(1);
+
+ //the V0
+ AliAODMCParticle* mcPartDaughterV0=0;
+ if (fGenLcOption==kCountLambdapi) {
+ mcPartDaughterV0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0lc)); // strange baryon (0)
+ AliInfo(Form(" Case Lc->Lambda+pi: (V0) daughter0_Lc=%d (%d)",daughter0lc,mcPartDaughterV0->GetPdgCode()));
+ }
+ else if (fGenLcOption==kCountK0Sp) {
+ AliAODMCParticle* mcPartDaughterK0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1lc)); // strange meson (1)
+ if (!mcPartDaughterK0) return bGenValues;
+ if (TMath::Abs(mcPartDaughterK0->GetPdgCode())!=311) return bGenValues;
+ Int_t daughterK0 = mcPartDaughterK0->GetDaughter(0);
+ mcPartDaughterV0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterK0));
+ if (TMath::Abs(mcPartDaughterV0->GetPdgCode())!=310) return bGenValues;
+ AliInfo(Form(" Case Lc->K0S+p: (V0) daughter1_Lc=%d (%d to %d)",daughter1lc,mcPartDaughterK0->GetPdgCode(),
+ mcPartDaughterV0->GetPdgCode()));
+ }
+
+ //the bachelor
+ AliAODMCParticle* mcPartDaughterBachelor=0;
+ if (fGenLcOption==kCountLambdapi) {
+ mcPartDaughterBachelor = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1lc)); // meson (1)
+ AliInfo(Form(" Case Lc->Lambda+pi: (bachelor) daughter1_Lc=%d (%d)",daughter1lc,mcPartDaughterBachelor->GetPdgCode()));
+ }
+ if (fGenLcOption==kCountK0Sp) {
+ mcPartDaughterBachelor = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0lc)); // baryon (0)
+ AliInfo(Form(" Case Lc->K0S+p: (bachelor) daughter0_Lc=%d (%d)",daughter0lc,mcPartDaughterBachelor->GetPdgCode()));
+ }
+
+ if (!mcPartDaughterV0 || !mcPartDaughterBachelor) return bGenValues;
+
+ if (fGenLcOption==kCountLambdapi) {
+ if ( !(mcPartDaughterV0->GetPdgCode()==3122 &&
+ mcPartDaughterBachelor->GetPdgCode()==211) )
+ AliInfo("It isn't a Lc->Lambda+pion candidate");
+ }
+ if (fGenLcOption==kCountK0Sp) {
+ if ( !(mcPartDaughterV0->GetPdgCode()==310 &&
+ mcPartDaughterBachelor->GetPdgCode()==2212) )
+ AliInfo("It isn't a Lc->K0S+proton candidate");
+ }
+
+ Double_t vtx1[3] = {0,0,0}; // primary vertex
+ fmcPartCandidate->XvYvZv(vtx1); // cm
+
+ // getting vertex from daughters
+ Double_t vtx1daughter0[3] = {0,0,0}; // secondary vertex from daughter 0
+ Double_t vtx1daughter1[3] = {0,0,0}; // secondary vertex from daughter 1
+ mcPartDaughterBachelor->XvYvZv(vtx1daughter0); //cm
+ mcPartDaughterV0->XvYvZv(vtx1daughter1); //cm
+ if (TMath::Abs(vtx1daughter0[0]-vtx1daughter1[0])>1E-5 ||
+ TMath::Abs(vtx1daughter0[1]-vtx1daughter1[1])>1E-5 ||
+ TMath::Abs(vtx1daughter0[2]-vtx1daughter1[2])>1E-5) {
+ AliError("Daughters have different secondary vertex, skipping the track");
+ return bGenValues;
+ }
+
+ // getting the momentum from the daughters
+ Double_t px1[2] = {mcPartDaughterBachelor->Px(), mcPartDaughterV0->Px()};
+ Double_t py1[2] = {mcPartDaughterBachelor->Py(), mcPartDaughterV0->Py()};
+ Double_t pz1[2] = {mcPartDaughterBachelor->Pz(), mcPartDaughterV0->Pz()};
+
+ Int_t nprongs = 2;
+ Short_t charge = mcPartDaughterBachelor->Charge();
+ Double_t d0[2] = {0.,0.};
+ AliAODRecoDecayHF* decayLc = new AliAODRecoDecayHF(vtx1,vtx1daughter0,nprongs,charge,px1,py1,pz1,d0);
+ Double_t cosPAwrtPrimVtxLc = decayLc->CosPointingAngle();
+
+ //V0 daughters
+ Int_t daughter0 = mcPartDaughterV0->GetDaughter(0);
+ Int_t daughter1 = mcPartDaughterV0->GetDaughter(1);
+ AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0)); //V0daughter0
+ AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1)); //V0daughter1
+
+ if (!mcPartDaughter0 || !mcPartDaughter1) return kFALSE;
+
+ // getting vertex from daughters
+ Double_t vtx2daughter0[3] = {0,0,0}; // secondary vertex from daughter 0
+ Double_t vtx2daughter1[3] = {0,0,0}; // secondary vertex from daughter 1
+ mcPartDaughter0->XvYvZv(vtx2daughter0); //cm
+ mcPartDaughter1->XvYvZv(vtx2daughter1); //cm
+ if (TMath::Abs(vtx2daughter0[0]-vtx2daughter1[0])>1E-5 ||
+ TMath::Abs(vtx2daughter0[1]-vtx2daughter1[1])>1E-5 ||
+ TMath::Abs(vtx2daughter0[2]-vtx2daughter1[2])>1E-5) {
+ AliError("Daughters have different secondary vertex, skipping the track");
+ return bGenValues;
+ }
+
+ // always instantiate the AliAODRecoDecay with the positive daughter first, the negative second
+ AliAODMCParticle* positiveDaugh = mcPartDaughter0;
+ AliAODMCParticle* negativeDaugh = mcPartDaughter1;
+ if (mcPartDaughter0->GetPdgCode()<0 && mcPartDaughter1->GetPdgCode()>0){
+ // inverting in case the positive daughter is the second one
+ positiveDaugh = mcPartDaughter1;
+ negativeDaugh = mcPartDaughter0;
+ }
+ // getting the momentum from the daughters
+ Double_t px[2] = {positiveDaugh->Px(), negativeDaugh->Px()};
+ Double_t py[2] = {positiveDaugh->Py(), negativeDaugh->Py()};
+ Double_t pz[2] = {positiveDaugh->Pz(), negativeDaugh->Pz()};
+
+ nprongs = 2;
+ charge = 0;
+ AliAODRecoDecayHF* decay = new AliAODRecoDecayHF(vtx1,vtx2daughter0,nprongs,charge,px,py,pz,d0);
+ Double_t cosPAwrtPrimVtxV0 = decay->CosPointingAngle();
+
+ //ct
+ Double_t cTV0 = 0.;
+ if (fGenLcOption==kCountK0Sp) {
+ cTV0 = decay->Ct(310); // by default wrt Primary Vtx
+ } else if (fGenLcOption==kCountLambdapi) {
+ cTV0 = decay->Ct(3122); // by default wrt Primary Vtx
+ }
+
+ Double_t cTLc = Ctau(fmcPartCandidate); // by default wrt Primary Vtx
+
+ // get the bachelor pT
+ Double_t pTbach = 0.;
+
+ if (TMath::Abs(fmcPartCandidate->GetPdgCode()) == 4122)
+ pTbach = mcPartDaughterBachelor->Pt();
+
+ Double_t invMass = 0.;
+ if (fGenLcOption==kCountK0Sp) {
+ invMass = decay->InvMass2Prongs(0,1,211,211);
+ } else if (fGenLcOption==kCountLambdapi) {
+ if (fmcPartCandidate->GetPdgCode() == 4122)
+ invMass = decay->InvMass2Prongs(0,1,2212,211);
+ else if (fmcPartCandidate->GetPdgCode() ==-4122)
+ invMass = decay->InvMass2Prongs(0,1,211,2212);
+ }
+
+ switch (fConfiguration){
+ case AliCFTaskVertexingHF::kSnail:
+ vectorMC[0] = fmcPartCandidate->Pt();
+ vectorMC[1] = fmcPartCandidate->Y() ;
+ vectorMC[2] = fmcPartCandidate->Phi();
+ vectorMC[3] = cosPAwrtPrimVtxV0;
+ vectorMC[4] = 0; // dummy value x MC, onTheFlyStatus
+ vectorMC[5] = fCentValue; // reconstructed centrality
+ vectorMC[6] = 1; // dummy value x MC, fFake
+ vectorMC[7] = fMultiplicity; // reconstructed multiplicity
+
+ vectorMC[8] = pTbach;
+ vectorMC[9] = positiveDaugh->Pt();
+ vectorMC[10] = negativeDaugh->Pt();
+ vectorMC[11] = invMass;
+ vectorMC[12] = 0; // dummy value x MC, V0 DCA
+ vectorMC[13] = cTV0*1.E4; // in micron
+ vectorMC[14] = cTLc*1.E4; // in micron
+ vectorMC[15] = cosPAwrtPrimVtxLc;
+ break;
+ case AliCFTaskVertexingHF::kCheetah:
+ vectorMC[0] = fmcPartCandidate->Pt();
+ vectorMC[1] = fmcPartCandidate->Y() ;
+ vectorMC[2] = fmcPartCandidate->Phi();
+ vectorMC[3] = cosPAwrtPrimVtxV0;
+ vectorMC[4] = 0; // dummy value x MC, onTheFlyStatus
+ vectorMC[5] = fCentValue; // reconstructed centrality
+ vectorMC[6] = 1; // dummy value x MC, fFake
+ vectorMC[7] = fMultiplicity; // reconstructed multiplicity
+ break;
+ }
+
+ delete decay;
+ delete decayLc;
+
+ bGenValues = kTRUE;
+ return bGenValues;
+
+}
+//____________________________________________
+Bool_t AliCFVertexingHFLctoV0bachelor::GetRecoValuesFromCandidate(Double_t *vectorReco) const
+{
+ // read the variables for the container
+
+ Bool_t bFillRecoValues = kFALSE;
+
+ //Get the Lc and the V0 from Lc
+ AliAODRecoCascadeHF* lcV0bachelor = (AliAODRecoCascadeHF*)fRecoCandidate;
+
+ if ( !(lcV0bachelor->Getv0()) ) return bFillRecoValues;
+
+ AliAODVertex *vtx0 = (AliAODVertex*)lcV0bachelor->GetPrimaryVtx();
+ if (vtx0) AliDebug(1,"lcV0bachelor has primary vtx");
+
+ AliAODTrack* bachelor = (AliAODTrack*)lcV0bachelor->GetBachelor();
+ AliAODv0* v0toDaughters = (AliAODv0*)lcV0bachelor->Getv0();
+ Bool_t onTheFlyStatus = v0toDaughters->GetOnFlyStatus();
+ AliAODTrack* v0positiveTrack = (AliAODTrack*)lcV0bachelor->Getv0PositiveTrack();
+ AliAODTrack* v0negativeTrack = (AliAODTrack*)lcV0bachelor->Getv0NegativeTrack();
+
+ Double_t pt = lcV0bachelor->Pt();
+ Double_t rapidity = lcV0bachelor->Y(4122);
+ Double_t invMassV0 = 0.;
+ Double_t primVtxPos[3]; vtx0->GetXYZ(primVtxPos);
+ Double_t cosPAwrtPrimVtxV0 = v0toDaughters->CosPointingAngle(primVtxPos);
+ Double_t cTV0 = 0.;
+
+ Double_t pTbachelor = bachelor->Pt();
+ Double_t pTV0pos = v0positiveTrack->Pt();
+ Double_t pTV0neg = v0negativeTrack->Pt();
+ Double_t phi = lcV0bachelor->Phi();
+ Double_t dcaV0 = v0toDaughters->GetDCA();
+ Double_t cTLc = lcV0bachelor->Ct(4122); // wrt PrimVtx
+ //Double_t dcaLc = lcV0bachelor->GetDCA();
+ Double_t cosPointingAngleLc = lcV0bachelor->CosPointingAngle();
+
+ if (fGenLcOption==kCountK0Sp) {
+ cTV0 = v0toDaughters->Ct(310,primVtxPos);
+ } else if (fGenLcOption==kCountLambdapi) {
+ cTV0 = v0toDaughters->Ct(3122,primVtxPos);
+ }
+
+ Short_t bachelorCharge = bachelor->Charge();
+ if (fGenLcOption==kCountLambdapi) {
+ if (bachelorCharge==1) {
+ invMassV0 = v0toDaughters->MassLambda();
+ } else if (bachelorCharge==-1) {
+ invMassV0 = v0toDaughters->MassAntiLambda();
+ }
+
+ } else if (fGenLcOption==kCountK0Sp) {
+ invMassV0 = v0toDaughters->MassK0Short();
+ }
+
+ switch (fConfiguration){
+ case AliCFTaskVertexingHF::kSnail:
+ vectorReco[0] = pt;
+ vectorReco[1] = rapidity;
+ vectorReco[2] = phi;
+ vectorReco[3] = cosPAwrtPrimVtxV0;
+ vectorReco[4] = onTheFlyStatus;
+ vectorReco[5] = fCentValue;
+ vectorReco[6] = fFake; // whether the reconstructed candidate was a fake (fFake = 0) or not (fFake = 2)
+ vectorReco[7] = fMultiplicity;
+
+ vectorReco[8] = pTbachelor;
+ vectorReco[9] = pTV0pos;
+ vectorReco[10] = pTV0neg;
+ vectorReco[11] = invMassV0;
+ vectorReco[12] = dcaV0;
+ vectorReco[13] = cTV0*1.E4; // in micron
+ vectorReco[14] = cTLc*1.E4; // in micron
+ vectorReco[15] = cosPointingAngleLc;
+
+ break;
+ case AliCFTaskVertexingHF::kCheetah:
+ vectorReco[0] = pt;
+ vectorReco[1] = rapidity;
+ vectorReco[2] = phi;
+ vectorReco[3] = cosPAwrtPrimVtxV0;
+ vectorReco[4] = onTheFlyStatus;
+ vectorReco[5] = fCentValue;
+ vectorReco[6] = fFake; // whether the reconstructed candidate was a fake (fFake = 0) or not (fFake = 2)
+ vectorReco[7] = fMultiplicity;
+ break;
+ }
+
+ bFillRecoValues = kTRUE;
+
+ return bFillRecoValues;
+}
+
+//_____________________________________________________________
+Bool_t AliCFVertexingHFLctoV0bachelor::CheckMCChannelDecay() const
+{
+ // check the required decay channel
+
+ Bool_t checkCD = kFALSE;
+
+ if (fmcPartCandidate->GetNDaughters()!=2) return checkCD;
+
+ Int_t daughter0 = fmcPartCandidate->GetDaughter(0);
+ Int_t daughter1 = fmcPartCandidate->GetDaughter(1);
+ AliAODMCParticle* mcPartDaughter0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0));
+ AliAODMCParticle* mcPartDaughter1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1));
+
+ if (!mcPartDaughter0 || !mcPartDaughter1) {
+ AliDebug (2,"Problems in the MC Daughters\n");
+ return checkCD;
+ }
+
+ // Lc -> Lambda + pion AND cc
+ if (fGenLcOption==kCountLambdapi) {
+
+ if (!(TMath::Abs(mcPartDaughter0->GetPdgCode())==3122 &&
+ TMath::Abs(mcPartDaughter1->GetPdgCode())==211) &&
+ !(TMath::Abs(mcPartDaughter0->GetPdgCode())==211 &&
+ TMath::Abs(mcPartDaughter1->GetPdgCode())==3122)) {
+ AliDebug(2, "The Lc MC doesn't decay in Lambda+pion (or cc), skipping!!");
+ return checkCD;
+ } else if (TMath::Abs(mcPartDaughter0->GetPdgCode())==3122) {
+ if (mcPartDaughter0->GetNDaughters()!=2) return checkCD;
+ Int_t daughter0D0 = mcPartDaughter0->GetDaughter(0);
+ AliAODMCParticle* mcPartDaughter0D0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0D0));
+ Int_t daughter0D1 = mcPartDaughter0->GetDaughter(1);
+ AliAODMCParticle* mcPartDaughter0D1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter0D1));
+ if (!(TMath::Abs(mcPartDaughter0D0->GetPdgCode())==211 &&
+ TMath::Abs(mcPartDaughter0D1->GetPdgCode())==2212) &&
+ !(TMath::Abs(mcPartDaughter0D0->GetPdgCode())==2212 &&
+ TMath::Abs(mcPartDaughter0D1->GetPdgCode())==211)) {
+ AliDebug(2, "The Lambda MC doesn't decay in pi+proton (or cc), skipping!!");
+ return checkCD;
+ }
+ } else if (TMath::Abs(mcPartDaughter1->GetPdgCode())==3122) {
+ if (mcPartDaughter1->GetNDaughters()!=2) return checkCD;
+ Int_t daughter1D0 = mcPartDaughter1->GetDaughter(0);
+ AliAODMCParticle* mcPartDaughter1D0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1D0));
+ Int_t daughter1D1 = mcPartDaughter1->GetDaughter(1);
+ AliAODMCParticle* mcPartDaughter1D1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter1D1));
+ if (!(TMath::Abs(mcPartDaughter1D0->GetPdgCode())==211 &&
+ TMath::Abs(mcPartDaughter1D1->GetPdgCode())==2212) &&
+ !(TMath::Abs(mcPartDaughter1D0->GetPdgCode())==2212 &&
+ TMath::Abs(mcPartDaughter1D1->GetPdgCode())==211)) {
+ AliDebug(2, "The Lambda MC doesn't decay in pi+proton (or cc), skipping!!");
+ return checkCD;
+ }
+ }
+
+ } else if (fGenLcOption==kCountK0Sp) {
+ // Lc -> K0bar + proton AND cc
+
+ if (!(TMath::Abs(mcPartDaughter0->GetPdgCode())==311 &&
+ TMath::Abs(mcPartDaughter1->GetPdgCode())==2212) &&
+ !(TMath::Abs(mcPartDaughter0->GetPdgCode())==2212 &&
+ TMath::Abs(mcPartDaughter1->GetPdgCode())==311)) {
+ AliDebug(2, "The Lc MC doesn't decay in K0+proton (or cc), skipping!!");
+ return checkCD;
+ }
+
+ if (TMath::Abs(mcPartDaughter0->GetPdgCode())==311) {
+ Int_t daughter = mcPartDaughter0->GetDaughter(0);
+ AliAODMCParticle* mcPartDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter));
+ if (!(TMath::Abs(mcPartDaughter->GetPdgCode())==310)) {
+ AliDebug(2, "The K0 (or K0bar) MC doesn't go in K0S, skipping!!");
+ return checkCD;
+ }
+ if (mcPartDaughter->GetNDaughters()!=2) return checkCD;
+ Int_t daughterD0 = mcPartDaughter->GetDaughter(0);
+ AliAODMCParticle* mcPartDaughterD0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD0));
+ Int_t daughterD1 = mcPartDaughter->GetDaughter(1);
+ AliAODMCParticle* mcPartDaughterD1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD1));
+ if (!(TMath::Abs(mcPartDaughterD0->GetPdgCode())==211 &&
+ TMath::Abs(mcPartDaughterD1->GetPdgCode())==211)) {
+ AliDebug(2, "The K0S MC doesn't decay in pi+pi, skipping!!");
+ return checkCD;
+ }
+
+ }
+
+ if (TMath::Abs(mcPartDaughter1->GetPdgCode())==311) {
+ Int_t daughter = mcPartDaughter1->GetDaughter(0);
+ AliAODMCParticle* mcPartDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughter));
+ if (!(TMath::Abs(mcPartDaughter->GetPdgCode())==310)) {
+ AliDebug(2, "The K0 (or K0bar) MC doesn't go in K0S, skipping!!");
+ return checkCD;
+ }
+ if (mcPartDaughter->GetNDaughters()!=2) return checkCD;
+ Int_t daughterD0 = mcPartDaughter->GetDaughter(0);
+ AliAODMCParticle* mcPartDaughterD0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD0));
+ Int_t daughterD1 = mcPartDaughter->GetDaughter(1);
+ AliAODMCParticle* mcPartDaughterD1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daughterD1));
+ if (! ( TMath::Abs(mcPartDaughterD0->GetPdgCode())==211 &&
+ TMath::Abs(mcPartDaughterD1->GetPdgCode())==211 ) ) {
+ AliDebug(2, "The K0S MC doesn't decay in pi+pi, skipping!!");
+ return checkCD;
+ }
+
+ }
+
+ }
+
+ checkCD = kTRUE;
+ return checkCD;
+
+}
+
+//___________________________________________________________
+
+void AliCFVertexingHFLctoV0bachelor::SetPtAccCut(Float_t* ptAccCut)
+{
+ //
+ // setting the pt cut to be used in the Acceptance steps (MC+Reco)
+ //
+
+ AliInfo("The 1st element of the pt cut array will correspond to the cut applied to the bachelor - please check that it is correct");
+ if (fProngs>0){
+ for (Int_t iP=0; iP<fProngs; iP++){
+ fPtAccCut[iP]=ptAccCut[iP];
+ }
+ }
+ return;
+}
+
+//___________________________________________________________
+
+void AliCFVertexingHFLctoV0bachelor::SetEtaAccCut(Float_t* etaAccCut)
+{
+ //
+ // setting the eta cut to be used in the Acceptance steps (MC+Reco)
+ //
+
+ AliInfo("The 1st element of the pt cut array will correspond to the cut applied to the bachelor - please check that it is correct");
+ if (fProngs>0){
+ for (Int_t iP=0; iP<fProngs; iP++){
+ fEtaAccCut[iP]=etaAccCut[iP];
+ }
+ }
+ return;
+}
+
+//___________________________________________________________
+
+void AliCFVertexingHFLctoV0bachelor::SetAccCut(Float_t* ptAccCut, Float_t* etaAccCut)
+{
+ //
+ // setting the pt and eta cut to be used in the Acceptance steps (MC+Reco)
+ //
+
+ AliInfo("The 1st element of the pt cut array will correspond to the cut applied to the bachelor - please check that it is correct");
+ if (fProngs>0){
+ for (Int_t iP=0; iP<fProngs; iP++){
+ fPtAccCut[iP]=ptAccCut[iP];
+ fEtaAccCut[iP]=etaAccCut[iP];
+ }
+ }
+ return;
+}
+
+//___________________________________________________________
+
+void AliCFVertexingHFLctoV0bachelor::SetAccCut()
+{
+ //
+ // setting the pt and eta cut to be used in the Acceptance steps (MC+Reco)
+ //
+
+ AliAODMCParticle* mcPartDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fLabelArray[0])); // bachelor
+ if (!mcPartDaughter) return;
+ Int_t mother = mcPartDaughter->GetMother();
+ AliAODMCParticle* mcMother = dynamic_cast<AliAODMCParticle*>(fmcArray->At(mother));
+ if (!mcMother) return;
+ /*
+ if (fGenLcOption==kCountK0Sp) {
+ if ( !(TMath::Abs(mcPartDaughter->GetPdgCode())== 2212) )
+ AliFatal(Form("Apparently the proton bachelor is not in the first position (%d <- %d), causing a crash!!",
+ mcPartDaughter->GetPdgCode(),mcMother->GetPdgCode()));
+ }
+ else if (fGenLcOption==kCountLambdapi) {
+ if ( !(TMath::Abs(mcPartDaughter->GetPdgCode())== 211) )
+ AliFatal(Form("Apparently the pion bachelor is not in the first position (%d <- %d), causing a crash!!",
+ mcPartDaughter->GetPdgCode(),mcMother->GetPdgCode()));
+ }
+ */
+ if (fProngs>0) {
+ fPtAccCut[0]=0.3; // bachelor
+ fEtaAccCut[0]=0.9; // bachelor
+ for (Int_t iP=1; iP<fProngs; iP++){
+ fPtAccCut[iP]=0.1;
+ fEtaAccCut[iP]=0.9;
+ }
+ }
+
+ return;
+
+}
+
+//_____________________________________________________________
+Double_t AliCFVertexingHFLctoV0bachelor::GetEtaProng(Int_t iProng) const
+{
+ //
+ // getting eta of the prong - overload the mother class method
+ //
+
+ if (fRecoCandidate) {
+
+ AliAODRecoCascadeHF* lcV0bachelor = (AliAODRecoCascadeHF*)fRecoCandidate;
+
+ Double_t etaProng =-9999;
+ if (iProng==0) etaProng = lcV0bachelor->GetBachelor()->Eta();
+ else if (iProng==1) etaProng = lcV0bachelor->Getv0PositiveTrack()->Eta();
+ else if (iProng==2) etaProng = lcV0bachelor->Getv0NegativeTrack()->Eta();
+
+ return etaProng;
+
+ }
+ return 999999;
+}
+
+//_____________________________________________________________
+
+Double_t AliCFVertexingHFLctoV0bachelor::GetPtProng(Int_t iProng) const
+{
+ //
+ // getting pt of the prong
+ //
+
+ Double_t ptProng=-9999.;
+
+ if (!fRecoCandidate) return ptProng;
+
+ AliAODRecoCascadeHF* lcV0bachelor = (AliAODRecoCascadeHF*)fRecoCandidate;
+
+ if (iProng==0) ptProng = lcV0bachelor->GetBachelor()->Pt();
+ else if (iProng==1) ptProng = lcV0bachelor->Getv0PositiveTrack()->Pt();
+ else if (iProng==2) ptProng = lcV0bachelor->Getv0NegativeTrack()->Pt();
+
+ return ptProng;
+
+}
+
+//_____________________________________________________________
+
+Double_t AliCFVertexingHFLctoV0bachelor::Ctau(AliAODMCParticle *mcPartCandidate)
+{
+
+ Double_t cTau = 999999.;
+
+ Double_t vtx1[3] = {0,0,0}; // primary vertex
+ Bool_t hasProdVertex = mcPartCandidate->XvYvZv(vtx1); // cm
+
+ AliAODMCParticle *mcPartDaughter0 = (AliAODMCParticle*)fmcArray->At(mcPartCandidate->GetDaughter(0));
+ AliAODMCParticle *mcPartDaughter1 = (AliAODMCParticle*)fmcArray->At(mcPartCandidate->GetDaughter(1));
+ if (!mcPartDaughter0 && !mcPartDaughter1) return cTau;
+ Double_t vtx1daughter[3] = {0,0,0}; // secondary vertex
+ if (mcPartDaughter0)
+ hasProdVertex = hasProdVertex || mcPartDaughter0->XvYvZv(vtx1daughter); //cm
+ if (mcPartDaughter1)
+ hasProdVertex = hasProdVertex || mcPartDaughter1->XvYvZv(vtx1daughter); //cm
+
+ if (!hasProdVertex) return cTau;
+
+ Double_t decayLength = 0.;
+ for (Int_t ii=0; ii<3; ii++) decayLength += (vtx1daughter[ii]-vtx1[ii])*(vtx1daughter[ii]-vtx1[ii]);
+ decayLength = TMath::Sqrt(decayLength);
+
+ cTau = decayLength * mcPartCandidate->M()/mcPartCandidate->P();
+ AliDebug(2,Form(" cTau(4122)=%f",cTau));
+ return cTau;
+
+}
--- /dev/null
+//DEFINITION OF A FEW CONSTANTS
+const Double_t ptmin = 0.0;
+const Double_t ptmax = 6.0;
+const Double_t ymin = -1.2 ;
+const Double_t ymax = 1.2 ;
+const Double_t cosPAV0min = 0.95;
+const Double_t cosPAV0max = 1.02;
+const Int_t onFlymin = 0;
+const Int_t onFlymax = 2;
+const Float_t centmin = 0.;
+//const Float_t centmin_0_10 = 0.;
+//const Float_t centmax_0_10 = 10.;
+//const Float_t centmin_10_60 = 10.;
+//const Float_t centmax_10_60 = 60.;
+//const Float_t centmin_60_100 = 60.;
+//const Float_t centmax_60_100 = 100.;
+const Float_t centmax = 100.;
+const Int_t fakemin = 0;
+const Int_t fakemax = 3;
+const Float_t multmin = 0;
+//const Float_t multmin_0_20 = 0;
+//const Float_t multmax_0_20 = 20;
+//const Float_t multmin_20_50 = 20;
+//const Float_t multmax_20_50 = 50;
+//const Float_t multmin_50_102 = 50;
+//const Float_t multmax_50_102 = 102;
+const Float_t multmax = 102;
+
+const Double_t ptBachmin = 0.0;
+const Double_t ptBachmax = 5.0;
+const Double_t ptV0posmin = 0.0;
+const Double_t ptV0posmax = 5.0;
+const Double_t ptV0negmin = 0.0;
+const Double_t ptV0negmax = 5.0;
+const Double_t dcaV0min = 0.0; // nSigma
+const Double_t dcaV0max = 1.5; // nSigma
+const Double_t cTV0min = 0.0; // micron
+const Double_t cTV0max = 300; // micron
+const Double_t cTmin = 0.0; // micron
+const Double_t cTmax = 300; // micron
+const Float_t cosPAmin =-1.02;
+const Float_t cosPAmax = 1.02;
+
+const Double_t etamin = -0.9;
+const Double_t etamax = 0.9;
+//const Double_t zmin = -15.;
+//const Double_t zmax = 15.;
+
+
+//----------------------------------------------------
+
+AliCFTaskVertexingHF *AddTaskCFVertexingHFLctoV0bachelor(const char* cutFile = "./LctoV0bachelorCuts.root",
+ Int_t configuration = AliCFTaskVertexingHF::kCheetah, Bool_t isKeepDfromB = kTRUE,
+ Bool_t isKeepDfromBOnly = kFALSE, Int_t pdgCode = 4122, Char_t isSign = 0, TString usercomment = "username")
+{
+
+
+ printf("Adding CF task using cuts from file %s\n",cutFile);
+ if (configuration == AliCFTaskVertexingHF::kSnail){
+ printf("The configuration is set to be SLOW --> all the variables will be used to fill the CF\n");
+ }
+ else if (configuration == AliCFTaskVertexingHF::kCheetah){
+ printf("The configuration is set to be FAST --> using only pt, y, ct, phi, zvtx, centrality, fake, multiplicity to fill the CF\n");
+ }
+ else{
+ printf("The configuration is not defined! returning\n");
+ return;
+ }
+
+ gSystem->Sleep(2000);
+
+ // isSign = 0 --> K0S only
+ // isSign = 1 --> Lambda only
+ // isSign = 2 --> LambdaBar only
+ // isSign = 3 --> K0S and Lambda and LambdaBar
+
+ TString expected;
+ if (isSign == 0 && pdgCode < 0){
+ AliError(Form("Error setting PDG code (%d) and sign (0 --> K0S only): they are not compatible, returning"));
+ return 0x0;
+ }
+ else if (isSign == 1 && pdgCode < 0){
+ AliError(Form("Error setting PDG code (%d) and sign (1 --> Lambda only): they are not compatible, returning"));
+ return 0x0;
+ }
+ else if (isSign == 2 && pdgCode > 0){
+ AliError(Form("Error setting PDG code (%d) and sign (2 --> LambdaBar only): they are not compatible, returning"));
+ return 0x0;
+ }
+ else if (isSign > 3 || isSign < 0){
+ AliError(Form("Sign not valid (%d, possible values are 0, 1, 2, 3), returning"));
+ return 0x0;
+ }
+
+ TFile* fileCuts = TFile::Open(cutFile);
+ AliRDHFCuts *cutsLctoV0 = (AliRDHFCutsLctoV0*)fileCuts->Get("LctoV0AnalysisCuts");
+
+ // check that the fKeepD0fromB flag is set to true when the fKeepD0fromBOnly flag is true
+ // for now the binning is the same than for all D's
+ if (isKeepDfromBOnly) isKeepDfromB = true;
+
+ Double_t massV0min = 0.47;
+ Double_t massV0max = 1.14;
+ if (isSign==0) {
+ massV0min = 0.47 ;
+ massV0max = 0.53 ;
+ } else if (isSign==1 || isSign==2) {
+ massV0min = 1.09;
+ massV0max = 1.14;
+ }
+
+ const Double_t phimin = 0.0;
+ const Double_t phimax = 2.*TMath::Pi();
+
+ const Int_t nbinpt = 6; //bins in pt from 0 to 6 GeV
+ const Int_t nbiny = 24; //bins in y
+ Int_t nbininvMassV0 = 60; //bins in invMassV0
+ if (isSign==3) nbininvMassV0=134; //bins in invMassV0
+ const Int_t nbinpointingV0 = 12; //bins in cosPointingAngleV0
+ const Int_t nbinonFly = 2; //bins in onFlyStatus x V0
+
+ const Int_t nbincent = 18; //bins in centrality (total number)
+ //const Int_t nbincent_0_10 = 4; //bins in centrality between 0 and 10
+ //const Int_t nbincent_10_60 = 10; //bins in centrality between 10 and 60
+ //const Int_t nbincent_60_100 = 4; //bins in centrality between 60 and 100
+
+ const Int_t nbinfake = 3; //bins in fake
+
+ const Int_t nbinmult = 48; //bins in multiplicity (total number)
+ //const Int_t nbinmult_0_20 = 20; //bins in multiplicity between 0 and 20
+ //const Int_t nbinmult_20_50 = 15; //bins in multiplicity between 20 and 50
+ //const Int_t nbinmult_50_102 = 13; //bins in multiplicity between 50 and 102
+
+
+ const Int_t nbinptBach = 50; //bins in pt from 0 to 5 GeV
+ const Int_t nbinptV0pos = 50; //bins in pt from 0 to 5 GeV
+ const Int_t nbinptV0neg = 50; //bins in pt from 0 to 5 GeV
+ const Int_t nbinphi = 18; //bins in Phi
+ const Int_t nbindcaV0 = 15; //bins in dcaV0
+ const Int_t nbincTV0 = 15; //bins in cTV0
+ const Int_t nbincT = 15; //bins in cT
+ const Int_t nbinpointing = 51; //bins in cosPointingAngle
+
+ //the sensitive variables, their indices
+
+ // variables' indices
+ const UInt_t ipT = 0;
+ const UInt_t iy = 1;
+ const UInt_t iphi = 2;
+ const UInt_t icosPAxV0 = 3;
+ const UInt_t ionFly = 4;
+ const UInt_t imult = 5;
+ const UInt_t ifake = 6;
+ const UInt_t icent = 7;
+
+ const UInt_t ipTbach = 8;
+ const UInt_t ipTposV0 = 9;
+ const UInt_t ipTnegV0 = 10;
+ const UInt_t iinvMassV0= 11;
+ const UInt_t idcaV0 = 12;
+ const UInt_t icTv0 = 13;
+ const UInt_t icT = 14;
+ const UInt_t icosPA = 15;
+
+ //Setting the bins: pt, ptPi, and ptK are considered seprately because for them you can either define the binning by hand, or using the cuts file
+
+ //arrays for the number of bins in each dimension
+
+ //if ( configuration ==AliCFTaskVertexingHF::kSnail)
+ const Int_t nvarTot = 16 ; //number of variables on the grid:pt, y, cosThetaStar, pTpi, pTk, cT, dca, d0pi, d0K, d0xd0, cosPointingAngle, phi, z, centrality, fake, cosPointingAngleXY, normDecayLengthXY, multiplicity
+ //if ( configuration ==AliCFTaskVertexingHF::kCheetah)
+ //const Int_t nvarTot = 8 ; //number of variables on the grid:pt, y, cosThetaStar, pTpi, pTk, cT, dca, d0pi, d0K, d0xd0, cosPointingAngle, phi, z, centrality, fake, cosPointingAngleXY, normDecayLengthXY, multiplicity
+
+ Int_t iBin[nvarTot];
+
+ //OPTION 1: defining the pt, ptPi, ptK bins by hand...
+ iBin[ipT]=nbinpt;
+ iBin[iy]=nbiny;
+ iBin[iphi]=nbinphi;
+ iBin[icosPAxV0]=nbinpointingV0;
+ iBin[ionFly]=nbinonFly;
+ iBin[imult]=nbinmult;
+ iBin[ifake]=nbinfake;
+ iBin[icent]=nbincent;
+
+ iBin[ipTbach]=nbinptBach;
+ iBin[ipTposV0]=nbinptV0pos;
+ iBin[ipTnegV0]=nbinptV0neg;
+ iBin[iinvMassV0]=nbininvMassV0;
+ iBin[idcaV0]=nbindcaV0;
+ iBin[icTv0]=nbincTV0;
+ iBin[icT]=nbincT;
+ iBin[icosPA]=nbinpointing;
+
+ // values for bin lower bounds
+
+ // pt
+ Double_t *binLimpT=new Double_t[iBin[0]+1];
+ for(Int_t i=0; i<=iBin[0]; i++) binLimpT[i]=(Double_t)ptmin + (ptmax-ptmin)/iBin[0]*(Double_t)i ;
+
+ // y
+ Double_t *binLimy=new Double_t[iBin[1]+1];
+ for(Int_t i=0; i<=iBin[1]; i++) binLimy[i]=(Double_t)ymin + (ymax-ymin)/iBin[1]*(Double_t)i ;
+
+ // phi
+ Double_t *binLimphi=new Double_t[iBin[2]+1];
+ for(Int_t i=0; i<=iBin[2]; i++) binLimphi[i]=(Double_t)phimin + (phimax-phimin)/iBin[2]*(Double_t)i ;
+
+ // cosPointingAngleV0
+ Double_t *binLimcosPAV0=new Double_t[iBin[3]+1];
+ for(Int_t i=0; i<=iBin[3]; i++) binLimcosPAV0[i]=(Double_t)cosPAV0min + (cosPAV0max-cosPAV0min)/iBin[3]*(Double_t)i ;
+
+ // onTheFlyV0
+ Double_t *binLimonFlyV0=new Double_t[iBin[4]+1];
+ for(Int_t i=0; i<=iBin[4]; i++) binLimonFlyV0[i]=(Double_t)onFlymin + (onFlymax-onFlymin)/iBin[4]*(Double_t)i ;
+
+ // centrality
+ Double_t *binLimcent=new Double_t[iBin[5]+1];
+ for(Int_t i=0; i<=iBin[5]; i++) binLimcent[i]=(Double_t)centmin + (centmax-centmin)/iBin[5]*(Double_t)i ;
+ /*
+ Double_t *binLimcent_0_10=new Double_t[iBin[icent]+1];
+ for(Int_t i=0; i<=nbincent_0_10; i++) binLimcent[i]=(Double_t)centmin_0_10 + (centmax_0_10-centmin_0_10)/nbincent_0_10*(Double_t)i ;
+ if (binLimcent[nbincent_0_10] != centmin_10_60) {
+ Error("AliCFHeavyFlavourTaskMultiVarMultiStep","Calculated bin lim for cent - 1st range - differs from expected!\n");
+ }
+
+ Double_t *binLimcent_10_60=new Double_t[iBin[icent]+1];
+ for(Int_t i=0; i<=nbincent_10_60; i++) binLimcent[i+nbincent_0_10]=(Double_t)centmin_10_60 + (centmax_10_60-centmin_10_60)/nbincent_10_60*(Double_t)i ;
+ if (binLimcent[nbincent_0_10+nbincent_10_60] != centmin_60_100) {
+ Error("AliCFHeavyFlavourTaskMultiVarMultiStep","Calculated bin lim for cent - 2st range - differs from expected!\n");
+ }
+
+ Double_t *binLimcent_60_100=new Double_t[iBin[icent]+1];
+ for(Int_t i=0; i<=nbincent_60_100; i++) binLimcent[i+nbincent_10_60]=(Double_t)centmin_60_100 + (centmax_60_100-centmin_60_100)/nbincent_60_100*(Double_t)i ;
+ */
+
+ // fake
+ Double_t *binLimfake=new Double_t[iBin[6]+1];
+ for(Int_t i=0; i<=iBin[6]; i++) binLimfake[i]=(Double_t)fakemin + (fakemax-fakemin)/iBin[6] * (Double_t)i;
+
+ // multiplicity
+ Double_t *binLimmult=new Double_t[iBin[7]+1];
+ for(Int_t i=0; i<=iBin[7]; i++) binLimmult[i]=(Double_t)multmin + (multmax-multmin)/iBin[7]*(Double_t)i ;
+ /*
+ for(Int_t i=0; i<=nbinmult_0_20; i++) binLimmult[i]=(Double_t)multmin_0_20 + (multmax_0_20-multmin_0_20)/nbinmult_0_20*(Double_t)i ;
+ if (binLimmult[nbinmult_0_20] != multmin_20_50) {
+ Error("AliCFHeavyFlavourTaskMultiVarMultiStep","Calculated bin lim for mult - 1st range - differs from expected!\n");
+ }
+ for(Int_t i=0; i<=nbinmult_20_50; i++) binLimmult[i+nbinmult_0_20]=(Double_t)multmin_20_50 + (multmax_20_50-multmin_20_50)/nbinmult_20_50*(Double_t)i ;
+ if (binLimmult[nbinmult_0_20+nbinmult_20_50] != multmin_50_102) {
+ Error("AliCFHeavyFlavourTaskMultiVarMultiStep","Calculated bin lim for mult - 2nd range - differs from expected!\n");
+ }
+ for(Int_t i=0; i<=nbinmult_50_102; i++) binLimmult[i+nbinmult_0_20+nbinmult_20_50]=(Double_t)multmin_50_102 + (multmax_50_102-multmin_50_102)/nbinmult_50_102*(Double_t)i ;
+ */
+
+
+ // ptBach
+ Double_t *binLimpTbach=new Double_t[iBin[8]+1];
+ for(Int_t i=0; i<=iBin[8]; i++) binLimpTbach[i]=(Double_t)ptBachmin + (ptBachmax-ptBachmin)/iBin[8]*(Double_t)i ;
+
+ // ptV0pos
+ Double_t *binLimpTV0pos=new Double_t[iBin[9]+1];
+ for(Int_t i=0; i<=iBin[9]; i++) binLimpTV0pos[i]=(Double_t)ptV0posmin + (ptV0posmax-ptV0posmin)/iBin[9]*(Double_t)i ;
+
+ // ptV0neg
+ Double_t *binLimpTV0neg=new Double_t[iBin[10]+1];
+ for(Int_t i=0; i<=iBin[10]; i++) binLimpTV0neg[i]=(Double_t)ptV0negmin + (ptV0negmax-ptV0negmin)/iBin[10]*(Double_t)i ;
+
+ // invMassV0
+ Double_t *binLimInvMassV0=new Double_t[iBin[11]+1];
+ for(Int_t i=0; i<=iBin[11]; i++) binLimInvMassV0[i]=(Double_t)massV0min + (massV0max-massV0min)/iBin[11]*(Double_t)i ;
+
+ // dcaV0
+ Double_t *binLimdcaV0=new Double_t[iBin[12]+1];
+ for(Int_t i=0; i<=iBin[12]; i++) binLimdcaV0[i]=(Double_t)dcaV0min + (dcaV0max-dcaV0min)/iBin[12]*(Double_t)i ;
+
+ // cTV0
+ Double_t *binLimcTV0=new Double_t[iBin[13]+1];
+ for(Int_t i=0; i<=iBin[13]; i++) binLimcTV0[i]=(Double_t)cTV0min + (cTV0max-cTV0min)/iBin[13]*(Double_t)i ;
+
+ // cT
+ Double_t *binLimcT=new Double_t[iBin[14]+1];
+ for(Int_t i=0; i<=iBin[14]; i++) binLimcT[i]=(Double_t)cTmin + (cTmax-cTmin)/iBin[14]*(Double_t)i ;
+
+ // cosPointingAngle
+ Double_t *binLimcosPA=new Double_t[iBin[15]+1];
+ for(Int_t i=0; i<=iBin[15]; i++) binLimcosPA[i]=(Double_t)cosPAmin + (cosPAmax-cosPAmin)/iBin[15]*(Double_t)i ;
+
+
+
+ // z Primary Vertex
+ //for(Int_t i=0; i<=nbinzvtx; i++) binLimzvtx[i]=(Double_t)zmin + (zmax-zmin) /nbinzvtx*(Double_t)i ;
+
+
+ //OPTION 2: ...or from the cuts file
+ /*
+ const Int_t nbinpt = cutsLctoV0->GetNPtBins(); // bins in pT
+ iBin[ipT]=nbinpt;
+ iBin[ipTpi]=nbinpt;
+ iBin[ipTk]=nbinpt;
+ Double_t *binLimpT=new Double_t[iBin[ipT]+1];
+ Double_t *binLimpTpi=new Double_t[iBin[ipTpi]+1];
+ Double_t *binLimpTk=new Double_t[iBin[ipTk]+1];
+ // values for bin lower bounds
+ Float_t* floatbinLimpT = cutsLctoV0->GetPtBinLimits();
+ for (Int_t ibin0 = 0 ; ibin0<iBin[ipT]+1; ibin0++){
+ binLimpT[ibin0] = (Double_t)floatbinLimpT[ibin0];
+ binLimpTpi[ibin0] = (Double_t)floatbinLimpT[ibin0];
+ binLimpTk[ibin0] = (Double_t)floatbinLimpT[ibin0];
+ }
+ for(Int_t i=0; i<=nbinpt; i++) printf("binLimpT[%d]=%f\n",i,binLimpT[i]);
+
+ printf("pT: nbin (from cuts file) = %d\n",nbinpt);
+
+ // defining now the binning for the other variables:
+
+ iBin[iy]=nbiny;
+ iBin[icosThetaStar]=nbincosThetaStar;
+ iBin[icT]=nbincT;
+ iBin[idca]=nbindca;
+ iBin[id0xd0]=nbind0xd0;
+ iBin[ipointing]=nbinpointing;
+ iBin[iphi]=nbinphi;
+ iBin[izvtx]=nbinzvtx;
+ iBin[icent]=nbincent;
+ iBin[ifake]=nbinfake;
+ iBin[ipointingXY]=nbinpointingXY;
+ iBin[inormDecayLXY]=nbinnormDecayLXY;
+ iBin[imult]=nbinmult;
+
+ //arrays for lower bounds :
+ Double_t *binLimy=new Double_t[iBin[iy]+1];
+ Double_t *binLimcosThetaStar=new Double_t[iBin[icosThetaStar]+1];
+ Double_t *binLimcT=new Double_t[iBin[icT]+1];
+ Double_t *binLimdca=new Double_t[iBin[idca]+1];
+ Double_t *binLimd0xd0=new Double_t[iBin[id0xd0]+1];
+ Double_t *binLimpointing=new Double_t[iBin[ipointing]+1];
+ Double_t *binLimphi=new Double_t[iBin[iphi]+1];
+ Double_t *binLimzvtx=new Double_t[iBin[izvtx]+1];
+ Double_t *binLimcent=new Double_t[iBin[icent]+1];
+ Double_t *binLimfake=new Double_t[iBin[ifake]+1];
+ Double_t *binLimpointingXY=new Double_t[iBin[ipointingXY]+1];
+ Double_t *binLimnormDecayLXY=new Double_t[iBin[inormDecayLXY]+1];
+ Double_t *binLimmult=new Double_t[iBin[imult]+1];
+
+ // y
+ for(Int_t i=0; i<=nbiny; i++) binLimy[i]=(Double_t)ymin + (ymax-ymin) /nbiny*(Double_t)i ;
+
+ // cosThetaStar
+ for(Int_t i=0; i<=nbincosThetaStar; i++) binLimcosThetaStar[i]=(Double_t)cosminTS + (cosmaxTS-cosminTS) /nbincosThetaStar*(Double_t)i ;
+
+ // cT
+ for(Int_t i=0; i<=nbincT; i++) binLimcT[i]=(Double_t)cTmin + (cTmax-cTmin) /nbincT*(Double_t)i ;
+
+ // dca
+ for(Int_t i=0; i<=nbindca; i++) binLimdca[i]=(Double_t)dcamin + (dcamax-dcamin) /nbindca*(Double_t)i ;
+
+ // d0xd0
+ for(Int_t i=0; i<=nbind0xd0; i++) binLimd0xd0[i]=(Double_t)d0xd0min + (d0xd0max-d0xd0min) /nbind0xd0*(Double_t)i ;
+
+ // cosPointingAngle
+ for(Int_t i=0; i<=nbinpointing; i++) binLimpointing[i]=(Double_t)cosmin + (cosmax-cosmin) /nbinpointing*(Double_t)i ;
+
+ // Phi
+ for(Int_t i=0; i<=nbinphi; i++) binLimphi[i]=(Double_t)phimin + (phimax-phimin) /nbinphi*(Double_t)i ;
+
+ // z Primary Vertex
+ for(Int_t i=0; i<=nbinzvtx; i++) {
+ binLimzvtx[i]=(Double_t)zmin + (zmax-zmin) /nbinzvtx*(Double_t)i ;
+ }
+
+ // centrality
+ for(Int_t i=0; i<=nbincent_0_10; i++) binLimcent[i]=(Double_t)centmin_0_10 + (centmax_0_10-centmin_0_10)/nbincent_0_10*(Double_t)i ;
+ if (binLimcent[nbincent_0_10] != centmin_10_60) {
+ Error("AliCFHeavyFlavourTaskMultiVarMultiStep","Calculated bin lim for cent - 1st range - differs from expected!\n");
+ }
+ for(Int_t i=0; i<=nbincent_10_60; i++) binLimcent[i+nbincent_0_10]=(Double_t)centmin_10_60 + (centmax_10_60-centmin_10_60)/nbincent_10_60*(Double_t)i ;
+ if (binLimcent[nbincent_0_10+nbincent_10_60] != centmin_60_100) {
+ Error("AliCFHeavyFlavourTaskMultiVarMultiStep","Calculated bin lim for cent - 2st range - differs from expected!\n");
+ }
+ for(Int_t i=0; i<=nbincent_60_100; i++) binLimcent[i+nbincent_10_60]=(Double_t)centmin_60_100 + (centmax_60_100-centmin_60_100)/nbincent_60_100*(Double_t)i ;
+
+ // fake
+ for(Int_t i=0; i<=nbinfake; i++) {
+ binLimfake[i]=(Double_t)fakemin + (fakemax-fakemin)/nbinfake * (Double_t)i;
+ }
+
+ // multiplicity
+ for(Int_t i=0; i<=nbinmult_0_20; i++) binLimmult[i]=(Double_t)multmin_0_20 + (multmax_0_20-multmin_0_20)/nbinmult_0_20*(Double_t)i ;
+ if (binLimmult[nbinmult_0_20] != multmin_20_50) {
+ Error("AliCFHeavyFlavourTaskMultiVarMultiStep","Calculated bin lim for mult - 1st range - differs from expected!\n");
+ }
+ for(Int_t i=0; i<=nbinmult_20_50; i++) binLimmult[i+nbinmult_0_20]=(Double_t)multmin_20_50 + (multmax_20_50-multmin_20_50)/nbinmult_20_50*(Double_t)i ;
+ if (binLimmult[nbinmult_0_20+nbinmult_20_50] != multmin_50_102) {
+ Error("AliCFHeavyFlavourTaskMultiVarMultiStep","Calculated bin lim for mult - 2nd range - differs from expected!\n");
+ }
+ for(Int_t i=0; i<=nbinmult_50_102; i++) binLimmult[i+nbinmult_0_20+nbinmult_20_50]=(Double_t)multmin_50_102 + (multmax_50_102-multmin_50_102)/nbinmult_50_102*(Double_t)i ;
+ */
+
+
+
+ //one "container" for MC
+ TString nameContainer="";
+ if (!isKeepDfromB) {
+ nameContainer="CFHFccontainer0_CommonFramework_"+usercomment;
+ }
+ else if (isKeepDfromBOnly) {
+ nameContainer="CFHFccontainer0LcfromB_CommonFramework_"+usercomment;
+ }
+ else {
+ nameContainer="CFHFccontainer0allLc_CommonFramework_"+usercomment;
+ }
+
+ //Setting up the container grid...
+
+ //CONTAINER DEFINITION
+ Info("AliCFTaskVertexingHF","SETUP CONTAINER");
+ UInt_t nstep = 10; //number of selection steps: MC with limited acceptance, MC, Acceptance, Vertex, Refit, Reco (no cuts), RecoAcceptance, RecoITSClusters (RecoAcceptance included), RecoPPR (RecoAcceptance+RecoITSCluster included), RecoPID
+
+ AliCFContainer* container;
+ if (configuration == AliCFTaskVertexingHF::kSnail) {
+ container = new AliCFContainer(nameContainer,"container for tracks",nstep,nvarTot,iBin);
+ }
+ else if (configuration == AliCFTaskVertexingHF::kCheetah) {
+ container = new AliCFContainer(nameContainer,"container for tracks",nstep,8,iBin);
+ }
+
+ //setting the bin limits
+ container -> SetBinLimits(0,binLimpT);
+ container -> SetBinLimits(1,binLimy);
+ container -> SetBinLimits(2,binLimphi);
+ container -> SetBinLimits(3,binLimcosPAV0);
+ container -> SetBinLimits(4,binLimonFlyV0);
+ container -> SetBinLimits(5,binLimcent);
+ container -> SetBinLimits(6,binLimfake);
+ container -> SetBinLimits(7,binLimmult);
+
+ container -> SetVarTitle(0,"pt");
+ container -> SetVarTitle(1,"y");
+ container -> SetVarTitle(2,"phi");
+ container -> SetVarTitle(3,"cosPA -V0-");
+ container -> SetVarTitle(4,"onFlyV0");
+ container -> SetVarTitle(5,"centrality");
+ container -> SetVarTitle(6,"fake");
+ container -> SetVarTitle(7,"multiplicity");
+
+ if (configuration == AliCFTaskVertexingHF::kSnail) {
+ container -> SetBinLimits(8,binLimpTbach);
+ container -> SetBinLimits(9,binLimpTV0pos);
+ container -> SetBinLimits(10,binLimpTV0neg);
+ container -> SetBinLimits(11,binLimInvMassV0);
+ container -> SetBinLimits(12,binLimdcaV0);
+ container -> SetBinLimits(13,binLimcTV0);
+ container -> SetBinLimits(14,binLimcT);
+ container -> SetBinLimits(15,binLimcosPA);
+
+ container -> SetVarTitle(8,"ptBachelor");
+ container -> SetVarTitle(9,"ptV0pos");
+ container -> SetVarTitle(10,"ptV0neg");
+ container -> SetVarTitle(11,"mV0");
+ container -> SetVarTitle(12,"DCA -V0-");
+ container -> SetVarTitle(13,"c#tau -V0-");
+ container -> SetVarTitle(14,"c#tau");
+ container -> SetVarTitle(15,"cosPA");
+ }
+
+ container -> SetStepTitle(0, "MCLimAcc");
+ container -> SetStepTitle(1, "MC");
+ container -> SetStepTitle(2, "MCAcc");
+ container -> SetStepTitle(3, "RecoVertex");
+ container -> SetStepTitle(4, "RecoRefit");
+ container -> SetStepTitle(5, "Reco");
+ container -> SetStepTitle(6, "RecoAcc");
+ container -> SetStepTitle(7, "RecoITSCluster");
+ container -> SetStepTitle(8, "RecoCuts");
+ container -> SetStepTitle(9, "RecoPID");
+
+ //return container;
+
+ //CREATE THE CUTS -----------------------------------------------
+
+ // Gen-Level kinematic cuts
+ AliCFTrackKineCuts *mcKineCuts = new AliCFTrackKineCuts("mcKineCuts","MC-level kinematic cuts");
+
+ //Particle-Level cuts:
+ AliCFParticleGenCuts* mcGenCuts = new AliCFParticleGenCuts("mcGenCuts","MC particle generation cuts");
+ Bool_t useAbsolute = kTRUE;
+ if (isSign != 3 && isSign!=0) {
+ useAbsolute = kFALSE;
+ }
+ mcGenCuts->SetRequirePdgCode(pdgCode, useAbsolute); // kTRUE set in order to include Lc-
+ mcGenCuts->SetAODMC(1); //special flag for reading MC in AOD tree (important)
+
+ // Acceptance cuts:
+ AliCFAcceptanceCuts* accCuts = new AliCFAcceptanceCuts("accCuts", "Acceptance cuts");
+ AliCFTrackKineCuts * kineAccCuts = new AliCFTrackKineCuts("kineAccCuts","Kine-Acceptance cuts");
+ kineAccCuts->SetPtRange(ptmin,ptmax);
+ kineAccCuts->SetEtaRange(etamin,etamax);
+
+ // Rec-Level kinematic cuts
+ AliCFTrackKineCuts *recKineCuts = new AliCFTrackKineCuts("recKineCuts","rec-level kine cuts");
+
+ AliCFTrackQualityCuts *recQualityCuts = new AliCFTrackQualityCuts("recQualityCuts","rec-level quality cuts");
+
+ AliCFTrackIsPrimaryCuts *recIsPrimaryCuts = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts","rec-level isPrimary cuts");
+
+ printf("CREATE MC KINE CUTS\n");
+ TObjArray* mcList = new TObjArray(0) ;
+ mcList->AddLast(mcKineCuts);
+ mcList->AddLast(mcGenCuts);
+
+ printf("CREATE ACCEPTANCE CUTS\n");
+ TObjArray* accList = new TObjArray(0) ;
+ accList->AddLast(kineAccCuts);
+
+ printf("CREATE RECONSTRUCTION CUTS\n");
+ TObjArray* recList = new TObjArray(0) ; // not used!!
+ recList->AddLast(recKineCuts);
+ recList->AddLast(recQualityCuts);
+ recList->AddLast(recIsPrimaryCuts);
+
+ TObjArray* emptyList = new TObjArray(0);
+
+ //CREATE THE INTERFACE TO CORRECTION FRAMEWORK USED IN THE TASK
+ printf("CREATE INTERFACE AND CUTS\n");
+ AliCFManager* man = new AliCFManager() ;
+ man->SetParticleContainer(container);
+ man->SetParticleCutsList(0 , mcList); // MC, Limited Acceptance
+ man->SetParticleCutsList(1 , mcList); // MC
+ man->SetParticleCutsList(2 , accList); // Acceptance
+ man->SetParticleCutsList(3 , emptyList); // Vertex
+ man->SetParticleCutsList(4 , emptyList); // Refit
+ man->SetParticleCutsList(5 , emptyList); // AOD
+ man->SetParticleCutsList(6 , emptyList); // AOD in Acceptance
+ man->SetParticleCutsList(7 , emptyList); // AOD with required n. of ITS clusters
+ man->SetParticleCutsList(8 , emptyList); // AOD Reco (PPR cuts implemented in Task)
+ man->SetParticleCutsList(9 , emptyList); // AOD Reco PID
+
+ // Get the pointer to the existing analysis manager via the static access method.
+ //==============================================================================
+ AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+ if (!mgr) {
+ ::Error("AddTaskCompareHF", "No analysis manager to connect to.");
+ return NULL;
+ }
+ //CREATE THE TASK
+ printf("CREATE TASK\n");
+
+ // create the task
+ AliCFTaskVertexingHF *task = new AliCFTaskVertexingHF("AliCFTaskVertexingHF",cutsLctoV0);
+ task->SetConfiguration(configuration);
+ task->SetFillFromGenerated(kFALSE);
+ task->SetCFManager(man); //here is set the CF manager
+ task->SetDecayChannel(22);//kLctoV0bachelor
+ task->SetUseWeight(kFALSE);
+ task->SetSign(isSign);
+ task->SetCentralitySelection(kFALSE);
+ task->SetFakeSelection(0);
+ task->SetRejectCandidateIfNotFromQuark(kTRUE); // put to false if you want to keep HIJING D0!!
+ task->SetUseMCVertex(kFALSE); // put to true if you want to do studies on pp
+
+ if (isKeepDfromB && !isKeepDfromBOnly) task->SetDselection(2);
+ if (isKeepDfromB && isKeepDfromBOnly) task->SetDselection(1);
+
+ TF1* funcWeight = 0x0;
+ if (task->GetUseWeight()) {
+ funcWeight = (TF1*)fileCuts->Get("funcWeight");
+ if (funcWeight == 0x0){
+ Printf("FONLL Weights will be used");
+ }
+ else {
+ task->SetWeightFunction(funcWeight);
+ Printf("User-defined Weights will be used. The function being:");
+ task->GetWeightFunction(funcWeight)->Print();
+ }
+ }
+
+ Printf("***************** CONTAINER SETTINGS *****************");
+ Printf("decay channel = %d",(Int_t)task->GetDecayChannel());
+ Printf("FillFromGenerated = %d",(Int_t)task->GetFillFromGenerated());
+ Printf("Dselection = %d",(Int_t)task->GetDselection());
+ Printf("UseWeight = %d",(Int_t)task->GetUseWeight());
+ if (task->GetUseWeight()) {
+ Printf("User-defined Weight function:");
+ task->GetWeightFunction(funcWeight)->Print();
+ }
+ else {
+ Printf("FONLL will be used for the weights");
+ }
+ Printf("Sign = %d",(Int_t)task->GetSign());
+ Printf("Centrality selection = %d",(Int_t)task->GetCentralitySelection());
+ Printf("Fake selection = %d",(Int_t)task->GetFakeSelection());
+ Printf("RejectCandidateIfNotFromQuark selection = %d",(Int_t)task->GetRejectCandidateIfNotFromQuark());
+ Printf("UseMCVertex selection = %d",(Int_t)task->GetUseMCVertex());
+ Printf("***************END CONTAINER SETTINGS *****************\n");
+
+ //-----------------------------------------------------------//
+ // create correlation matrix for unfolding - only eta-pt //
+ //-----------------------------------------------------------//
+
+ Bool_t AcceptanceUnf = kTRUE; // unfold at acceptance level, otherwise PPR
+
+ Int_t thnDim[4];
+
+ //first half : reconstructed
+ //second half : MC
+
+ thnDim[0] = iBin[0];
+ thnDim[2] = iBin[0];
+ thnDim[1] = iBin[1];
+ thnDim[3] = iBin[1];
+
+ TString nameCorr="";
+ if (!isKeepDfromB) {
+ nameCorr="CFHFcorr0_CommonFramework_"+usercomment;
+ }
+ else if (isKeepDfromBOnly) {
+ nameCorr= "CFHFcorr0KeepDfromBOnly_CommonFramework_"+usercomment;
+ }
+ else {
+ nameCorr="CFHFcorr0allLc_CommonFramework_"+usercomment;
+ }
+
+ THnSparseD* correlation = new THnSparseD(nameCorr,"THnSparse with correlations",4,thnDim);
+ Double_t** binEdges = new Double_t[2];
+
+ // set bin limits
+
+ binEdges[0]= binLimpT;
+ binEdges[1]= binLimy;
+
+ correlation->SetBinEdges(0,binEdges[0]);
+ correlation->SetBinEdges(2,binEdges[0]);
+
+ correlation->SetBinEdges(1,binEdges[1]);
+ correlation->SetBinEdges(3,binEdges[1]);
+
+ correlation->Sumw2();
+
+ // correlation matrix ready
+ //------------------------------------------------//
+
+ task->SetCorrelationMatrix(correlation); // correlation matrix for unfolding
+
+ // Create and connect containers for input/output
+
+ // ------ input data ------
+ AliAnalysisDataContainer *cinput0 = mgr->GetCommonInputContainer();
+
+ // ----- output data -----
+
+ TString outputfile = AliAnalysisManager::GetCommonFileName();
+ TString output1name="", output2name="", output3name="",output4name="";
+ output2name=nameContainer;
+ output3name=nameCorr;
+ if (!isKeepDfromB) {
+ outputfile += ":PWG3_D2H_CFtaskLctoK0Sp_CommonFramework_"+usercomment;
+ output1name="CFHFchist0_CommonFramework_"+usercomment;
+ output4name= "Cuts_CommonFramework_"+usercomment;
+ }
+ else if (isKeepDfromBOnly) {
+ outputfile += ":PWG3_D2H_CFtaskLctoK0SpKeepDfromBOnly_CommonFramework_"+usercomment;
+ output1name="CFHFchist0DfromB_CommonFramework_"+usercomment;
+ output4name= "Cuts_CommonFramework_DfromB_"+usercomment;
+ }
+ else {
+ outputfile += ":PWG3_D2H_CFtaskLctoK0SpKeepDfromB_CommonFramework_"+usercomment;
+ output1name="CFHFchist0allLc_CommonFramework_"+usercomment;
+ output4name= "Cuts_CommonFramework_allLc_"+usercomment;
+ }
+
+ //now comes user's output objects :
+ // output TH1I for event counting
+ AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(output1name, TH1I::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data());
+ // output Correction Framework Container (for acceptance & efficiency calculations)
+ AliAnalysisDataContainer *coutput2 = mgr->CreateContainer(output2name, AliCFContainer::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data());
+ // Unfolding - correlation matrix
+ AliAnalysisDataContainer *coutput3 = mgr->CreateContainer(output3name, THnSparseD::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data());
+ // cuts
+ AliAnalysisDataContainer *coutput4 = mgr->CreateContainer(output4name, AliRDHFCuts::Class(),AliAnalysisManager::kOutputContainer, outputfile.Data());
+
+ mgr->AddTask(task);
+
+ mgr->ConnectInput(task,0,mgr->GetCommonInputContainer());
+ mgr->ConnectOutput(task,1,coutput1);
+ mgr->ConnectOutput(task,2,coutput2);
+ mgr->ConnectOutput(task,3,coutput3);
+ mgr->ConnectOutput(task,4,coutput4);
+ return task;
+
+}