Lines getting the matched track moved to a method in AliCalorimeterUtils. Lines copie...
[u/mrichter/AliRoot.git] / PWG2 / FEMTOSCOPY / macros / Train / TwoTrackQA / ConfigFemtoAnalysis.C
1 /*********************************************************************
2  *                                                                   *
3  * ConfigFemtoAnalysis.C - configuration macro for the femtoscopic   *
4  * analysis, meant as a QA process for two-particle effects          *
5  *                                                                   *
6  * Author: Adam Kisiel (Adam.Kisiel@cern.ch)                         *
7  *                                                                   *
8  *********************************************************************/
9
10 #if !defined(__CINT__) || defined(__MAKECINT_)
11 #include "AliFemtoManager.h"
12 #include "AliFemtoEventReaderESDChain.h"
13 #include "AliFemtoEventReaderESDChainKine.h"
14 #include "AliFemtoSimpleAnalysis.h"
15 #include "AliFemtoBasicEventCut.h"
16 #include "AliFemtoESDTrackCut.h"
17 #include "AliFemtoCorrFctn.h"
18 #include "AliFemtoCutMonitorParticleYPt.h"
19 #include "AliFemtoCutMonitorParticleVertPos.h"
20 #include "AliFemtoCutMonitorParticleMomRes.h"
21 #include "AliFemtoCutMonitorEventMult.h"
22 #include "AliFemtoCutMonitorEventVertex.h"
23 #include "AliFemtoShareQualityTPCEntranceSepPairCut.h"
24 #include "AliFemtoQinvCorrFctn.h"
25 #include "AliFemtoShareQualityCorrFctn.h"
26 #include "AliFemtoTPCInnerCorrFctn.h"
27 #include "AliFemtoVertexMultAnalysis.h"
28 #include "AliFemtoCorrFctn3DSpherical.h"
29 #include "AliFemtoChi2CorrFctn.h"
30 #include "AliFemtoCorrFctnTPCNcls.h"
31 #include "AliFemtoModelBPLCMSCorrFctn.h"
32 #include "AliFemtoModelCorrFctn3DSpherical.h"
33 #include "AliFemtoModelGausLCMSFreezeOutGenerator.h"
34 #include "AliFemtoModelGausRinvFreezeOutGenerator.h"
35 #include "AliFemtoModelManager.h"
36 #include "AliFemtoModelWeightGeneratorBasic.h"
37 #include "AliFemtoModelWeightGeneratorLednicky.h"
38 #include "AliFemtoCorrFctnDirectYlm.h"
39 #include "AliFemtoModelCorrFctnDirectYlm.h"
40 #include "AliFemtoModelCorrFctnSource.h"
41 #include "AliFemtoCutMonitorParticlePtPDG.h"
42 #include "AliFemtoKTPairCut.h"
43 #endif
44
45 //________________________________________________________________________
46 AliFemtoManager* ConfigFemtoAnalysis() {
47
48   double PionMass = 0.13956995;
49   double KaonMass = 0.493677;
50   
51 //   AliFemtoEventReaderESDChainKine* Reader=new AliFemtoEventReaderESDChainKine();
52 //   Reader->SetConstrained(true);
53 //   Reader->SetUseTPCOnly(false);
54
55   AliFemtoEventReaderESDChain* Reader=new AliFemtoEventReaderESDChain();
56   Reader->SetConstrained(true);
57   Reader->SetUseTPCOnly(false);
58
59   AliFemtoManager* Manager=new AliFemtoManager();
60   Manager->SetEventReader(Reader);
61         
62   // Switches for QA analyses
63   int runPiPlusStandard = 1;
64   int runPiMinusStandard = 1;
65   int runPositiveQA = 1;
66   int runNegativeQA = 1;
67   int runPositiveTPCQA = 1;
68   int runNegativeTPCQA = 1;
69
70   // *** First QA task - standard HBT analysis with all the cut on ***
71
72   // *** Begin pion-pion (positive) analysis ***
73   if (runPiPlusStandard) {
74     AliFemtoVertexMultAnalysis *anpipstd = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
75     anpipstd->SetNumEventsToMix(10);
76     anpipstd->SetMinSizePartCollection(2);
77
78     AliFemtoBasicEventCut* mecpipstd = new AliFemtoBasicEventCut();
79     mecpipstd->SetEventMult(2,100000);
80     mecpipstd->SetVertZPos(-1000,1000);
81         
82     AliFemtoESDTrackCut* dtcpipstd = new AliFemtoESDTrackCut();
83     dtcpipstd->SetPidProbPion(0.2,1.001);
84     dtcpipstd->SetPidProbMuon(0.0,1.0);
85     dtcpipstd->SetPidProbKaon(0.0,1.0);
86     dtcpipstd->SetPidProbProton(0.0,1.0);
87     dtcpipstd->SetCharge(1.0);
88     dtcpipstd->SetPt(0.15,0.5);
89     dtcpipstd->SetMass(PionMass);
90     // Track quality cuts
91     dtcpipstd->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
92     //  dtcpipstd->SetStatus(AliESDtrack::kTPCrefit);
93     dtcpipstd->SetminTPCncls(50);
94     dtcpipstd->SetRemoveKinks(kTRUE);
95     dtcpipstd->SetLabel(kFALSE);
96     dtcpipstd->SetMaxITSChiNdof(2.5);
97     dtcpipstd->SetMaxTPCChiNdof(3.0);
98     dtcpipstd->SetMaxImpactXY(3.0);
99     dtcpipstd->SetMaxImpactZ(3.0);
100
101     AliFemtoCutMonitorParticleYPt *cutPassYPtpipstd = new AliFemtoCutMonitorParticleYPt("cutPasspipstd", 0.13957);
102     AliFemtoCutMonitorParticleYPt *cutFailYPtpipstd = new AliFemtoCutMonitorParticleYPt("cutFailpipstd", 0.13957);
103     dtcpipstd->AddCutMonitor(cutPassYPtpipstd, cutFailYPtpipstd);
104
105     AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpipstd = new AliFemtoShareQualityTPCEntranceSepPairCut();
106     sqpcpipstd->SetShareQualityMax(0.0);
107     sqpcpipstd->SetShareFractionMax(0.02);
108     sqpcpipstd->SetRemoveSameLabel(kFALSE);
109     sqpcpipstd->SetTPCEntranceSepMinimum(2.0);
110
111     anpipstd->SetEventCut(mecpipstd);
112     anpipstd->SetFirstParticleCut(dtcpipstd);
113     anpipstd->SetSecondParticleCut(dtcpipstd);
114     anpipstd->SetPairCut(sqpcpipstd);
115     
116     AliFemtoShareQualityCorrFctn *csqqinvpipstd= new AliFemtoShareQualityCorrFctn("sqqinvcfpipstd",40,0.0,0.4);
117     anpipstd->AddCorrFctn(csqqinvpipstd);
118
119     AliFemtoCorrFctnDirectYlm *cylmpipstd = new AliFemtoCorrFctnDirectYlm("cylmpipstd",3,60,0.0,0.3,1);
120     anpipstd->AddCorrFctn(cylmpipstd);
121     
122     AliFemtoQinvCorrFctn *cqinvpipstd = new AliFemtoQinvCorrFctn("qinvcfpipstd", 100,0.0,1.0);
123     anpipstd->AddCorrFctn(cqinvpipstd);
124
125     AliFemtoChi2CorrFctn *cchiqinvpipstd= new AliFemtoChi2CorrFctn("chicfpipstd",40,0.0,0.4);
126     anpipstd->AddCorrFctn(cchiqinvpipstd);
127
128     AliFemtoCorrFctnTPCNcls *cqtpcnclspipstd = new AliFemtoCorrFctnTPCNcls("cqtpcnclspipstd",40,0.0,0.4);
129     anpipstd->AddCorrFctn(cqtpcnclspipstd);
130
131     Manager->AddAnalysis(anpipstd);     
132   }
133   // *** End pion-pion (positive) analysis
134
135   // *** Begin pion-pion (negative) analysis ***
136   if (runPiMinusStandard) {
137     AliFemtoVertexMultAnalysis *anpimstd = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
138     anpimstd->SetNumEventsToMix(10);
139     anpimstd->SetMinSizePartCollection(2);
140
141     AliFemtoBasicEventCut* mecpimstd = new AliFemtoBasicEventCut();
142     mecpimstd->SetEventMult(2,100000);
143     mecpimstd->SetVertZPos(-1000,1000);
144         
145     AliFemtoESDTrackCut* dtcpimstd = new AliFemtoESDTrackCut();
146     dtcpimstd->SetPidProbPion(0.2,1.001);
147     dtcpimstd->SetPidProbMuon(0.0,1.0);
148     dtcpimstd->SetPidProbKaon(0.0,1.0);
149     dtcpimstd->SetPidProbProton(0.0,1.0);
150     dtcpimstd->SetCharge(-1.0);
151     dtcpimstd->SetPt(0.15,0.5);
152     dtcpimstd->SetMass(PionMass);
153     // Track quality cuts
154     dtcpimstd->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
155     //  dtcpimstd->SetStatus(AliESDtrack::kTPCrefit);
156     dtcpimstd->SetminTPCncls(50);
157     dtcpimstd->SetRemoveKinks(kTRUE);
158     dtcpimstd->SetLabel(kFALSE);
159     dtcpimstd->SetMaxITSChiNdof(2.5);
160     dtcpimstd->SetMaxTPCChiNdof(3.0);
161     dtcpimstd->SetMaxImpactXY(3.0);
162     dtcpimstd->SetMaxImpactZ(3.0);
163
164     AliFemtoCutMonitorParticleYPt *cutPassYPtpimstd = new AliFemtoCutMonitorParticleYPt("cutPasspimstd", 0.13957);
165     AliFemtoCutMonitorParticleYPt *cutFailYPtpimstd = new AliFemtoCutMonitorParticleYPt("cutFailpimstd", 0.13957);
166     dtcpimstd->AddCutMonitor(cutPassYPtpimstd, cutFailYPtpimstd);
167
168     AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpimstd = new AliFemtoShareQualityTPCEntranceSepPairCut();
169     sqpcpimstd->SetShareQualityMax(0.0);
170     sqpcpimstd->SetShareFractionMax(0.02);
171     sqpcpimstd->SetRemoveSameLabel(kFALSE);
172     sqpcpimstd->SetTPCEntranceSepMinimum(2.0);
173
174     anpimstd->SetEventCut(mecpimstd);
175     anpimstd->SetFirstParticleCut(dtcpimstd);
176     anpimstd->SetSecondParticleCut(dtcpimstd);
177     anpimstd->SetPairCut(sqpcpimstd);
178     
179     AliFemtoShareQualityCorrFctn *csqqinvpimstd= new AliFemtoShareQualityCorrFctn("sqqinvcfpimstd",40,0.0,0.4);
180     anpimstd->AddCorrFctn(csqqinvpimstd);
181
182     AliFemtoCorrFctnDirectYlm *cylmpimstd = new AliFemtoCorrFctnDirectYlm("cylmpimstd",3,60,0.0,0.3,1);
183     anpimstd->AddCorrFctn(cylmpimstd);
184     
185     AliFemtoQinvCorrFctn *cqinvpimstd = new AliFemtoQinvCorrFctn("qinvcfpimstd", 100,0.0,1.0);
186     anpimstd->AddCorrFctn(cqinvpimstd);
187
188     AliFemtoChi2CorrFctn *cchiqinvpimstd= new AliFemtoChi2CorrFctn("chicfpimstd",40,0.0,0.4);
189     anpimstd->AddCorrFctn(cchiqinvpimstd);
190
191     AliFemtoCorrFctnTPCNcls *cqtpcnclspimstd = new AliFemtoCorrFctnTPCNcls("cqtpcnclspimstd",40,0.0,0.4);
192     anpimstd->AddCorrFctn(cqtpcnclspimstd);
193
194     Manager->AddAnalysis(anpimstd);     
195   }
196   // *** End pion-pion (negative) analysis
197
198   // *** Second QA task - HBT analysis with all pair cuts off ***
199   // *** Begin pion-pion (positive) analysis ***
200   if (runPositiveQA) {
201     AliFemtoVertexMultAnalysis *anpipnct = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
202     anpipnct->SetNumEventsToMix(10);
203     anpipnct->SetMinSizePartCollection(2);
204
205     AliFemtoBasicEventCut* mecpipnct = new AliFemtoBasicEventCut();
206     mecpipnct->SetEventMult(2,100000);
207     mecpipnct->SetVertZPos(-1000,1000);
208         
209     AliFemtoESDTrackCut* dtcpipnct = new AliFemtoESDTrackCut();
210     dtcpipnct->SetPidProbPion(0.0,1.001);
211     dtcpipnct->SetPidProbMuon(0.0,1.0);
212     dtcpipnct->SetPidProbKaon(0.0,1.0);
213     dtcpipnct->SetPidProbProton(0.0,1.0);
214     dtcpipnct->SetCharge(1.0);
215     dtcpipnct->SetPt(0.1,1.0);
216     dtcpipnct->SetMass(PionMass);
217     // Track quality cuts
218     dtcpipnct->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
219     //  dtcpipnct->SetStatus(AliESDtrack::kTPCrefit);
220     dtcpipnct->SetminTPCncls(50);
221     dtcpipnct->SetRemoveKinks(kTRUE);
222     dtcpipnct->SetLabel(kFALSE);
223     dtcpipnct->SetMaxITSChiNdof(6.0);
224     dtcpipnct->SetMaxTPCChiNdof(6.0);
225     dtcpipnct->SetMaxImpactXY(3.0);
226     dtcpipnct->SetMaxImpactZ(3.0);
227
228     AliFemtoCutMonitorParticleYPt *cutPassYPtpipnct = new AliFemtoCutMonitorParticleYPt("cutPasspipnct", 0.13957);
229     AliFemtoCutMonitorParticleYPt *cutFailYPtpipnct = new AliFemtoCutMonitorParticleYPt("cutFailpipnct", 0.13957);
230     dtcpipnct->AddCutMonitor(cutPassYPtpipnct, cutFailYPtpipnct);
231
232     AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpipnct = new AliFemtoShareQualityTPCEntranceSepPairCut();
233     sqpcpipnct->SetShareQualityMax(1.0);
234     sqpcpipnct->SetShareFractionMax(1.0);
235     sqpcpipnct->SetRemoveSameLabel(kFALSE);
236     sqpcpipnct->SetTPCEntranceSepMinimum(0.0);
237
238     anpipnct->SetEventCut(mecpipnct);
239     anpipnct->SetFirstParticleCut(dtcpipnct);
240     anpipnct->SetSecondParticleCut(dtcpipnct);
241     anpipnct->SetPairCut(sqpcpipnct);
242     
243     AliFemtoShareQualityCorrFctn *csqqinvpipnct= new AliFemtoShareQualityCorrFctn("sqqinvcfpipnct",40,0.0,0.4);
244     anpipnct->AddCorrFctn(csqqinvpipnct);
245
246     AliFemtoCorrFctnDirectYlm *cylmpipnct = new AliFemtoCorrFctnDirectYlm("cylmpipnct",3,60,0.0,0.3,1);
247     anpipnct->AddCorrFctn(cylmpipnct);
248     
249     AliFemtoQinvCorrFctn *cqinvpipnct = new AliFemtoQinvCorrFctn("qinvcfpipnct", 100,0.0,1.0);
250     anpipnct->AddCorrFctn(cqinvpipnct);
251
252     AliFemtoChi2CorrFctn *cchiqinvpipnct= new AliFemtoChi2CorrFctn("chicfpipnct",40,0.0,0.4);
253     anpipnct->AddCorrFctn(cchiqinvpipnct);
254
255     AliFemtoCorrFctnTPCNcls *cqtpcnclspipnct = new AliFemtoCorrFctnTPCNcls("cqtpcnclspipnct",40,0.0,0.4);
256     anpipnct->AddCorrFctn(cqtpcnclspipnct);
257
258     Manager->AddAnalysis(anpipnct);     
259   }
260   // *** End pion-pion (positive) analysis
261
262   // *** Begin pion-pion (negative) analysis ***
263   if (runNegativeQA) {
264     AliFemtoVertexMultAnalysis *anpimnct = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
265     anpimnct->SetNumEventsToMix(10);
266     anpimnct->SetMinSizePartCollection(2);
267
268     AliFemtoBasicEventCut* mecpimnct = new AliFemtoBasicEventCut();
269     mecpimnct->SetEventMult(2,100000);
270     mecpimnct->SetVertZPos(-1000,1000);
271         
272     AliFemtoESDTrackCut* dtcpimnct = new AliFemtoESDTrackCut();
273     dtcpimnct->SetPidProbPion(0.0,1.001);
274     dtcpimnct->SetPidProbMuon(0.0,1.0);
275     dtcpimnct->SetPidProbKaon(0.0,1.0);
276     dtcpimnct->SetPidProbProton(0.0,1.0);
277     dtcpimnct->SetCharge(-1.0);
278     dtcpimnct->SetPt(0.1,1.0);
279     dtcpimnct->SetMass(PionMass);
280     // Track quality cuts
281     dtcpimnct->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
282     //  dtcpimnct->SetStatus(AliESDtrack::kTPCrefit);
283     dtcpimnct->SetminTPCncls(50);
284     dtcpimnct->SetRemoveKinks(kTRUE);
285     dtcpimnct->SetLabel(kFALSE);
286     dtcpimnct->SetMaxITSChiNdof(6.0);
287     dtcpimnct->SetMaxTPCChiNdof(6.0);
288     dtcpimnct->SetMaxImpactXY(3.0);
289     dtcpimnct->SetMaxImpactZ(3.0);
290
291     AliFemtoCutMonitorParticleYPt *cutPassYPtpimnct = new AliFemtoCutMonitorParticleYPt("cutPasspimnct", 0.13957);
292     AliFemtoCutMonitorParticleYPt *cutFailYPtpimnct = new AliFemtoCutMonitorParticleYPt("cutFailpimnct", 0.13957);
293     dtcpimnct->AddCutMonitor(cutPassYPtpimnct, cutFailYPtpimnct);
294
295     AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpimnct = new AliFemtoShareQualityTPCEntranceSepPairCut();
296     sqpcpimnct->SetShareQualityMax(1.0);
297     sqpcpimnct->SetShareFractionMax(1.0);
298     sqpcpimnct->SetRemoveSameLabel(kFALSE);
299     sqpcpimnct->SetTPCEntranceSepMinimum(0.0);
300
301     anpimnct->SetEventCut(mecpimnct);
302     anpimnct->SetFirstParticleCut(dtcpimnct);
303     anpimnct->SetSecondParticleCut(dtcpimnct);
304     anpimnct->SetPairCut(sqpcpimnct);
305     
306     AliFemtoShareQualityCorrFctn *csqqinvpimnct= new AliFemtoShareQualityCorrFctn("sqqinvcfpimnct",40,0.0,0.4);
307     anpimnct->AddCorrFctn(csqqinvpimnct);
308
309     AliFemtoCorrFctnDirectYlm *cylmpimnct = new AliFemtoCorrFctnDirectYlm("cylmpimnct",3,60,0.0,0.3,1);
310     anpimnct->AddCorrFctn(cylmpimnct);
311     
312     AliFemtoQinvCorrFctn *cqinvpimnct = new AliFemtoQinvCorrFctn("qinvcfpimnct", 100,0.0,1.0);
313     anpimnct->AddCorrFctn(cqinvpimnct);
314
315     AliFemtoChi2CorrFctn *cchiqinvpimnct= new AliFemtoChi2CorrFctn("chicfpimnct",40,0.0,0.4);
316     anpimnct->AddCorrFctn(cchiqinvpimnct);
317
318     AliFemtoCorrFctnTPCNcls *cqtpcnclspimnct = new AliFemtoCorrFctnTPCNcls("cqtpcnclspimnct",40,0.0,0.4);
319     anpimnct->AddCorrFctn(cqtpcnclspimnct);
320
321     Manager->AddAnalysis(anpimnct);     
322   }
323   // *** End pion-pion (negative) analysis
324
325   // *** Third QA task - HBT analysis with all pair cuts off, TPC only ***
326   // *** Begin pion-pion (positive) analysis ***
327   if (runPositiveTPCQA) {
328     AliFemtoVertexMultAnalysis *anpiptpc = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
329     anpiptpc->SetNumEventsToMix(10);
330     anpiptpc->SetMinSizePartCollection(2);
331
332     AliFemtoBasicEventCut* mecpiptpc = new AliFemtoBasicEventCut();
333     mecpiptpc->SetEventMult(2,100000);
334     mecpiptpc->SetVertZPos(-1000,1000);
335         
336     AliFemtoESDTrackCut* dtcpiptpc = new AliFemtoESDTrackCut();
337     dtcpiptpc->SetPidProbPion(0.0,1.001);
338     dtcpiptpc->SetPidProbMuon(0.0,1.0);
339     dtcpiptpc->SetPidProbKaon(0.0,1.0);
340     dtcpiptpc->SetPidProbProton(0.0,1.0);
341     dtcpiptpc->SetCharge(1.0);
342     dtcpiptpc->SetPt(0.1,1.0);
343     dtcpiptpc->SetMass(PionMass);
344     // Track quality cuts
345     //dtcpiptpc->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
346     dtcpiptpc->SetStatus(AliESDtrack::kTPCrefit);
347     dtcpiptpc->SetminTPCncls(50);
348     dtcpiptpc->SetRemoveKinks(kTRUE);
349     dtcpiptpc->SetLabel(kFALSE);
350     dtcpiptpc->SetMaxITSChiNdof(6.0);
351     dtcpiptpc->SetMaxTPCChiNdof(6.0);
352     dtcpiptpc->SetMaxImpactXY(3.0);
353     dtcpiptpc->SetMaxImpactZ(3.0);
354
355     AliFemtoCutMonitorParticleYPt *cutPassYPtpiptpc = new AliFemtoCutMonitorParticleYPt("cutPasspiptpc", 0.13957);
356     AliFemtoCutMonitorParticleYPt *cutFailYPtpiptpc = new AliFemtoCutMonitorParticleYPt("cutFailpiptpc", 0.13957);
357     dtcpiptpc->AddCutMonitor(cutPassYPtpiptpc, cutFailYPtpiptpc);
358
359     AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpiptpc = new AliFemtoShareQualityTPCEntranceSepPairCut();
360     sqpcpiptpc->SetShareQualityMax(1.0);
361     sqpcpiptpc->SetShareFractionMax(1.0);
362     sqpcpiptpc->SetRemoveSameLabel(kFALSE);
363     sqpcpiptpc->SetTPCEntranceSepMinimum(0.0);
364
365     anpiptpc->SetEventCut(mecpiptpc);
366     anpiptpc->SetFirstParticleCut(dtcpiptpc);
367     anpiptpc->SetSecondParticleCut(dtcpiptpc);
368     anpiptpc->SetPairCut(sqpcpiptpc);
369     
370     AliFemtoShareQualityCorrFctn *csqqinvpiptpc= new AliFemtoShareQualityCorrFctn("sqqinvcfpiptpc",40,0.0,0.4);
371     anpiptpc->AddCorrFctn(csqqinvpiptpc);
372
373     AliFemtoCorrFctnDirectYlm *cylmpiptpc = new AliFemtoCorrFctnDirectYlm("cylmpiptpc",3,60,0.0,0.3,1);
374     anpiptpc->AddCorrFctn(cylmpiptpc);
375     
376     AliFemtoQinvCorrFctn *cqinvpiptpc = new AliFemtoQinvCorrFctn("qinvcfpiptpc", 100,0.0,1.0);
377     anpiptpc->AddCorrFctn(cqinvpiptpc);
378
379     AliFemtoChi2CorrFctn *cchiqinvpiptpc= new AliFemtoChi2CorrFctn("chicfpiptpc",40,0.0,0.4);
380     anpiptpc->AddCorrFctn(cchiqinvpiptpc);
381
382     AliFemtoCorrFctnTPCNcls *cqtpcnclspiptpc = new AliFemtoCorrFctnTPCNcls("cqtpcnclspiptpc",40,0.0,0.4);
383     anpiptpc->AddCorrFctn(cqtpcnclspiptpc);
384
385     Manager->AddAnalysis(anpiptpc);     
386   }
387   // *** End pion-pion (positive) analysis
388
389   // *** Begin pion-pion (negative) analysis ***
390   if (runNegativeTPCQA) {
391     AliFemtoVertexMultAnalysis *anpimtpc = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
392     anpimtpc->SetNumEventsToMix(10);
393     anpimtpc->SetMinSizePartCollection(2);
394
395     AliFemtoBasicEventCut* mecpimtpc = new AliFemtoBasicEventCut();
396     mecpimtpc->SetEventMult(2,100000);
397     mecpimtpc->SetVertZPos(-1000,1000);
398         
399     AliFemtoESDTrackCut* dtcpimtpc = new AliFemtoESDTrackCut();
400     dtcpimtpc->SetPidProbPion(0.0,1.001);
401     dtcpimtpc->SetPidProbMuon(0.0,1.0);
402     dtcpimtpc->SetPidProbKaon(0.0,1.0);
403     dtcpimtpc->SetPidProbProton(0.0,1.0);
404     dtcpimtpc->SetCharge(-1.0);
405     dtcpimtpc->SetPt(0.1,1.0);
406     dtcpimtpc->SetMass(PionMass);
407     // Track quality cuts
408     //dtcpimtpc->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
409     dtcpimtpc->SetStatus(AliESDtrack::kTPCrefit);
410     dtcpimtpc->SetminTPCncls(50);
411     dtcpimtpc->SetRemoveKinks(kTRUE);
412     dtcpimtpc->SetLabel(kFALSE);
413     dtcpimtpc->SetMaxITSChiNdof(6.0);
414     dtcpimtpc->SetMaxTPCChiNdof(6.0);
415     dtcpimtpc->SetMaxImpactXY(3.0);
416     dtcpimtpc->SetMaxImpactZ(3.0);
417
418     AliFemtoCutMonitorParticleYPt *cutPassYPtpimtpc = new AliFemtoCutMonitorParticleYPt("cutPasspimtpc", 0.13957);
419     AliFemtoCutMonitorParticleYPt *cutFailYPtpimtpc = new AliFemtoCutMonitorParticleYPt("cutFailpimtpc", 0.13957);
420     dtcpimtpc->AddCutMonitor(cutPassYPtpimtpc, cutFailYPtpimtpc);
421
422     AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpimtpc = new AliFemtoShareQualityTPCEntranceSepPairCut();
423     sqpcpimtpc->SetShareQualityMax(1.0);
424     sqpcpimtpc->SetShareFractionMax(1.0);
425     sqpcpimtpc->SetRemoveSameLabel(kFALSE);
426     sqpcpimtpc->SetTPCEntranceSepMinimum(0.0);
427
428     anpimtpc->SetEventCut(mecpimtpc);
429     anpimtpc->SetFirstParticleCut(dtcpimtpc);
430     anpimtpc->SetSecondParticleCut(dtcpimtpc);
431     anpimtpc->SetPairCut(sqpcpimtpc);
432     
433     AliFemtoShareQualityCorrFctn *csqqinvpimtpc= new AliFemtoShareQualityCorrFctn("sqqinvcfpimtpc",40,0.0,0.4);
434     anpimtpc->AddCorrFctn(csqqinvpimtpc);
435
436     AliFemtoCorrFctnDirectYlm *cylmpimtpc = new AliFemtoCorrFctnDirectYlm("cylmpimtpc",3,60,0.0,0.3,1);
437     anpimtpc->AddCorrFctn(cylmpimtpc);
438     
439     AliFemtoQinvCorrFctn *cqinvpimtpc = new AliFemtoQinvCorrFctn("qinvcfpimtpc", 100,0.0,1.0);
440     anpimtpc->AddCorrFctn(cqinvpimtpc);
441
442     AliFemtoChi2CorrFctn *cchiqinvpimtpc= new AliFemtoChi2CorrFctn("chicfpimtpc",40,0.0,0.4);
443     anpimtpc->AddCorrFctn(cchiqinvpimtpc);
444
445     AliFemtoCorrFctnTPCNcls *cqtpcnclspimtpc = new AliFemtoCorrFctnTPCNcls("cqtpcnclspimtpc",40,0.0,0.4);
446     anpimtpc->AddCorrFctn(cqtpcnclspimtpc);
447
448     Manager->AddAnalysis(anpimtpc);     
449   }
450   // *** End pion-pion (negative) analysis
451
452
453
454   return Manager;
455 }                         
456