[u/mrichter/AliRoot.git] / PWG4 / JetTasks / AliPWG4CosmicCandidates.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

//-----------------------------------------------------------------------
// Analysis task looking for cosmic candidates embedded in a p-p event
// Task checks if particles are back-to-back in eta and phi
//             
//
// Author : Marta Verweij - UU - marta.verweij@cern.ch
//-----------------------------------------------------------------------

#include "TVector3.h"
#include <iostream>
#include "TH1.h"
#include "TH2.h"
#include "TH3.h"
#include "TList.h"
#include "TChain.h"

#include "AliAnalysisTask.h"
#include "AliAnalysisManager.h"

#include "AliESDEvent.h"
#include "AliESDInputHandler.h"

#include "AliESDtrack.h"
#include "AliESDtrackCuts.h"
#include "AliExternalTrackParam.h"

#include "AliLog.h"

#include "AliPWG4CosmicCandidates.h"

//using namespace std; //required for resolving the 'cout' symbol
using namespace std;

ClassImp(AliPWG4CosmicCandidates)

//________________________________________________________________________
AliPWG4CosmicCandidates::AliPWG4CosmicCandidates()
: AliAnalysisTaskSE(),
  fTrackCuts(0), 
  fPtMin(5.),
  fMaxCosmicAngle(0.002),
  fNEventAll(0),
  fNEventSel(0),
  fPtSignedCosmicCandidates(0),
  fDeltaPtCosmicCandidates(0),
  fDeltaPhiSumEta(0),
  fDCAZCosmicCandidates(0),
  fDCARCosmicCandidates(0),
  fTheta(0),
  fThetaZoom(0),
  fThetaPt1Pt2(0),
  fThetaPt1Pt2Signed(0),
  fDeltaPhiSumEtaPt1(0),
  fDeltaPhiSumEtaPt2(0),
  fThetaDCAZ1DCAZ2(0),
  fRisol(0),
  fRisolTheta(0),
  fHistListCosmics(0)
{
  //
  // Default constructor
  //
}

//________________________________________________________________________
AliPWG4CosmicCandidates::AliPWG4CosmicCandidates(const char *name)
  : AliAnalysisTaskSE(name),
    fTrackCuts(0), 
    fPtMin(5.),
    fMaxCosmicAngle(0.002),
    fNEventAll(0),
    fNEventSel(0),
    fPtSignedCosmicCandidates(0),
    fDeltaPtCosmicCandidates(0),
    fDeltaPhiSumEta(0),
    fDCAZCosmicCandidates(0),
    fDCARCosmicCandidates(0),
    fTheta(0),
    fThetaZoom(0),
    fThetaPt1Pt2(0),
    fThetaPt1Pt2Signed(0),
    fDeltaPhiSumEtaPt1(0),
    fDeltaPhiSumEtaPt2(0),
    fThetaDCAZ1DCAZ2(0),
    fRisol(0),
    fRisolTheta(0),
    fHistListCosmics(0)
{
  // Constructor. Initialization of Inputs and Outputs

  // Define input and output slots here
  // Input slot #0 works with a TChain
  //  DefineInput(0, TChain::Class());
  // Output slot #0 id reserved by the base class for AOD
  // Output slot #1 writes into a TList
  DefineOutput(1, TList::Class());
  // Output slot #2 writes into a AliESDtrackCuts
  DefineOutput(2, AliESDtrackCuts::Class());
}
//________________________________________________________________________
AliPWG4CosmicCandidates::AliPWG4CosmicCandidates(const AliPWG4CosmicCandidates &res)
  : AliAnalysisTaskSE(res),
    fTrackCuts(0), 
    fPtMin(5.),
    fMaxCosmicAngle(0.002),
    fNEventAll(0),
    fNEventSel(0),
    fPtSignedCosmicCandidates(0),
    fDeltaPtCosmicCandidates(0),
    fDeltaPhiSumEta(0),
    fDCAZCosmicCandidates(0),
    fDCARCosmicCandidates(0),
    fTheta(0),
    fThetaZoom(0),
    fThetaPt1Pt2(0),
    fThetaPt1Pt2Signed(0),
    fDeltaPhiSumEtaPt1(0),
    fDeltaPhiSumEtaPt2(0),
    fThetaDCAZ1DCAZ2(0),
    fRisol(0),
    fRisolTheta(0),
    fHistListCosmics(0)
{
  //
  // Dummy copy constructor
  //
}

//________________________________________________________________________
AliPWG4CosmicCandidates& AliPWG4CosmicCandidates::operator=(const AliPWG4CosmicCandidates& /*trclass*/)
{
  //
  // Dummy assignment operator
  //
  return *this;
}

//________________________________________________________________________
void AliPWG4CosmicCandidates::LocalInit()
{
  //
  // Only called once at beginning
  //
  PostData(2,fTrackCuts);
}

//________________________________________________________________________
void AliPWG4CosmicCandidates::UserCreateOutputObjects()
{
  //
  // Create output objects
  // Called once
  AliDebug(2,Form(">> AliPWG4CosmicCandidates::UserCreateOutputObjects \n")); 

  Bool_t oldStatus = TH1::AddDirectoryStatus();
  TH1::AddDirectory(kFALSE); 
  
  OpenFile(1);
  fHistListCosmics = new TList();

  fNEventAll = new TH1F("fNEventAll","NEventAll",1,-0.5,0.5);
  fHistListCosmics->Add(fNEventAll);
  fNEventSel = new TH1F("fNEventSel","NEvent Selected for analysis",1,-0.5,0.5);
  fHistListCosmics->Add(fNEventSel);

  Float_t fgkPtMin=0.;
  Float_t fgkPtMax=100.;
  Int_t fgkNPtBins= (int)(fgkPtMax-fgkPtMin);

  Int_t fgkNPhiBins=18;
  Float_t kMinPhi = -0.5*TMath::Pi();
  Float_t kMaxPhi = 3./2.*TMath::Pi();

  Int_t fgkNThetaBins=fgkNPhiBins*8;
  Float_t kMinTheta = -0.5*TMath::Pi();
  Float_t kMaxTheta = 3./2.*TMath::Pi();

  Int_t fgkNDCARBins=80;
  Float_t fgkDCARMin = -0.2;
  Float_t fgkDCARMax = 0.2;
  Float_t *binsDCAR=new Float_t[fgkNDCARBins+1];
  for(Int_t i=0; i<=fgkNDCARBins; i++) binsDCAR[i]=(Float_t)fgkDCARMin + (fgkDCARMax-fgkDCARMin)/fgkNDCARBins*(Float_t)i ;

  Int_t fgkNDCAZBins=80;
  Float_t fgkDCAZMin = -2.;
  Float_t fgkDCAZMax = 2.;
  Float_t *binsDCAZ=new Float_t[fgkNDCAZBins+1];
  for(Int_t i=0; i<=fgkNDCAZBins; i++) binsDCAZ[i]=(Float_t)fgkDCAZMin + (fgkDCAZMax-fgkDCAZMin)/fgkNDCAZBins*(Float_t)i ;

  fPtSignedCosmicCandidates = new TH1F("fPtSignedCosmicCandidates","fPtSignedCosmicCandidates",2*(int)(fgkPtMax-fgkPtMin), -1.*fgkPtMax, fgkPtMax);
  fHistListCosmics->Add(fPtSignedCosmicCandidates);  

  fDeltaPtCosmicCandidates = new TH1F("fDeltaPtCosmicCandidates","fDeltaPtCosmicCandidates",fgkNPtBins, -50., 50.);
  fHistListCosmics->Add(fDeltaPtCosmicCandidates);  

  fDeltaPhiSumEta = new TH2F("fDeltaPhiSumEta","fDeltaPhiSumEta",fgkNPhiBins*4,kMinPhi,kMaxPhi,80, -2.,2.);
  fHistListCosmics->Add(fDeltaPhiSumEta);  

  fDCAZCosmicCandidates = new TH2F("fDCAZCosmicCandidates","fDCAZCosmicCandidates",fgkNDCAZBins,binsDCAZ,fgkNDCAZBins,binsDCAZ);
  fHistListCosmics->Add(fDCAZCosmicCandidates);

  fDCARCosmicCandidates = new TH2F("fDCARCosmicCandidates","fDCARCosmicCandidates",fgkNDCARBins,binsDCAR,fgkNDCARBins,binsDCAR);
  fHistListCosmics->Add(fDCARCosmicCandidates);

  fTheta = new TH1F("fTheta","fTheta",fgkNThetaBins,kMinTheta,kMaxTheta);
  fHistListCosmics->Add(fTheta);

  fThetaZoom = new TH1F("fThetaZoom","fThetaZoom",100,TMath::Pi()-1.,TMath::Pi()+1.);
  fHistListCosmics->Add(fThetaZoom);

  fThetaPt1Pt2 = new TH3F("fThetaPt1Pt2","fThetaPt1Pt2",fgkNThetaBins,kMinTheta,kMaxTheta,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
  fHistListCosmics->Add(fThetaPt1Pt2);

  fThetaPt1Pt2Signed = new TH3F("fThetaPt1Pt2Signed","fThetaPt1Pt2Signed",fgkNThetaBins,kMinTheta,kMaxTheta,4*(int)(fgkPtMax-fgkPtMin),-1.*fgkPtMax,fgkPtMax,4*(int)(fgkPtMax-fgkPtMin),-1.*fgkPtMax,fgkPtMax);
  fHistListCosmics->Add(fThetaPt1Pt2Signed);

  fDeltaPhiSumEtaPt1 = new TH3F("fDeltaPhiSumEtaPt1","fDeltaPhiSumEtaPt1",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
  fHistListCosmics->Add(fDeltaPhiSumEtaPt1);

  fDeltaPhiSumEtaPt2 = new TH3F("fDeltaPhiSumEtaPt2","fDeltaPhiSumEtaPt2",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
  fHistListCosmics->Add(fDeltaPhiSumEtaPt2);

  fThetaDCAZ1DCAZ2 = new TH3F("fThetaDCAZ1DCAZ2","fThetaDCAZ1DCAZ2",fgkNThetaBins,kMinTheta,kMaxTheta,fgkNDCAZBins,-2.,2.,fgkNDCAZBins,-2.,2.);
  fHistListCosmics->Add(fThetaDCAZ1DCAZ2);

  fRisol = new TH1F("fRisol","fRisol",100,0.,10.);
  fHistListCosmics->Add(fRisol);

  fRisolTheta = new TH2F("fRisolTheta","fRisolTheta",100,0.,10.,fgkNThetaBins,kMinTheta,kMaxTheta);
  fHistListCosmics->Add(fRisolTheta);

  TH1::AddDirectory(oldStatus);   

  if(binsDCAR) delete [] binsDCAR;
  if(binsDCAZ) delete [] binsDCAZ;

}

//________________________________________________________________________
void AliPWG4CosmicCandidates::UserExec(Option_t *) 
{
  // Main loop
  // Called for each event
  //

  // All events without selection
  fNEventAll->Fill(0.);
   
  if (!fInputEvent) {
    AliDebug(2,Form("ERROR: fESD not available"));
    cout << "ERROR: fESD not available" << endl;
    PostData(1, fHistListCosmics);
    return;
  }

  const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
  // Need vertex cut
  TString vtxName(vtx->GetName());
  if(!vtx || vtx->GetNContributors() < 2 || (vtxName.Contains("TPCVertex")) ) {
    // Post output data
    PostData(1, fHistListCosmics);
    return;
  }

  //  AliDebug(2,Form("Vertex title %s, status %d, nCont %d\n",vtx->GetTitle(), vtx->GetStatus(), vtx->GetNContributors()));
  double primVtx[3];
  vtx->GetXYZ(primVtx);
  if(TMath::Sqrt(primVtx[0]*primVtx[0] + primVtx[1]*primVtx[1])>1. || TMath::Abs(primVtx[2]>10.)){
    // Post output data
    PostData(1, fHistListCosmics);
    return;
  }
  if(!fInputEvent->GetNumberOfTracks() || fInputEvent->GetNumberOfTracks()<2){ 
    // Post output data
    PostData(1, fHistListCosmics);
    return;
  }
  Int_t nTracks = fInputEvent->GetNumberOfTracks();

  if(!fTrackCuts) {
   // Post output data
    PostData(1, fHistListCosmics);
    return;
  }

  fNEventSel->Fill(0.);

  Float_t dcaR[2] = {0.,0.};
  Float_t dcaZ[2] = {0.,0.};

  for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) {

    AliESDtrack* track1 = (AliESDtrack*)fInputEvent->GetTrack(iTrack1);
    if (!track1)  continue;
    if(!(fTrackCuts->AcceptTrack(track1))) { continue; }
    if(track1->Pt()<fPtMin) continue;
    //Start 2nd track loop to look for correlations
    for (Int_t iTrack2 = iTrack1+1; iTrack2 < nTracks; iTrack2++) {
      AliESDtrack *track2 = (AliESDtrack*)fInputEvent->GetTrack(iTrack2);
      if(!track2) continue;
      if(!(fTrackCuts->AcceptTrack(track2))) { continue; }
         
      //Check if back-to-back
      Double_t mom1[3],mom2[3];
      track1->GetPxPyPz(mom1);
      track2->GetPxPyPz(mom2);
      //     Double_t cosTheta = (mom1[0]*mom2[0]+mom1[1]*mom2[1]+mom1[2]*mom2[2])/( TMath::Sqrt(mom1[0]*mom1[0]+mom1[1]*mom1[1]+mom1[2]*mom1[2])*TMath::Sqrt(mom2[0]*mom2[0]+mom2[1]*mom2[1]+mom2[2]*mom2[2]) );
      TVector3 momv1(mom1[0],mom1[1],mom1[2]);
      TVector3 momv2(mom2[0],mom2[1],mom2[2]);
      //Double_t theta = momv1.Angle(momv2);
      Double_t theta = momv1.Phi()-momv2.Phi();
      if(theta<-0.5*TMath::Pi()) theta+=2.*TMath::Pi();
    
      fDeltaPtCosmicCandidates->Fill(track1->Pt()-track2->Pt());
      Float_t deltaPhi = track1->Phi()-track2->Phi();
      if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
      fDeltaPhiSumEta->Fill(deltaPhi,track1->Eta()+track2->Eta());

      track1->GetImpactParameters(dcaR[0],dcaZ[0]);
      track2->GetImpactParameters(dcaR[1],dcaZ[1]);

      if(track2->Pt()<0.5) continue;
      Double_t rIsol = TMath::Sqrt( deltaPhi*deltaPhi+(track1->Eta()-track2->Eta())*(track1->Eta()-track2->Eta()) );
      fRisol->Fill(rIsol); //Fill R histogram
      if(track2->Pt()<fPtMin) continue;

	fTheta->Fill(theta);
	fThetaZoom->Fill(theta);
	fThetaPt1Pt2->Fill(theta,track1->Pt(),track2->Pt());
	fThetaPt1Pt2Signed->Fill(theta,track1->GetSign()*track1->Pt(),track2->GetSign()*track2->Pt());
	fDeltaPhiSumEtaPt1->Fill(deltaPhi,track1->Eta()+track2->Eta(),track1->Pt());
	fDeltaPhiSumEtaPt2->Fill(deltaPhi,track1->Eta()+track2->Eta(),track2->Pt());
	fThetaDCAZ1DCAZ2->Fill(theta,dcaZ[0],dcaZ[1]);
	fRisolTheta->Fill(rIsol,theta);
	if(TMath::Abs(TMath::Pi()-theta)<fMaxCosmicAngle) {
	  fDCAZCosmicCandidates->Fill(dcaZ[0],dcaZ[1]);
	  fDCARCosmicCandidates->Fill(dcaR[0],dcaR[1]);
	  fPtSignedCosmicCandidates->Fill(track1->GetSign()*track1->Pt());
	}

    } // track2 loop
    
  } // track1 loop 

  PostData(1,fHistListCosmics);
}      

//________________________________________________________________________
void AliPWG4CosmicCandidates::Terminate(Option_t *) 
{
  //
  // Called once at the end of the query
  //
  
}
Commit	Line	Data
1f329128	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
	16	//-----------------------------------------------------------------------
9b58297c	17	// Analysis task looking for cosmic candidates embedded in a p-p event
1f329128	18	// Task checks if particles are back-to-back in eta and phi
	19	//
	20	//
	21	// Author : Marta Verweij - UU - marta.verweij@cern.ch
	22	//-----------------------------------------------------------------------
	23
	24	#include "TVector3.h"
	25	#include <iostream>
	26	#include "TH1.h"
	27	#include "TH2.h"
	28	#include "TH3.h"
	29	#include "TList.h"
	30	#include "TChain.h"
	31
	32	#include "AliAnalysisTask.h"
	33	#include "AliAnalysisManager.h"
	34
	35	#include "AliESDEvent.h"
	36	#include "AliESDInputHandler.h"
	37
	38	#include "AliESDtrack.h"
	39	#include "AliESDtrackCuts.h"
	40	#include "AliExternalTrackParam.h"
	41
	42	#include "AliLog.h"
	43
	44	#include "AliPWG4CosmicCandidates.h"
	45
ffc1a68e	46	//using namespace std; //required for resolving the 'cout' symbol
1f329128	47	using namespace std;
1f329128	48
ffc1a68e	49	ClassImp(AliPWG4CosmicCandidates)
1f329128	50
	51	//________________________________________________________________________
	52	AliPWG4CosmicCandidates::AliPWG4CosmicCandidates()
ffc1a68e	53	: AliAnalysisTaskSE(),
	54	fTrackCuts(0),
	55	fPtMin(5.),
	56	fMaxCosmicAngle(0.002),
	57	fNEventAll(0),
	58	fNEventSel(0),
	59	fPtSignedCosmicCandidates(0),
	60	fDeltaPtCosmicCandidates(0),
	61	fDeltaPhiSumEta(0),
	62	fDCAZCosmicCandidates(0),
	63	fDCARCosmicCandidates(0),
	64	fTheta(0),
	65	fThetaZoom(0),
	66	fThetaPt1Pt2(0),
	67	fThetaPt1Pt2Signed(0),
	68	fDeltaPhiSumEtaPt1(0),
	69	fDeltaPhiSumEtaPt2(0),
	70	fThetaDCAZ1DCAZ2(0),
	71	fRisol(0),
	72	fRisolTheta(0),
	73	fHistListCosmics(0)
1f329128	74	{
ffc1a68e	75	//
1f329128	76	// Default constructor
ffc1a68e	77	//
1f329128	78	}
	79
	80	//________________________________________________________________________
	81	AliPWG4CosmicCandidates::AliPWG4CosmicCandidates(const char *name)
	82	: AliAnalysisTaskSE(name),
	83	fTrackCuts(0),
	84	fPtMin(5.),
	85	fMaxCosmicAngle(0.002),
ffc1a68e	86	fNEventAll(0),
ffc1a68e	87	fNEventSel(0),
1f329128	88	fPtSignedCosmicCandidates(0),
	89	fDeltaPtCosmicCandidates(0),
	90	fDeltaPhiSumEta(0),
	91	fDCAZCosmicCandidates(0),
	92	fDCARCosmicCandidates(0),
	93	fTheta(0),
	94	fThetaZoom(0),
	95	fThetaPt1Pt2(0),
	96	fThetaPt1Pt2Signed(0),
	97	fDeltaPhiSumEtaPt1(0),
	98	fDeltaPhiSumEtaPt2(0),
	99	fThetaDCAZ1DCAZ2(0),
	100	fRisol(0),
	101	fRisolTheta(0),
	102	fHistListCosmics(0)
	103	{
	104	// Constructor. Initialization of Inputs and Outputs
	105
	106	// Define input and output slots here
	107	// Input slot #0 works with a TChain
	108	// DefineInput(0, TChain::Class());
	109	// Output slot #0 id reserved by the base class for AOD
	110	// Output slot #1 writes into a TList
	111	DefineOutput(1, TList::Class());
	112	// Output slot #2 writes into a AliESDtrackCuts
	113	DefineOutput(2, AliESDtrackCuts::Class());
	114	}
	115	//________________________________________________________________________
	116	AliPWG4CosmicCandidates::AliPWG4CosmicCandidates(const AliPWG4CosmicCandidates &res)
	117	: AliAnalysisTaskSE(res),
	118	fTrackCuts(0),
	119	fPtMin(5.),
	120	fMaxCosmicAngle(0.002),
ffc1a68e	121	fNEventAll(0),
ffc1a68e	122	fNEventSel(0),
1f329128	123	fPtSignedCosmicCandidates(0),
	124	fDeltaPtCosmicCandidates(0),
	125	fDeltaPhiSumEta(0),
	126	fDCAZCosmicCandidates(0),
	127	fDCARCosmicCandidates(0),
	128	fTheta(0),
	129	fThetaZoom(0),
	130	fThetaPt1Pt2(0),
	131	fThetaPt1Pt2Signed(0),
	132	fDeltaPhiSumEtaPt1(0),
	133	fDeltaPhiSumEtaPt2(0),
	134	fThetaDCAZ1DCAZ2(0),
	135	fRisol(0),
	136	fRisolTheta(0),
	137	fHistListCosmics(0)
	138	{
9b58297c	139	//
	140	// Dummy copy constructor
	141	//
1f329128	142	}
	143
	144	//________________________________________________________________________
	145	AliPWG4CosmicCandidates& AliPWG4CosmicCandidates::operator=(const AliPWG4CosmicCandidates& /trclass/)
	146	{
9b58297c	147	//
	148	// Dummy assignment operator
	149	//
	150	return *this;
1f329128	151	}
	152
	153	//________________________________________________________________________
	154	void AliPWG4CosmicCandidates::LocalInit()
	155	{
	156	//
	157	// Only called once at beginning
	158	//
	159	PostData(2,fTrackCuts);
	160	}
	161
	162	//________________________________________________________________________
	163	void AliPWG4CosmicCandidates::UserCreateOutputObjects()
	164	{
9b58297c	165	//
9b58297c	166	// Create output objects
1f329128	167	// Called once
	168	AliDebug(2,Form(">> AliPWG4CosmicCandidates::UserCreateOutputObjects \n"));
	169
	170	Bool_t oldStatus = TH1::AddDirectoryStatus();
	171	TH1::AddDirectory(kFALSE);
	172
	173	OpenFile(1);
	174	fHistListCosmics = new TList();
	175
	176	fNEventAll = new TH1F("fNEventAll","NEventAll",1,-0.5,0.5);
	177	fHistListCosmics->Add(fNEventAll);
	178	fNEventSel = new TH1F("fNEventSel","NEvent Selected for analysis",1,-0.5,0.5);
	179	fHistListCosmics->Add(fNEventSel);
	180
	181	Float_t fgkPtMin=0.;
	182	Float_t fgkPtMax=100.;
	183	Int_t fgkNPtBins= (int)(fgkPtMax-fgkPtMin);
	184
	185	Int_t fgkNPhiBins=18;
	186	Float_t kMinPhi = -0.5*TMath::Pi();
	187	Float_t kMaxPhi = 3./2.*TMath::Pi();
	188
	189	Int_t fgkNThetaBins=fgkNPhiBins*8;
	190	Float_t kMinTheta = -0.5*TMath::Pi();
	191	Float_t kMaxTheta = 3./2.*TMath::Pi();
	192
	193	Int_t fgkNDCARBins=80;
	194	Float_t fgkDCARMin = -0.2;
	195	Float_t fgkDCARMax = 0.2;
60829498	196	Float_t *binsDCAR=new Float_t[fgkNDCARBins+1];
60829498	197	for(Int_t i=0; i<=fgkNDCARBins; i++) binsDCAR[i]=(Float_t)fgkDCARMin + (fgkDCARMax-fgkDCARMin)/fgkNDCARBins*(Float_t)i ;
1f329128	198
	199	Int_t fgkNDCAZBins=80;
	200	Float_t fgkDCAZMin = -2.;
	201	Float_t fgkDCAZMax = 2.;
60829498	202	Float_t *binsDCAZ=new Float_t[fgkNDCAZBins+1];
60829498	203	for(Int_t i=0; i<=fgkNDCAZBins; i++) binsDCAZ[i]=(Float_t)fgkDCAZMin + (fgkDCAZMax-fgkDCAZMin)/fgkNDCAZBins*(Float_t)i ;
1f329128	204
	205	fPtSignedCosmicCandidates = new TH1F("fPtSignedCosmicCandidates","fPtSignedCosmicCandidates",2(int)(fgkPtMax-fgkPtMin), -1.fgkPtMax, fgkPtMax);
	206	fHistListCosmics->Add(fPtSignedCosmicCandidates);
	207
	208	fDeltaPtCosmicCandidates = new TH1F("fDeltaPtCosmicCandidates","fDeltaPtCosmicCandidates",fgkNPtBins, -50., 50.);
	209	fHistListCosmics->Add(fDeltaPtCosmicCandidates);
	210
	211	fDeltaPhiSumEta = new TH2F("fDeltaPhiSumEta","fDeltaPhiSumEta",fgkNPhiBins*4,kMinPhi,kMaxPhi,80, -2.,2.);
	212	fHistListCosmics->Add(fDeltaPhiSumEta);
	213
	214	fDCAZCosmicCandidates = new TH2F("fDCAZCosmicCandidates","fDCAZCosmicCandidates",fgkNDCAZBins,binsDCAZ,fgkNDCAZBins,binsDCAZ);
	215	fHistListCosmics->Add(fDCAZCosmicCandidates);
	216
	217	fDCARCosmicCandidates = new TH2F("fDCARCosmicCandidates","fDCARCosmicCandidates",fgkNDCARBins,binsDCAR,fgkNDCARBins,binsDCAR);
	218	fHistListCosmics->Add(fDCARCosmicCandidates);
	219
	220	fTheta = new TH1F("fTheta","fTheta",fgkNThetaBins,kMinTheta,kMaxTheta);
	221	fHistListCosmics->Add(fTheta);
	222
	223	fThetaZoom = new TH1F("fThetaZoom","fThetaZoom",100,TMath::Pi()-1.,TMath::Pi()+1.);
	224	fHistListCosmics->Add(fThetaZoom);
	225
	226	fThetaPt1Pt2 = new TH3F("fThetaPt1Pt2","fThetaPt1Pt2",fgkNThetaBins,kMinTheta,kMaxTheta,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
	227	fHistListCosmics->Add(fThetaPt1Pt2);
	228
	229	fThetaPt1Pt2Signed = new TH3F("fThetaPt1Pt2Signed","fThetaPt1Pt2Signed",fgkNThetaBins,kMinTheta,kMaxTheta,4(int)(fgkPtMax-fgkPtMin),-1.fgkPtMax,fgkPtMax,4(int)(fgkPtMax-fgkPtMin),-1.fgkPtMax,fgkPtMax);
	230	fHistListCosmics->Add(fThetaPt1Pt2Signed);
	231
	232	fDeltaPhiSumEtaPt1 = new TH3F("fDeltaPhiSumEtaPt1","fDeltaPhiSumEtaPt1",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
	233	fHistListCosmics->Add(fDeltaPhiSumEtaPt1);
	234
	235	fDeltaPhiSumEtaPt2 = new TH3F("fDeltaPhiSumEtaPt2","fDeltaPhiSumEtaPt2",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
	236	fHistListCosmics->Add(fDeltaPhiSumEtaPt2);
	237
	238	fThetaDCAZ1DCAZ2 = new TH3F("fThetaDCAZ1DCAZ2","fThetaDCAZ1DCAZ2",fgkNThetaBins,kMinTheta,kMaxTheta,fgkNDCAZBins,-2.,2.,fgkNDCAZBins,-2.,2.);
	239	fHistListCosmics->Add(fThetaDCAZ1DCAZ2);
	240
	241	fRisol = new TH1F("fRisol","fRisol",100,0.,10.);
	242	fHistListCosmics->Add(fRisol);
	243
	244	fRisolTheta = new TH2F("fRisolTheta","fRisolTheta",100,0.,10.,fgkNThetaBins,kMinTheta,kMaxTheta);
	245	fHistListCosmics->Add(fRisolTheta);
	246
	247	TH1::AddDirectory(oldStatus);
	248
	249	if(binsDCAR) delete [] binsDCAR;
	250	if(binsDCAZ) delete [] binsDCAZ;
	251
	252	}
	253
	254	//________________________________________________________________________
	255	void AliPWG4CosmicCandidates::UserExec(Option_t *)
	256	{
9b58297c	257	// Main loop
	258	// Called for each event
	259	//
1f329128	260
	261	// All events without selection
	262	fNEventAll->Fill(0.);
	263
	264	if (!fInputEvent) {
	265	AliDebug(2,Form("ERROR: fESD not available"));
	266	cout << "ERROR: fESD not available" << endl;
	267	PostData(1, fHistListCosmics);
	268	return;
	269	}
	270
	271	const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
	272	// Need vertex cut
75946d5d	273	TString vtxName(vtx->GetName());
75946d5d	274	if(!vtx \|\| vtx->GetNContributors() < 2 \|\| (vtxName.Contains("TPCVertex")) ) {
1f329128	275	// Post output data
	276	PostData(1, fHistListCosmics);
	277	return;
	278	}
	279
	280	// AliDebug(2,Form("Vertex title %s, status %d, nCont %d\n",vtx->GetTitle(), vtx->GetStatus(), vtx->GetNContributors()));
	281	double primVtx[3];
	282	vtx->GetXYZ(primVtx);
	283	if(TMath::Sqrt(primVtx[0]primVtx[0] + primVtx[1]primVtx[1])>1. \|\| TMath::Abs(primVtx[2]>10.)){
	284	// Post output data
	285	PostData(1, fHistListCosmics);
	286	return;
	287	}
	288	if(!fInputEvent->GetNumberOfTracks() \|\| fInputEvent->GetNumberOfTracks()<2){
	289	// Post output data
	290	PostData(1, fHistListCosmics);
	291	return;
	292	}
	293	Int_t nTracks = fInputEvent->GetNumberOfTracks();
	294
	295	if(!fTrackCuts) {
	296	// Post output data
	297	PostData(1, fHistListCosmics);
	298	return;
	299	}
	300
9b58297c	301	fNEventSel->Fill(0.);
1f329128	302
	303	Float_t dcaR[2] = {0.,0.};
	304	Float_t dcaZ[2] = {0.,0.};
	305
1f329128	306	for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) {
	307
	308	AliESDtrack* track1 = (AliESDtrack*)fInputEvent->GetTrack(iTrack1);
	309	if (!track1) continue;
	310	if(!(fTrackCuts->AcceptTrack(track1))) { continue; }
	311	if(track1->Pt()<fPtMin) continue;
	312	//Start 2nd track loop to look for correlations
	313	for (Int_t iTrack2 = iTrack1+1; iTrack2 < nTracks; iTrack2++) {
	314	AliESDtrack track2 = (AliESDtrack)fInputEvent->GetTrack(iTrack2);
	315	if(!track2) continue;
	316	if(!(fTrackCuts->AcceptTrack(track2))) { continue; }
	317
	318	//Check if back-to-back
	319	Double_t mom1[3],mom2[3];
	320	track1->GetPxPyPz(mom1);
	321	track2->GetPxPyPz(mom2);
	322	// Double_t cosTheta = (mom1[0]mom2[0]+mom1[1]mom2[1]+mom1[2]mom2[2])/( TMath::Sqrt(mom1[0]mom1[0]+mom1[1]mom1[1]+mom1[2]mom1[2])TMath::Sqrt(mom2[0]mom2[0]+mom2[1]mom2[1]+mom2[2]mom2[2]) );
	323	TVector3 momv1(mom1[0],mom1[1],mom1[2]);
	324	TVector3 momv2(mom2[0],mom2[1],mom2[2]);
	325	//Double_t theta = momv1.Angle(momv2);
	326	Double_t theta = momv1.Phi()-momv2.Phi();
	327	if(theta<-0.5TMath::Pi()) theta+=2.TMath::Pi();
	328
	329	fDeltaPtCosmicCandidates->Fill(track1->Pt()-track2->Pt());
	330	Float_t deltaPhi = track1->Phi()-track2->Phi();
	331	if(deltaPhi<-0.5TMath::Pi()) deltaPhi+=2.TMath::Pi();
	332	fDeltaPhiSumEta->Fill(deltaPhi,track1->Eta()+track2->Eta());
	333
	334	track1->GetImpactParameters(dcaR[0],dcaZ[0]);
	335	track2->GetImpactParameters(dcaR[1],dcaZ[1]);
	336
	337	if(track2->Pt()<0.5) continue;
	338	Double_t rIsol = TMath::Sqrt( deltaPhideltaPhi+(track1->Eta()-track2->Eta())(track1->Eta()-track2->Eta()) );
	339	fRisol->Fill(rIsol); //Fill R histogram
	340	if(track2->Pt()<fPtMin) continue;
	341
	342	fTheta->Fill(theta);
	343	fThetaZoom->Fill(theta);
	344	fThetaPt1Pt2->Fill(theta,track1->Pt(),track2->Pt());
	345	fThetaPt1Pt2Signed->Fill(theta,track1->GetSign()track1->Pt(),track2->GetSign()track2->Pt());
	346	fDeltaPhiSumEtaPt1->Fill(deltaPhi,track1->Eta()+track2->Eta(),track1->Pt());
	347	fDeltaPhiSumEtaPt2->Fill(deltaPhi,track1->Eta()+track2->Eta(),track2->Pt());
	348	fThetaDCAZ1DCAZ2->Fill(theta,dcaZ[0],dcaZ[1]);
	349	fRisolTheta->Fill(rIsol,theta);
	350	if(TMath::Abs(TMath::Pi()-theta)<fMaxCosmicAngle) {
	351	fDCAZCosmicCandidates->Fill(dcaZ[0],dcaZ[1]);
	352	fDCARCosmicCandidates->Fill(dcaR[0],dcaR[1]);
	353	fPtSignedCosmicCandidates->Fill(track1->GetSign()*track1->Pt());
	354	}
	355
	356	} // track2 loop
	357
	358	} // track1 loop
	359
	360	PostData(1,fHistListCosmics);
	361	}
	362
	363	//________________________________________________________________________
	364	void AliPWG4CosmicCandidates::Terminate(Option_t *)
	365	{
9b58297c	366	//
1f329128	367	// Called once at the end of the query
9b58297c	368	//
1f329128	369
1f329128	370	}