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
  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),DEBUG(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),fNFitLow(0),fNFitHigh(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),tPhi(NULL),tEta(NULL),tPt(NULL),tCharge(NULL),tEff(NULL),tPtAssoc3(NULL),tNPtAssoc3(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   // tPhi=NULL;
  85   //  tEta=NULL;
  86   // tPt=NULL;
  87
  88   //fPi=3.1415926535898;
  89
  90
  91   }
  92 //--------------------------------------
  93 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, char **TrigIDArray){
  94   fMinClustersTPC=MinClustersTPC;
  95   fMinClusterRatio=MinClusterRatio;
  96   fMaxTPCchi2=MaxTPCchi2;
  97   fMinClustersITS=MinClustersITS;
  98   fEtaCut=EtaCut;
  99   fTrigEtaCut=TrigEtaCut;
 100   fNearPhiCut=NearPhiCut;
 101   fXECut=XECut;
 102   fMaxDCA=MaxDCA;
 103   fMaxDCAXY=MaxDCAXY;
 104   fMaxDCAZ=MaxDCAZ;
 105   fDCA2D=DCA2D;
 106   fTPCRefit=TPCRefit;
 107   fITSRefit=ITSRefit;
 108   fSPDCut=SPDCut;
 109   fMinPtAssoc=MinPtAssoc;
 110   fMaxPtAssoc=MaxPtAssoc;
 111   fVzCut=VzCut;
 112   fNIDs=NIDs;
 113    fTrigIDArray=(char**)TrigIDArray;
 114
 115    //  Printf("TPC%d R%2.1f Chi%2.1f ITS%d Eta%2.1f Near%2.1f NearX%2.1f DCA%2.1f DCAR%2.1f DCAZ%2.1f DCAM%d TPCFit%d ITSFit%d SPD%d MinPt%2.1f MaxPt%2.1f Vz%2.1f Trig_%s",fMinClustersTPC,fMinClusterRatio,fMaxTPCchi2,fMinClustersITS,fEtaCut,fNearPhiCut,fXECut,fMaxDCA,fMaxDCAXY,fMaxDCAZ,fDCA2D,fTPCRefit,fITSRefit,fSPDCut,fMinPtAssoc,fMaxPtAssoc,fVzCut,fTrigIDArray[0]);
 116 }
 117 //--------------------------------------------------------
 118 void AliAnalysisTaskDiHadron::SetOptions(Int_t EfficiencyCorr, Int_t fDEBUG,  Int_t MCHistos){
 119   fEfficiencyCorr=EfficiencyCorr;
 120   DEBUG=fDEBUG;
 121   fMCHistos=MCHistos;
 122   //Printf("Eff%d DEBUG%d MCHistos%d",fEfficiencyCorr,DEBUG,fMCHistos);
 123 }
 124 //------------------------------------------------------
 125 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){
 126   fnBinPhi=nBinPhi;
 127   fnBinEta=nBinEta;
 128   fnBinPhiEtaPhi=nBinPhiEtaPhi;
 129   fnBinPhiEtaEta=nBinPhiEtaEta;
 130   fnBinPhi3=nBinPhi3;
 131   fnBinEta3=nBinEta3;
 132   fdPhiMin=dPhiMin;
 133   fdPhiMax=dPhiMax;
 134   fNTPtBins=NTPtBins;
 135   fNMix=NMixBins;
 136   fNCentBins=NCentBins;
 137   fNAPtBins=NAPtBins;
 138   fNAPt3Bins=NAPt3Bins;
 139   fNVertexBins=NVertexBins;
 140   fNXEBins=NXEBins;
 141   fPtTrigArray=(Float_t*)PtTrigArray;
 142   fPtAssocArray=(Float_t*)PtAssocArray;
 143   fPtAssoc3Array1=(Float_t*)PtAssoc3Array1;
 144   fPtAssoc3Array2=(Float_t*)PtAssoc3Array2;
 145   fCentArrayMin=(Int_t*)CentArrayMin;
 146   fCentArrayMax=(Int_t*)CentArrayMax;
 147   fXEArray=(Float_t*)XEArray;
 148  for(int i=0;i<=fNVertexBins;i++)fVertexArray[i]=(2.*i/fNVertexBins-1)*fVzCut;
 149  //Printf("BinPhi%d BinEta%d BinPhiEta%d %d Bin3Phi%d Bin3Eta%d VertexArray%2.1f %2.1f",fnBinPhi,fnBinEta,fnBinPhiEtaPhi,fnBinPhiEtaEta,fnBinPhi3,fnBinEta3,fVertexArray[0],fVertexArray[fNVertexBins]);
 150 }
 151 //-------------------------------------------------------
 152 void AliAnalysisTaskDiHadron::SetEfficiencies(Float_t EffFitPt, TF1 *FitLow, TF1 *FitHigh, Int_t NFitLowParam, Int_t NFitHighParam, Float_t *FitLowParam, Float_t *FitHighParam){
 153   fEffFitPt=EffFitPt;
 154   fFitLow=(TF1*)FitLow;
 155   fFitHigh=(TF1*)FitHigh;
 156   fNFitLowParam=NFitLowParam;
 157   fNFitHighParam=NFitHighParam;
 158   // Printf("NCent%d NFitLow%d NFitHigh%d",NCent,fNFitLowParam,fNFitHighParam);
 159   //  Printf("FitLowParm%2.4f",FitLowParam[0]);
 160   fFitLowParam=(Float_t*)FitLowParam;
 161   fFitHighParam=(Float_t*)FitHighParam;
 162   //Printf("fEffFitPt%2.1f FitParam%2.1f",fEffFitPt,fFitLowParam[0]);
 163 }
 164
 165 //-----------------------------------------------------------
 166 void AliAnalysisTaskDiHadron::ConnectInputData(Option_t *){
 167   //Connect to ESD
 168   if(DEBUG)Printf("Connecting");
 169    TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
 170    if (!tree&&DEBUG) {Printf("ERROR: Could not read chain from input slot 0");}
 171   else {
 172     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
 173     if (!esdH&&DEBUG) 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(DEBUG)Printf("Connected");
 181 }
 182
 183 //---------------------------------------------------------
 184 void AliAnalysisTaskDiHadron::CreateOutputObjects(){
 185   if(DEBUG)Printf("Output");
 186   //Creates the histograms and list
 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   char *cmc1[2]={"","_MC"};
 194   char *cmc2[2]={""," MC"};
 195   char *sign1[3]={"","_LS","_ULS"};
 196   char *sign2[3]={""," Like-Sign"," Unlike-Sign"};
 197   char *sign31[4]={"","_LS","_ULT","_ULA"};
 198   char *sign32[4]={""," Like-Sign"," Trigger-Diff"," Assoc-Diff"};
 199   Float_t EtaEdge=2*fEtaCut;
 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(DEBUG)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(fPtTrigArray[i]==fPtAssocArray[j])lptbins=j;
 335         }
 336         //if(DEBUG)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(DEBUG)Printf("OutPut Created");
 630 }//CreateOutputObjects
 631
 632 Int_t AliAnalysisTaskDiHadron::CheckVertex(AliESDEvent *rESD){
 633   Int_t rGood=-1;
 634   Float_t Vtx[3];
 635   Vtx[0]=rESD->GetPrimaryVertex()->GetX();
 636   Vtx[1]=rESD->GetPrimaryVertex()->GetY();
 637   Vtx[2]=rESD->GetPrimaryVertex()->GetZ();
 638   if((Vtx[0]*Vtx[0]+Vtx[1]*Vtx[1])<9) rGood=0; //vertex out of beam pipe
 639   if(fabs(Vtx[2])<fVzCut)rGood=0;//Vertex Z cut
 640   if(DEBUG)Printf("VtxZ %f",Vtx[2]);
 641   for(int i=0;i<fNVertexBins;i++){
 642     if(Vtx[2]>fVertexArray[i]&&Vtx[2]<=fVertexArray[i+1]&&rGood==0)rGood=i;
 643   }
 644   return rGood;
 645 }
 646
 647 Int_t AliAnalysisTaskDiHadron::CheckTrigger(AliESDEvent *rESD){
 648   Int_t rGood=0;
 649   TString trigID=rESD->GetFiredTriggerClasses();
 650   for(int i=0;i<fNIDs;i++){
 651     if(trigID.Contains(fTrigIDArray[i])) rGood=1;
 652   }
 653     return rGood;
 654 }
 655
 656 Int_t AliAnalysisTaskDiHadron::TrackCuts(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){
 657     //fills arrays with all of the tracks passing cuts
 658   rGoodTracks[0]=0;
 659   Int_t Lead=0;
 660   Float_t LeadPt=0;
 661   Int_t rTrack=fESD->GetNumberOfTracks();
 662   Float_t sPt, sEta, sPhi, sChi, sb[2], sbCov[3];
 663   Int_t sNcls, sNclsF, sITScls;
 664   Short_t sCharge;
 665   for(int iTrack=0;iTrack<rTrack;iTrack++){
 666     AliESDtrack *ESDtrack=rESD->GetTrack(iTrack);
 667     const AliExternalTrackParam *ConTrack = ESDtrack->GetConstrainedParam();
 668     if(!ConTrack)continue;
 669     sPt=ConTrack->Pt();
 670     //if(DEBUG)Printf("Pt%f",rPt);
 671     sEta=ConTrack->Eta();
 672     sPhi=ConTrack->Phi();
 673     sCharge=ConTrack->Charge();
 674     if(sPhi<fdPhiMin)sPhi+=2*fPi;
 675     if(sPhi>fdPhiMax)sPhi-=2*fPi;
 676     if(sPt<fMinPtAssoc||sPt>fMaxPtAssoc)continue;//set Pt range
 677     if(fabs(sEta)>fEtaCut)continue;//set Eta Range
 678     if(!sCharge)continue;
 679     sNcls=ESDtrack->GetTPCNcls();
 680     //if(DEBUG)Printf("NCLS%d",sNcls);
 681     if(sNcls<fMinClustersTPC)continue;
 682     sNclsF=ESDtrack->GetTPCNclsF();
 683     if((1.0*sNcls/sNclsF)<fMinClusterRatio)continue;//Clusters fit/ Possible
 684     sChi=(ESDtrack->GetTPCchi2())/sNcls;
 685     if(sChi>fMaxTPCchi2)continue;
 686     sITScls=ESDtrack->GetNcls(0);
 687     if(sITScls<fMinClustersITS)continue;
 688     ESDtrack->GetImpactParameters(sb,sbCov);
 689     if(!fDCA2D&&(sb[0]*sb[0]+sb[1]*sb[1])>(fMaxDCA*fMaxDCA))continue;//DCA cut
 690     if(fDCA2D==1&&(sb[0]*sb[0]/fMaxDCAXY/fMaxDCAXY+sb[1]*sb[1]/fMaxDCAZ/fMaxDCAZ)>1)continue;
 691     if(fDCA2D==2&&(0.35+0.42*std::pow(double(sPt),-0.9))<(sb[0]*sb[0]))continue;
 692     if(ESDtrack->GetKinkIndex(0)>0)continue;//removes kinked tracks
 693     if(!ESDtrack->GetStatus()&AliESDtrack::kTPCrefit&&fTPCRefit)continue;//refit in TPC
 694     if((fITSRefit==1||(fITSRefit==2&&sPt>5))&&!ESDtrack->GetStatus()&AliESDtrack::kITSrefit)continue;//refit of its tracks either for none,all, or >5 GeV/c
 695     if(fSPDCut&&!ESDtrack->HasPointOnITSLayer(0)&&!ESDtrack->HasPointOnITSLayer(1))continue;
 696     rPt[rGoodTracks[0]]=sPt;
 697     rEta[rGoodTracks[0]]=sEta;
 698     rPhi[rGoodTracks[0]]=sPhi;
 699     rCharge[rGoodTracks[0]]=sCharge;
 700     if(fEfficiencyCorr){
 701     if(sPt<fEffFitPt)rEff[rGoodTracks[0]]=1./fFitLow->Eval(sPt);
 702     else rEff[rGoodTracks[0]]=1./fFitHigh->Eval(sPt);
 703     }
 704     else rEff[rGoodTracks[0]]=1;
 705     if(sPt>LeadPt)Lead=rGoodTracks[0];
 706     //rPtAssoc3[rGoodTracks[0]]=new Int_t [10];
 707     rNPtAssoc3[rGoodTracks[0]]=0;
 708     for(int apt3=0;apt3<fNAPt3Bins;apt3++){
 709       if(sPt<fPtAssoc3Array2[apt3]&&sPt>=fPtAssoc3Array1[apt3]){
 710         rPtAssoc3[rGoodTracks[0]][rNPtAssoc3[rGoodTracks[0]]]=apt3;
 711         rNPtAssoc3[rGoodTracks[0]]++;
 712       }
 713     }
 714
 715     rGoodTracks[0]++;
 716
 717   }
 718   return Lead;
 719 }
 720
 721 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){//Fills Arrays of MC particles
 722   rGoodTracks[1]=0;
 723   AliStack *rStack=rMC->Stack();
 724   Int_t rTrack=rStack->GetNtrack();
 725   Float_t sPt, sEta, sPhi;
 726   Short_t sCharge;
 727   Int_t Lead=0;
 728   Float_t LeadPt=0;
 729   for(int iTrack=0;iTrack<rTrack;iTrack++){
 730     TParticle *rParticle=rStack->Particle(iTrack);
 731     sPt=rParticle->Pt();
 732     //if(DEBUG)Printf("MCPt%f",rPt);
 733     sEta=rParticle->Eta();
 734     sPhi=rParticle->Phi();
 735     sCharge=rMC->GetTrack(iTrack)->Charge();
 736     if(sPhi<fdPhiMin)sPhi+=2*fPi;
 737     if(sPhi>fdPhiMax)sPhi-=2*fPi;
 738     if(sPt<fMinPtAssoc||sPt>fMaxPtAssoc)continue;//set Pt range
 739     if(fabs(sEta)>fEtaCut)continue;//set Eta Range
 740     if(!rStack->IsPhysicalPrimary(iTrack))continue;//primary particles only
 741     if(!sCharge)continue;//only charged particles kept
 742     rPt[rGoodTracks[1]]=sPt;
 743     rEta[rGoodTracks[1]]=sEta;
 744     rPhi[rGoodTracks[1]]=sPhi;
 745     rCharge[rGoodTracks[1]]=sCharge;
 746     rEff[rGoodTracks[1]]=1;
 747     if(sPt>LeadPt)Lead=rGoodTracks[1];
 748     rNPtAssoc3[rGoodTracks[1]]=0;
 749     for(int apt3=0;apt3<fNAPt3Bins;apt3++){
 750       if(sPt<fPtAssoc3Array2[apt3]&&sPt>=fPtAssoc3Array1[apt3]){
 751         rPtAssoc3[rGoodTracks[1]][rNPtAssoc3[rGoodTracks[1]]]=apt3;
 752         rNPtAssoc3[rGoodTracks[1]]++;
 753       }
 754     }
 755     rGoodTracks[1]++;
 756   }
 757   return Lead;
 758 }
 759 //------------------------------------------------------------
 760 void AliAnalysisTaskDiHadron::Exec(Option_t *)
 761 {
 762   if(DEBUG)Printf("Exec");
 763
 764   const int NTPtBins=fNTPtBins;
 765   const int NCentBins=fNCentBins;
 766   for(int ievent=0;ievent<=1;ievent++){
 767
 768   if(!fESD&&ievent==0){
 769     if(DEBUG)Printf("Error: fESD not found");
 770     break;
 771   }
 772   if(!fMC&&ievent==1){
 773     break;
 774   }
 775   if(ievent==1&&!fMCHistos)break;//break out if MC event and we don't have fill of those set
 776   //Secondary check
 777   if(ievent==0){
 778     if(fESD->GetNumberOfTracks()<=0){
 779       if(DEBUG)Printf("Error: no tracks");
 780       break;
 781     }
 782   }
 783   //The previous check doesn't seem to work as a fMC is bad not NULL
 784   if(ievent==1){
 785     if(fMC->GetNumberOfTracks()<=0){
 786       if(DEBUG)Printf("<=0 MCTracks");
 787       break;
 788     }
 789   }
 790
 791   //Check for Trigger only on real data
 792   if(!fMC){
 793     if(!CheckTrigger(fESD)) break;
 794   }
 795   //I'll only cut on the reconstructed vertex since these are the events that will be used
 796   int VertexBin;
 797   VertexBin=CheckVertex(fESD);
 798   //else VertexBin=CheckVertex(fMC);
 799   if(VertexBin<0)break;
 800
 801   Int_t nGoodTracks[2]={0,0}, nTriggers[NTPtBins][NCentBins][2];
 802   Int_t nTrack;
 803   if(!ievent)nTrack=fESD->GetNumberOfTracks();
 804   else nTrack=fMC->Stack()->GetNtrack();
 805
 806   Float_t tdPhi, tdEta, tXE;
 807   Float_t tdPhi2, tdEta2;
 808   tPhi=new Float_t [nTrack];
 809   tEta=new Float_t [nTrack];
 810   tPt=new Float_t [nTrack];
 811   tCharge=new Short_t [nTrack];
 812   tEff=new Float_t [nTrack];
 813   tPtAssoc3=new Int_t *[nTrack];
 814   for(int i=0;i<nTrack;i++){
 815     tPtAssoc3[i]=new Int_t [10];
 816   }
 817   tNPtAssoc3=new Int_t [nTrack];
 818   Short_t Sign;
 819
 820   //will need to do something exta for the effieiency once it comes from embedding
 821   for(int i=0;i<fNTPtBins;i++){
 822     for(int c=0;c<fNCentBins;c++){
 823     nTriggers[i][c][ievent]=0;
 824     }
 825   }
 826   Int_t tMult=fESD->GetMultiplicity()->GetNumberOfTracklets();//I think this is the correct multiplicity to use
 827   if(DEBUG)Printf("Mult%d",tMult);
 828
 829   //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
 830   Int_t MultArray[4]={0,0,0,0};
 831   Int_t MaxArray=0;
 832   for(int imult=0;imult<fNCentBins;imult++){
 833     if(tMult>=fCentArrayMin[imult]&&tMult<fCentArrayMax[imult]){
 834       MultArray[MaxArray]=imult;
 835       MaxArray++;
 836     }
 837   }
 838   if(DEBUG)Printf("MaxArray%d",MaxArray);
 839   //Set Efficiency for the centrality bin (lowest bin used in array if multiple overlap)
 840   for(int ipar=0;ipar<fNFitLowParam;ipar++){
 841     fFitLow->SetParameter(ipar,fFitLowParam[MultArray[0]*fNCentBins+ipar]);
 842   }
 843    for(int ipar=0;ipar<fNFitHighParam;ipar++){
 844     fFitHigh->SetParameter(ipar,fFitHighParam[MultArray[0]*fNCentBins+ipar]);
 845   }
 846   fHistMult[ievent]->Fill(tMult);
 847   for(int c=0;c<MaxArray;c++){fHistNEvents[MultArray[c]][ievent]->Fill(0);}//count the number of events used
 848   Int_t LeadPart;
 849
 850   //returns arrays filled up to nGoodTracks with tracks passing cuts
 851   if(!ievent)LeadPart=TrackCuts(fESD,tPt,tEta,tPhi,tCharge,tEff,tPtAssoc3,tNPtAssoc3,nGoodTracks);
 852   else LeadPart=TrackCutsMC(fMC,tPt,tEta,tPhi,tCharge,tEff,tPtAssoc3,tNPtAssoc3,nGoodTracks);
 853   int NearEta=0,NearXE=0;
 854   int NearEta2=0;
 855   if(DEBUG)Printf("Track Loop");
 856   for(int iTrack=0;iTrack<nGoodTracks[ievent];iTrack++){
 857     if(DEBUG)Printf("Track%d Pt%f",iTrack,tPt[iTrack]);
 858     //if(tPhi[iTrack]<fdPhiMin)tPhi[iTrack]+=2*fPi;
 859     //if(tPhi[iTrack]>fdPhiMax)tPhi[iTrack]-=2*fPi;
 860     for(int c=0;c<MaxArray;c++){
 861       fHistPt[MultArray[c]][ievent]->Fill(tPt[iTrack],tEff[iTrack]);
 862     fHistPtEff[MultArray[c]][ievent]->Fill(tPt[iTrack]);
 863     fHistPhi[MultArray[c]][ievent]->Fill(tPhi[iTrack],tPt[iTrack],tEff[iTrack]);
 864     fHistPhiPt[MultArray[c]][ievent]->Fill(tPhi[iTrack],tPt[iTrack],tPt[iTrack]*tEff[iTrack]);
 865     fHistEta[MultArray[c]][ievent]->Fill(tEta[iTrack],tPt[iTrack],tEff[iTrack]);
 866     fHistEtaPt[MultArray[c]][ievent]->Fill(tEta[iTrack],tPt[iTrack],tPt[iTrack]*tEff[iTrack]);
 867     fHistPhiEta[MultArray[c]][ievent]->Fill(tPhi[iTrack],tEta[iTrack],tPt[iTrack],tEff[iTrack]);
 868     fHistPhiEtaPt[MultArray[c]][ievent]->Fill(tPhi[iTrack],tEta[iTrack],tPt[iTrack],tPt[iTrack]*tEff[iTrack]);
 869     }
 870     for(int i=0;i<fNTPtBins;i++){
 871       if(tPt[iTrack]>fPtTrigArray[i]&&tPt[iTrack]<fPtTrigArray[i+1]&&fabs(tEta[iTrack])<fTrigEtaCut){
 872         if(DEBUG)Printf("In %fpt%f",fPtTrigArray[i],fPtTrigArray[i+1]);
 873         fHistMultTrig[i][ievent]->Fill(tMult);
 874         for(int c=0;c<MaxArray;c++){
 875           nTriggers[i][MultArray[c]][ievent]++;
 876           fHistNTrigger[MultArray[c]][ievent]->Fill(i);
 877           fHistNTriggerPt[MultArray[c]][ievent]->Fill(i,tPt[iTrack]);
 878         }
 879         if(DEBUG)Printf("Assiciated Particle Loop");
 880         for(int iTrack2=0;iTrack2<nGoodTracks[ievent];iTrack2++){
 881           if(iTrack==iTrack2) continue;
 882           if(tPt[iTrack2]>tPt[iTrack])continue;
 883           tdPhi=tPhi[iTrack]-tPhi[iTrack2];
 884           if(tdPhi<-fPi)tdPhi+=2*fPi;
 885           if(tdPhi>fPi)tdPhi-=2*fPi;
 886           if(fabs(tdPhi)<fNearPhiCut)NearEta=1;
 887           else NearEta=0;
 888           if(fabs(tdPhi)<fXECut)NearXE=1;
 889           else NearXE=0;
 890           if(fabs(tdPhi)<(fPi/2))tdEta=tEta[iTrack]-tEta[iTrack2];
 891           else tdEta=tEta[iTrack]+tEta[iTrack2];
 892           if(tdPhi<fdPhiMin)tdPhi+=2*fPi;
 893           if(tdPhi>fdPhiMax)tdPhi-=2*fPi;
 894           if((tCharge[iTrack]<0&&tCharge[iTrack2]<0)||(tCharge[iTrack]>0&&tCharge[iTrack2]>0))Sign=1;
 895           else Sign=2;
 896           if(DEBUG) Printf("dPhi %f  dEta %f",tdPhi,tdEta);
 897           for(int c=0;c<MaxArray;c++){//loop over multiplicity bins
 898             fHistPtTrig[i][MultArray[c]][ievent]->Fill(tPt[iTrack2],tEff[iTrack2]);
 899             fHistPhiTrig[i][MultArray[c]][ievent]->Fill(tPhi[iTrack2],tPt[iTrack2],tEff[iTrack2]);
 900             fHistPhiTrigPt[i][MultArray[c]][ievent]->Fill(tPhi[iTrack2],tPt[iTrack2],tPt[iTrack2]*tEff[iTrack2]);
 901             fHistEtaTrig[i][MultArray[c]][ievent]->Fill(tEta[iTrack2],tPt[iTrack2],tEff[iTrack2]);
 902             fHistEtaTrigPt[i][MultArray[c]][ievent]->Fill(tEta[iTrack2],tPt[iTrack2],tPt[iTrack2]*tEff[iTrack2]);
 903
 904             fHistPhiEtaTrig[i][MultArray[c]][ievent]->Fill(tPhi[iTrack2],tEta[iTrack2],tPt[iTrack2],tEff[iTrack2]);
 905             fHistPhiEtaTrigPt[i][MultArray[c]][ievent]->Fill(tPhi[iTrack2],tEta[iTrack2],tPt[iTrack2],tPt[iTrack2]*tEff[iTrack2]);
 906             fHistDeltaPhi[i][MultArray[c]][0][ievent]->Fill(tdPhi,tPt[iTrack2],tEff[iTrack2]);
 907             fHistDeltaPhiPt[i][MultArray[c]][0][ievent]->Fill(tdPhi,tPt[iTrack2],tPt[iTrack2]*tEff[iTrack2]);
 908             fHistDeltaPhi[i][MultArray[c]][Sign][ievent]->Fill(tdPhi,tPt[iTrack2],tEff[iTrack2]);
 909             fHistDeltaPhiPt[i][MultArray[c]][Sign][ievent]->Fill(tdPhi,tPt[iTrack2],tPt[iTrack2]*tEff[iTrack2]);
 910
 911             if(NearEta){
 912             fHistDeltaEtaN[i][MultArray[c]][0][ievent]->Fill(tdEta,tPt[iTrack2],tEff[iTrack2]);
 913             fHistDeltaEtaNPt[i][MultArray[c]][0][ievent]->Fill(tdEta,tPt[iTrack2],tPt[iTrack2]*tEff[iTrack2]);
 914             fHistDeltaEtaN[i][MultArray[c]][Sign][ievent]->Fill(tdEta,tPt[iTrack2],tEff[iTrack2]);
 915             fHistDeltaEtaNPt[i][MultArray[c]][Sign][ievent]->Fill(tdEta,tPt[iTrack2],tPt[iTrack2]*tEff[iTrack2]);
 916             }
 917             else{
 918             fHistDeltaEtaA[i][MultArray[c]][0][ievent]->Fill(tdEta,tPt[iTrack2],tEff[iTrack2]);
 919             fHistDeltaEtaAPt[i][MultArray[c]][0][ievent]->Fill(tdEta,tPt[iTrack2],tPt[iTrack2]*tEff[iTrack2]);
 920             fHistDeltaEtaA[i][MultArray[c]][Sign][ievent]->Fill(tdEta,tPt[iTrack2],tEff[iTrack2]);
 921             fHistDeltaEtaAPt[i][MultArray[c]][Sign][ievent]->Fill(tdEta,tPt[iTrack2],tPt[iTrack2]*tEff[iTrack2]);
 922             }
 923             fHistDeltaPhiEta[i][MultArray[c]][ievent]->Fill(tdPhi,tdEta,tPt[iTrack2],tEff[iTrack2]);
 924             fHistDeltaPhiEtaPt[i][MultArray[c]][ievent]->Fill(tdPhi,tdEta,tPt[iTrack2],tPt[iTrack2]*tEff[iTrack2]);
 925
 926             //only fill these if trigger particle is the leading particle
 927             if(iTrack==LeadPart){
 928               if(NearXE){
 929                 tXE=tPt[iTrack2]*cos(tdPhi)/tPt[iTrack];
 930                 fHistXEN[i][MultArray[c]][ievent]->Fill(tXE,tEff[iTrack2]);
 931               }
 932               else{
 933                 tXE=tPt[iTrack2]*cos(tdPhi+fPi)/tPt[iTrack];
 934                 fHistXEA[i][MultArray[c]][ievent]->Fill(tXE,tEff[iTrack2]);
 935               }
 936             }
 937
 938           }//Centrality loop (c)
 939
 940           //3-particle Correlations
 941           for(int iTrack3=0;iTrack3<nGoodTracks[ievent];iTrack3++){
 942             if(iTrack2==iTrack3)continue;
 943             if(tPt[iTrack3]>tPt[iTrack])continue;
 944             tdPhi2=tPhi[iTrack]-tPhi[iTrack3];
 945             if(tdPhi2<-fPi)tdPhi2+=2*fPi;
 946             if(tdPhi2>fPi)tdPhi2-=2*fPi;
 947             if(fabs(tdPhi2)<fNearPhiCut&&NearEta==1)NearEta2=1;
 948             else NearEta2=0;
 949             //if(fabs(tdPhi)<fXECut)NearXE=1;
 950             //else NearXE=0;
 951             if(fabs(tdPhi2)<(fPi/2))tdEta2=tEta[iTrack]-tEta[iTrack3];
 952             else tdEta2=tEta[iTrack]+tEta[iTrack3];
 953             if(tdPhi2<fdPhiMin)tdPhi2+=2*fPi;
 954             if(tdPhi2>fdPhiMax)tdPhi2-=2*fPi;
 955             // if((tCharge[iTrack]<0&&tCharge[iTrack2]<0)||(tCharge[iTrack]>0&&tCharge[iTrack2]>0))Sign=1;
 956             if((tCharge[iTrack]<0&&tCharge[iTrack2]<0&&tCharge[iTrack3]<0)||(tCharge[iTrack]>0&&tCharge[iTrack2]>0&&tCharge[iTrack3]>0))Sign=1;
 957             else if((tCharge[iTrack3]<0&&tCharge[iTrack2]<0)||(tCharge[iTrack3]>0&&tCharge[iTrack2]>0))Sign=2;
 958             else Sign=3;
 959             for(int e=0;e<tNPtAssoc3[iTrack2];e++){//check associated pT bin
 960               for(int f=0;f<tNPtAssoc3[iTrack3];f++){
 961                 if(tPtAssoc3[iTrack2][e]==tPtAssoc3[iTrack3][f]){
 962                   for(int c=0;c<MaxArray;c++){//loop over multiplicity bins
 963                     fHistDeltaPhiPhi[i][tPtAssoc3[iTrack2][e]][MultArray[c]][0][ievent]->Fill(tdPhi,tdPhi2,tEff[iTrack2]*tEff[iTrack3]);
 964                     fHistDeltaPhiPhi[i][tPtAssoc3[iTrack2][e]][MultArray[c]][Sign][ievent]->Fill(tdPhi2,tdPhi,tEff[iTrack2]*tEff[iTrack3]);
 965
 966
 967                     if(NearEta2){
 968                       fHistDeltaEtaEta[i][tPtAssoc3[iTrack2][e]][MultArray[c]][0][ievent]->Fill(tdEta,tdEta2,tEff[iTrack2]*tEff[iTrack3]);
 969                       fHistDeltaEtaEta[i][tPtAssoc3[iTrack2][e]][MultArray[c]][Sign][ievent]->Fill(tdEta,tdEta2,tEff[iTrack2]*tEff[iTrack3]);
 970                     }
 971                   }//multiplicity loop (c)
 972                 }
 973               }
 974             }//track checking loops
 975           }//iTrack3
 976         }//iTrack2 (associated track loop)
 977
 978         if(DEBUG)Printf("Mixed Event Loop");
 979         for(int c=0;c<MaxArray;c++){
 980           int d=MultArray[c];//Centrality bin we are in
 981           if(fMixEnd[d][VertexBin][ievent]>=0){//check if there are any mixed events for this bin
 982             for(int imix=0;imix<=fMixEnd[d][VertexBin][ievent];imix++){//loop over the stored mixed events
 983               fHistNMix[d][ievent]->Fill(i);
 984               for(int iTrack2=0;iTrack2<fMixTrack[imix][d][VertexBin][ievent];iTrack2++){
 985                 if(tPt[iTrack]<fMPt[imix][d][VertexBin][ievent][iTrack2])continue;
 986                 tdPhi=tPhi[iTrack]-fMPhi[imix][d][VertexBin][ievent][iTrack2];
 987                 if(tdPhi<-fPi)tdPhi+=2*fPi;
 988                 if(tdPhi>fPi)tdPhi-=2*fPi;
 989                 if(fabs(tdPhi)<fNearPhiCut)NearEta=1;
 990                 else NearEta=0;
 991                 if(fabs(tdPhi)<fXECut)NearXE=1;
 992                 else NearXE=0;
 993                 if(fabs(tdPhi)<(fPi/2))tdEta=tEta[iTrack]-fMEta[imix][d][VertexBin][ievent][iTrack2];
 994                 else tdEta=tEta[iTrack]+fMEta[imix][d][VertexBin][ievent][iTrack2];
 995                 if(tdPhi<fdPhiMin)tdPhi+=2*fPi;
 996                 if(tdPhi>fdPhiMax)tdPhi-=2*fPi;
 997                 if((tCharge[iTrack]<0&&fMCharge[imix][d][VertexBin][ievent][iTrack2]<0)||(tCharge[iTrack]>0&&fMCharge[imix][d][VertexBin][ievent][iTrack2]>0))Sign=1;
 998                 else Sign=2;
 999
1000                 fHistDeltaPhiMix[i][d][0][ievent]->Fill(tdPhi,fMPt[imix][d][VertexBin][ievent][iTrack2],fMEff[imix][d][VertexBin][ievent][iTrack2]);
1001                 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]);
1002                 fHistDeltaPhiMix[i][d][Sign][ievent]->Fill(tdPhi,fMPt[imix][d][VertexBin][ievent][iTrack2],fMEff[imix][d][VertexBin][ievent][iTrack2]);
1003                 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]);
1004                 if(NearEta){
1005                   fHistDeltaEtaNMix[i][d][0][ievent]->Fill(tdEta,fMPt[imix][d][VertexBin][ievent][iTrack2],fMEff[imix][d][VertexBin][ievent][iTrack2]);
1006                 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]);
1007                 fHistDeltaEtaNMix[i][d][Sign][ievent]->Fill(tdEta,fMPt[imix][d][VertexBin][ievent][iTrack2],fMEff[imix][d][VertexBin][ievent][iTrack2]);
1008                 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]);
1009                 }
1010                 else{
1011                   fHistDeltaEtaAMix[i][d][0][ievent]->Fill(tdEta,fMPt[imix][d][VertexBin][ievent][iTrack2],fMEff[imix][d][VertexBin][ievent][iTrack2]);
1012                 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]);
1013                 fHistDeltaEtaAMix[i][d][Sign][ievent]->Fill(tdEta,fMPt[imix][d][VertexBin][ievent][iTrack2],fMEff[imix][d][VertexBin][ievent][iTrack2]);
1014                 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]);
1015                 }
1016
1017                 fHistDeltaPhiEtaMix[i][d][ievent]->Fill(tdPhi,tdEta,fMPt[imix][d][VertexBin][ievent][iTrack2],fMEff[imix][d][VertexBin][ievent][iTrack2]);
1018                 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]);
1019
1020                 if(iTrack==LeadPart){
1021                   if(NearXE){
1022                     tXE=fMPt[imix][d][VertexBin][ievent][iTrack2]*cos(tdPhi)/tPt[iTrack];
1023                     fHistXENMix[i][d][ievent]->Fill(tXE,fMEff[imix][d][VertexBin][ievent][iTrack2]);
1024                   }
1025                   else{
1026                     tXE=fMPt[imix][d][VertexBin][ievent][iTrack2]*cos(tdPhi+fPi)/tPt[iTrack];
1027                     fHistXEAMix[i][MultArray[c]][ievent]->Fill(tXE,fMEff[imix][d][VertexBin][ievent][iTrack2]);
1028                   }
1029                 }
1030                 //3-particle correlation soft-soft term (both associated from the same event)
1031                 for(int iTrack3=0;iTrack3<fMixTrack[imix][d][VertexBin][ievent];iTrack3++){
1032                   if(iTrack3==iTrack2)continue;
1033                   if(tPt[iTrack]<fMPt[imix][d][VertexBin][ievent][iTrack3])continue;
1034                 tdPhi2=tPhi[iTrack]-fMPhi[imix][d][VertexBin][ievent][iTrack3];
1035                 if(tdPhi2<-fPi)tdPhi2+=2*fPi;
1036                 if(tdPhi2>fPi)tdPhi2-=2*fPi;
1037                 if(fabs(tdPhi2)<fNearPhiCut&&NearEta)NearEta2=1;
1038                 else NearEta2=0;
1039                 if(fabs(tdPhi2)<(fPi/2))tdEta2=tEta[iTrack]-fMEta[imix][d][VertexBin][ievent][iTrack3];
1040                 else tdEta2=tEta[iTrack]+fMEta[imix][d][VertexBin][ievent][iTrack3];
1041                 if(tdPhi2<fdPhiMin)tdPhi2+=2*fPi;
1042                 if(tdPhi2>fdPhiMax)tdPhi2-=2*fPi;
1043                 //if((tCharge[iTrack]<0&&fMCharge[imix][d][VertexBin][ievent][iTrack2]<0)||(tCharge[iTrack]>0&&fMCharge[imix][d][VertexBin][ievent][iTrack2]>0))Sign=1;
1044                 //else Sign=2;
1045                 if((tCharge[iTrack]<0&&fMCharge[imix][d][VertexBin][ievent][iTrack2]<0&&fMCharge[imix][d][VertexBin][ievent][iTrack3]<0)||(tCharge[iTrack]>0&&fMCharge[imix][d][VertexBin][ievent][iTrack2]>0&&fMCharge[imix][d][VertexBin][ievent][iTrack3]>0))Sign=1;
1046                 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;
1047                 else Sign=3;
1048                 for(int e=0;e<fMNPtAssoc3[imix][d][VertexBin][ievent][iTrack2];e++){//check associated pT bin
1049                   for(int f=0;f<fMNPtAssoc3[imix][d][VertexBin][ievent][iTrack3];f++){
1050                     if(fMPtAssoc3[imix][d][VertexBin][ievent][e][iTrack2]==fMPtAssoc3[imix][d][VertexBin][ievent][f][iTrack3]){
1051                       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]);
1052                       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]);
1053
1054                       if(NearEta2){
1055                         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]);
1056                         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]);
1057                       }//near-side
1058                     }
1059                   }
1060                 }//associated pt bin
1061                 }//iTrack3
1062
1063                 //3-particle mixed event (associated from different events)
1064                 //for(int imix2=0;imix2<=fMixEnd[d][VertexBin][ievent];imix2++){//loop over the stored mixed events
1065                 //if(imix2==imix)continue;
1066                   int imix2=imix+1;
1067                   if(imix2>=fMixEnd[d][VertexBin][ievent])imix2=0;
1068                   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)))
1069                   for(int iTrack3=0;iTrack3<fMixTrack[imix2][d][VertexBin][ievent];iTrack3++){
1070                     if(tPt[iTrack]<fMPt[imix2][d][VertexBin][ievent][iTrack3])continue;
1071                     tdPhi2=tPhi[iTrack]-fMPhi[imix2][d][VertexBin][ievent][iTrack3];
1072                     if(tdPhi2<-fPi)tdPhi2+=2*fPi;
1073                     if(tdPhi2>fPi)tdPhi2-=2*fPi;
1074                     if(fabs(tdPhi2)<fNearPhiCut&&NearEta)NearEta2=1;
1075                     else NearEta2=0;
1076                     if(fabs(tdPhi2)<(fPi/2))tdEta2=tEta[iTrack]-fMEta[imix2][d][VertexBin][ievent][iTrack3];
1077                     else tdEta2=tEta[iTrack]+fMEta[imix2][d][VertexBin][ievent][iTrack3];
1078                     if(tdPhi2<fdPhiMin)tdPhi2+=2*fPi;
1079                     if(tdPhi2>fdPhiMax)tdPhi2-=2*fPi;
1080                     //if((tCharge[iTrack]<0&&fMCharge[imix2][d][VertexBin][ievent][iTrack2]<0)||(tCharge[iTrack]>0&&fMCharge[imix2][d][VertexBin][ievent][iTrack2]>0))Sign=1;
1081                     //else Sign=2;
1082                     if((tCharge[iTrack]<0&&fMCharge[imix][d][VertexBin][ievent][iTrack2]<0&&fMCharge[imix2][d][VertexBin][ievent][iTrack3]<0)||(tCharge[iTrack]>0&&fMCharge[imix][d][VertexBin][ievent][iTrack2]>0&&fMCharge[imix2][d][VertexBin][ievent][iTrack3]>0))Sign=1;
1083                     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;
1084                     else Sign=3;
1085                     for(int e=0;e<fMNPtAssoc3[imix][d][VertexBin][ievent][iTrack2];e++){//check associated pT bin
1086                       for(int f=0;f<fMNPtAssoc3[imix2][d][VertexBin][ievent][iTrack3];f++){
1087                         if(fMPtAssoc3[imix][d][VertexBin][ievent][e][iTrack2]==fMPtAssoc3[imix2][d][VertexBin][ievent][f][iTrack3]){
1088                           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]);
1089                           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]);
1090                           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
1091                           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]);
1092                       if(NearEta2){
1093                         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]);
1094                         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]);
1095                         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]);
1096                         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]);
1097                       }//near-side
1098                         }
1099                       }
1100                     }//associated pt bin
1101                   }//iTrack3
1102
1103               }//iTrack2
1104               }//imix
1105             }//fMixEnd
1106           }//Centrality bins (c)
1107         }//pt trig cuts
1108       }//i Pt Trig
1109     }//itrack
1110
1111   //now store this event for mixing (using these dynamic arrays instead of the other fix to save memory)
1112     if(DEBUG)Printf("Store Event For Mixing");
1113   for(int c=0;c<MaxArray;c++){//loops over centrality bins
1114     int d=MultArray[c];//too many nested arrays looked confusing d=which centrality bin
1115     if(fMixEnd[d][VertexBin][ievent]<(fNMix-1))fMixEnd[d][VertexBin][ievent]++;
1116     if(fMixPointer[d][VertexBin][ievent]<(fNMix-1))fMixPointer[d][VertexBin][ievent]++;
1117     else fMixPointer[d][VertexBin][ievent]=0;
1118     int e=fMixPointer[d][VertexBin][ievent];//nested arrays (e is event number in pool)
1119     delete [] fMPt[e][d][VertexBin][ievent];
1120     delete [] fMPhi[e][d][VertexBin][ievent];
1121     delete [] fMEta[e][d][VertexBin][ievent];
1122     delete [] fMCharge[e][d][VertexBin][ievent];
1123     delete [] fMEff[e][d][VertexBin][ievent];
1124     delete [] fMNPtAssoc3[e][d][VertexBin][ievent];
1125     for(int jj=0;jj<10;jj++){
1126     delete [] fMPtAssoc3[e][d][VertexBin][ievent][jj];
1127     }
1128     fMPt[e][d][VertexBin][ievent]=new Float_t [nGoodTracks[ievent]];
1129     fMPhi[e][d][VertexBin][ievent]=new Float_t [nGoodTracks[ievent]];
1130     fMEta[e][d][VertexBin][ievent]=new Float_t [nGoodTracks[ievent]];
1131     fMCharge[e][d][VertexBin][ievent]=new Short_t [nGoodTracks[ievent]];
1132     fMEff[e][d][VertexBin][ievent]=new Float_t [nGoodTracks[ievent]];
1133     fMixTrack[e][d][VertexBin][ievent]=nGoodTracks[ievent];
1134     fMNPtAssoc3[e][d][VertexBin][ievent]=new Short_t [nGoodTracks[ievent]];
1135     for(int jj=0;jj<10;jj++){
1136     fMPtAssoc3[e][d][VertexBin][ievent][jj]=new Short_t [nGoodTracks[ievent]];
1137     }
1138
1139     for(int iTrack=0;iTrack<nGoodTracks[ievent];iTrack++){
1140       fMPt[e][d][VertexBin][ievent][iTrack]=tPt[iTrack];
1141       fMPhi[e][d][VertexBin][ievent][iTrack]=tPhi[iTrack];
1142       fMEta[e][d][VertexBin][ievent][iTrack]=tEta[iTrack];
1143       fMCharge[e][d][VertexBin][ievent][iTrack]=tCharge[iTrack];
1144       fMEff[e][d][VertexBin][ievent][iTrack]=tEff[iTrack];
1145       fMNPtAssoc3[e][d][VertexBin][ievent][iTrack]=tNPtAssoc3[iTrack];
1146       // fMPtAssoc3[e][d][VertexBin][ievent][iTrack]=new Int_t [tNPtAssoc3[iTrack]];
1147       for(int jj=0;jj<tNPtAssoc3[iTrack];jj++){
1148         //if(DEBUG) Printf("%d",tPtAssoc3[iTrack][jj]);
1149         fMPtAssoc3[e][d][VertexBin][ievent][jj][iTrack]=tPtAssoc3[iTrack][jj];
1150       }
1151
1152     }//iTracks
1153   }//Centrality (c)
1154   }//ievent
1155   //track=0;
1156   //track2=0;
1157
1158   PostData(0, fOutput);
1159   //get rid of these arrays from memory
1160   delete [] tPhi;
1161   delete [] tEta;
1162   delete [] tPt;
1163   delete [] tCharge;
1164   delete [] tEff;
1165   delete [] tNPtAssoc3;
1166   delete [] tPtAssoc3;
1167   tPhi=NULL;
1168   tEta=NULL;
1169   tPt=NULL;
1170   tCharge=NULL;
1171   tEff=NULL;
1172   tNPtAssoc3=NULL;
1173   tPtAssoc3=NULL;
1174
1175 }//Exec
1176
1177 //---------------------------------------------------
1178 void AliAnalysisTaskDiHadron::Terminate(Option_t *){
1179   for(int ii=0;ii<fNMix;ii++){
1180     for(int cc=0;cc<fNCentBins;cc++){
1181       for(int vtx=0;vtx<fNVertexBins;vtx++){
1182         for(int jj=0;jj<2;jj++){
1183           delete [] fMPt[ii][cc][vtx][jj];
1184           delete [] fMPhi[ii][cc][vtx][jj];
1185           delete [] fMEta[ii][cc][vtx][jj];
1186           delete [] fMCharge[ii][cc][vtx][jj];
1187           delete [] fMEff[ii][cc][vtx][jj];
1188           delete [] fMNPtAssoc3[ii][cc][vtx][jj];
1189           for(int qq=0;qq<10;qq++){
1190             delete [] fMPtAssoc3[ii][cc][vtx][jj][qq];
1191             fMPtAssoc3[ii][cc][vtx][jj][qq]=NULL;
1192           }
1193           fMPt[ii][cc][vtx][jj]=NULL;
1194           fMPhi[ii][cc][vtx][jj]=NULL;
1195           fMEta[ii][cc][vtx][jj]=NULL;
1196           fMCharge[ii][cc][vtx][jj]=NULL;
1197           fMEff[ii][cc][vtx][jj]=NULL;
1198           fMNPtAssoc3[ii][cc][vtx][jj]=NULL;
1199
1200         }
1201       }
1202     }
1203   }
1204   Printf("Terminate AliAnalysisTaskDiHadron");
1205 }