[u/mrichter/AliRoot.git] / PWGPP / ITS / AliTrackMatchingTPCITSCosmics.cxx

/**************************************************************************
 * Copyright(c) 1998-2009, 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.                  *
 **************************************************************************/

/////////////////////////////////////////////////////////////
//
// AliAnalysisTask to study the TPC-ITS track matching
//
// Author: A.Dainese, andrea.dainese@pd.infn.it
/////////////////////////////////////////////////////////////

#include <TTree.h>
#include <TChain.h>
#include <TNtuple.h>
#include <TList.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TVector3.h>
#include <TGeoManager.h>

#include "AliLog.h"
#include "AliGeomManager.h"
#include "AliTrackPointArray.h"
#include "AliESDInputHandler.h"
#include "AliESDVertex.h"
#include "AliESDtrack.h"
#include "AliESDEvent.h"
#include "AliESDfriend.h"
#include "AliAnalysisTask.h"
#include "AliAnalysisManager.h"
#include "AliExternalTrackParam.h"
#include "AliTrackMatchingTPCITSCosmics.h"


ClassImp(AliTrackMatchingTPCITSCosmics)


//________________________________________________________________________
AliTrackMatchingTPCITSCosmics::AliTrackMatchingTPCITSCosmics(const char *name):
AliAnalysisTask(name,"task"),
fOnlySPDFO(kFALSE),
fReadHLTESD(kFALSE),
fGeometryFileName("geometry.root"),
fESD(0),
fList(0),
fHistEvCount(0),
fntTrks(0),
fntPairs(0),
fntITSPairs(0)
{
  // Constructor

  // Define input and output slots here
  // Input slot #0 works with a TChain
  DefineInput(0, TChain::Class());
  // Output slot #0 writes into a TList
  DefineOutput(0,TList::Class());  //My private output
}

//________________________________________________________________________
AliTrackMatchingTPCITSCosmics::~AliTrackMatchingTPCITSCosmics()
{
  // Destructor
  if (fList) {
    delete fList;
    fList = 0;
  }
  if (fHistEvCount) {
    delete fHistEvCount;
    fHistEvCount = 0;
  }
  if (fntTrks) {
    delete fntTrks;
    fntTrks = 0;
  }
  if (fntPairs) {
    delete fntPairs;
    fntPairs = 0;
  }
  if (fntITSPairs) {
    delete fntITSPairs;
    fntITSPairs = 0;
  }
}  

//________________________________________________________________________
void AliTrackMatchingTPCITSCosmics::ConnectInputData(Option_t *) 
{
  // Connect ESD
  // Called once

  TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
  if(!tree) {
    printf("ERROR: Could not read chain from input slot 0\n");
  } else {
    AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
    
    if(!esdH) {
      printf("ERROR: Could not get ESDInputHandler\n");
    } else {
      fESD = esdH->GetEvent();
    }
  }
  
  return;
}

//________________________________________________________________________
void AliTrackMatchingTPCITSCosmics::Init()
{
  // Initialization

  return;
}

//________________________________________________________________________
void AliTrackMatchingTPCITSCosmics::CreateOutputObjects()
{
  // Create the output container
  //

  // load the geometry  
  if(!gGeoManager) {    
    AliGeomManager::LoadGeometry(fGeometryFileName.Data());
    if(!gGeoManager) { 
      printf("AliTrackMatchingTPCITSCosmics::CreateOutputObjects(): no geometry loaded \n");
      return;
    }
  }

  // Several histograms are more conveniently managed in a TList
  fList = new TList();
  fList->SetOwner();

  fHistEvCount = new TH1F("fHistEvCount","0: all, 1: SPDFO, 2: vtx, 3: TPCtrks, 4: ITStrks, 5: HLT, 6: HLT&TPCtrks, 7: HLT&ITStrks",11,-0.5,10.5);
  fList->Add(fHistEvCount);

  fntTrks = new TNtuple("fntTrks","TPC-ITS matching","ptTPC:nClsTPC:nClsITSSA:nClsITSMI:SSD1:SSD2:phi:z:dx:dy:dz:drphi:dphi:dtgl");
  fList->Add(fntTrks);
  fntPairs = new TNtuple("fntPairs","pairs at vertex","nClsITSSAin:nClsITSSAout:nClsITSMIin:nClsITSMIout:dxyITSSA:dzITSSA:dxyITSMI:dzITSMI:ptITSSAin:ptITSSAout:ptITSMIin:ptITSMIout:sigmad0MIin:sigmad0MIout");
  fList->Add(fntPairs);
  fntITSPairs = new TNtuple("fntITSPairs","pairs at vertex","nClsITSSAin:nClsITSSAout:dxyITSSA:dzITSSA:ptITSSAin:ptITSSAout:d0mu");
  fList->Add(fntITSPairs);

  return;
}

//________________________________________________________________________
void AliTrackMatchingTPCITSCosmics::Exec(Option_t */*option*/)
{
  // Execute analysis for current event:
  // write ITS AliTrackPoints for selected tracks to fspTree
  
  // check the geometry  
  if(!gGeoManager) { 
    printf("AliTrackMatchingTPCITSCosmics::Exec(): no geometry loaded \n");
    return;
  }

  if(!fESD) {
    printf("AliTrackMatchingTPCITSCosmics::Exec(): no ESD \n");
    return;
  } 


  fHistEvCount->Fill(0);
  PostData(0,fList);

  // check if event is triggered by HLT
  Bool_t hltTrigg=kFALSE;
  if(fReadHLTESD) {
    AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
    if(!esdH) {
      printf("AliTrackMatchingTPCITSCosmics::Exec(): no ESD \n");
      return;
    } 
    AliESDEvent *hltESD = esdH->GetHLTEvent();
    if(!hltESD) {
      printf("AliTrackMatchingTPCITSCosmics::Exec(): no HLT ESD \n");
      return;
    } 
    if(hltESD->IsHLTTriggerFired()) {fHistEvCount->Fill(5);hltTrigg=kTRUE;}
  }

  
  TString triggeredClass = fESD->GetFiredTriggerClasses(); 
  if(triggeredClass.Contains("C0SCO-ABCE-NOPF-CENT")) {
    fHistEvCount->Fill(1);
  } else {
    if(fOnlySPDFO) {  PostData(0,fList); return; }
  }
  
  Int_t ntracks = fESD->GetNumberOfTracks();
  printf("CONTR: %d\n",fESD->GetVertex()->GetNContributors());
  if(fESD->GetVertex()->GetNContributors()>=0) {
    fHistEvCount->Fill(2);
  } 
  
  Int_t nTrksITSSA=0,nTrksTPC=0,nTrksTPCITS=0;
  Int_t idxITSSA[10000],idxTPC[10000];
  Int_t idxITSSAin=-1,idxITSSAout=-1;
  Int_t idxTPCin=-1,idxTPCout=-1;
  for(Int_t itr=0; itr<ntracks; itr++) {
    AliESDtrack *t = fESD->GetTrack(itr);
    Int_t nClsTPC = t->GetNcls(1);
    Int_t nClsITS = t->GetNcls(0);
    if(nClsITS>=2 && nClsTPC==0) { // ITS SA
      idxITSSA[nTrksITSSA]=itr;
      nTrksITSSA++;
      //printf("%d\n",itr);
    }
    if(nClsTPC>=50) {              // TPC
      idxTPC[nTrksTPC]=itr;
      nTrksTPC++;
    }
    if(nClsTPC>=50 && nClsITS>=2) {// TPC+ITS
      nTrksTPCITS++;
      /*printf("  --> TPC+ITS: %d + %d (%d,%d,%d,%d,%d,%d) pt = %f\n",nClsTPC,nClsITS,
	(Int_t)(t->HasPointOnITSLayer(0)),
	(Int_t)(t->HasPointOnITSLayer(1)),
	(Int_t)(t->HasPointOnITSLayer(2)),
	(Int_t)(t->HasPointOnITSLayer(3)),
	(Int_t)(t->HasPointOnITSLayer(4)),
	(Int_t)(t->HasPointOnITSLayer(5)),
	t->Pt());*/
    }
  }
  
  //printf("nTrksTPC %d  nTrksITSSA %d\n",nTrksTPC,nTrksITSSA);
  if(nTrksITSSA>0) fHistEvCount->Fill(4);
  if(hltTrigg && nTrksITSSA>0) fHistEvCount->Fill(7);
  if(nTrksTPC>0)   fHistEvCount->Fill(3);
  if(hltTrigg && nTrksTPC>0)   fHistEvCount->Fill(6);
  if(nTrksTPCITS>0)   fHistEvCount->Fill(8);
  if(hltTrigg && nTrksTPCITS>0)   fHistEvCount->Fill(9);
  
  if(nTrksTPC>2 || nTrksTPC==0 || nTrksITSSA>2 || nTrksITSSA==0) { PostData(0,fList); return; }
  
  for(Int_t itr=0; itr<nTrksTPC; itr++) {
    AliESDtrack *t = fESD->GetTrack(idxTPC[itr]);
    if(t->Py()>0) idxTPCin=idxTPC[itr];
    if(t->Py()<0) idxTPCout=idxTPC[itr];
  }
  
  
  for(Int_t itr=0; itr<nTrksITSSA; itr++) {
    AliESDtrack *t = fESD->GetTrack(idxITSSA[itr]);
    if(t->Py()>0) idxITSSAin=idxITSSA[itr];
    if(t->Py()<0) idxITSSAout=idxITSSA[itr];
  }
  
  Double_t xyzITS[3],xyzTPC[3];

  // analysis for inward track
  if(idxITSSAin>-1 && idxTPCin>-1) { 
    //printf(" %d\n",idxITSSAin);
    AliESDtrack *tITS = fESD->GetTrack(idxITSSAin);
    const AliExternalTrackParam *outerITS = tITS->GetOuterParam();
    xyzITS[0]=-1000.;xyzITS[1]=-1000.;xyzITS[2]=-1000.;
    if(outerITS) outerITS->GetXYZAt(50.,fESD->GetMagneticField(),xyzITS);
    AliESDtrack *tTPC = fESD->GetTrack(idxTPCin);
    AliExternalTrackParam *innerTPC = new AliExternalTrackParam(*(tTPC->GetInnerParam()));
    innerTPC->Rotate(tITS->GetAlpha());
    innerTPC->GetXYZAt(50.,fESD->GetMagneticField(),xyzTPC);
    fntTrks->Fill(innerTPC->Pt(),
		 tTPC->GetNcls(1),
		 tITS->GetNcls(0),
		 tTPC->GetNcls(0),
		 tITS->HasPointOnITSLayer(4),
		 tITS->HasPointOnITSLayer(5),
		 TMath::ATan2(xyzITS[1],xyzITS[0]),
		 xyzITS[2],
		 xyzITS[0]-xyzTPC[0],
		 xyzITS[1]-xyzTPC[1],
		 xyzITS[2]-xyzTPC[2],
		 0.5*(TMath::Sqrt(xyzITS[0]*xyzITS[0]+xyzITS[1]*xyzITS[1])+TMath::Sqrt(xyzTPC[0]*xyzTPC[0]+xyzTPC[1]*xyzTPC[1]))*(TMath::ATan2(xyzITS[1],xyzITS[0])-TMath::ATan2(xyzTPC[1],xyzTPC[0])),
		 TMath::ATan2(tITS->Py(),tITS->Px())-TMath::ATan2(innerTPC->Py(),innerTPC->Px()),
		 tITS->Pz()/tITS->Pt()-innerTPC->Pz()/innerTPC->Pt());
    delete innerTPC; innerTPC=NULL;
  }

  // analysis for outward track
  if(idxITSSAout>-1 && idxTPCout>-1) { 
    //printf(" %d\n",idxITSSAout);
    AliESDtrack *tITS = fESD->GetTrack(idxITSSAout);
    const AliExternalTrackParam *outerITS = tITS->GetOuterParam();
    xyzITS[0]=-1000.;xyzITS[1]=-1000.;xyzITS[2]=-1000.;
    if(outerITS) outerITS->GetXYZAt(50.,fESD->GetMagneticField(),xyzITS);
    AliESDtrack *tTPC = fESD->GetTrack(idxTPCout);
    AliExternalTrackParam *innerTPC = new AliExternalTrackParam(*(tTPC->GetInnerParam()));
    innerTPC->Rotate(tITS->GetAlpha());
    innerTPC->GetXYZAt(50.,fESD->GetMagneticField(),xyzTPC);
    fntTrks->Fill(innerTPC->Pt(),
		 tTPC->GetNcls(1),
		 tITS->GetNcls(0),
		 tTPC->GetNcls(0),
		 tITS->HasPointOnITSLayer(4),
		 tITS->HasPointOnITSLayer(5),
		 TMath::ATan2(xyzITS[1],xyzITS[0]),
		 xyzITS[2],
		 xyzITS[0]-xyzTPC[0],
		 xyzITS[1]-xyzTPC[1],
		 xyzITS[2]-xyzTPC[2],
		 0.5*(TMath::Sqrt(xyzITS[0]*xyzITS[0]+xyzITS[1]*xyzITS[1])+TMath::Sqrt(xyzTPC[0]*xyzTPC[0]+xyzTPC[1]*xyzTPC[1]))*(TMath::ATan2(xyzITS[1],xyzITS[0])-TMath::ATan2(xyzTPC[1],xyzTPC[0])),
		 TMath::ATan2(tITS->Py(),tITS->Px())-TMath::ATan2(innerTPC->Py(),innerTPC->Px()),
		 tITS->Pz()/tITS->Pt()-innerTPC->Pz()/innerTPC->Pt()
		 );
    delete innerTPC; innerTPC=NULL;
  }

  // inward-outward track ITS only
  if(idxITSSAin>-1 && idxITSSAout>-1) {
    AliESDtrack *tITSin = fESD->GetTrack(idxITSSAin);
    AliESDtrack *tITSout = fESD->GetTrack(idxITSSAout);
    if(tITSin->HasPointOnITSLayer(0) && 
       tITSin->HasPointOnITSLayer(1) &&
       tITSout->HasPointOnITSLayer(0) &&
       tITSout->HasPointOnITSLayer(1)) {
      Double_t alpha = TMath::ATan2(tITSin->Py(),tITSin->Px());
      tITSin->Propagate(alpha,0.,fESD->GetMagneticField());
      tITSout->Propagate(alpha,0.,fESD->GetMagneticField());
      Float_t d0[2],z0[2];
      d0[0] = tITSin->GetY();
      z0[0] = tITSin->GetZ();
      d0[1] = tITSout->GetY();
      z0[1] = tITSout->GetZ();
      Float_t dxyITS = -(d0[0]-d0[1]);
      Float_t dzITS = z0[0]-z0[1];
      fntITSPairs->Fill(tITSin->GetNcls(0),tITSout->GetNcls(0),dxyITS,dzITS,tITSin->Pt(),tITSout->Pt(),TMath::Abs(d0[0]));
    }
  }

  // inward-outward track TPC+ITS
  if(idxITSSAin>-1 && idxTPCin>-1 && idxITSSAout>-1 && idxTPCout>-1) {
    AliESDtrack *tTPCin = fESD->GetTrack(idxTPCin);
    AliESDtrack *tITSin = fESD->GetTrack(idxITSSAin);
    AliESDtrack *tTPCout = fESD->GetTrack(idxTPCout);
    AliESDtrack *tITSout = fESD->GetTrack(idxITSSAout);
    if(tTPCin->HasPointOnITSLayer(0) && 
       tTPCin->HasPointOnITSLayer(1) &&
       tITSin->HasPointOnITSLayer(0) && 
       tITSin->HasPointOnITSLayer(1) && 
       tTPCout->HasPointOnITSLayer(0) &&
       tTPCout->HasPointOnITSLayer(1) &&
       tITSout->HasPointOnITSLayer(0) &&
       tITSout->HasPointOnITSLayer(1)) {
      Double_t alpha = TMath::ATan2(tTPCin->Py(),tTPCin->Px());
      tTPCin->Propagate(alpha,0.,fESD->GetMagneticField());
      tTPCout->Propagate(alpha,0.,fESD->GetMagneticField());
      Float_t d0[2],z0[2];
      d0[0] = tTPCin->GetY();
      z0[0] = tTPCin->GetZ();
      d0[1] = tTPCout->GetY();
      z0[1] = tTPCout->GetZ();
      Float_t dxyTPC = -(d0[0]-d0[1]);
      Float_t dzTPC = z0[0]-z0[1];
      alpha = TMath::ATan2(tITSin->Py(),tITSin->Px());
      tITSin->Propagate(alpha,0.,fESD->GetMagneticField());
      tITSout->Propagate(alpha,0.,fESD->GetMagneticField());
      d0[0] = tITSin->GetY();
      z0[0] = tITSin->GetZ();
      d0[1] = tITSout->GetY();
      z0[1] = tITSout->GetZ();
      Float_t dxyITS = -(d0[0]-d0[1]);
      Float_t dzITS = z0[0]-z0[1];
      fntPairs->Fill(tITSin->GetNcls(0),tITSout->GetNcls(0),tTPCin->GetNcls(0),tTPCout->GetNcls(0),dxyITS,dzITS,dxyTPC,dzTPC,tITSin->Pt(),tITSout->Pt(),tTPCin->Pt(),tTPCout->Pt(),tTPCin->GetSigmaY2(),tTPCout->GetSigmaY2());
    }
  }

  PostData(0,fList);

  return;
}

//________________________________________________________________________
void AliTrackMatchingTPCITSCosmics::Terminate(Option_t */*option*/)
{
  // Terminate analysis
  //
  AliDebug(2,"AliTrackMatchingTPCITSCosmics: Terminate() \n");

  fList = dynamic_cast<TList*> (GetOutputData(0));
  if (!fList) {     
    printf("ERROR: fList not available\n");
    return;
  }

  return;
}
Commit	Line	Data
dae0a357	1	/**************************************************************************
	2	* Copyright(c) 1998-2009, 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	/////////////////////////////////////////////////////////////
	17	//
	18	// AliAnalysisTask to study the TPC-ITS track matching
	19	//
	20	// Author: A.Dainese, andrea.dainese@pd.infn.it
	21	/////////////////////////////////////////////////////////////
	22
	23	#include <TTree.h>
	24	#include <TChain.h>
	25	#include <TNtuple.h>
	26	#include <TList.h>
	27	#include <TH1F.h>
	28	#include <TH2F.h>
	29	#include <TVector3.h>
	30	#include <TGeoManager.h>
	31
	32	#include "AliLog.h"
	33	#include "AliGeomManager.h"
	34	#include "AliTrackPointArray.h"
	35	#include "AliESDInputHandler.h"
	36	#include "AliESDVertex.h"
	37	#include "AliESDtrack.h"
	38	#include "AliESDEvent.h"
	39	#include "AliESDfriend.h"
	40	#include "AliAnalysisTask.h"
	41	#include "AliAnalysisManager.h"
	42	#include "AliExternalTrackParam.h"
	43	#include "AliTrackMatchingTPCITSCosmics.h"
	44
	45
	46	ClassImp(AliTrackMatchingTPCITSCosmics)
	47
	48
	49	//________________________________________________________________________
	50	AliTrackMatchingTPCITSCosmics::AliTrackMatchingTPCITSCosmics(const char *name):
	51	AliAnalysisTask(name,"task"),
77ec32bb	52	fOnlySPDFO(kFALSE),
dae0a357	53	fReadHLTESD(kFALSE),
	54	fGeometryFileName("geometry.root"),
	55	fESD(0),
	56	fList(0),
	57	fHistEvCount(0),
	58	fntTrks(0),
	59	fntPairs(0),
	60	fntITSPairs(0)
	61	{
	62	// Constructor
	63
	64	// Define input and output slots here
	65	// Input slot #0 works with a TChain
	66	DefineInput(0, TChain::Class());
	67	// Output slot #0 writes into a TList
	68	DefineOutput(0,TList::Class()); //My private output
	69	}
	70
	71	//________________________________________________________________________
	72	AliTrackMatchingTPCITSCosmics::~AliTrackMatchingTPCITSCosmics()
	73	{
	74	// Destructor
	75	if (fList) {
	76	delete fList;
	77	fList = 0;
	78	}
	79	if (fHistEvCount) {
	80	delete fHistEvCount;
	81	fHistEvCount = 0;
	82	}
	83	if (fntTrks) {
	84	delete fntTrks;
	85	fntTrks = 0;
	86	}
	87	if (fntPairs) {
	88	delete fntPairs;
	89	fntPairs = 0;
	90	}
	91	if (fntITSPairs) {
	92	delete fntITSPairs;
	93	fntITSPairs = 0;
	94	}
	95	}
	96
	97	//________________________________________________________________________
	98	void AliTrackMatchingTPCITSCosmics::ConnectInputData(Option_t *)
	99	{
	100	// Connect ESD
	101	// Called once
	102
	103	TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
	104	if(!tree) {
	105	printf("ERROR: Could not read chain from input slot 0\n");
	106	} else {
	107	AliESDInputHandler esdH = dynamic_cast<AliESDInputHandler> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
	108
	109	if(!esdH) {
	110	printf("ERROR: Could not get ESDInputHandler\n");
	111	} else {
77ec32bb	112	fESD = esdH->GetEvent();
dae0a357	113	}
	114	}
	115
	116	return;
	117	}
	118
	119	//________________________________________________________________________
	120	void AliTrackMatchingTPCITSCosmics::Init()
	121	{
	122	// Initialization
	123
	124	return;
	125	}
	126
	127	//________________________________________________________________________
	128	void AliTrackMatchingTPCITSCosmics::CreateOutputObjects()
	129	{
	130	// Create the output container
	131	//
	132
c1605e49	133	// load the geometry
	134	if(!gGeoManager) {
	135	AliGeomManager::LoadGeometry(fGeometryFileName.Data());
	136	if(!gGeoManager) {
	137	printf("AliTrackMatchingTPCITSCosmics::CreateOutputObjects(): no geometry loaded \n");
	138	return;
	139	}
	140	}
	141
dae0a357	142	// Several histograms are more conveniently managed in a TList
	143	fList = new TList();
	144	fList->SetOwner();
	145
77ec32bb	146	fHistEvCount = new TH1F("fHistEvCount","0: all, 1: SPDFO, 2: vtx, 3: TPCtrks, 4: ITStrks, 5: HLT, 6: HLT&TPCtrks, 7: HLT&ITStrks",11,-0.5,10.5);
dae0a357	147	fList->Add(fHistEvCount);
	148
	149	fntTrks = new TNtuple("fntTrks","TPC-ITS matching","ptTPC:nClsTPC:nClsITSSA:nClsITSMI:SSD1:SSD2:phi:z:dx:dy:dz:drphi:dphi:dtgl");
	150	fList->Add(fntTrks);
	151	fntPairs = new TNtuple("fntPairs","pairs at vertex","nClsITSSAin:nClsITSSAout:nClsITSMIin:nClsITSMIout:dxyITSSA:dzITSSA:dxyITSMI:dzITSMI:ptITSSAin:ptITSSAout:ptITSMIin:ptITSMIout:sigmad0MIin:sigmad0MIout");
	152	fList->Add(fntPairs);
	153	fntITSPairs = new TNtuple("fntITSPairs","pairs at vertex","nClsITSSAin:nClsITSSAout:dxyITSSA:dzITSSA:ptITSSAin:ptITSSAout:d0mu");
	154	fList->Add(fntITSPairs);
	155
	156	return;
	157	}
	158
	159	//________________________________________________________________________
	160	void AliTrackMatchingTPCITSCosmics::Exec(Option_t /option*/)
	161	{
	162	// Execute analysis for current event:
	163	// write ITS AliTrackPoints for selected tracks to fspTree
	164
c1605e49	165	// check the geometry
	166	if(!gGeoManager) {
	167	printf("AliTrackMatchingTPCITSCosmics::Exec(): no geometry loaded \n");
	168	return;
dae0a357	169	}
	170
	171	if(!fESD) {
	172	printf("AliTrackMatchingTPCITSCosmics::Exec(): no ESD \n");
	173	return;
	174	}
	175
	176
	177	fHistEvCount->Fill(0);
	178	PostData(0,fList);
	179
77ec32bb	180	// check if event is triggered by HLT
	181	Bool_t hltTrigg=kFALSE;
	182	if(fReadHLTESD) {
	183	AliESDInputHandler esdH = dynamic_cast<AliESDInputHandler> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
ca6660ea	184	if(!esdH) {
	185	printf("AliTrackMatchingTPCITSCosmics::Exec(): no ESD \n");
	186	return;
	187	}
77ec32bb	188	AliESDEvent *hltESD = esdH->GetHLTEvent();
	189	if(!hltESD) {
	190	printf("AliTrackMatchingTPCITSCosmics::Exec(): no HLT ESD \n");
	191	return;
	192	}
	193	if(hltESD->IsHLTTriggerFired()) {fHistEvCount->Fill(5);hltTrigg=kTRUE;}
	194	}
dae0a357	195
77ec32bb	196
dae0a357	197	TString triggeredClass = fESD->GetFiredTriggerClasses();
	198	if(triggeredClass.Contains("C0SCO-ABCE-NOPF-CENT")) {
	199	fHistEvCount->Fill(1);
	200	} else {
	201	if(fOnlySPDFO) { PostData(0,fList); return; }
	202	}
77ec32bb	203
dae0a357	204	Int_t ntracks = fESD->GetNumberOfTracks();
77ec32bb	205	printf("CONTR: %d\n",fESD->GetVertex()->GetNContributors());
dae0a357	206	if(fESD->GetVertex()->GetNContributors()>=0) {
	207	fHistEvCount->Fill(2);
	208	}
	209
	210	Int_t nTrksITSSA=0,nTrksTPC=0,nTrksTPCITS=0;
	211	Int_t idxITSSA[10000],idxTPC[10000];
	212	Int_t idxITSSAin=-1,idxITSSAout=-1;
	213	Int_t idxTPCin=-1,idxTPCout=-1;
	214	for(Int_t itr=0; itr<ntracks; itr++) {
	215	AliESDtrack *t = fESD->GetTrack(itr);
	216	Int_t nClsTPC = t->GetNcls(1);
	217	Int_t nClsITS = t->GetNcls(0);
	218	if(nClsITS>=2 && nClsTPC==0) { // ITS SA
	219	idxITSSA[nTrksITSSA]=itr;
	220	nTrksITSSA++;
	221	//printf("%d\n",itr);
	222	}
	223	if(nClsTPC>=50) { // TPC
	224	idxTPC[nTrksTPC]=itr;
	225	nTrksTPC++;
	226	}
	227	if(nClsTPC>=50 && nClsITS>=2) {// TPC+ITS
	228	nTrksTPCITS++;
	229	/*printf(" --> TPC+ITS: %d + %d (%d,%d,%d,%d,%d,%d) pt = %f\n",nClsTPC,nClsITS,
	230	(Int_t)(t->HasPointOnITSLayer(0)),
	231	(Int_t)(t->HasPointOnITSLayer(1)),
	232	(Int_t)(t->HasPointOnITSLayer(2)),
	233	(Int_t)(t->HasPointOnITSLayer(3)),
	234	(Int_t)(t->HasPointOnITSLayer(4)),
	235	(Int_t)(t->HasPointOnITSLayer(5)),
	236	t->Pt());*/
	237	}
	238	}
	239
	240	//printf("nTrksTPC %d nTrksITSSA %d\n",nTrksTPC,nTrksITSSA);
77ec32bb	241	if(nTrksITSSA>0) fHistEvCount->Fill(4);
	242	if(hltTrigg && nTrksITSSA>0) fHistEvCount->Fill(7);
	243	if(nTrksTPC>0) fHistEvCount->Fill(3);
	244	if(hltTrigg && nTrksTPC>0) fHistEvCount->Fill(6);
	245	if(nTrksTPCITS>0) fHistEvCount->Fill(8);
	246	if(hltTrigg && nTrksTPCITS>0) fHistEvCount->Fill(9);
dae0a357	247
	248	if(nTrksTPC>2 \|\| nTrksTPC==0 \|\| nTrksITSSA>2 \|\| nTrksITSSA==0) { PostData(0,fList); return; }
	249
	250	for(Int_t itr=0; itr<nTrksTPC; itr++) {
	251	AliESDtrack *t = fESD->GetTrack(idxTPC[itr]);
	252	if(t->Py()>0) idxTPCin=idxTPC[itr];
	253	if(t->Py()<0) idxTPCout=idxTPC[itr];
	254	}
	255
	256
	257	for(Int_t itr=0; itr<nTrksITSSA; itr++) {
	258	AliESDtrack *t = fESD->GetTrack(idxITSSA[itr]);
	259	if(t->Py()>0) idxITSSAin=idxITSSA[itr];
	260	if(t->Py()<0) idxITSSAout=idxITSSA[itr];
	261	}
	262
	263	Double_t xyzITS[3],xyzTPC[3];
	264
	265	// analysis for inward track
	266	if(idxITSSAin>-1 && idxTPCin>-1) {
	267	//printf(" %d\n",idxITSSAin);
	268	AliESDtrack *tITS = fESD->GetTrack(idxITSSAin);
	269	const AliExternalTrackParam *outerITS = tITS->GetOuterParam();
	270	xyzITS[0]=-1000.;xyzITS[1]=-1000.;xyzITS[2]=-1000.;
	271	if(outerITS) outerITS->GetXYZAt(50.,fESD->GetMagneticField(),xyzITS);
	272	AliESDtrack *tTPC = fESD->GetTrack(idxTPCin);
	273	AliExternalTrackParam innerTPC = new AliExternalTrackParam((tTPC->GetInnerParam()));
	274	innerTPC->Rotate(tITS->GetAlpha());
	275	innerTPC->GetXYZAt(50.,fESD->GetMagneticField(),xyzTPC);
	276	fntTrks->Fill(innerTPC->Pt(),
	277	tTPC->GetNcls(1),
	278	tITS->GetNcls(0),
	279	tTPC->GetNcls(0),
	280	tITS->HasPointOnITSLayer(4),
	281	tITS->HasPointOnITSLayer(5),
	282	TMath::ATan2(xyzITS[1],xyzITS[0]),
	283	xyzITS[2],
	284	xyzITS[0]-xyzTPC[0],
	285	xyzITS[1]-xyzTPC[1],
	286	xyzITS[2]-xyzTPC[2],
	287	0.5(TMath::Sqrt(xyzITS[0]xyzITS[0]+xyzITS[1]xyzITS[1])+TMath::Sqrt(xyzTPC[0]xyzTPC[0]+xyzTPC[1]xyzTPC[1]))(TMath::ATan2(xyzITS[1],xyzITS[0])-TMath::ATan2(xyzTPC[1],xyzTPC[0])),
	288	TMath::ATan2(tITS->Py(),tITS->Px())-TMath::ATan2(innerTPC->Py(),innerTPC->Px()),
	289	tITS->Pz()/tITS->Pt()-innerTPC->Pz()/innerTPC->Pt());
	290	delete innerTPC; innerTPC=NULL;
	291	}
	292
	293	// analysis for outward track
	294	if(idxITSSAout>-1 && idxTPCout>-1) {
	295	//printf(" %d\n",idxITSSAout);
	296	AliESDtrack *tITS = fESD->GetTrack(idxITSSAout);
	297	const AliExternalTrackParam *outerITS = tITS->GetOuterParam();
	298	xyzITS[0]=-1000.;xyzITS[1]=-1000.;xyzITS[2]=-1000.;
	299	if(outerITS) outerITS->GetXYZAt(50.,fESD->GetMagneticField(),xyzITS);
	300	AliESDtrack *tTPC = fESD->GetTrack(idxTPCout);
	301	AliExternalTrackParam innerTPC = new AliExternalTrackParam((tTPC->GetInnerParam()));
	302	innerTPC->Rotate(tITS->GetAlpha());
	303	innerTPC->GetXYZAt(50.,fESD->GetMagneticField(),xyzTPC);
	304	fntTrks->Fill(innerTPC->Pt(),
	305	tTPC->GetNcls(1),
	306	tITS->GetNcls(0),
	307	tTPC->GetNcls(0),
	308	tITS->HasPointOnITSLayer(4),
	309	tITS->HasPointOnITSLayer(5),
	310	TMath::ATan2(xyzITS[1],xyzITS[0]),
311	xyzITS[2],
312	xyzITS[0]-xyzTPC[0],
313	xyzITS[1]-xyzTPC[1],
314	xyzITS[2]-xyzTPC[2],
315	0.5(TMath::Sqrt(xyzITS[0]xyzITS[0]+xyzITS[1]xyzITS[1])+TMath::Sqrt(xyzTPC[0]xyzTPC[0]+xyzTPC[1]xyzTPC[1]))(TMath::ATan2(xyzITS[1],xyzITS[0])-TMath::ATan2(xyzTPC[1],xyzTPC[0])),
316	TMath::ATan2(tITS->Py(),tITS->Px())-TMath::ATan2(innerTPC->Py(),innerTPC->Px()),
317	tITS->Pz()/tITS->Pt()-innerTPC->Pz()/innerTPC->Pt()
318	);
319	delete innerTPC; innerTPC=NULL;
320	}
321
322	// inward-outward track ITS only
323	if(idxITSSAin>-1 && idxITSSAout>-1) {
324	AliESDtrack *tITSin = fESD->GetTrack(idxITSSAin);
325	AliESDtrack *tITSout = fESD->GetTrack(idxITSSAout);
326	if(tITSin->HasPointOnITSLayer(0) &&
327	tITSin->HasPointOnITSLayer(1) &&
328	tITSout->HasPointOnITSLayer(0) &&
329	tITSout->HasPointOnITSLayer(1)) {
330	Double_t alpha = TMath::ATan2(tITSin->Py(),tITSin->Px());
331	tITSin->Propagate(alpha,0.,fESD->GetMagneticField());
332	tITSout->Propagate(alpha,0.,fESD->GetMagneticField());
333	Float_t d0[2],z0[2];
334	d0[0] = tITSin->GetY();
335	z0[0] = tITSin->GetZ();
336	d0[1] = tITSout->GetY();
337	z0[1] = tITSout->GetZ();
338	Float_t dxyITS = -(d0[0]-d0[1]);
339	Float_t dzITS = z0[0]-z0[1];
340	fntITSPairs->Fill(tITSin->GetNcls(0),tITSout->GetNcls(0),dxyITS,dzITS,tITSin->Pt(),tITSout->Pt(),TMath::Abs(d0[0]));
341	}
342	}
343
344	// inward-outward track TPC+ITS
345	if(idxITSSAin>-1 && idxTPCin>-1 && idxITSSAout>-1 && idxTPCout>-1) {
346	AliESDtrack *tTPCin = fESD->GetTrack(idxTPCin);
347	AliESDtrack *tITSin = fESD->GetTrack(idxITSSAin);
348	AliESDtrack *tTPCout = fESD->GetTrack(idxTPCout);
349	AliESDtrack *tITSout = fESD->GetTrack(idxITSSAout);
350	if(tTPCin->HasPointOnITSLayer(0) &&
351	tTPCin->HasPointOnITSLayer(1) &&
352	tITSin->HasPointOnITSLayer(0) &&
353	tITSin->HasPointOnITSLayer(1) &&
354	tTPCout->HasPointOnITSLayer(0) &&
355	tTPCout->HasPointOnITSLayer(1) &&
356	tITSout->HasPointOnITSLayer(0) &&
357	tITSout->HasPointOnITSLayer(1)) {
358	Double_t alpha = TMath::ATan2(tTPCin->Py(),tTPCin->Px());
359	tTPCin->Propagate(alpha,0.,fESD->GetMagneticField());
360	tTPCout->Propagate(alpha,0.,fESD->GetMagneticField());
361	Float_t d0[2],z0[2];
362	d0[0] = tTPCin->GetY();
363	z0[0] = tTPCin->GetZ();
364	d0[1] = tTPCout->GetY();
365	z0[1] = tTPCout->GetZ();
366	Float_t dxyTPC = -(d0[0]-d0[1]);
367	Float_t dzTPC = z0[0]-z0[1];
368	alpha = TMath::ATan2(tITSin->Py(),tITSin->Px());
369	tITSin->Propagate(alpha,0.,fESD->GetMagneticField());
370	tITSout->Propagate(alpha,0.,fESD->GetMagneticField());
371	d0[0] = tITSin->GetY();
372	z0[0] = tITSin->GetZ();
373	d0[1] = tITSout->GetY();
374	z0[1] = tITSout->GetZ();
375	Float_t dxyITS = -(d0[0]-d0[1]);
376	Float_t dzITS = z0[0]-z0[1];
377	fntPairs->Fill(tITSin->GetNcls(0),tITSout->GetNcls(0),tTPCin->GetNcls(0),tTPCout->GetNcls(0),dxyITS,dzITS,dxyTPC,dzTPC,tITSin->Pt(),tITSout->Pt(),tTPCin->Pt(),tTPCout->Pt(),tTPCin->GetSigmaY2(),tTPCout->GetSigmaY2());
378	}
379	}
380
381	PostData(0,fList);
382
383	return;
384	}
385
386	//________________________________________________________________________
387	void AliTrackMatchingTPCITSCosmics::Terminate(Option_t /option*/)
388	{
389	// Terminate analysis
390	//
391	AliDebug(2,"AliTrackMatchingTPCITSCosmics: Terminate() \n");
392
393	fList = dynamic_cast<TList*> (GetOutputData(0));
394	if (!fList) {
395	printf("ERROR: fList not available\n");
396	return;
397	}
398
399	return;
400	}
401