[u/mrichter/AliRoot.git] / PWG1 / 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
  //

  // 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
  
  // load the geometry  
  if(!gGeoManager) {    
    AliGeomManager::LoadGeometry(fGeometryFileName.Data());
    if(!gGeoManager) { 
      printf("AliAlignmentDataFilterITS::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());
    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
	133	// Several histograms are more conveniently managed in a TList
	134	fList = new TList();
	135	fList->SetOwner();
	136
77ec32bb	137	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	138	fList->Add(fHistEvCount);
	139
	140	fntTrks = new TNtuple("fntTrks","TPC-ITS matching","ptTPC:nClsTPC:nClsITSSA:nClsITSMI:SSD1:SSD2:phi:z:dx:dy:dz:drphi:dphi:dtgl");
	141	fList->Add(fntTrks);
	142	fntPairs = new TNtuple("fntPairs","pairs at vertex","nClsITSSAin:nClsITSSAout:nClsITSMIin:nClsITSMIout:dxyITSSA:dzITSSA:dxyITSMI:dzITSMI:ptITSSAin:ptITSSAout:ptITSMIin:ptITSMIout:sigmad0MIin:sigmad0MIout");
	143	fList->Add(fntPairs);
	144	fntITSPairs = new TNtuple("fntITSPairs","pairs at vertex","nClsITSSAin:nClsITSSAout:dxyITSSA:dzITSSA:ptITSSAin:ptITSSAout:d0mu");
	145	fList->Add(fntITSPairs);
	146
	147	return;
	148	}
	149
	150	//________________________________________________________________________
	151	void AliTrackMatchingTPCITSCosmics::Exec(Option_t /option*/)
	152	{
	153	// Execute analysis for current event:
	154	// write ITS AliTrackPoints for selected tracks to fspTree
	155
	156	// load the geometry
	157	if(!gGeoManager) {
	158	AliGeomManager::LoadGeometry(fGeometryFileName.Data());
	159	if(!gGeoManager) {
	160	printf("AliAlignmentDataFilterITS::Exec(): no geometry loaded \n");
	161	return;
	162	}
	163	}
	164
	165	if(!fESD) {
	166	printf("AliTrackMatchingTPCITSCosmics::Exec(): no ESD \n");
	167	return;
	168	}
	169
	170
	171	fHistEvCount->Fill(0);
	172	PostData(0,fList);
	173
77ec32bb	174	// check if event is triggered by HLT
	175	Bool_t hltTrigg=kFALSE;
	176	if(fReadHLTESD) {
	177	AliESDInputHandler esdH = dynamic_cast<AliESDInputHandler> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
	178	AliESDEvent *hltESD = esdH->GetHLTEvent();
	179	if(!hltESD) {
	180	printf("AliTrackMatchingTPCITSCosmics::Exec(): no HLT ESD \n");
	181	return;
	182	}
	183	if(hltESD->IsHLTTriggerFired()) {fHistEvCount->Fill(5);hltTrigg=kTRUE;}
	184	}
dae0a357	185
77ec32bb	186
dae0a357	187	TString triggeredClass = fESD->GetFiredTriggerClasses();
	188	if(triggeredClass.Contains("C0SCO-ABCE-NOPF-CENT")) {
	189	fHistEvCount->Fill(1);
	190	} else {
	191	if(fOnlySPDFO) { PostData(0,fList); return; }
	192	}
77ec32bb	193
dae0a357	194	Int_t ntracks = fESD->GetNumberOfTracks();
77ec32bb	195	printf("CONTR: %d\n",fESD->GetVertex()->GetNContributors());
dae0a357	196	if(fESD->GetVertex()->GetNContributors()>=0) {
	197	fHistEvCount->Fill(2);
	198	}
	199
	200	Int_t nTrksITSSA=0,nTrksTPC=0,nTrksTPCITS=0;
	201	Int_t idxITSSA[10000],idxTPC[10000];
	202	Int_t idxITSSAin=-1,idxITSSAout=-1;
	203	Int_t idxTPCin=-1,idxTPCout=-1;
	204	for(Int_t itr=0; itr<ntracks; itr++) {
	205	AliESDtrack *t = fESD->GetTrack(itr);
	206	Int_t nClsTPC = t->GetNcls(1);
	207	Int_t nClsITS = t->GetNcls(0);
	208	if(nClsITS>=2 && nClsTPC==0) { // ITS SA
	209	idxITSSA[nTrksITSSA]=itr;
	210	nTrksITSSA++;
	211	//printf("%d\n",itr);
	212	}
	213	if(nClsTPC>=50) { // TPC
	214	idxTPC[nTrksTPC]=itr;
	215	nTrksTPC++;
	216	}
	217	if(nClsTPC>=50 && nClsITS>=2) {// TPC+ITS
	218	nTrksTPCITS++;
	219	/*printf(" --> TPC+ITS: %d + %d (%d,%d,%d,%d,%d,%d) pt = %f\n",nClsTPC,nClsITS,
	220	(Int_t)(t->HasPointOnITSLayer(0)),
	221	(Int_t)(t->HasPointOnITSLayer(1)),
	222	(Int_t)(t->HasPointOnITSLayer(2)),
	223	(Int_t)(t->HasPointOnITSLayer(3)),
	224	(Int_t)(t->HasPointOnITSLayer(4)),
	225	(Int_t)(t->HasPointOnITSLayer(5)),
	226	t->Pt());*/
	227	}
	228	}
	229
	230	//printf("nTrksTPC %d nTrksITSSA %d\n",nTrksTPC,nTrksITSSA);
77ec32bb	231	if(nTrksITSSA>0) fHistEvCount->Fill(4);
	232	if(hltTrigg && nTrksITSSA>0) fHistEvCount->Fill(7);
	233	if(nTrksTPC>0) fHistEvCount->Fill(3);
	234	if(hltTrigg && nTrksTPC>0) fHistEvCount->Fill(6);
	235	if(nTrksTPCITS>0) fHistEvCount->Fill(8);
	236	if(hltTrigg && nTrksTPCITS>0) fHistEvCount->Fill(9);
dae0a357	237
	238	if(nTrksTPC>2 \|\| nTrksTPC==0 \|\| nTrksITSSA>2 \|\| nTrksITSSA==0) { PostData(0,fList); return; }
	239
	240	for(Int_t itr=0; itr<nTrksTPC; itr++) {
	241	AliESDtrack *t = fESD->GetTrack(idxTPC[itr]);
	242	if(t->Py()>0) idxTPCin=idxTPC[itr];
	243	if(t->Py()<0) idxTPCout=idxTPC[itr];
	244	}
	245
	246
	247	for(Int_t itr=0; itr<nTrksITSSA; itr++) {
	248	AliESDtrack *t = fESD->GetTrack(idxITSSA[itr]);
	249	if(t->Py()>0) idxITSSAin=idxITSSA[itr];
	250	if(t->Py()<0) idxITSSAout=idxITSSA[itr];
	251	}
	252
	253	Double_t xyzITS[3],xyzTPC[3];
	254
	255	// analysis for inward track
	256	if(idxITSSAin>-1 && idxTPCin>-1) {
	257	//printf(" %d\n",idxITSSAin);
	258	AliESDtrack *tITS = fESD->GetTrack(idxITSSAin);
	259	const AliExternalTrackParam *outerITS = tITS->GetOuterParam();
	260	xyzITS[0]=-1000.;xyzITS[1]=-1000.;xyzITS[2]=-1000.;
	261	if(outerITS) outerITS->GetXYZAt(50.,fESD->GetMagneticField(),xyzITS);
	262	AliESDtrack *tTPC = fESD->GetTrack(idxTPCin);
	263	AliExternalTrackParam innerTPC = new AliExternalTrackParam((tTPC->GetInnerParam()));
	264	innerTPC->Rotate(tITS->GetAlpha());
	265	innerTPC->GetXYZAt(50.,fESD->GetMagneticField(),xyzTPC);
	266	fntTrks->Fill(innerTPC->Pt(),
	267	tTPC->GetNcls(1),
	268	tITS->GetNcls(0),
	269	tTPC->GetNcls(0),
	270	tITS->HasPointOnITSLayer(4),
	271	tITS->HasPointOnITSLayer(5),
	272	TMath::ATan2(xyzITS[1],xyzITS[0]),
	273	xyzITS[2],
	274	xyzITS[0]-xyzTPC[0],
	275	xyzITS[1]-xyzTPC[1],
	276	xyzITS[2]-xyzTPC[2],
	277	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])),
	278	TMath::ATan2(tITS->Py(),tITS->Px())-TMath::ATan2(innerTPC->Py(),innerTPC->Px()),
	279	tITS->Pz()/tITS->Pt()-innerTPC->Pz()/innerTPC->Pt());
	280	delete innerTPC; innerTPC=NULL;
	281	}
	282
	283	// analysis for outward track
	284	if(idxITSSAout>-1 && idxTPCout>-1) {
	285	//printf(" %d\n",idxITSSAout);
	286	AliESDtrack *tITS = fESD->GetTrack(idxITSSAout);
	287	const AliExternalTrackParam *outerITS = tITS->GetOuterParam();
	288	xyzITS[0]=-1000.;xyzITS[1]=-1000.;xyzITS[2]=-1000.;
	289	if(outerITS) outerITS->GetXYZAt(50.,fESD->GetMagneticField(),xyzITS);
	290	AliESDtrack *tTPC = fESD->GetTrack(idxTPCout);
	291	AliExternalTrackParam innerTPC = new AliExternalTrackParam((tTPC->GetInnerParam()));
	292	innerTPC->Rotate(tITS->GetAlpha());
	293	innerTPC->GetXYZAt(50.,fESD->GetMagneticField(),xyzTPC);
	294	fntTrks->Fill(innerTPC->Pt(),
	295	tTPC->GetNcls(1),
	296	tITS->GetNcls(0),
	297	tTPC->GetNcls(0),
	298	tITS->HasPointOnITSLayer(4),
	299	tITS->HasPointOnITSLayer(5),
	300	TMath::ATan2(xyzITS[1],xyzITS[0]),
301	xyzITS[2],
302	xyzITS[0]-xyzTPC[0],
303	xyzITS[1]-xyzTPC[1],
304	xyzITS[2]-xyzTPC[2],
305	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])),
306	TMath::ATan2(tITS->Py(),tITS->Px())-TMath::ATan2(innerTPC->Py(),innerTPC->Px()),
307	tITS->Pz()/tITS->Pt()-innerTPC->Pz()/innerTPC->Pt()
308	);
309	delete innerTPC; innerTPC=NULL;
310	}
311
312	// inward-outward track ITS only
313	if(idxITSSAin>-1 && idxITSSAout>-1) {
314	AliESDtrack *tITSin = fESD->GetTrack(idxITSSAin);
315	AliESDtrack *tITSout = fESD->GetTrack(idxITSSAout);
316	if(tITSin->HasPointOnITSLayer(0) &&
317	tITSin->HasPointOnITSLayer(1) &&
318	tITSout->HasPointOnITSLayer(0) &&
319	tITSout->HasPointOnITSLayer(1)) {
320	Double_t alpha = TMath::ATan2(tITSin->Py(),tITSin->Px());
321	tITSin->Propagate(alpha,0.,fESD->GetMagneticField());
322	tITSout->Propagate(alpha,0.,fESD->GetMagneticField());
323	Float_t d0[2],z0[2];
324	d0[0] = tITSin->GetY();
325	z0[0] = tITSin->GetZ();
326	d0[1] = tITSout->GetY();
327	z0[1] = tITSout->GetZ();
328	Float_t dxyITS = -(d0[0]-d0[1]);
329	Float_t dzITS = z0[0]-z0[1];
330	fntITSPairs->Fill(tITSin->GetNcls(0),tITSout->GetNcls(0),dxyITS,dzITS,tITSin->Pt(),tITSout->Pt(),TMath::Abs(d0[0]));
331	}
332	}
333
334	// inward-outward track TPC+ITS
335	if(idxITSSAin>-1 && idxTPCin>-1 && idxITSSAout>-1 && idxTPCout>-1) {
336	AliESDtrack *tTPCin = fESD->GetTrack(idxTPCin);
337	AliESDtrack *tITSin = fESD->GetTrack(idxITSSAin);
338	AliESDtrack *tTPCout = fESD->GetTrack(idxTPCout);
339	AliESDtrack *tITSout = fESD->GetTrack(idxITSSAout);
340	if(tTPCin->HasPointOnITSLayer(0) &&
341	tTPCin->HasPointOnITSLayer(1) &&
342	tITSin->HasPointOnITSLayer(0) &&
343	tITSin->HasPointOnITSLayer(1) &&
344	tTPCout->HasPointOnITSLayer(0) &&
345	tTPCout->HasPointOnITSLayer(1) &&
346	tITSout->HasPointOnITSLayer(0) &&
347	tITSout->HasPointOnITSLayer(1)) {
348	Double_t alpha = TMath::ATan2(tTPCin->Py(),tTPCin->Px());
349	tTPCin->Propagate(alpha,0.,fESD->GetMagneticField());
350	tTPCout->Propagate(alpha,0.,fESD->GetMagneticField());
351	Float_t d0[2],z0[2];
352	d0[0] = tTPCin->GetY();
353	z0[0] = tTPCin->GetZ();
354	d0[1] = tTPCout->GetY();
355	z0[1] = tTPCout->GetZ();
356	Float_t dxyTPC = -(d0[0]-d0[1]);
357	Float_t dzTPC = z0[0]-z0[1];
358	alpha = TMath::ATan2(tITSin->Py(),tITSin->Px());
359	tITSin->Propagate(alpha,0.,fESD->GetMagneticField());
360	tITSout->Propagate(alpha,0.,fESD->GetMagneticField());
361	d0[0] = tITSin->GetY();
362	z0[0] = tITSin->GetZ();
363	d0[1] = tITSout->GetY();
364	z0[1] = tITSout->GetZ();
365	Float_t dxyITS = -(d0[0]-d0[1]);
366	Float_t dzITS = z0[0]-z0[1];
367	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());
368	}
369	}
370
371	PostData(0,fList);
372
373	return;
374	}
375
376	//________________________________________________________________________
377	void AliTrackMatchingTPCITSCosmics::Terminate(Option_t /option*/)
378	{
379	// Terminate analysis
380	//
381	AliDebug(2,"AliTrackMatchingTPCITSCosmics: Terminate() \n");
382
383	fList = dynamic_cast<TList*> (GetOutputData(0));
384	if (!fList) {
385	printf("ERROR: fList not available\n");
386	return;
387	}
388
389	return;
390	}
391