[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
  //

  // 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());
    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());
	184	AliESDEvent *hltESD = esdH->GetHLTEvent();
	185	if(!hltESD) {
	186	printf("AliTrackMatchingTPCITSCosmics::Exec(): no HLT ESD \n");
	187	return;
	188	}
	189	if(hltESD->IsHLTTriggerFired()) {fHistEvCount->Fill(5);hltTrigg=kTRUE;}
	190	}
dae0a357	191
77ec32bb	192
dae0a357	193	TString triggeredClass = fESD->GetFiredTriggerClasses();
	194	if(triggeredClass.Contains("C0SCO-ABCE-NOPF-CENT")) {
	195	fHistEvCount->Fill(1);
	196	} else {
	197	if(fOnlySPDFO) { PostData(0,fList); return; }
	198	}
77ec32bb	199
dae0a357	200	Int_t ntracks = fESD->GetNumberOfTracks();
77ec32bb	201	printf("CONTR: %d\n",fESD->GetVertex()->GetNContributors());
dae0a357	202	if(fESD->GetVertex()->GetNContributors()>=0) {
	203	fHistEvCount->Fill(2);
	204	}
	205
	206	Int_t nTrksITSSA=0,nTrksTPC=0,nTrksTPCITS=0;
	207	Int_t idxITSSA[10000],idxTPC[10000];
	208	Int_t idxITSSAin=-1,idxITSSAout=-1;
	209	Int_t idxTPCin=-1,idxTPCout=-1;
	210	for(Int_t itr=0; itr<ntracks; itr++) {
	211	AliESDtrack *t = fESD->GetTrack(itr);
	212	Int_t nClsTPC = t->GetNcls(1);
	213	Int_t nClsITS = t->GetNcls(0);
	214	if(nClsITS>=2 && nClsTPC==0) { // ITS SA
	215	idxITSSA[nTrksITSSA]=itr;
	216	nTrksITSSA++;
	217	//printf("%d\n",itr);
	218	}
	219	if(nClsTPC>=50) { // TPC
	220	idxTPC[nTrksTPC]=itr;
	221	nTrksTPC++;
	222	}
	223	if(nClsTPC>=50 && nClsITS>=2) {// TPC+ITS
	224	nTrksTPCITS++;
	225	/*printf(" --> TPC+ITS: %d + %d (%d,%d,%d,%d,%d,%d) pt = %f\n",nClsTPC,nClsITS,
	226	(Int_t)(t->HasPointOnITSLayer(0)),
	227	(Int_t)(t->HasPointOnITSLayer(1)),
	228	(Int_t)(t->HasPointOnITSLayer(2)),
	229	(Int_t)(t->HasPointOnITSLayer(3)),
	230	(Int_t)(t->HasPointOnITSLayer(4)),
	231	(Int_t)(t->HasPointOnITSLayer(5)),
	232	t->Pt());*/
	233	}
	234	}
	235
	236	//printf("nTrksTPC %d nTrksITSSA %d\n",nTrksTPC,nTrksITSSA);
77ec32bb	237	if(nTrksITSSA>0) fHistEvCount->Fill(4);
	238	if(hltTrigg && nTrksITSSA>0) fHistEvCount->Fill(7);
	239	if(nTrksTPC>0) fHistEvCount->Fill(3);
	240	if(hltTrigg && nTrksTPC>0) fHistEvCount->Fill(6);
	241	if(nTrksTPCITS>0) fHistEvCount->Fill(8);
	242	if(hltTrigg && nTrksTPCITS>0) fHistEvCount->Fill(9);
dae0a357	243
	244	if(nTrksTPC>2 \|\| nTrksTPC==0 \|\| nTrksITSSA>2 \|\| nTrksITSSA==0) { PostData(0,fList); return; }
	245
	246	for(Int_t itr=0; itr<nTrksTPC; itr++) {
	247	AliESDtrack *t = fESD->GetTrack(idxTPC[itr]);
	248	if(t->Py()>0) idxTPCin=idxTPC[itr];
	249	if(t->Py()<0) idxTPCout=idxTPC[itr];
	250	}
	251
	252
	253	for(Int_t itr=0; itr<nTrksITSSA; itr++) {
	254	AliESDtrack *t = fESD->GetTrack(idxITSSA[itr]);
	255	if(t->Py()>0) idxITSSAin=idxITSSA[itr];
	256	if(t->Py()<0) idxITSSAout=idxITSSA[itr];
	257	}
	258
	259	Double_t xyzITS[3],xyzTPC[3];
	260
	261	// analysis for inward track
	262	if(idxITSSAin>-1 && idxTPCin>-1) {
	263	//printf(" %d\n",idxITSSAin);
	264	AliESDtrack *tITS = fESD->GetTrack(idxITSSAin);
	265	const AliExternalTrackParam *outerITS = tITS->GetOuterParam();
	266	xyzITS[0]=-1000.;xyzITS[1]=-1000.;xyzITS[2]=-1000.;
	267	if(outerITS) outerITS->GetXYZAt(50.,fESD->GetMagneticField(),xyzITS);
	268	AliESDtrack *tTPC = fESD->GetTrack(idxTPCin);
	269	AliExternalTrackParam innerTPC = new AliExternalTrackParam((tTPC->GetInnerParam()));
	270	innerTPC->Rotate(tITS->GetAlpha());
	271	innerTPC->GetXYZAt(50.,fESD->GetMagneticField(),xyzTPC);
	272	fntTrks->Fill(innerTPC->Pt(),
	273	tTPC->GetNcls(1),
	274	tITS->GetNcls(0),
	275	tTPC->GetNcls(0),
	276	tITS->HasPointOnITSLayer(4),
	277	tITS->HasPointOnITSLayer(5),
	278	TMath::ATan2(xyzITS[1],xyzITS[0]),
	279	xyzITS[2],
	280	xyzITS[0]-xyzTPC[0],
	281	xyzITS[1]-xyzTPC[1],
	282	xyzITS[2]-xyzTPC[2],
	283	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])),
	284	TMath::ATan2(tITS->Py(),tITS->Px())-TMath::ATan2(innerTPC->Py(),innerTPC->Px()),
	285	tITS->Pz()/tITS->Pt()-innerTPC->Pz()/innerTPC->Pt());
	286	delete innerTPC; innerTPC=NULL;
	287	}
	288
	289	// analysis for outward track
	290	if(idxITSSAout>-1 && idxTPCout>-1) {
	291	//printf(" %d\n",idxITSSAout);
	292	AliESDtrack *tITS = fESD->GetTrack(idxITSSAout);
	293	const AliExternalTrackParam *outerITS = tITS->GetOuterParam();
	294	xyzITS[0]=-1000.;xyzITS[1]=-1000.;xyzITS[2]=-1000.;
	295	if(outerITS) outerITS->GetXYZAt(50.,fESD->GetMagneticField(),xyzITS);
	296	AliESDtrack *tTPC = fESD->GetTrack(idxTPCout);
	297	AliExternalTrackParam innerTPC = new AliExternalTrackParam((tTPC->GetInnerParam()));
	298	innerTPC->Rotate(tITS->GetAlpha());
	299	innerTPC->GetXYZAt(50.,fESD->GetMagneticField(),xyzTPC);
	300	fntTrks->Fill(innerTPC->Pt(),
	301	tTPC->GetNcls(1),
	302	tITS->GetNcls(0),
	303	tTPC->GetNcls(0),
	304	tITS->HasPointOnITSLayer(4),
	305	tITS->HasPointOnITSLayer(5),
	306	TMath::ATan2(xyzITS[1],xyzITS[0]),
307	xyzITS[2],
308	xyzITS[0]-xyzTPC[0],
309	xyzITS[1]-xyzTPC[1],
310	xyzITS[2]-xyzTPC[2],
311	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])),
312	TMath::ATan2(tITS->Py(),tITS->Px())-TMath::ATan2(innerTPC->Py(),innerTPC->Px()),
313	tITS->Pz()/tITS->Pt()-innerTPC->Pz()/innerTPC->Pt()
314	);
315	delete innerTPC; innerTPC=NULL;
316	}
317
318	// inward-outward track ITS only
319	if(idxITSSAin>-1 && idxITSSAout>-1) {
320	AliESDtrack *tITSin = fESD->GetTrack(idxITSSAin);
321	AliESDtrack *tITSout = fESD->GetTrack(idxITSSAout);
322	if(tITSin->HasPointOnITSLayer(0) &&
323	tITSin->HasPointOnITSLayer(1) &&
324	tITSout->HasPointOnITSLayer(0) &&
325	tITSout->HasPointOnITSLayer(1)) {
326	Double_t alpha = TMath::ATan2(tITSin->Py(),tITSin->Px());
327	tITSin->Propagate(alpha,0.,fESD->GetMagneticField());
328	tITSout->Propagate(alpha,0.,fESD->GetMagneticField());
329	Float_t d0[2],z0[2];
330	d0[0] = tITSin->GetY();
331	z0[0] = tITSin->GetZ();
332	d0[1] = tITSout->GetY();
333	z0[1] = tITSout->GetZ();
334	Float_t dxyITS = -(d0[0]-d0[1]);
335	Float_t dzITS = z0[0]-z0[1];
336	fntITSPairs->Fill(tITSin->GetNcls(0),tITSout->GetNcls(0),dxyITS,dzITS,tITSin->Pt(),tITSout->Pt(),TMath::Abs(d0[0]));
337	}
338	}
339
340	// inward-outward track TPC+ITS
341	if(idxITSSAin>-1 && idxTPCin>-1 && idxITSSAout>-1 && idxTPCout>-1) {
342	AliESDtrack *tTPCin = fESD->GetTrack(idxTPCin);
343	AliESDtrack *tITSin = fESD->GetTrack(idxITSSAin);
344	AliESDtrack *tTPCout = fESD->GetTrack(idxTPCout);
345	AliESDtrack *tITSout = fESD->GetTrack(idxITSSAout);
346	if(tTPCin->HasPointOnITSLayer(0) &&
347	tTPCin->HasPointOnITSLayer(1) &&
348	tITSin->HasPointOnITSLayer(0) &&
349	tITSin->HasPointOnITSLayer(1) &&
350	tTPCout->HasPointOnITSLayer(0) &&
351	tTPCout->HasPointOnITSLayer(1) &&
352	tITSout->HasPointOnITSLayer(0) &&
353	tITSout->HasPointOnITSLayer(1)) {
354	Double_t alpha = TMath::ATan2(tTPCin->Py(),tTPCin->Px());
355	tTPCin->Propagate(alpha,0.,fESD->GetMagneticField());
356	tTPCout->Propagate(alpha,0.,fESD->GetMagneticField());
357	Float_t d0[2],z0[2];
358	d0[0] = tTPCin->GetY();
359	z0[0] = tTPCin->GetZ();
360	d0[1] = tTPCout->GetY();
361	z0[1] = tTPCout->GetZ();
362	Float_t dxyTPC = -(d0[0]-d0[1]);
363	Float_t dzTPC = z0[0]-z0[1];
364	alpha = TMath::ATan2(tITSin->Py(),tITSin->Px());
365	tITSin->Propagate(alpha,0.,fESD->GetMagneticField());
366	tITSout->Propagate(alpha,0.,fESD->GetMagneticField());
367	d0[0] = tITSin->GetY();
368	z0[0] = tITSin->GetZ();
369	d0[1] = tITSout->GetY();
370	z0[1] = tITSout->GetZ();
371	Float_t dxyITS = -(d0[0]-d0[1]);
372	Float_t dzITS = z0[0]-z0[1];
373	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());
374	}
375	}
376
377	PostData(0,fList);
378
379	return;
380	}
381
382	//________________________________________________________________________
383	void AliTrackMatchingTPCITSCosmics::Terminate(Option_t /option*/)
384	{
385	// Terminate analysis
386	//
387	AliDebug(2,"AliTrackMatchingTPCITSCosmics: Terminate() \n");
388
389	fList = dynamic_cast<TList*> (GetOutputData(0));
390	if (!fList) {
391	printf("ERROR: fList not available\n");
392	return;
393	}
394
395	return;
396	}
397