04fe7352976c8f23190c370764d038950419fe72
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / macros / AddTaskSignificance.C
1 AliAnalysisTaskSESignificance *AddTaskSignificance(TString filename="cuts4SignifMaximDplus.root",Int_t decCh=0,Bool_t readMC=kFALSE,Int_t flagOPartAntiPart=0,Int_t nofsteps=8,AliAnalysisTaskSESignificance::FeedDownEnum fromcb=AliAnalysisTaskSESignificance::kBoth)
2 {
3   //
4   // Test macro for the AliAnalysisTaskSE for D meson candidates
5   // Invariant mass histogram and
6   // association with MC truth (using MC info in AOD)
7   //  R. Bala, bala@to.infn.it
8   // C. Bianchin, cbianchi@pd.infn.it
9   // Get the pointer to the existing analysis manager via the static access method.
10   //============================================================================
11   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
12   if (!mgr) {
13     ::Error("AddTaskSignificance", "No analysis manager to connect to.");
14     return NULL;
15   }
16
17   TFile* filecuts=new TFile(filename.Data());
18   if(!filecuts->IsOpen()){
19     cout<<"Input file not found: exit"<<endl;
20     return;
21   }
22   
23   AliRDHFCuts *analysiscuts=0x0;
24   TString suffix="";
25   TString suffix2="";
26   if(flagOPartAntiPart==1) suffix2="P"; //P=particle, A=antiparticle
27   if(flagOPartAntiPart==-1) suffix2="A";
28   if(fromcb==AliAnalysisTaskSESignificance::kCharmOnly) suffix2+="prompt";
29   if(fromcb==AliAnalysisTaskSESignificance::kBeautyOnly) suffix2+="feeddown";
30
31   TString cutsobjname="loosercuts";
32   //Analysis cuts
33   switch (decCh){
34   case 0:
35     analysiscuts = (AliRDHFCutsDplustoKpipi*)filecuts->Get(cutsobjname);
36     suffix=Form("Dplus%s",suffix2.Data());
37     break;
38   case 1:
39     analysiscuts = (AliRDHFCutsD0toKpi*)filecuts->Get(cutsobjname);
40     suffix=Form("D0%s",suffix2.Data());
41     break;
42   case 2:
43     analysiscuts = (AliRDHFCutsDstartoKpipi*)filecuts->Get(cutsobjname);
44     suffix=Form("Dstar%s",suffix2.Data());
45     break;
46   case 3:
47     analysiscuts = (AliRDHFCutsDstoKKpi*)filecuts->Get(cutsobjname);
48     suffix=Form("Ds%s",suffix2.Data());
49     break;
50   case 4:
51     analysiscuts = (AliRDHFCutsD0toKpipipi*)filecuts->Get(cutsobjname);
52     suffix=Form("D04%s",suffix2.Data());
53     break;
54   case 5:
55     analysiscuts = (AliRDHFCutsLctopKpi*)filecuts->Get(cutsobjname);
56     suffix=Form("Lc%s",suffix2.Data());
57     break;
58   }
59   //ptbins
60
61   if(!analysiscuts){
62     cout<<"Specific AliRDHFCuts not found"<<endl;
63     return;
64   }
65   
66   const Int_t nptbins=analysiscuts->GetNPtBins();
67   Float_t* ptbins=analysiscuts->GetPtBinLimits();
68
69   // Analysis task    
70   const Int_t npars=analysiscuts->GetNVarsForOpt();//numbers of var for opt
71   Bool_t* varsforopt=analysiscuts->GetVarsForOpt();
72   cout<<"pt bins= "<<nptbins<<"  varsforopt= "<<npars<<"  nvars= "<<analysiscuts->GetNVars()<<endl;
73   Int_t* nofcells=new Int_t[npars];//={4,4,4,4};
74   Float_t** looses;
75   looses=new Float_t*[npars];
76   Float_t** tights;
77   tights=new Float_t*[npars];
78   TString axisTitle[npars];
79   TString parname="";
80
81   for (Int_t ivop=0;ivop<npars;ivop++){
82     looses[ivop]=new Float_t[nptbins];
83     tights[ivop]=new Float_t[nptbins];
84     nofcells[ivop]=nofsteps;
85   }
86
87   Int_t count=0;
88   
89   for (Int_t ip=0;ip<nptbins;ip++){
90     for(Int_t iv=0;iv<analysiscuts->GetNVars();iv++){
91       if(varsforopt[iv]){
92         looses[count][ip]=analysiscuts->GetCutValue(iv,ip);
93         axisTitle[count]=(analysiscuts->GetVarNames())[iv];
94         parname=Form("par%dptbin%d",count,ip);
95         TParameter<float>* par=(TParameter<float>*)filecuts->Get(parname.Data());
96         tights[count][ip]=par->GetVal();
97         count++;
98       }
99     }
100     count=0;
101   }
102   
103   //creation TList of AliMultiDimVector
104   TList *listMDV=new TList();
105   listMDV->SetOwner();
106   listMDV->SetName("listMDV");
107
108   for (Int_t ip=0;ip<nptbins;ip++){
109     Float_t *loosescut=new Float_t[npars];
110     Float_t *tightscut=new Float_t[npars];
111     for(Int_t i=0;i<npars;i++){
112       loosescut[i]=looses[i][ip];
113       tightscut[i]=tights[i][ip];
114     }
115     Float_t ptbincut[2]={ptbins[ip],ptbins[ip+1]};
116     TString mdvname=Form("multiDimVectorPtBin%d",ip);
117     AliMultiDimVector *mv=new AliMultiDimVector(mdvname.Data(),"MultiDimVector",1,ptbincut,npars,nofcells,loosescut,tightscut,axisTitle);
118     listMDV->Add(mv);
119   }
120   
121   AliAnalysisTaskSESignificance *sigTask = new AliAnalysisTaskSESignificance("SignificanceAnalysis",listMDV,analysiscuts,decCh,AliRDHFCuts::kAll);//AliRDHFCuts::kCandidate
122   sigTask->SetReadMC(readMC);
123   //sigTask->SetDoLikeSign(kTRUE);
124   sigTask->SetBFeedDown(fromcb);
125   sigTask->SetDebugLevel(3);
126   sigTask->SetFillWithPartAntiPartBoth(flagOPartAntiPart);
127   mgr->AddTask(sigTask);
128
129   TString contname=Form("cinputSig%s",suffix.Data());
130   // Create containers for input/output
131   AliAnalysisDataContainer *cinputSig = mgr->CreateContainer(contname.Data(),TChain::Class(),AliAnalysisManager::kInputContainer);
132   TString outputfile = AliAnalysisManager::GetCommonFileName();
133   TString outputhistos = outputfile += ":PWG3_D2H_Significance"; 
134   contname=Form("coutputSig%s",suffix.Data());
135   AliAnalysisDataContainer *coutputSig = mgr->CreateContainer(contname.Data(),TList::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data());
136   contname=Form("coutputmv%s",suffix.Data());
137   AliAnalysisDataContainer *coutputmv = mgr->CreateContainer(contname.Data(),TList::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data());
138   contname=Form("cloosecuts%s",suffix.Data());
139   AliAnalysisDataContainer *coutputcuts = mgr->CreateContainer(contname.Data(),AliRDHFCuts::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data());
140
141
142   mgr->ConnectInput(sigTask,0,mgr->GetCommonInputContainer());
143   
144   mgr->ConnectOutput(sigTask,1,coutputSig);
145   
146   mgr->ConnectOutput(sigTask,2,coutputmv);
147
148   mgr->ConnectOutput(sigTask,3,coutputcuts);
149  
150   return sigTask;
151 }
152