b129735616ed10d01bfbcfc239bc179d81f64fac
[u/mrichter/AliRoot.git] / PWG4 / JCORRAN / AliAnalysisTaskDiHadron.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15 //2- and 3-particle trigger particle correlation analysis
16 //Author: Jason Glyndwr Ulery, ulery@uni-frankfurt.de
17 //version: 3.4,  last revised: 2010/08/15
18
19 #include "Riostream.h"
20 #include "TChain.h"
21 #include "TTree.h"
22 #include "TH1F.h"
23 #include "TH2F.h"
24 #include "TH3F.h"
25 #include "TFormula.h"
26 #include "TF1.h"
27 #include "TF2.h"
28 #include "TF3.h"
29 #include "TVector3.h"
30 #include "TMath.h"
31
32
33 #include "AliAnalysisTask.h"
34 #include "AliAnalysisManager.h"
35
36 #include "AliESDEvent.h"
37 #include "AliESDInputHandler.h"
38 #include "AliMCEvent.h"
39 #include "AliMCEventHandler.h"
40 #include "AliMCParticle.h"
41 #include "AliStack.h"
42 #include "AliESDVertex.h"
43 #include "AliMultiplicity.h"
44 #include "TParticle.h"
45
46
47 //#include "AliHeader.h"
48 //#include "AliGenEventHeader.h"
49
50
51
52 #include "AliAnalysisTaskDiHadron.h"
53
54
55 ClassImp(AliAnalysisTaskDiHadron)
56
57 //----------------------------------------
58 AliAnalysisTaskDiHadron::AliAnalysisTaskDiHadron(const char *name):
59 AliAnalysisTask(name,""), fESD(0), fMC(0), fOutput(0),fMinClustersTPC(0),fMinClusterRatio(0),fMaxTPCchi2(0),fMinClustersITS(0),fEtaCut(0),fTrigEtaCut(0),fNearPhiCut(0),fXECut(0),fMaxDCA(0),fMaxDCAXY(0),fMaxDCAZ(0),fDCA2D(0),fTPCRefit(0),fITSRefit(0),fSPDCut(0),fMinPtAssoc(0),fMaxPtAssoc(0),fVzCut(0),fEfficiencyCorr(0),fDEBUG(0),fnBinPhi(0),fnBinEta(0),fnBinPhiEtaPhi(0),fnBinPhiEtaEta(0),fnBinPhi3(0),fnBinEta3(0),fPi(3.1415926535898),fdPhiMin(0),fdPhiMax(0),fNTPtBins(0),fNMix(0),fNCentBins(0),fNAPtBins(0),fNAPt3Bins(0),fNVertexBins(0),fNXEBins(0),fNIDs(0),fEffFitPt(0),fNFitLowParam(0),fNFitHighParam(0),fMCHistos(0),fFitLow(NULL),fFitHigh(NULL),fFitLowParam(NULL),fFitHighParam(NULL),fPtTrigArray(NULL),fPtAssocArray(NULL),fPtAssoc3Array1(NULL),fPtAssoc3Array2(NULL),fCentArrayMin(NULL),fCentArrayMax(NULL),fXEArray(NULL),fTrigIDArray(NULL),ftPhi(NULL),ftEta(NULL),ftPt(NULL),ftCharge(NULL),ftEff(NULL),ftPtAssoc3(NULL),ftNPtAssoc3(NULL)
60   
61   {
62  
63   //IO Slots
64   DefineInput(0, TChain::Class());
65   DefineOutput(0,TList::Class());
66   
67   
68   for(int c=0;c<fNCentBins;c++){
69     for(int v=0;v<fNVertexBins;v++){
70       for(int jmc=0;jmc<2;jmc++){
71         fMixPointer[c][v][jmc]=-1;
72         fMixEnd[c][v][jmc]=-1;
73         for(int ievts=0;ievts<fNMix;ievts++){
74           fMPt[ievts][c][v][jmc]=NULL;
75           fMPhi[ievts][c][v][jmc]=NULL;
76           fMEta[ievts][c][v][jmc]=NULL;
77           for(int dd=0;dd<10;dd++)fMPtAssoc3[ievts][c][v][jmc][dd]=NULL;
78           fMNPtAssoc3[ievts][c][v][jmc]=NULL;
79           fMixTrack[ievts][c][v][jmc]=0;
80         }
81       }
82     }
83   }
84  
85   }
86 //--------------------------------------
87 void AliAnalysisTaskDiHadron::SetCuts(Int_t MinClustersTPC,  Float_t MinClusterRatio, Float_t MaxTPCchi2, Int_t MinClustersITS, Float_t EtaCut, Float_t TrigEtaCut, Float_t NearPhiCut, Float_t XECut, Float_t MaxDCA, Float_t MaxDCAXY, Float_t MaxDCAZ, Int_t DCA2D, Int_t TPCRefit, Int_t ITSRefit, Int_t SPDCut, Float_t MinPtAssoc, Float_t MaxPtAssoc, Float_t VzCut, Int_t NIDs, const char * TrigIDArray){
88 //Sets the varibles for track and event cuts
89   fMinClustersTPC=MinClustersTPC;
90   fMinClusterRatio=MinClusterRatio;
91   fMaxTPCchi2=MaxTPCchi2;
92   fMinClustersITS=MinClustersITS;
93   fEtaCut=EtaCut;
94   fTrigEtaCut=TrigEtaCut;
95   fNearPhiCut=NearPhiCut;
96   fXECut=XECut;
97   fMaxDCA=MaxDCA;
98   fMaxDCAXY=MaxDCAXY;
99   fMaxDCAZ=MaxDCAZ;
100   fDCA2D=DCA2D;
101   fTPCRefit=TPCRefit;
102   fITSRefit=ITSRefit;
103   fSPDCut=SPDCut;
104   fMinPtAssoc=MinPtAssoc;
105   fMaxPtAssoc=MaxPtAssoc;
106   fVzCut=VzCut;
107   fNIDs=NIDs;
108   fTrigIDArray=(char*)TrigIDArray;
109 }
110 //--------------------------------------------------------
111 void AliAnalysisTaskDiHadron::SetOptions(Int_t EfficiencyCorr, Int_t ffDEBUG,  Int_t MCHistos){
112 //Sets some options
113   fEfficiencyCorr=EfficiencyCorr;
114   fDEBUG=ffDEBUG;
115   fMCHistos=MCHistos;
116 }
117 //------------------------------------------------------
118 void AliAnalysisTaskDiHadron::SetBins(Int_t nBinPhi, Int_t nBinEta, Int_t nBinPhiEtaPhi, Int_t nBinPhiEtaEta, Int_t nBinPhi3, Int_t nBinEta3,Float_t dPhiMin, Float_t dPhiMax, Int_t NTPtBins, Int_t NMixBins, Int_t NCentBins,Int_t NAPtBins, Int_t NAPt3Bins, Int_t NVertexBins, Int_t NXEBins,Float_t *PtTrigArray, Float_t *PtAssocArray,Float_t *PtAssoc3Array1, Float_t *PtAssoc3Array2, Int_t *CentArrayMin, Int_t *CentArrayMax, Float_t *XEArray){
119 //sets up the histogram binning
120   fnBinPhi=nBinPhi;
121   fnBinEta=nBinEta;
122   fnBinPhiEtaPhi=nBinPhiEtaPhi;
123   fnBinPhiEtaEta=nBinPhiEtaEta;
124   fnBinPhi3=nBinPhi3;
125   fnBinEta3=nBinEta3;
126   fdPhiMin=dPhiMin;
127   fdPhiMax=dPhiMax;
128   fNTPtBins=NTPtBins;
129   fNMix=NMixBins;
130   fNCentBins=NCentBins;
131   fNAPtBins=NAPtBins;
132   fNAPt3Bins=NAPt3Bins;
133   fNVertexBins=NVertexBins;
134   fNXEBins=NXEBins;
135   fPtTrigArray=new Float_t [fNTPtBins];
136   for(int i=0;i<fNTPtBins;i++)fPtTrigArray[i]=PtTrigArray[i];
137   fPtAssocArray=new Float_t [fNAPtBins];
138   for(int i=0;i<fNAPtBins;i++)fPtAssocArray[i]=PtAssocArray[i];
139   fPtAssoc3Array1=new Float_t [fNAPt3Bins];
140   for(int i=0;i<fNAPt3Bins;i++)fPtAssoc3Array1[i]=PtAssoc3Array1[i];
141   fPtAssoc3Array2=new Float_t [fNAPt3Bins];
142   for(int i=0;i<fNAPt3Bins;i++)fPtAssoc3Array2[i]=PtAssoc3Array2[i];
143   fCentArrayMin=new Int_t [fNCentBins];
144   for(int i=0;i<NCentBins;i++)fCentArrayMin[i]=CentArrayMin[i];
145   fCentArrayMax=new Int_t [fNCentBins];
146   for(int i=0;i<NCentBins;i++)fCentArrayMax[i]=CentArrayMax[i];
147   fXEArray=new Float_t [fNXEBins];
148   for(int i=0;i<fNXEBins;i++)fXEArray[i]=XEArray[i];
149  for(int i=0;i<=fNVertexBins;i++)fVertexArray[i]=(2.*i/fNVertexBins-1)*fVzCut;
150 }
151 //-------------------------------------------------------
152 void AliAnalysisTaskDiHadron::SetEfficiencies(Float_t EffFitPt, const TF1 *FitLow, const TF1 *FitHigh, Int_t NFitLowParam, Int_t NFitHighParam, Float_t *FitLowParam, Float_t *FitHighParam){
153 //Sets up the efficiency corrections
154   fEffFitPt=EffFitPt;
155   fFitLow=(TF1*)FitLow;
156   fFitHigh=(TF1*)FitHigh;
157   fNFitLowParam=NFitLowParam;
158   fNFitHighParam=NFitHighParam;
159   fFitLowParam=new Float_t [fNFitLowParam*fNCentBins];
160   for(int i=0;i<fNFitLowParam*fNCentBins;i++)fFitLowParam[i]=FitLowParam[i];
161   fFitHighParam=new Float_t [fNFitHighParam*fNCentBins];
162   for(int i=0;i<fNFitHighParam*fNCentBins;i++)fFitHighParam[i]=FitHighParam[i];
163 }
164
165 //-----------------------------------------------------------
166 void AliAnalysisTaskDiHadron::ConnectInputData(Option_t *){
167   //Connect to ESD
168   if(fDEBUG)Printf("Connecting");
169    TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
170    if (!tree&&fDEBUG) {Printf("ERROR: Could not read chain from input slot 0");} 
171   else {
172     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
173     if (!esdH&&fDEBUG) Printf("ERROR: Could not get ESDInputHandler");
174     else fESD = esdH->GetEvent();
175     
176     //MC Data handler (so one can calcualte eff)
177     AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler*>((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
178     if(mcH)fMC=mcH->MCEvent();
179   }
180   if(fDEBUG)Printf("Connected");
181 }
182   
183 //---------------------------------------------------------  
184 void AliAnalysisTaskDiHadron::CreateOutputObjects(){
185   //Creates the histograms and list
186   if(fDEBUG)Printf("Output");
187   fOutput=new TList();
188   fOutput->SetName(GetName());
189   char histname[100];
190   char histtitle[200];
191   int nptbins=fNAPtBins;
192   int lptbins=0;
193   const char *cmc1[2]={"","_MC"};
194   const char *cmc2[2]={""," MC"};
195   const char *sign1[3]={"","_LS","_ULS"};
196   const char *sign2[3]={""," Like-Sign"," Unlike-Sign"};
197   const char *sign31[4]={"","_LS","_ULT","_ULA"};
198   const char *sign32[4]={""," Like-Sign"," Trigger-Diff"," Assoc-Diff"};
199   Float_t etaEdge=fEtaCut+fTrigEtaCut;
200   Float_t phiArray[fnBinPhi+1];
201   Float_t etaArray[fnBinEta+1];
202   Float_t phiEtaArrayPhi[fnBinPhiEtaPhi+1];
203   Float_t phiEtaArrayEta[fnBinPhiEtaEta+1];
204   for(int iphi=0;iphi<=fnBinPhi;iphi++){
205     phiArray[iphi]=fdPhiMin+iphi*2*fPi/fnBinPhi;
206   }
207   for(int ieta=0;ieta<=fnBinEta;ieta++){
208     etaArray[ieta]=-etaEdge+ieta*2*etaEdge/fnBinEta;
209   }
210   for(int iphi=0;iphi<=fnBinPhiEtaPhi;iphi++){
211     phiEtaArrayPhi[iphi]=fdPhiMin+iphi*2*fPi/fnBinPhiEtaPhi;
212   }
213   for(int ieta=0;ieta<=fnBinPhiEtaEta;ieta++){
214     phiEtaArrayEta[ieta]=-etaEdge+ieta*2*etaEdge/fnBinPhiEtaEta;
215   }
216   for(int imc=0;imc<=1;imc++){//MC loop
217     if(imc==1&&!fMCHistos) continue;
218     //Create the histograms
219     sprintf(histname,"fHistMult%s",cmc1[imc]);
220     sprintf(histtitle,"Multiplicity%s",cmc2[imc]);
221     fHistMult[imc]=new TH1F(histname,histtitle,2000,-0.5,1999.5);
222     fHistMult[imc]->Sumw2();
223     fHistMult[imc]->GetXaxis()->SetTitle("Number of tracks");
224     fHistMult[imc]->GetYaxis()->SetTitle("Counts");
225     fOutput->Add(fHistMult[imc]);
226     
227     for(int imult=0;imult<fNCentBins;imult++){//loop for multiplicity bins
228       
229       //Histograms that are independent of the trigger
230       sprintf(histname,"fHistPt_C%d%s",imult,cmc1[imc]);
231       sprintf(histtitle,"P_{T} Distribution of Tracks %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
232       fHistPt[imult][imc]=new TH1F(histname,histtitle,nptbins,fPtAssocArray);
233       fHistPt[imult][imc]->Sumw2();
234       fHistPt[imult][imc]->GetXaxis()->SetTitle("p_{T}");
235       fHistPt[imult][imc]->GetYaxis()->SetTitle("Counts");
236       fOutput->Add(fHistPt[imult][imc]);
237
238  //Histograms that are independent of the trigger
239       sprintf(histname,"fHistPtEff_C%d%s",imult,cmc1[imc]);
240       sprintf(histtitle,"P_{T} Distribution of Tracks %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
241       fHistPtEff[imult][imc]=new TH1F(histname,histtitle,1000,0,100);
242       fHistPtEff[imult][imc]->Sumw2();
243       fHistPtEff[imult][imc]->GetXaxis()->SetTitle("p_{T}");
244       fHistPtEff[imult][imc]->GetYaxis()->SetTitle("Counts");
245       fOutput->Add(fHistPtEff[imult][imc]);
246       
247       sprintf(histname,"fHistPhi_C%d%s",imult,cmc1[imc]);
248       sprintf(histtitle,"#phi Distribution of Tracks %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
249       fHistPhi[imult][imc]=new TH2F(histname,histtitle,fnBinPhi,phiArray,nptbins,fPtAssocArray);
250       fHistPhi[imult][imc]->Sumw2();
251       fHistPhi[imult][imc]->GetXaxis()->SetTitle("#phi");
252       fHistPhi[imult][imc]->GetYaxis()->SetTitle("P_{T}");
253       fOutput->Add(fHistPhi[imult][imc]);
254       
255       sprintf(histname,"fHistPhiPt_C%d%s",imult,cmc1[imc]);
256       sprintf(histtitle,"P_{T} weighted #phi Distribution of tracks %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
257       fHistPhiPt[imult][imc]=new TH2F(histname,histtitle,fnBinPhi,phiArray,nptbins,fPtAssocArray);
258       fHistPhiPt[imult][imc]->Sumw2();
259       fHistPhiPt[imult][imc]->GetXaxis()->SetTitle("#phi");
260       fHistPhiPt[imult][imc]->GetYaxis()->SetTitle("P_{T}");
261       fOutput->Add(fHistPhiPt[imult][imc]);
262       
263       sprintf(histname,"fHistEta_C%d%s",imult,cmc1[imc]);
264       sprintf(histtitle,"#eta Distribution of Tracks %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
265       fHistEta[imult][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,nptbins,fPtAssocArray);
266       fHistEta[imult][imc]->Sumw2();
267       fHistEta[imult][imc]->GetXaxis()->SetTitle("#eta");
268       fHistEta[imult][imc]->GetYaxis()->SetTitle("P_{T}");
269       fOutput->Add(fHistEta[imult][imc]);
270       
271       sprintf(histname,"fHistEtaPt_C%d%s",imult,cmc1[imc]);
272       sprintf(histtitle,"P_{T} weighted #eta Distribution of tracks %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
273       fHistEtaPt[imult][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,nptbins,fPtAssocArray);
274       fHistEtaPt[imult][imc]->Sumw2();
275       fHistEtaPt[imult][imc]->GetXaxis()->SetTitle("#eta");
276       fHistEtaPt[imult][imc]->GetYaxis()->SetTitle("P_{T}");
277       fOutput->Add(fHistEtaPt[imult][imc]);
278       
279       sprintf(histname,"fHistNEvents_C%d%s",imult,cmc1[imc]);
280       sprintf(histtitle,"Number of Events and Number Passing Cuts %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
281       fHistNEvents[imult][imc]=new TH1F(histname,histtitle,2,-0.5,1.5);
282       fHistNEvents[imult][imc]->Sumw2();
283       fHistNEvents[imult][imc]->GetXaxis()->SetTitle("Events,Passing Cuts");
284       fHistNEvents[imult][imc]->GetYaxis()->SetTitle("Number of Events");
285       fOutput->Add(fHistNEvents[imult][imc]);
286       
287       sprintf(histname,"fHistNTrigger_C%d%s",imult,cmc1[imc]);
288       sprintf(histtitle,"Number of Triggers %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
289       fHistNTrigger[imult][imc]=new TH1F(histname,histtitle,fNTPtBins,-0.5,fNTPtBins-0.5);
290       fHistNTrigger[imult][imc]->Sumw2();
291       fHistNTrigger[imult][imc]->GetXaxis()->SetTitle("Trigger Number");
292       fHistNTrigger[imult][imc]->GetYaxis()->SetTitle("Number of Triggers");
293       fOutput->Add(fHistNTrigger[imult][imc]);
294       
295       sprintf(histname,"fHistNTriggerPt_C%d%s",imult,cmc1[imc]);
296       sprintf(histtitle,"P_{T} Weighted Number of Triggers %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
297       fHistNTriggerPt[imult][imc]=new TH1F(histname,histtitle,fNTPtBins,-0.5,fNTPtBins-0.5);
298       fHistNTriggerPt[imult][imc]->Sumw2();
299       fHistNTriggerPt[imult][imc]->GetXaxis()->SetTitle("Trigger Number");
300       fHistNTriggerPt[imult][imc]->GetYaxis()->SetTitle("Number of Triggers");
301       fOutput->Add(fHistNTriggerPt[imult][imc]);
302       
303       sprintf(histname,"fHistNMix_C%d%s",imult,cmc1[imc]);
304       sprintf(histtitle,"Number of Mixed Events %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
305       fHistNMix[imult][imc]=new TH1F(histname,histtitle,fNTPtBins,-0.5,fNTPtBins-0.5);
306       fHistNMix[imult][imc]->Sumw2();
307       fHistNMix[imult][imc]->GetXaxis()->SetTitle("Trigger Number");
308       fHistNMix[imult][imc]->GetYaxis()->SetTitle("Number of Mixed Events");
309       fOutput->Add(fHistNMix[imult][imc]);
310       
311       sprintf(histname,"fHistPhiEta_C%d%s",imult,cmc1[imc]);
312       sprintf(histtitle,"#phi-#eta distribution of tracks %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
313       fHistPhiEta[imult][imc]=new TH3F(histname, histtitle,fnBinPhiEtaPhi,phiEtaArrayPhi,fnBinPhiEtaEta,phiEtaArrayEta,nptbins,fPtAssocArray);
314       fHistPhiEta[imult][imc]->Sumw2();
315       fHistPhiEta[imult][imc]->GetXaxis()->SetTitle("#phi");
316       fHistPhiEta[imult][imc]->GetYaxis()->SetTitle("#eta");
317       fHistPhiEta[imult][imc]->GetZaxis()->SetTitle("p_{T}");
318       fOutput->Add(fHistPhiEta[imult][imc]);
319       
320       sprintf(histname,"fHistPhiEtaPt_C%d%s",imult,cmc1[imc]);
321       sprintf(histtitle,"Pt Weighted #phi-#eta distribution of tracks %dMult%d%s",fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
322       fHistPhiEtaPt[imult][imc]=new TH3F(histname, histtitle,fnBinPhiEtaPhi,phiEtaArrayPhi,fnBinPhiEtaEta,phiEtaArrayEta,nptbins,fPtAssocArray);
323       fHistPhiEtaPt[imult][imc]->Sumw2();
324       fHistPhiEtaPt[imult][imc]->GetXaxis()->SetTitle("#phi");
325       fHistPhiEtaPt[imult][imc]->GetYaxis()->SetTitle("#eta");
326       fHistPhiEtaPt[imult][imc]->GetZaxis()->SetTitle("p_{T}");
327       fOutput->Add(fHistPhiEtaPt[imult][imc]);
328       
329       //if(fDEBUG)Printf("OutPut2");
330       //Histograms with a trigger dependence
331       
332       for(int i=0;i<fNTPtBins;i++){
333         for(int j=1;j<fNAPtBins;j++){
334           if(fabs(fPtTrigArray[i]-fPtAssocArray[j])<1E-5)lptbins=j;
335         }
336         //if(fDEBUG)Printf("Loop: %d Pt %3.2f",i,fPtTrigArray[i]/fPtBinWidth);
337         
338         //Ones with no centrality binning
339         if(imult==0){
340           sprintf(histname,"fHistMultTrig_P%d%s",i,cmc1[imc]);
341           sprintf(histtitle,"Distrubition of number of tracks in triggered events with %3.1f<p_{T}^{Trig}<%3.1f%s",fPtTrigArray[i],fPtTrigArray[i+1],cmc2[imc]);
342           fHistMultTrig[i][imc]=new TH1F(histname,histtitle,2000,0,2000);
343           fHistMultTrig[i][imc]->Sumw2();
344           fHistMultTrig[i][imc]->GetXaxis()->SetTitle("Number of Tracks");
345           fHistMultTrig[i][imc]->GetYaxis()->SetTitle("Counts");
346           fOutput->Add(fHistMultTrig[i][imc]);
347         }
348         sprintf(histname,"fHistPtTrig_P%d_C%d%s",i,imult,cmc1[imc]);
349         sprintf(histtitle,"P_{T} distribution in triggered events with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
350         fHistPtTrig[i][imult][imc]=new TH1F(histname,histtitle,nptbins,fPtAssocArray);
351         fHistPtTrig[i][imult][imc]->Sumw2();
352         fHistPtTrig[i][imult][imc]->GetXaxis()->SetTitle("p_{T}");
353         fHistPtTrig[i][imult][imc]->GetYaxis()->SetTitle("Counts");
354         fOutput->Add(fHistPtTrig[i][imult][imc]);
355         
356         sprintf(histname,"fHistPhiTrig_P%d_C%d%s",i,imult,cmc1[imc]);
357         sprintf(histtitle,"Phi Distribution of triggered events with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
358         fHistPhiTrig[i][imult][imc]=new TH2F(histname,histtitle,fnBinPhi,phiArray,lptbins,fPtAssocArray);
359         fHistPhiTrig[i][imult][imc]->Sumw2();
360         fHistPhiTrig[i][imult][imc]->GetXaxis()->SetTitle("#phi");
361         fHistPhiTrig[i][imult][imc]->GetYaxis()->SetTitle("p_{T}");
362         fOutput->Add(fHistPhiTrig[i][imult][imc]);
363         
364         sprintf(histname,"fHistPhiTrigPt_P%d_C%d%s",i,imult,cmc1[imc]);
365         sprintf(histtitle,"P_{T} Weighted Phi Distribution of triggered events with %3.1f<p_{T}^{Trig}<%3.1f%s",fPtTrigArray[i],fPtTrigArray[i+1],cmc2[imc]);
366         fHistPhiTrigPt[i][imult][imc]=new TH2F(histname,histtitle,fnBinPhi,phiArray,lptbins,fPtAssocArray);
367         fHistPhiTrigPt[i][imult][imc]->Sumw2();
368         fHistPhiTrigPt[i][imult][imc]->GetXaxis()->SetTitle("#phi");
369         fHistPhiTrigPt[i][imult][imc]->GetYaxis()->SetTitle("p_{T}");
370         fOutput->Add(fHistPhiTrigPt[i][imult][imc]);
371         
372         sprintf(histname,"fHistEtaTrig_P%d_C%d%s",i,imult,cmc1[imc]);
373         sprintf(histtitle,"Eta Distribution of triggered events with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
374         fHistEtaTrig[i][imult][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,lptbins,fPtAssocArray);
375         fHistEtaTrig[i][imult][imc]->Sumw2();
376         fHistEtaTrig[i][imult][imc]->GetXaxis()->SetTitle("#eta");
377         fHistEtaTrig[i][imult][imc]->GetYaxis()->SetTitle("p_{T}");
378         fOutput->Add(fHistEtaTrig[i][imult][imc]);
379         
380         sprintf(histname,"fHistEtaTrigPt_P%d_C%d%s",i,imult,cmc1[imc]);
381         sprintf(histtitle,"P_{T} Weighted Eta Distribution of triggered events with %3.1f<p_{T}^{Trig}<%3.1f%s",fPtTrigArray[i],fPtTrigArray[i+1],cmc2[imc]);
382         fHistEtaTrigPt[i][imult][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,lptbins,fPtAssocArray);
383         fHistEtaTrigPt[i][imult][imc]->Sumw2();
384         fHistEtaTrigPt[i][imult][imc]->GetXaxis()->SetTitle("#eta");
385         fHistEtaTrigPt[i][imult][imc]->GetYaxis()->SetTitle("p_{T}");
386         fOutput->Add(fHistEtaTrigPt[i][imult][imc]);
387         
388         sprintf(histname,"fHistPhiEtaTrig_P%d_C%d%s",i,imult,cmc1[imc]);
389         sprintf(histtitle,"#phi-#eta distribution in triggered events %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
390         fHistPhiEtaTrig[i][imult][imc]=new TH3F(histname,histtitle,fnBinPhiEtaPhi,phiEtaArrayPhi,fnBinPhiEtaEta,phiEtaArrayEta,lptbins,fPtAssocArray);
391         fHistPhiEtaTrig[i][imult][imc]->Sumw2();
392         fHistPhiEtaTrig[i][imult][imc]->GetXaxis()->SetTitle("#phi");
393         fHistPhiEtaTrig[i][imult][imc]->GetYaxis()->SetTitle("#eta");
394         fHistPhiEtaTrig[i][imult][imc]->GetZaxis()->SetTitle("p_{T}");
395         fOutput->Add(fHistPhiEtaTrig[i][imult][imc]);
396
397         sprintf(histname,"fHistXEN_P%d_C%d%s",i,imult,cmc1[imc]);
398         sprintf(histtitle,"Near-Side X_{E} distribution for %3.1f<p_{T}^{Lead}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
399         fHistXEN[i][imult][imc]=new TH1F(histname,histtitle,fNXEBins,fXEArray);
400         fHistXEN[i][imult][imc]->Sumw2();
401         fHistXEN[i][imult][imc]->GetXaxis()->SetTitle("X_{E}");
402         fOutput->Add(fHistXEN[i][imult][imc]);
403
404         sprintf(histname,"fHistXENMixed_P%d_C%d%s",i,imult,cmc1[imc]);
405         sprintf(histtitle,"Mixed Near-Side X_{E} distribution for %3.1f<p_{T}^{Lead}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
406         fHistXENMix[i][imult][imc]=new TH1F(histname,histtitle,fNXEBins,fXEArray);
407         fHistXENMix[i][imult][imc]->Sumw2();
408         fHistXENMix[i][imult][imc]->GetXaxis()->SetTitle("X_{E}");
409         fOutput->Add(fHistXENMix[i][imult][imc]);
410         
411         sprintf(histname,"fHistXEA_P%d_C%d%s",i,imult,cmc1[imc]);
412         sprintf(histtitle,"Away-Side X_{E} distribution for %3.1f<p_{T}^{Lead}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
413         fHistXEA[i][imult][imc]=new TH1F(histname,histtitle,fNXEBins,fXEArray);
414         fHistXEA[i][imult][imc]->Sumw2();
415         fHistXEA[i][imult][imc]->GetXaxis()->SetTitle("X_{E}");
416         fOutput->Add(fHistXEA[i][imult][imc]);
417
418         sprintf(histname,"fHistXEAMixed_P%d_C%d%s",i,imult,cmc1[imc]);
419         sprintf(histtitle,"Mixed Away-Side X_{E} distribution for %3.1f<p_{T}^{Lead}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
420         fHistXEAMix[i][imult][imc]=new TH1F(histname,histtitle,fNXEBins,fXEArray);
421         fHistXEAMix[i][imult][imc]->Sumw2();
422         fHistXEAMix[i][imult][imc]->GetXaxis()->SetTitle("X_{E}");
423         fOutput->Add(fHistXEAMix[i][imult][imc]);
424
425         //signloop
426         for(int isign=0;isign<3;isign++){
427           sprintf(histname,"fHistDeltaPhi_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
428           sprintf(histtitle,"#Delta#phi Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
429           fHistDeltaPhi[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinPhi,phiArray,lptbins,fPtAssocArray);
430           fHistDeltaPhi[i][imult][isign][imc]->Sumw2();
431           fHistDeltaPhi[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#phi");
432           fHistDeltaPhi[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
433           fOutput->Add(fHistDeltaPhi[i][imult][isign][imc]);
434           
435           sprintf(histname,"fHistDeltaPhiPt_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
436           sprintf(histtitle,"P_{T} Weighted #Delta#phi Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
437           fHistDeltaPhiPt[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinPhi,phiArray,lptbins,fPtAssocArray);
438           fHistDeltaPhiPt[i][imult][isign][imc]->Sumw2();
439           fHistDeltaPhiPt[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#phi");
440           fHistDeltaPhiPt[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
441           fOutput->Add(fHistDeltaPhiPt[i][imult][isign][imc]);
442           
443           sprintf(histname,"fHistDeltaPhiMix_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
444           sprintf(histtitle,"#Delta#phi Mixed Event Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
445           fHistDeltaPhiMix[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinPhi,phiArray,lptbins,fPtAssocArray);
446           fHistDeltaPhiMix[i][imult][isign][imc]->Sumw2();
447           fHistDeltaPhiMix[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#phi");
448           fHistDeltaPhiMix[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
449           fOutput->Add(fHistDeltaPhiMix[i][imult][isign][imc]);
450           
451           sprintf(histname,"fHistDeltaPhiMixPt_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
452           sprintf(histtitle,"P_{T} Weighted #Delta#phi Mixed Event Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
453           fHistDeltaPhiMixPt[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinPhi,phiArray,lptbins,fPtAssocArray);
454           fHistDeltaPhiMixPt[i][imult][isign][imc]->Sumw2();
455           fHistDeltaPhiMixPt[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#phi");
456           fHistDeltaPhiMixPt[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
457           fOutput->Add(fHistDeltaPhiMixPt[i][imult][isign][imc]);
458           
459           //etaNear
460           sprintf(histname,"fHistDeltaEtaN_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
461           sprintf(histtitle,"Near-Side #Delta#eta Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
462           fHistDeltaEtaN[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,lptbins,fPtAssocArray);
463           fHistDeltaEtaN[i][imult][isign][imc]->Sumw2();
464           fHistDeltaEtaN[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#eta");
465           fHistDeltaEtaN[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
466           fOutput->Add(fHistDeltaEtaN[i][imult][isign][imc]);
467           
468           sprintf(histname,"fHistDeltaEtaNPt_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
469           sprintf(histtitle,"Near-Side P_{T} Weighted #Delta#eta Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
470           fHistDeltaEtaNPt[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,lptbins,fPtAssocArray);
471           fHistDeltaEtaNPt[i][imult][isign][imc]->Sumw2();
472           fHistDeltaEtaNPt[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#eta");
473           fHistDeltaEtaNPt[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
474           fOutput->Add(fHistDeltaEtaNPt[i][imult][isign][imc]);
475           
476           sprintf(histname,"fHistDeltaEtaNMix_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
477           sprintf(histtitle,"Near-Side #Delta#eta Mixed Event Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
478           fHistDeltaEtaNMix[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,lptbins,fPtAssocArray);
479           fHistDeltaEtaNMix[i][imult][isign][imc]->Sumw2();
480           fHistDeltaEtaNMix[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#eta");
481           fHistDeltaEtaNMix[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
482           fOutput->Add(fHistDeltaEtaNMix[i][imult][isign][imc]);
483           
484           sprintf(histname,"fHistDeltaEtaNMixPt_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
485           sprintf(histtitle,"Near-Side P_{T} Weighted #Delta#eta Mixed Event Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
486           fHistDeltaEtaNMixPt[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,lptbins,fPtAssocArray);
487           fHistDeltaEtaNMixPt[i][imult][isign][imc]->Sumw2();
488           fHistDeltaEtaNMixPt[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#eta");
489           fHistDeltaEtaNMixPt[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
490           fOutput->Add(fHistDeltaEtaNMixPt[i][imult][isign][imc]);
491           
492           //Away Eta
493           sprintf(histname,"fHistDeltaEtaA_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
494           sprintf(histtitle,"Away-Side #Delta#eta Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
495           fHistDeltaEtaA[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,lptbins,fPtAssocArray);
496           fHistDeltaEtaA[i][imult][isign][imc]->Sumw2();
497           fHistDeltaEtaA[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#eta");
498           fHistDeltaEtaA[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
499           fOutput->Add(fHistDeltaEtaA[i][imult][isign][imc]);
500           
501           sprintf(histname,"fHistDeltaEtaAPt_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
502           sprintf(histtitle,"Away-Side P_{T} Weighted #Delta#eta Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
503           fHistDeltaEtaAPt[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,lptbins,fPtAssocArray);
504           fHistDeltaEtaAPt[i][imult][isign][imc]->Sumw2();
505           fHistDeltaEtaAPt[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#eta");
506           fHistDeltaEtaAPt[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
507           fOutput->Add(fHistDeltaEtaAPt[i][imult][isign][imc]);
508           
509           sprintf(histname,"fHistDeltaEtaAMix_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
510           sprintf(histtitle,"Away-Side #Delta#eta Mixed Event Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
511           fHistDeltaEtaAMix[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,lptbins,fPtAssocArray);
512           fHistDeltaEtaAMix[i][imult][isign][imc]->Sumw2();
513           fHistDeltaEtaAMix[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#eta");
514           fHistDeltaEtaAMix[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
515           fOutput->Add(fHistDeltaEtaAMix[i][imult][isign][imc]);
516           
517           sprintf(histname,"fHistDeltaEtaAMixPt_P%d_C%d%s%s",i,imult,sign1[isign],cmc1[imc]);
518           sprintf(histtitle,"Away-Side P_{T} Weighted #Delta#eta Mixed Event Distribution with %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],sign2[isign],cmc2[imc]);
519           fHistDeltaEtaAMixPt[i][imult][isign][imc]=new TH2F(histname,histtitle,fnBinEta,etaArray,lptbins,fPtAssocArray);
520           fHistDeltaEtaAMixPt[i][imult][isign][imc]->Sumw2();
521           fHistDeltaEtaAMixPt[i][imult][isign][imc]->GetXaxis()->SetTitle("#Delta#eta");
522           fHistDeltaEtaAMixPt[i][imult][isign][imc]->GetYaxis()->SetTitle("p_{T}");
523           fOutput->Add(fHistDeltaEtaAMixPt[i][imult][isign][imc]);
524
525
526       //====
527         }//end isignloop
528       sprintf(histname,"fHistDeltaPhiEta_P%d_C%d%s",i,imult,cmc1[imc]);
529       sprintf(histtitle,"#Delta#phi-#Delta#eta %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
530       fHistDeltaPhiEta[i][imult][imc]=new TH3F(histname,histtitle,fnBinPhiEtaPhi,phiEtaArrayPhi,fnBinPhiEtaEta,phiEtaArrayEta,lptbins,fPtAssocArray);
531       fHistDeltaPhiEta[i][imult][imc]->Sumw2();
532       fHistDeltaPhiEta[i][imult][imc]->GetXaxis()->SetTitle("#phi");
533       fHistDeltaPhiEta[i][imult][imc]->GetYaxis()->SetTitle("#eta");
534       fHistDeltaPhiEta[i][imult][imc]->GetZaxis()->SetTitle("p_{T}");
535       fOutput->Add(fHistDeltaPhiEta[i][imult][imc]);
536       
537       sprintf(histname,"fHistDeltaPhiEtaMix_P%d_C%d%s",i,imult,cmc1[imc]);
538       sprintf(histtitle,"#Delta#phi-#Delta#eta from Mixed Events %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
539       fHistDeltaPhiEtaMix[i][imult][imc]=new TH3F(histname,histtitle,fnBinPhiEtaPhi,phiEtaArrayPhi,fnBinPhiEtaEta,phiEtaArrayEta,lptbins,fPtAssocArray);
540       fHistDeltaPhiEtaMix[i][imult][imc]->Sumw2();
541       fHistDeltaPhiEtaMix[i][imult][imc]->GetXaxis()->SetTitle("#phi");
542       fHistDeltaPhiEtaMix[i][imult][imc]->GetYaxis()->SetTitle("#eta");
543       fHistDeltaPhiEtaMix[i][imult][imc]->GetZaxis()->SetTitle("p_{T}");
544       fOutput->Add(fHistDeltaPhiEtaMix[i][imult][imc]);
545       
546       sprintf(histname,"fHistPhiEtaTrigPt_P%d_C%d%s",i,imult,cmc1[imc]);
547       sprintf(histtitle,"P_{T}-Weighted #phi-#eta distribution in triggered events %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
548       fHistPhiEtaTrigPt[i][imult][imc]=new TH3F(histname,histtitle,fnBinPhiEtaPhi,phiEtaArrayPhi,fnBinPhiEtaEta,phiEtaArrayEta,lptbins,fPtAssocArray);
549       fHistPhiEtaTrigPt[i][imult][imc]->Sumw2();
550       fHistPhiEtaTrigPt[i][imult][imc]->GetXaxis()->SetTitle("#phi");
551       fHistPhiEtaTrigPt[i][imult][imc]->GetYaxis()->SetTitle("#eta");
552       fHistPhiEtaTrigPt[i][imult][imc]->GetZaxis()->SetTitle("p_{T}");
553       fOutput->Add(fHistPhiEtaTrigPt[i][imult][imc]);
554     
555       sprintf(histname,"fHistDeltaPhiEtaPt_P%d_C%d%s",i,imult,cmc1[imc]);
556       sprintf(histtitle,"P_{T}-Weighted #Delta#phi-#Delta#eta %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
557       fHistDeltaPhiEtaPt[i][imult][imc]=new TH3F(histname,histtitle,fnBinPhiEtaPhi,phiEtaArrayPhi,fnBinPhiEtaEta,phiEtaArrayEta,lptbins,fPtAssocArray);
558       fHistDeltaPhiEtaPt[i][imult][imc]->Sumw2();
559       fHistDeltaPhiEtaPt[i][imult][imc]->GetXaxis()->SetTitle("#phi");
560       fHistDeltaPhiEtaPt[i][imult][imc]->GetYaxis()->SetTitle("#eta");
561       fHistDeltaPhiEtaPt[i][imult][imc]->GetZaxis()->SetTitle("p_{T}");
562       fOutput->Add(fHistDeltaPhiEtaPt[i][imult][imc]);
563       
564       sprintf(histname,"fHistDeltaPhiEtaMixPt_P%d_C%d%s",i,imult,cmc1[imc]);
565       sprintf(histtitle,"P_{T}-Weighted #Delta#phi-#Delta#eta from Mixed Events %3.1f<p_{T}^{Trig}<%3.1f %dMult%d%s",fPtTrigArray[i],fPtTrigArray[i+1],fCentArrayMin[imult],fCentArrayMax[imult],cmc2[imc]);
566       fHistDeltaPhiEtaMixPt[i][imult][imc]=new TH3F(histname,histtitle,fnBinPhiEtaPhi,phiEtaArrayPhi,fnBinPhiEtaEta,phiEtaArrayEta,lptbins,fPtAssocArray);
567       fHistDeltaPhiEtaMixPt[i][imult][imc]->Sumw2();
568       fHistDeltaPhiEtaMixPt[i][imult][imc]->GetXaxis()->SetTitle("#phi");
569       fHistDeltaPhiEtaMixPt[i][imult][imc]->GetYaxis()->SetTitle("#eta");
570       fHistDeltaPhiEtaMixPt[i][imult][imc]->GetZaxis()->SetTitle("p_{T}");
571       fOutput->Add(fHistDeltaPhiEtaMixPt[i][imult][imc]);
572
573       //Three-Particle Histograms
574       for(int ipt=0;ipt<fNAPt3Bins;ipt++){
575         for(int iSign=0;iSign<4;iSign++){
576           sprintf(histname,"fHistDeltaPhiPhi_P%dp%d_C%d%s%s",i,ipt,imult,cmc1[imc],sign31[iSign]);
577           sprintf(histtitle,"Raw #Delta#phi-#Delta#phi %3.1f<p_{T}^{Trig}<%3.1f %3.2f<p_{T}^{Assoc}<%3.2f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fPtAssoc3Array1[ipt],fPtAssoc3Array2[ipt],fCentArrayMin[imult],fCentArrayMax[imult],sign32[iSign],cmc2[imc]);
578           fHistDeltaPhiPhi[i][ipt][imult][iSign][imc]=new TH2F(histname,histtitle,fnBinPhi3,fdPhiMin,fdPhiMax,fnBinPhi3,fdPhiMin,fdPhiMax);
579           fHistDeltaPhiPhi[i][ipt][imult][iSign][imc]->Sumw2();
580           fHistDeltaPhiPhi[i][ipt][imult][iSign][imc]->GetXaxis()->SetTitle("#Delta#phi_{1}");
581           fHistDeltaPhiPhi[i][ipt][imult][iSign][imc]->GetYaxis()->SetTitle("#Delta#phi_{2}");
582           fOutput->Add(fHistDeltaPhiPhi[i][ipt][imult][iSign][imc]);
583           
584           sprintf(histname,"fHistDeltaPhiPhiMix_P%dp%d_C%d%s%s",i,ipt,imult,cmc1[imc],sign31[iSign]);
585           sprintf(histtitle,"Mixed #Delta#phi-#Delta#phi %3.1f<p_{T}^{Trig}<%3.1f %3.2f<p_{T}^{Assoc}<%3.2f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fPtAssoc3Array1[ipt],fPtAssoc3Array2[ipt],fCentArrayMin[imult],fCentArrayMax[imult],sign32[iSign],cmc2[imc]);
586           fHistDeltaPhiPhiMix[i][ipt][imult][iSign][imc]=new TH2F(histname,histtitle,fnBinPhi3,fdPhiMin,fdPhiMax,fnBinPhi3,fdPhiMin,fdPhiMax);
587           fHistDeltaPhiPhiMix[i][ipt][imult][iSign][imc]->Sumw2();
588           fHistDeltaPhiPhiMix[i][ipt][imult][iSign][imc]->GetXaxis()->SetTitle("#Delta#phi_{1}");
589           fHistDeltaPhiPhiMix[i][ipt][imult][iSign][imc]->GetYaxis()->SetTitle("#Delta#phi_{2}");
590           fOutput->Add(fHistDeltaPhiPhiMix[i][ipt][imult][iSign][imc]);
591
592           sprintf(histname,"fHistDeltaPhiPhiSS_P%dp%d_C%d%s%s",i,ipt,imult,cmc1[imc],sign31[iSign]);
593           sprintf(histtitle,"Soft-Soft #Delta#phi-#Delta#phi %3.1f<p_{T}^{Trig}<%3.1f %3.2f<p_{T}^{Assoc}<%3.2f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fPtAssoc3Array1[ipt],fPtAssoc3Array2[ipt],fCentArrayMin[imult],fCentArrayMax[imult],sign32[iSign],cmc2[imc]);
594           fHistDeltaPhiPhiSS[i][ipt][imult][iSign][imc]=new TH2F(histname,histtitle,fnBinPhi3,fdPhiMin,fdPhiMax,fnBinPhi3,fdPhiMin,fdPhiMax);
595           fHistDeltaPhiPhiSS[i][ipt][imult][iSign][imc]->Sumw2();
596           fHistDeltaPhiPhiSS[i][ipt][imult][iSign][imc]->GetXaxis()->SetTitle("#Delta#phi_{1}");
597           fHistDeltaPhiPhiSS[i][ipt][imult][iSign][imc]->GetYaxis()->SetTitle("#Delta#phi_{2}");
598           fOutput->Add(fHistDeltaPhiPhiSS[i][ipt][imult][iSign][imc]);
599
600           sprintf(histname,"fHistDeltaEtaEta_P%dp%d_C%d%s%s",i,ipt,imult,cmc1[imc],sign31[iSign]);
601           sprintf(histtitle,"Raw #Delta#eta-#Delta#eta %3.1f<p_{T}^{Trig}<%3.1f %3.2f<p_{T}^{Assoc}<%3.2f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fPtAssoc3Array1[ipt],fPtAssoc3Array2[ipt],fCentArrayMin[imult],fCentArrayMax[imult],sign32[iSign],cmc2[imc]);
602           fHistDeltaEtaEta[i][ipt][imult][iSign][imc]=new TH2F(histname,histtitle,fnBinEta3,-etaEdge,etaEdge,fnBinEta3,-etaEdge,etaEdge);
603           fHistDeltaEtaEta[i][ipt][imult][iSign][imc]->Sumw2();
604           fHistDeltaEtaEta[i][ipt][imult][iSign][imc]->GetXaxis()->SetTitle("#Delta#eta_{1}");
605           fHistDeltaEtaEta[i][ipt][imult][iSign][imc]->GetYaxis()->SetTitle("#Delta#eta_{2}");
606           fOutput->Add(fHistDeltaEtaEta[i][ipt][imult][iSign][imc]);
607
608 sprintf(histname,"fHistDeltaEtaEtaMix_P%dp%d_C%d%s%s",i,ipt,imult,cmc1[imc],sign31[iSign]);
609           sprintf(histtitle,"Mixed #Delta#eta-#Delta#eta %3.1f<p_{T}^{Trig}<%3.1f %3.2f<p_{T}^{Assoc}<%3.2f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fPtAssoc3Array1[ipt],fPtAssoc3Array2[ipt],fCentArrayMin[imult],fCentArrayMax[imult],sign32[iSign],cmc2[imc]);
610           fHistDeltaEtaEtaMix[i][ipt][imult][iSign][imc]=new TH2F(histname,histtitle,fnBinEta3,-etaEdge,etaEdge,fnBinEta3,-etaEdge,etaEdge);
611           fHistDeltaEtaEtaMix[i][ipt][imult][iSign][imc]->Sumw2();
612           fHistDeltaEtaEtaMix[i][ipt][imult][iSign][imc]->GetXaxis()->SetTitle("#Delta#eta_{1}");
613           fHistDeltaEtaEtaMix[i][ipt][imult][iSign][imc]->GetYaxis()->SetTitle("#Delta#eta_{2}");
614           fOutput->Add(fHistDeltaEtaEtaMix[i][ipt][imult][iSign][imc]);
615
616 sprintf(histname,"fHistDeltaEtaEtaSS_P%dp%d_C%d%s%s",i,ipt,imult,cmc1[imc],sign31[iSign]);
617           sprintf(histtitle,"Soft-Soft #Delta#eta-#Delta#eta %3.1f<p_{T}^{Trig}<%3.1f %3.2f<p_{T}^{Assoc}<%3.2f %dMult%d%s%s",fPtTrigArray[i],fPtTrigArray[i+1],fPtAssoc3Array1[ipt],fPtAssoc3Array2[ipt],fCentArrayMin[imult],fCentArrayMax[imult],sign32[iSign],cmc2[imc]);
618           fHistDeltaEtaEtaSS[i][ipt][imult][iSign][imc]=new TH2F(histname,histtitle,fnBinEta3,-etaEdge,etaEdge,fnBinEta3,-etaEdge,etaEdge);
619           fHistDeltaEtaEtaSS[i][ipt][imult][iSign][imc]->Sumw2();
620           fHistDeltaEtaEtaSS[i][ipt][imult][iSign][imc]->GetXaxis()->SetTitle("#Delta#eta_{1}");
621           fHistDeltaEtaEtaSS[i][ipt][imult][iSign][imc]->GetYaxis()->SetTitle("#Delta#eta_{2}");
622           fOutput->Add(fHistDeltaEtaEtaSS[i][ipt][imult][iSign][imc]);
623
624         }//iSign
625       }//associated pt (ipt)
626       }//pt loop (i) 
627     }//centrality loop (imult)
628   }//imc
629   if(fDEBUG)Printf("OutPut Created");
630 }//CreateOutputObjects    
631
632 Int_t AliAnalysisTaskDiHadron::CheckVertex(const AliESDEvent *rESD){
633   //checks whether the vertex passes cuts
634   Int_t rGood=-1;
635   Float_t vtx[3];
636   vtx[0]=rESD->GetPrimaryVertex()->GetX();
637   vtx[1]=rESD->GetPrimaryVertex()->GetY();
638   vtx[2]=rESD->GetPrimaryVertex()->GetZ();
639   if((vtx[0]*vtx[0]+vtx[1]*vtx[1])<9) rGood=0; //vertex out of beam pipe
640   if(fabs(vtx[2])<fVzCut)rGood=0;//Vertex Z cut
641   if(fDEBUG)Printf("vtxZ %f",vtx[2]);
642   for(int i=0;i<fNVertexBins;i++){
643     if(vtx[2]>fVertexArray[i]&&vtx[2]<=fVertexArray[i+1]&&rGood==0)rGood=i;
644   }
645   return rGood;
646 }
647
648 Int_t AliAnalysisTaskDiHadron::CheckTrigger(const AliESDEvent *rESD){
649   //checks whether the trigger passes cuts
650   Int_t rGood=0;
651   TString trigID=rESD->GetFiredTriggerClasses();
652   int count=0;
653   char trigID2[50];
654   int stop=0;//in as a safety
655
656   for(int i=0;i<fNIDs;i++){
657     if(stop==1)continue;
658     for(int j=0;j<50;j++){
659       if(fTrigIDArray[count]==',')trigID2[j]='\0';
660       else if(fTrigIDArray[count]=='\0'){trigID2[j]='\0';stop=1;}
661       else trigID2[j]=fTrigIDArray[count];
662       count++;
663       if(trigID2[j]=='\0') break;
664       }
665       if(trigID.Contains(trigID2)) rGood=1;
666   }
667     return rGood;
668 }
669
670 Int_t AliAnalysisTaskDiHadron::TrackCuts(const AliESDEvent *rESD, Float_t *rPt, Float_t *rEta, Float_t *rPhi, Short_t *rCharge, Float_t *rEff, Int_t **rPtAssoc3, Int_t *rNPtAssoc3, Int_t *rGoodTracks){
671     //fills arrays with all of the tracks passing cuts
672   rGoodTracks[0]=0;
673   Int_t lead=0;
674   Float_t leadPt=0;
675   Int_t rTrack=fESD->GetNumberOfTracks();
676   Float_t sPt, sEta, sPhi, sChi, sb[2], sbCov[3];
677   Int_t sNcls, sNclsF, sITScls;
678   Short_t sCharge;
679   for(int iTrack=0;iTrack<rTrack;iTrack++){
680     AliESDtrack *eSDtrack=rESD->GetTrack(iTrack);
681     const AliExternalTrackParam *conTrack = eSDtrack->GetConstrainedParam();
682     if(!conTrack)continue;
683     sPt=conTrack->Pt();
684     //if(fDEBUG)Printf("Pt%f",rPt);
685     sEta=conTrack->Eta();
686     sPhi=conTrack->Phi();
687     sCharge=conTrack->Charge();
688     if(sPhi<fdPhiMin)sPhi+=2*fPi;
689     if(sPhi>fdPhiMax)sPhi-=2*fPi;
690     if(sPt<fMinPtAssoc||sPt>fMaxPtAssoc)continue;//set Pt range
691     if(fabs(sEta)>fEtaCut)continue;//set Eta Range
692     if(!sCharge)continue;
693     sNcls=eSDtrack->GetTPCNcls();
694     //if(fDEBUG)Printf("NCLS%d",sNcls);
695     if(sNcls<fMinClustersTPC)continue;
696     sNclsF=eSDtrack->GetTPCNclsF();
697     if((1.0*sNcls/sNclsF)<fMinClusterRatio)continue;//Clusters fit/ Possible
698     sChi=(eSDtrack->GetTPCchi2())/sNcls;
699     if(sChi>fMaxTPCchi2)continue;
700     sITScls=eSDtrack->GetNcls(0);
701     if(sITScls<fMinClustersITS)continue;
702     eSDtrack->GetImpactParameters(sb,sbCov);
703     if(!fDCA2D&&(sb[0]*sb[0]+sb[1]*sb[1])>(fMaxDCA*fMaxDCA))continue;//DCA cut
704     if(fDCA2D==1&&(sb[0]*sb[0]/fMaxDCAXY/fMaxDCAXY+sb[1]*sb[1]/fMaxDCAZ/fMaxDCAZ)>1)continue;
705     if(fDCA2D==2&&(0.35+0.42*std::pow(double(sPt),-0.9))<(sb[0]*sb[0]))continue;
706     if(eSDtrack->GetKinkIndex(0)>0)continue;//removes kinked tracks
707     if(!eSDtrack->GetStatus()&AliESDtrack::kTPCrefit&&fTPCRefit)continue;//refit in TPC
708     if((fITSRefit==1||(fITSRefit==2&&sPt>5))&&!eSDtrack->GetStatus()&AliESDtrack::kITSrefit)continue;//refit of its tracks either for none,all, or >5 GeV/c
709     if(fSPDCut&&!eSDtrack->HasPointOnITSLayer(0)&&!eSDtrack->HasPointOnITSLayer(1))continue;
710     rPt[rGoodTracks[0]]=sPt;
711     rEta[rGoodTracks[0]]=sEta;
712     rPhi[rGoodTracks[0]]=sPhi;
713     rCharge[rGoodTracks[0]]=sCharge;
714     if(fEfficiencyCorr){
715     if(sPt<fEffFitPt)rEff[rGoodTracks[0]]=1./fFitLow->Eval(sPt);
716     else rEff[rGoodTracks[0]]=1./fFitHigh->Eval(sPt);
717     }
718     else rEff[rGoodTracks[0]]=1;
719     if(sPt>leadPt)lead=rGoodTracks[0];
720     //rPtAssoc3[rGoodTracks[0]]=new Int_t [10];
721     rNPtAssoc3[rGoodTracks[0]]=0;
722     for(int apt3=0;apt3<fNAPt3Bins;apt3++){
723       if(sPt<fPtAssoc3Array2[apt3]&&sPt>=fPtAssoc3Array1[apt3]){
724         rPtAssoc3[rGoodTracks[0]][rNPtAssoc3[rGoodTracks[0]]]=apt3;
725         rNPtAssoc3[rGoodTracks[0]]++;
726       }
727     }
728
729     rGoodTracks[0]++;
730     
731   }
732   return lead;
733 }
734
735 Int_t AliAnalysisTaskDiHadron::TrackCutsMC(AliMCEvent *rMC, Float_t *rPt, Float_t *rEta, Float_t *rPhi, Short_t *rCharge, Float_t *rEff, Int_t **rPtAssoc3, Int_t *rNPtAssoc3, Int_t *rGoodTracks){
736 //Fills Arrays of MC particles
737   rGoodTracks[1]=0;
738   AliStack *rStack=rMC->Stack();
739   Int_t rTrack=rStack->GetNtrack();
740   Float_t sPt, sEta, sPhi;
741   Short_t sCharge;
742   Int_t lead=0;
743   Float_t leadPt=0;
744   for(int iTrack=0;iTrack<rTrack;iTrack++){
745     TParticle *rParticle=rStack->Particle(iTrack);
746     sPt=rParticle->Pt();
747     //if(fDEBUG)Printf("MCPt%f",rPt);
748     sEta=rParticle->Eta();
749     sPhi=rParticle->Phi();
750     sCharge=rMC->GetTrack(iTrack)->Charge();
751     if(sPhi<fdPhiMin)sPhi+=2*fPi;
752     if(sPhi>fdPhiMax)sPhi-=2*fPi;
753     if(sPt<fMinPtAssoc||sPt>fMaxPtAssoc)continue;//set Pt range
754     if(fabs(sEta)>fEtaCut)continue;//set Eta Range
755     if(!rStack->IsPhysicalPrimary(iTrack))continue;//primary particles only
756     if(!sCharge)continue;//only charged particles kept
757     rPt[rGoodTracks[1]]=sPt;
758     rEta[rGoodTracks[1]]=sEta;
759     rPhi[rGoodTracks[1]]=sPhi;
760     rCharge[rGoodTracks[1]]=sCharge;
761     rEff[rGoodTracks[1]]=1;
762     if(sPt>leadPt)lead=rGoodTracks[1];
763     rNPtAssoc3[rGoodTracks[1]]=0;
764     for(int apt3=0;apt3<fNAPt3Bins;apt3++){
765       if(sPt<fPtAssoc3Array2[apt3]&&sPt>=fPtAssoc3Array1[apt3]){
766         rPtAssoc3[rGoodTracks[1]][rNPtAssoc3[rGoodTracks[1]]]=apt3;
767         rNPtAssoc3[rGoodTracks[1]]++;
768       }
769     }
770     rGoodTracks[1]++;
771   }
772   return lead;
773 }
774 //------------------------------------------------------------
775 void AliAnalysisTaskDiHadron::Exec(Option_t *)
776
777   //Main executable
778   if(fDEBUG)Printf("Exec");
779
780   const int nTPtBins=fNTPtBins;
781   const int nCentBins=fNCentBins;
782   for(int ievent=0;ievent<=1;ievent++){
783
784   if(!fESD&&ievent==0){
785     if(fDEBUG)Printf("Error: fESD not found");
786     break;
787   }
788   if(!fMC&&ievent==1){
789     break;
790   }
791   if(ievent==1&&!fMCHistos)break;//break out if MC event and we don't have fill of those set
792   //Secondary check
793   if(ievent==0){
794     if(fESD->GetNumberOfTracks()<=0){
795       if(fDEBUG)Printf("Error: no tracks");
796       break;
797     }
798   }
799   //The previous check doesn't seem to work as a fMC is bad not NULL
800   if(ievent==1){
801     if(fMC->GetNumberOfTracks()<=0){
802       if(fDEBUG)Printf("<=0 MCTracks");
803       break;
804     }
805   }
806   
807   //Check for Trigger only on real data
808   if(!fMC){
809     if(!CheckTrigger(fESD)) break;
810   }
811   //I'll only cut on the reconstructed vertex since these are the events that will be used
812   int vertexBin;
813   vertexBin=CheckVertex(fESD);
814   //else vertexBin=CheckVertex(fMC);
815   if(vertexBin<0)break;
816
817   Int_t nGoodTracks[2]={0,0}, nTriggers[nTPtBins][nCentBins][2];
818   Int_t nTrack;
819   if(!ievent)nTrack=fESD->GetNumberOfTracks();
820   else nTrack=fMC->Stack()->GetNtrack();
821  
822   Float_t tdPhi, tdEta, tXE;
823   Float_t tdPhi2, tdEta2;
824   ftPhi=new Float_t [nTrack];
825   ftEta=new Float_t [nTrack];
826   ftPt=new Float_t [nTrack];
827   ftCharge=new Short_t [nTrack];
828   ftEff=new Float_t [nTrack];
829   ftPtAssoc3=new Int_t *[nTrack];
830   for(int i=0;i<nTrack;i++){
831     ftPtAssoc3[i]=new Int_t [10];
832   }
833   ftNPtAssoc3=new Int_t [nTrack];
834   Short_t sign;
835
836   //will need to do something exta for the effieiency once it comes from embedding
837   for(int i=0;i<fNTPtBins;i++){
838     for(int c=0;c<fNCentBins;c++){
839     nTriggers[i][c][ievent]=0;
840     }
841   }
842   Int_t tMult=fESD->GetMultiplicity()->GetNumberOfTracklets();//I think this is the correct multiplicity to use
843   if(fDEBUG)Printf("Mult%d",tMult);
844   
845   //Decide what multiplicy bins are filled with this event, note set to max of 4 as I didn't think more then 2 overlapping bins likely at one time, easliy changed
846   Int_t multArray[4]={0,0,0,0};
847   Int_t maxArray=0;
848   for(int imult=0;imult<fNCentBins;imult++){
849     if(tMult>=fCentArrayMin[imult]&&tMult<fCentArrayMax[imult]){
850       multArray[maxArray]=imult;
851       maxArray++;
852     }
853   }
854   if(fDEBUG)Printf("maxArray%d",maxArray);
855   //Set Efficiency for the centrality bin (lowest bin used in array if multiple overlap)
856   for(int ipar=0;ipar<fNFitLowParam;ipar++){
857     fFitLow->SetParameter(ipar,fFitLowParam[multArray[0]*fNCentBins+ipar]);
858   }
859    for(int ipar=0;ipar<fNFitHighParam;ipar++){
860     fFitHigh->SetParameter(ipar,fFitHighParam[multArray[0]*fNCentBins+ipar]);
861   }
862   fHistMult[ievent]->Fill(tMult);
863   for(int c=0;c<maxArray;c++){fHistNEvents[multArray[c]][ievent]->Fill(0);}//count the number of events used
864   Int_t leadPart;
865   
866   //returns arrays filled up to nGoodTracks with tracks passing cuts
867   if(!ievent)leadPart=TrackCuts(fESD,ftPt,ftEta,ftPhi,ftCharge,ftEff,ftPtAssoc3,ftNPtAssoc3,nGoodTracks);
868   else leadPart=TrackCutsMC(fMC,ftPt,ftEta,ftPhi,ftCharge,ftEff,ftPtAssoc3,ftNPtAssoc3,nGoodTracks);
869   int nearEta=0,NearXE=0;
870   int nearEta2=0;
871   if(fDEBUG)Printf("Track Loop");
872   for(int iTrack=0;iTrack<nGoodTracks[ievent];iTrack++){
873     if(fDEBUG)Printf("Track%d Pt%f",iTrack,ftPt[iTrack]);
874     //if(ftPhi[iTrack]<fdPhiMin)ftPhi[iTrack]+=2*fPi;
875     //if(ftPhi[iTrack]>fdPhiMax)ftPhi[iTrack]-=2*fPi;
876     for(int c=0;c<maxArray;c++){
877       fHistPt[multArray[c]][ievent]->Fill(ftPt[iTrack],ftEff[iTrack]);
878     fHistPtEff[multArray[c]][ievent]->Fill(ftPt[iTrack]);
879     fHistPhi[multArray[c]][ievent]->Fill(ftPhi[iTrack],ftPt[iTrack],ftEff[iTrack]);
880     fHistPhiPt[multArray[c]][ievent]->Fill(ftPhi[iTrack],ftPt[iTrack],ftPt[iTrack]*ftEff[iTrack]);
881     fHistEta[multArray[c]][ievent]->Fill(ftEta[iTrack],ftPt[iTrack],ftEff[iTrack]);
882     fHistEtaPt[multArray[c]][ievent]->Fill(ftEta[iTrack],ftPt[iTrack],ftPt[iTrack]*ftEff[iTrack]);
883     fHistPhiEta[multArray[c]][ievent]->Fill(ftPhi[iTrack],ftEta[iTrack],ftPt[iTrack],ftEff[iTrack]);
884     fHistPhiEtaPt[multArray[c]][ievent]->Fill(ftPhi[iTrack],ftEta[iTrack],ftPt[iTrack],ftPt[iTrack]*ftEff[iTrack]);
885     }
886     for(int i=0;i<fNTPtBins;i++){
887       if(ftPt[iTrack]>fPtTrigArray[i]&&ftPt[iTrack]<fPtTrigArray[i+1]&&fabs(ftEta[iTrack])<fTrigEtaCut){
888         if(fDEBUG)Printf("In %fpt%f",fPtTrigArray[i],fPtTrigArray[i+1]);
889         fHistMultTrig[i][ievent]->Fill(tMult);
890         for(int c=0;c<maxArray;c++){
891           nTriggers[i][multArray[c]][ievent]++;
892           fHistNTrigger[multArray[c]][ievent]->Fill(i);
893           fHistNTriggerPt[multArray[c]][ievent]->Fill(i,ftPt[iTrack]);
894         }
895         if(fDEBUG)Printf("Assiciated Particle Loop");
896         for(int iTrack2=0;iTrack2<nGoodTracks[ievent];iTrack2++){
897           if(iTrack==iTrack2) continue;
898           if(ftPt[iTrack2]>ftPt[iTrack])continue;
899           tdPhi=ftPhi[iTrack]-ftPhi[iTrack2];
900           if(tdPhi<-fPi)tdPhi+=2*fPi;
901           if(tdPhi>fPi)tdPhi-=2*fPi;
902           if(fabs(tdPhi)<fNearPhiCut)nearEta=1;
903           else nearEta=0;
904           if(fabs(tdPhi)<fXECut)NearXE=1;
905           else NearXE=0;
906           if(fabs(tdPhi)<(fPi/2))tdEta=ftEta[iTrack]-ftEta[iTrack2];
907           else tdEta=ftEta[iTrack]+ftEta[iTrack2];
908           if(tdPhi<fdPhiMin)tdPhi+=2*fPi;
909           if(tdPhi>fdPhiMax)tdPhi-=2*fPi;
910           if((ftCharge[iTrack]<0&&ftCharge[iTrack2]<0)||(ftCharge[iTrack]>0&&ftCharge[iTrack2]>0))sign=1;
911           else sign=2;
912           if(fDEBUG) Printf("dPhi %f  dEta %f",tdPhi,tdEta);
913           for(int c=0;c<maxArray;c++){//loop over multiplicity bins
914             fHistPtTrig[i][multArray[c]][ievent]->Fill(ftPt[iTrack2],ftEff[iTrack2]);
915             fHistPhiTrig[i][multArray[c]][ievent]->Fill(ftPhi[iTrack2],ftPt[iTrack2],ftEff[iTrack2]);
916             fHistPhiTrigPt[i][multArray[c]][ievent]->Fill(ftPhi[iTrack2],ftPt[iTrack2],ftPt[iTrack2]*ftEff[iTrack2]);
917             fHistEtaTrig[i][multArray[c]][ievent]->Fill(ftEta[iTrack2],ftPt[iTrack2],ftEff[iTrack2]);
918             fHistEtaTrigPt[i][multArray[c]][ievent]->Fill(ftEta[iTrack2],ftPt[iTrack2],ftPt[iTrack2]*ftEff[iTrack2]);
919
920             fHistPhiEtaTrig[i][multArray[c]][ievent]->Fill(ftPhi[iTrack2],ftEta[iTrack2],ftPt[iTrack2],ftEff[iTrack2]);
921             fHistPhiEtaTrigPt[i][multArray[c]][ievent]->Fill(ftPhi[iTrack2],ftEta[iTrack2],ftPt[iTrack2],ftPt[iTrack2]*ftEff[iTrack2]);
922             fHistDeltaPhi[i][multArray[c]][0][ievent]->Fill(tdPhi,ftPt[iTrack2],ftEff[iTrack2]);
923             fHistDeltaPhiPt[i][multArray[c]][0][ievent]->Fill(tdPhi,ftPt[iTrack2],ftPt[iTrack2]*ftEff[iTrack2]);
924             fHistDeltaPhi[i][multArray[c]][sign][ievent]->Fill(tdPhi,ftPt[iTrack2],ftEff[iTrack2]);
925             fHistDeltaPhiPt[i][multArray[c]][sign][ievent]->Fill(tdPhi,ftPt[iTrack2],ftPt[iTrack2]*ftEff[iTrack2]);
926
927             if(nearEta){
928             fHistDeltaEtaN[i][multArray[c]][0][ievent]->Fill(tdEta,ftPt[iTrack2],ftEff[iTrack2]);
929             fHistDeltaEtaNPt[i][multArray[c]][0][ievent]->Fill(tdEta,ftPt[iTrack2],ftPt[iTrack2]*ftEff[iTrack2]);
930             fHistDeltaEtaN[i][multArray[c]][sign][ievent]->Fill(tdEta,ftPt[iTrack2],ftEff[iTrack2]);
931             fHistDeltaEtaNPt[i][multArray[c]][sign][ievent]->Fill(tdEta,ftPt[iTrack2],ftPt[iTrack2]*ftEff[iTrack2]);
932             }
933             else{
934             fHistDeltaEtaA[i][multArray[c]][0][ievent]->Fill(tdEta,ftPt[iTrack2],ftEff[iTrack2]);
935             fHistDeltaEtaAPt[i][multArray[c]][0][ievent]->Fill(tdEta,ftPt[iTrack2],ftPt[iTrack2]*ftEff[iTrack2]);
936             fHistDeltaEtaA[i][multArray[c]][sign][ievent]->Fill(tdEta,ftPt[iTrack2],ftEff[iTrack2]);
937             fHistDeltaEtaAPt[i][multArray[c]][sign][ievent]->Fill(tdEta,ftPt[iTrack2],ftPt[iTrack2]*ftEff[iTrack2]);
938             }
939             fHistDeltaPhiEta[i][multArray[c]][ievent]->Fill(tdPhi,tdEta,ftPt[iTrack2],ftEff[iTrack2]);
940             fHistDeltaPhiEtaPt[i][multArray[c]][ievent]->Fill(tdPhi,tdEta,ftPt[iTrack2],ftPt[iTrack2]*ftEff[iTrack2]);
941             
942             //only fill these if trigger particle is the leading particle
943             if(iTrack==leadPart){
944               if(NearXE){
945                 tXE=ftPt[iTrack2]*cos(tdPhi)/ftPt[iTrack];
946                 fHistXEN[i][multArray[c]][ievent]->Fill(tXE,ftEff[iTrack2]);
947               }
948               else{
949                 tXE=ftPt[iTrack2]*cos(tdPhi+fPi)/ftPt[iTrack];
950                 fHistXEA[i][multArray[c]][ievent]->Fill(tXE,ftEff[iTrack2]);
951               }
952             }
953
954           }//Centrality loop (c)
955
956           //3-particle Correlations
957           for(int iTrack3=0;iTrack3<nGoodTracks[ievent];iTrack3++){
958             if(iTrack2==iTrack3)continue;
959             if(ftPt[iTrack3]>ftPt[iTrack])continue;
960             tdPhi2=ftPhi[iTrack]-ftPhi[iTrack3];
961             if(tdPhi2<-fPi)tdPhi2+=2*fPi;
962             if(tdPhi2>fPi)tdPhi2-=2*fPi;
963             if(fabs(tdPhi2)<fNearPhiCut&&nearEta==1)nearEta2=1;
964             else nearEta2=0;
965             //if(fabs(tdPhi)<fXECut)NearXE=1;
966             //else NearXE=0;
967             if(fabs(tdPhi2)<(fPi/2))tdEta2=ftEta[iTrack]-ftEta[iTrack3];
968             else tdEta2=ftEta[iTrack]+ftEta[iTrack3];
969             if(tdPhi2<fdPhiMin)tdPhi2+=2*fPi;
970             if(tdPhi2>fdPhiMax)tdPhi2-=2*fPi;
971             // if((ftCharge[iTrack]<0&&ftCharge[iTrack2]<0)||(ftCharge[iTrack]>0&&ftCharge[iTrack2]>0))sign=1;
972             if((ftCharge[iTrack]<0&&ftCharge[iTrack2]<0&&ftCharge[iTrack3]<0)||(ftCharge[iTrack]>0&&ftCharge[iTrack2]>0&&ftCharge[iTrack3]>0))sign=1;
973             else if((ftCharge[iTrack3]<0&&ftCharge[iTrack2]<0)||(ftCharge[iTrack3]>0&&ftCharge[iTrack2]>0))sign=2;
974             else sign=3;
975             for(int e=0;e<ftNPtAssoc3[iTrack2];e++){//check associated pT bin
976               for(int f=0;f<ftNPtAssoc3[iTrack3];f++){
977                 if(ftPtAssoc3[iTrack2][e]==ftPtAssoc3[iTrack3][f]){
978                   for(int c=0;c<maxArray;c++){//loop over multiplicity bins
979                     fHistDeltaPhiPhi[i][ftPtAssoc3[iTrack2][e]][multArray[c]][0][ievent]->Fill(tdPhi,tdPhi2,ftEff[iTrack2]*ftEff[iTrack3]);
980                     fHistDeltaPhiPhi[i][ftPtAssoc3[iTrack2][e]][multArray[c]][sign][ievent]->Fill(tdPhi2,tdPhi,ftEff[iTrack2]*ftEff[iTrack3]);
981                  
982
983                     if(nearEta2){
984                       fHistDeltaEtaEta[i][ftPtAssoc3[iTrack2][e]][multArray[c]][0][ievent]->Fill(tdEta,tdEta2,ftEff[iTrack2]*ftEff[iTrack3]);
985                       fHistDeltaEtaEta[i][ftPtAssoc3[iTrack2][e]][multArray[c]][sign][ievent]->Fill(tdEta,tdEta2,ftEff[iTrack2]*ftEff[iTrack3]);
986                     }
987                   }//multiplicity loop (c)
988                 }
989               }
990             }//track checking loops
991           }//iTrack3
992         }//iTrack2 (associated track loop)
993         
994         if(fDEBUG)Printf("Mixed Event Loop");
995         for(int c=0;c<maxArray;c++){
996           int d=multArray[c];//Centrality bin we are in
997           if(fMixEnd[d][vertexBin][ievent]>=0){//check if there are any mixed events for this bin
998             for(int imix=0;imix<=fMixEnd[d][vertexBin][ievent];imix++){//loop over the stored mixed events
999               fHistNMix[d][ievent]->Fill(i);
1000               for(int iTrack2=0;iTrack2<fMixTrack[imix][d][vertexBin][ievent];iTrack2++){
1001                 if(ftPt[iTrack]<fMPt[imix][d][vertexBin][ievent][iTrack2])continue;
1002                 tdPhi=ftPhi[iTrack]-fMPhi[imix][d][vertexBin][ievent][iTrack2];
1003                 if(tdPhi<-fPi)tdPhi+=2*fPi;
1004                 if(tdPhi>fPi)tdPhi-=2*fPi;
1005                 if(fabs(tdPhi)<fNearPhiCut)nearEta=1;
1006                 else nearEta=0;
1007                 if(fabs(tdPhi)<fXECut)NearXE=1;
1008                 else NearXE=0;
1009                 if(fabs(tdPhi)<(fPi/2))tdEta=ftEta[iTrack]-fMEta[imix][d][vertexBin][ievent][iTrack2];
1010                 else tdEta=ftEta[iTrack]+fMEta[imix][d][vertexBin][ievent][iTrack2];
1011                 if(tdPhi<fdPhiMin)tdPhi+=2*fPi; 
1012                 if(tdPhi>fdPhiMax)tdPhi-=2*fPi;
1013                 if((ftCharge[iTrack]<0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]<0)||(ftCharge[iTrack]>0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]>0))sign=1;
1014                 else sign=2;
1015
1016                 fHistDeltaPhiMix[i][d][0][ievent]->Fill(tdPhi,fMPt[imix][d][vertexBin][ievent][iTrack2],fMEff[imix][d][vertexBin][ievent][iTrack2]);
1017                 fHistDeltaPhiMixPt[i][d][0][ievent]->Fill(tdPhi,fMPt[imix][d][vertexBin][ievent][iTrack2],fMPt[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack2]);
1018                 fHistDeltaPhiMix[i][d][sign][ievent]->Fill(tdPhi,fMPt[imix][d][vertexBin][ievent][iTrack2],fMEff[imix][d][vertexBin][ievent][iTrack2]);
1019                 fHistDeltaPhiMixPt[i][d][sign][ievent]->Fill(tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMPt[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack2]);
1020                 if(nearEta){
1021                   fHistDeltaEtaNMix[i][d][0][ievent]->Fill(tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMEff[imix][d][vertexBin][ievent][iTrack2]);
1022                 fHistDeltaEtaNMixPt[i][d][0][ievent]->Fill(tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMPt[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack2]);
1023                 fHistDeltaEtaNMix[i][d][sign][ievent]->Fill(tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMEff[imix][d][vertexBin][ievent][iTrack2]);
1024                 fHistDeltaEtaNMixPt[i][d][sign][ievent]->Fill(tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMPt[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack2]);
1025                 }
1026                 else{
1027                   fHistDeltaEtaAMix[i][d][0][ievent]->Fill(tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMEff[imix][d][vertexBin][ievent][iTrack2]);
1028                 fHistDeltaEtaAMixPt[i][d][0][ievent]->Fill(tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMPt[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack2]);
1029                 fHistDeltaEtaAMix[i][d][sign][ievent]->Fill(tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMEff[imix][d][vertexBin][ievent][iTrack2]);
1030                 fHistDeltaEtaAMixPt[i][d][sign][ievent]->Fill(tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMPt[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack2]);
1031                 }       
1032
1033                 fHistDeltaPhiEtaMix[i][d][ievent]->Fill(tdPhi,tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMEff[imix][d][vertexBin][ievent][iTrack2]);
1034                 fHistDeltaPhiEtaMixPt[i][d][ievent]->Fill(tdPhi,tdEta,fMPt[imix][d][vertexBin][ievent][iTrack2],fMPt[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack2]);
1035
1036                 if(iTrack==leadPart){
1037                   if(NearXE){
1038                     tXE=fMPt[imix][d][vertexBin][ievent][iTrack2]*cos(tdPhi)/ftPt[iTrack];
1039                     fHistXENMix[i][d][ievent]->Fill(tXE,fMEff[imix][d][vertexBin][ievent][iTrack2]);
1040                   }
1041                   else{
1042                     tXE=fMPt[imix][d][vertexBin][ievent][iTrack2]*cos(tdPhi+fPi)/ftPt[iTrack];
1043                     fHistXEAMix[i][multArray[c]][ievent]->Fill(tXE,fMEff[imix][d][vertexBin][ievent][iTrack2]);
1044                   }
1045                 }
1046                 //3-particle correlation soft-soft term (both associated from the same event)
1047                 for(int iTrack3=0;iTrack3<fMixTrack[imix][d][vertexBin][ievent];iTrack3++){
1048                   if(iTrack3==iTrack2)continue;
1049                   if(ftPt[iTrack]<fMPt[imix][d][vertexBin][ievent][iTrack3])continue;
1050                 tdPhi2=ftPhi[iTrack]-fMPhi[imix][d][vertexBin][ievent][iTrack3];
1051                 if(tdPhi2<-fPi)tdPhi2+=2*fPi;
1052                 if(tdPhi2>fPi)tdPhi2-=2*fPi;
1053                 if(fabs(tdPhi2)<fNearPhiCut&&nearEta)nearEta2=1;
1054                 else nearEta2=0;
1055                 if(fabs(tdPhi2)<(fPi/2))tdEta2=ftEta[iTrack]-fMEta[imix][d][vertexBin][ievent][iTrack3];
1056                 else tdEta2=ftEta[iTrack]+fMEta[imix][d][vertexBin][ievent][iTrack3];
1057                 if(tdPhi2<fdPhiMin)tdPhi2+=2*fPi;       
1058                 if(tdPhi2>fdPhiMax)tdPhi2-=2*fPi;
1059                 //if((ftCharge[iTrack]<0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]<0)||(ftCharge[iTrack]>0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]>0))sign=1;
1060                 //else sign=2;
1061                 if((ftCharge[iTrack]<0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]<0&&fMCharge[imix][d][vertexBin][ievent][iTrack3]<0)||(ftCharge[iTrack]>0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]>0&&fMCharge[imix][d][vertexBin][ievent][iTrack3]>0))sign=1;
1062                 else if((fMCharge[imix][d][vertexBin][ievent][iTrack3]<0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]<0)||(fMCharge[imix][d][vertexBin][ievent][iTrack3]>0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]>0))sign=2;
1063                 else sign=3;
1064                 for(int e=0;e<fMNPtAssoc3[imix][d][vertexBin][ievent][iTrack2];e++){//check associated pT bin
1065                   for(int f=0;f<fMNPtAssoc3[imix][d][vertexBin][ievent][iTrack3];f++){
1066                     if(fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]==fMPtAssoc3[imix][d][vertexBin][ievent][f][iTrack3]){
1067                       fHistDeltaPhiPhiSS[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][0][ievent]->Fill(tdPhi,tdPhi2,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack3]);
1068                       fHistDeltaPhiPhiSS[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][sign][ievent]->Fill(tdPhi,tdPhi2,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack3]); 
1069
1070                       if(nearEta2){
1071                         fHistDeltaEtaEtaSS[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][0][ievent]->Fill(tdEta,tdEta2,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack3]);
1072                         fHistDeltaEtaEtaSS[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][sign][ievent]->Fill(tdEta,tdEta2,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix][d][vertexBin][ievent][iTrack3]);
1073                       }//near-side
1074                     }
1075                   }
1076                 }//associated pt bin
1077                 }//iTrack3
1078
1079                 //3-particle mixed event (associated from different events)
1080                 //for(int imix2=0;imix2<=fMixEnd[d][vertexBin][ievent];imix2++){//loop over the stored mixed events
1081                 //if(imix2==imix)continue;
1082                   int imix2=imix+1;
1083                   if(imix2>=fMixEnd[d][vertexBin][ievent])imix2=0;
1084                   if(imix2==imix)continue;//will kill it when there is only 1 mixed event (remember to scale by 1 less then the number of mixed events in others (number of mixed-mixed is 2*(mixed-1)))
1085                   for(int iTrack3=0;iTrack3<fMixTrack[imix2][d][vertexBin][ievent];iTrack3++){
1086                     if(ftPt[iTrack]<fMPt[imix2][d][vertexBin][ievent][iTrack3])continue;
1087                     tdPhi2=ftPhi[iTrack]-fMPhi[imix2][d][vertexBin][ievent][iTrack3];
1088                     if(tdPhi2<-fPi)tdPhi2+=2*fPi;
1089                     if(tdPhi2>fPi)tdPhi2-=2*fPi;
1090                     if(fabs(tdPhi2)<fNearPhiCut&&nearEta)nearEta2=1;
1091                     else nearEta2=0;
1092                     if(fabs(tdPhi2)<(fPi/2))tdEta2=ftEta[iTrack]-fMEta[imix2][d][vertexBin][ievent][iTrack3];
1093                     else tdEta2=ftEta[iTrack]+fMEta[imix2][d][vertexBin][ievent][iTrack3];
1094                     if(tdPhi2<fdPhiMin)tdPhi2+=2*fPi;   
1095                     if(tdPhi2>fdPhiMax)tdPhi2-=2*fPi;
1096                     //if((ftCharge[iTrack]<0&&fMCharge[imix2][d][vertexBin][ievent][iTrack2]<0)||(ftCharge[iTrack]>0&&fMCharge[imix2][d][vertexBin][ievent][iTrack2]>0))sign=1;
1097                     //else sign=2;
1098                     if((ftCharge[iTrack]<0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]<0&&fMCharge[imix2][d][vertexBin][ievent][iTrack3]<0)||(ftCharge[iTrack]>0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]>0&&fMCharge[imix2][d][vertexBin][ievent][iTrack3]>0))sign=1;
1099                     else if((fMCharge[imix2][d][vertexBin][ievent][iTrack3]<0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]<0)||(fMCharge[imix2][d][vertexBin][ievent][iTrack3]>0&&fMCharge[imix][d][vertexBin][ievent][iTrack2]>0))sign=2;
1100                     else sign=3;
1101                     for(int e=0;e<fMNPtAssoc3[imix][d][vertexBin][ievent][iTrack2];e++){//check associated pT bin
1102                       for(int f=0;f<fMNPtAssoc3[imix2][d][vertexBin][ievent][iTrack3];f++){
1103                         if(fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]==fMPtAssoc3[imix2][d][vertexBin][ievent][f][iTrack3]){
1104                           fHistDeltaPhiPhiMix[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][0][ievent]->Fill(tdPhi,tdPhi2,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix2][d][vertexBin][ievent][iTrack3]);
1105                           fHistDeltaPhiPhiMix[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][sign][ievent]->Fill(tdPhi,tdPhi2,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix2][d][vertexBin][ievent][iTrack3]); 
1106                           fHistDeltaPhiPhiMix[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][0][ievent]->Fill(tdPhi2,tdPhi,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix2][d][vertexBin][ievent][iTrack3]);//free factor of 2 in statistics
1107                           fHistDeltaPhiPhiMix[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][sign][ievent]->Fill(tdPhi2,tdPhi,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix2][d][vertexBin][ievent][iTrack3]); 
1108                       if(nearEta2){
1109                         fHistDeltaEtaEtaMix[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][0][ievent]->Fill(tdEta,tdEta2,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix2][d][vertexBin][ievent][iTrack3]);
1110                         fHistDeltaEtaEtaMix[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][sign][ievent]->Fill(tdEta,tdEta2,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix2][d][vertexBin][ievent][iTrack3]);
1111                         fHistDeltaEtaEtaMix[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][0][ievent]->Fill(tdEta2,tdEta,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix2][d][vertexBin][ievent][iTrack3]);
1112                         fHistDeltaEtaEtaMix[i][fMPtAssoc3[imix][d][vertexBin][ievent][e][iTrack2]][d][sign][ievent]->Fill(tdEta2,tdEta,fMEff[imix][d][vertexBin][ievent][iTrack2]*fMEff[imix2][d][vertexBin][ievent][iTrack3]);
1113                       }//near-side
1114                         }
1115                       }
1116                     }//associated pt bin
1117                   }//iTrack3
1118
1119               }//iTrack2
1120               }//imix
1121             }//fMixEnd
1122           }//Centrality bins (c)
1123         }//pt trig cuts
1124       }//i Pt Trig
1125     }//itrack    
1126   
1127   //now store this event for mixing (using these dynamic arrays instead of the other fix to save memory)
1128     if(fDEBUG)Printf("Store Event For Mixing");
1129   for(int c=0;c<maxArray;c++){//loops over centrality bins
1130     int d=multArray[c];//too many nested arrays looked confusing d=which centrality bin
1131     if(fMixEnd[d][vertexBin][ievent]<(fNMix-1))fMixEnd[d][vertexBin][ievent]++;
1132     if(fMixPointer[d][vertexBin][ievent]<(fNMix-1))fMixPointer[d][vertexBin][ievent]++;
1133     else fMixPointer[d][vertexBin][ievent]=0;
1134     int e=fMixPointer[d][vertexBin][ievent];//nested arrays (e is event number in pool)
1135     delete [] fMPt[e][d][vertexBin][ievent];
1136     delete [] fMPhi[e][d][vertexBin][ievent];
1137     delete [] fMEta[e][d][vertexBin][ievent];
1138     delete [] fMCharge[e][d][vertexBin][ievent];
1139     delete [] fMEff[e][d][vertexBin][ievent];
1140     delete [] fMNPtAssoc3[e][d][vertexBin][ievent];
1141     for(int jj=0;jj<10;jj++){
1142     delete [] fMPtAssoc3[e][d][vertexBin][ievent][jj];
1143     }
1144     fMPt[e][d][vertexBin][ievent]=new Float_t [nGoodTracks[ievent]];
1145     fMPhi[e][d][vertexBin][ievent]=new Float_t [nGoodTracks[ievent]];
1146     fMEta[e][d][vertexBin][ievent]=new Float_t [nGoodTracks[ievent]];
1147     fMCharge[e][d][vertexBin][ievent]=new Short_t [nGoodTracks[ievent]];
1148     fMEff[e][d][vertexBin][ievent]=new Float_t [nGoodTracks[ievent]];
1149     fMixTrack[e][d][vertexBin][ievent]=nGoodTracks[ievent];
1150     fMNPtAssoc3[e][d][vertexBin][ievent]=new Short_t [nGoodTracks[ievent]];
1151     for(int jj=0;jj<10;jj++){
1152     fMPtAssoc3[e][d][vertexBin][ievent][jj]=new Short_t [nGoodTracks[ievent]];
1153     }
1154
1155     for(int iTrack=0;iTrack<nGoodTracks[ievent];iTrack++){
1156       fMPt[e][d][vertexBin][ievent][iTrack]=ftPt[iTrack];
1157       fMPhi[e][d][vertexBin][ievent][iTrack]=ftPhi[iTrack];
1158       fMEta[e][d][vertexBin][ievent][iTrack]=ftEta[iTrack];
1159       fMCharge[e][d][vertexBin][ievent][iTrack]=ftCharge[iTrack];
1160       fMEff[e][d][vertexBin][ievent][iTrack]=ftEff[iTrack];
1161       fMNPtAssoc3[e][d][vertexBin][ievent][iTrack]=ftNPtAssoc3[iTrack];
1162       // fMPtAssoc3[e][d][vertexBin][ievent][iTrack]=new Int_t [ftNPtAssoc3[iTrack]];
1163       for(int jj=0;jj<ftNPtAssoc3[iTrack];jj++){
1164         //if(fDEBUG) Printf("%d",ftPtAssoc3[iTrack][jj]);
1165         fMPtAssoc3[e][d][vertexBin][ievent][jj][iTrack]=ftPtAssoc3[iTrack][jj];
1166       }
1167       
1168     }//iTracks
1169   }//Centrality (c)
1170   }//ievent
1171   //track=0;
1172   //track2=0;
1173
1174   PostData(0, fOutput);
1175   //get rid of these arrays from memory
1176   delete [] ftPhi;
1177   delete [] ftEta;
1178   delete [] ftPt;
1179   delete [] ftCharge;
1180   delete [] ftEff;
1181   delete [] ftNPtAssoc3;
1182   delete [] ftPtAssoc3;
1183   ftPhi=NULL;
1184   ftEta=NULL;
1185   ftPt=NULL;
1186   ftCharge=NULL;
1187   ftEff=NULL;
1188   ftNPtAssoc3=NULL;
1189   ftPtAssoc3=NULL;
1190
1191 }//Exec
1192
1193 //---------------------------------------------------
1194 void AliAnalysisTaskDiHadron::Terminate(Option_t *){
1195   //Terminates the code, frees up memory
1196   for(int ii=0;ii<fNMix;ii++){
1197     for(int cc=0;cc<fNCentBins;cc++){
1198       for(int vtx=0;vtx<fNVertexBins;vtx++){
1199         for(int jj=0;jj<2;jj++){
1200           delete [] fMPt[ii][cc][vtx][jj];
1201           delete [] fMPhi[ii][cc][vtx][jj];
1202           delete [] fMEta[ii][cc][vtx][jj];
1203           delete [] fMCharge[ii][cc][vtx][jj];
1204           delete [] fMEff[ii][cc][vtx][jj];
1205           delete [] fMNPtAssoc3[ii][cc][vtx][jj];
1206           for(int qq=0;qq<10;qq++){
1207             delete [] fMPtAssoc3[ii][cc][vtx][jj][qq];
1208             fMPtAssoc3[ii][cc][vtx][jj][qq]=NULL;
1209           }
1210           fMPt[ii][cc][vtx][jj]=NULL;
1211           fMPhi[ii][cc][vtx][jj]=NULL;
1212           fMEta[ii][cc][vtx][jj]=NULL;
1213           fMCharge[ii][cc][vtx][jj]=NULL;
1214           fMEff[ii][cc][vtx][jj]=NULL;
1215           fMNPtAssoc3[ii][cc][vtx][jj]=NULL;
1216
1217         }
1218       }
1219     }
1220   }
1221   delete [] fFitLowParam;
1222   delete [] fFitHighParam;
1223   delete [] fPtTrigArray;
1224   delete [] fPtAssocArray;
1225   delete [] fPtAssoc3Array1;
1226   delete [] fPtAssoc3Array2;
1227   delete [] fCentArrayMin;
1228   delete [] fCentArrayMax;
1229   delete [] fXEArray;
1230   fFitLowParam=NULL;
1231   fFitHighParam=NULL;
1232   fPtTrigArray=NULL;
1233   fPtAssocArray=NULL;
1234   fPtAssoc3Array1=NULL;
1235   fPtAssoc3Array2=NULL;
1236   fCentArrayMin=NULL;
1237   fCentArrayMax=NULL;
1238   Printf("Terminate AliAnalysisTaskDiHadron");
1239 }