[u/mrichter/AliRoot.git] / TPC / AliTPCTracking.C

////////////////////////////////////////////////////////////////////////
//
// AliTPCTracking.C 
//
// date: 22.08.2002
// author: Jiri Chudoba based on code of Jourij Belikov
// version: 1.0
// description: 
//      reconstructs of tracks in TPC inthe following steps:
//         TPC cluster finding
//         TPC track finding
// input parameters: 
//        Int_t nEvents      ... nr of events to process
//        Int_t firstEventNr ... first event number (starts from 0)
//        Char_t* fileNameHits ... name of file with hits
//        Char_t* fileNameDigits .. name of file with TPC digits
//        Char_t* fileNameClusters .. name of file with TPC clusters (output)
//        Char_t* fileNameTracks .. name of file with TPC tracks (output)
//
//        default file names correspond to pp production (2002-04)
//
// History:
//
//     03.03.2003 ... SetFieldFactor moved to AliTracker class and
//                    LoadTPCParam moved to AliTPC class
//                    TString replaced by Char_t*
//
//     20.11.2002 ... Changes due to a changed interface of AliTPCtracker. 
//                    Use Riostream.h instead of iostream.h
//
//     22.08.2002 ... first version
//
////////////////////////////////////////////////////////////////////////

#if !defined(__CINT__) || defined(__MAKECINT__)
#include "Riostream.h"
#include "TTree.h"
#include "TSystem.h"
#include "TArrayF.h"
#include "TPC/alles.h"
#include "TPC/AliTPCtracker.h"
#include "TPC/AliTPCclusterer.h"
#include "TPC/AliTPC.h"
#include "STEER/AliRun.h"
#include "STEER/AliHeader.h"
#include "STEER/AliGenEventHeader.h"
#include "STEER/AliMagF.h"
#include "STEER/AliTracker.h"

#endif

Int_t gDEBUG = 2;

Int_t TPCFindClusters(Int_t nEvents=1, Int_t firstEvent=0,
		      Char_t* fileNameDigits="rfio:galiceSDR.root", 
		      Char_t* fileNameClusters="tpc.clusters.root");
Int_t TPCFindClusters(Int_t nEvents, Int_t firstEvent,
		      TFile* fileDigits, TFile* fileClusters, 
		      AliTPCParam* paramTPC=0);
Int_t TPCFindTracks(Int_t nEvents=1, Int_t firstEvent=0,
		    Char_t* fileNameClusters="tpc.clusters.root",
		    Char_t* fileNameTracks="tpc.tracks.root");
Int_t TPCFindTracks(Int_t nEvents, Int_t firstEvent,
		    TFile* fileClusters, TFile* fileTracks,
		    AliTPCParam* paramTPC=0);

void FindVertex(Int_t iEvent, Double_t *vertex);
void PrintVertex(TArrayF &primaryVertex);

Int_t AliTPCTracking(Int_t nEvents=1, Int_t firstEvent=0,
		     Char_t* fileNameHits="rfio:galice.root",
		     Char_t* fileNameDigits="rfio:galiceSDR.root",
		     Char_t* fileNameClusters="tpc.clusters.root",
		     Char_t* fileNameTracks="tpc.tracks.root");

////////////////////////////////////////////////////////////////////////
Int_t AliTPCTracking( Int_t nEvents, Int_t firstEvent,
		      Char_t* fileNameHits,
		      Char_t* fileNameDigits,
		      Char_t* fileNameClusters,
		      Char_t* fileNameTracks) {

  AliTracker::SetFieldFactor(fileNameHits,kFALSE);

// ********** Find TPC clusters *********** //
  if (TPCFindClusters(nEvents,firstEvent,fileNameDigits,fileNameClusters)) {
    cerr<<"Failed to get TPC clusters: !\n";
    return 1;
  }      

// ********** Find TPC tracks *********** //
  if (TPCFindTracks(nEvents,firstEvent,fileNameClusters,fileNameTracks)) {
    cerr<<"Failed to get TPC tracks !\n";
    return 2;
  }

  return 0;
}

////////////////////////////////////////////////////////////////////////
Int_t TPCFindClusters(Int_t nEvents, Int_t firstEvent,
		      Char_t* fileNameDigits, Char_t* fileNameClusters) {
  
  Int_t rc;
  const Char_t *name="TPCFindClusters";
  if (gDEBUG>1) cout<<name<<" starts...\n";
  if (gDEBUG>1) gBenchmark->Start(name);
  TFile *fileClusters = TFile::Open(fileNameClusters,"recreate");
  TFile *fileDigits = TFile::Open(fileNameDigits);
  if (!fileDigits->IsOpen()) {
    cerr<<"Cannnot open "<<fileNameDigits<<" !\n"; 
    return 1;
  }
  if (!fileClusters->IsOpen()) {
    cerr<<"Cannnot open "<<fileNameClusters<<" !\n"; 
    return 1;
  }

  rc = TPCFindClusters(nEvents,firstEvent,fileDigits,fileClusters);

  fileDigits->Close();
  fileClusters->Close();
  delete fileDigits;
  delete fileClusters;
  if (gDEBUG>1) gBenchmark->Stop(name);
  if (gDEBUG>1) gBenchmark->Show(name);

  return rc;
}
////////////////////////////////////////////////////////////////////////
Int_t TPCFindClusters(Int_t nEvents, Int_t firstEvent,
		      TFile* fileDigits, TFile* fileClusters, 
		      AliTPCParam* paramTPC) {

  fileDigits->cd();
  if (!paramTPC) paramTPC = AliTPC::LoadTPCParam(fileDigits);
  if (!paramTPC) return 1;

  for (Int_t iEvent = firstEvent; iEvent < firstEvent+nEvents; iEvent++){
    if (gDEBUG > 2) cout<<"TPCFindClusters: event "<<iEvent<<endl;
    AliTPCclusterer::Digits2Clusters(paramTPC, fileClusters, iEvent);
  }
  return 0;
}
////////////////////////////////////////////////////////////////////////
Int_t TPCFindTracks(Int_t nEvents, Int_t firstEvent,
		    Char_t* fileNameClusters, Char_t* fileNameTracks) {

  Int_t rc = 0;
  const Char_t *name="TPCFindTracks";
  if (gDEBUG>1) cout<<name<<" starts"<<endl;
  if (gDEBUG>1) gBenchmark->Start(name);
  TFile *fileTracks = TFile::Open(fileNameTracks,"recreate");
  TFile *fileClusters =TFile::Open(fileNameClusters);

  rc = TPCFindTracks(nEvents, firstEvent, fileClusters, fileTracks);

  fileClusters->Close();
  fileTracks->Close();
  delete fileClusters;
  delete fileTracks;
  if (gDEBUG>1) gBenchmark->Stop(name);
  if (gDEBUG>1) gBenchmark->Show(name);
  return rc;

}
////////////////////////////////////////////////////////////////////////
Int_t TPCFindTracks(Int_t nEvents, Int_t firstEvent,
		    TFile *fileClusters, TFile * fileTracks,
		    AliTPCParam* paramTPC) {

  Int_t rc = 0;
  if (!paramTPC) paramTPC = AliTPC::LoadTPCParam(fileClusters);
  if (!paramTPC) return 1;

  for (Int_t iEvent = firstEvent; iEvent < firstEvent+nEvents; iEvent++){
    if (gDEBUG > 2) cout<<"TPCFindTracks: event "<<iEvent<<endl;
    AliTPCtracker *tracker = new AliTPCtracker(paramTPC);
    tracker->SetEventNumber(iEvent);
    Double_t vertex[3];
    FindVertex(iEvent,vertex);
    tracker->SetVertex(vertex);
    fileClusters->cd();
    rc = tracker->Clusters2Tracks(0,fileTracks);
    delete tracker;
  }
  return rc;
}

////////////////////////////////////////////////////////////////////////
void FindVertex(Int_t eventNr, Double_t *vertex) {

  vertex[0] = vertex[1] = vertex[2] = 0.;
  if (!gAlice) {
    cerr<<"gAlice was not found! Using vertex position (0,0,0).\n";
    return;
  }
  
  gAlice->GetEvent(eventNr);
  AliHeader *header = gAlice->GetHeader();
  if (!header) {
    cerr<<"header was not found!\n";
    return;
  } 
  AliGenEventHeader* genEventHeader = header->GenEventHeader();
  if (!genEventHeader) {
    cerr<<"AliGenEventHeader was not found!\n";
    return;
  } 

  TArrayF primaryVertex(3);
  genEventHeader->PrimaryVertex(primaryVertex);
  PrintVertex(primaryVertex);
  vertex[0] = static_cast<Double_t>(primaryVertex[0]);
  vertex[1] = static_cast<Double_t>(primaryVertex[1]);
  vertex[2] = static_cast<Double_t>(primaryVertex[2]);
//  delete header;
  delete genEventHeader;
  return;
     
}
////////////////////////////////////////////////////////////////////////
void PrintVertex(TArrayF &primaryVertex) 
{
  cout <<"Vertex: "
       <<primaryVertex[0]<<" "
       <<primaryVertex[1]<<" "
       <<primaryVertex[2]<<" "<<endl;
  exit;
} 
////////////////////////////////////////////////////////////////////////
Commit	Line	Data
171e422d	1	////////////////////////////////////////////////////////////////////////
	2	//
	3	// AliTPCTracking.C
	4	//
	5	// date: 22.08.2002
	6	// author: Jiri Chudoba based on code of Jourij Belikov
	7	// version: 1.0
	8	// description:
	9	// reconstructs of tracks in TPC inthe following steps:
	10	// TPC cluster finding
	11	// TPC track finding
	12	// input parameters:
	13	// Int_t nEvents ... nr of events to process
	14	// Int_t firstEventNr ... first event number (starts from 0)
2d813ef1	15	// Char_t* fileNameHits ... name of file with hits
	16	// Char_t* fileNameDigits .. name of file with TPC digits
	17	// Char_t* fileNameClusters .. name of file with TPC clusters (output)
	18	// Char_t* fileNameTracks .. name of file with TPC tracks (output)
171e422d	19	//
	20	// default file names correspond to pp production (2002-04)
	21	//
	22	// History:
	23	//
2d813ef1	24	// 03.03.2003 ... SetFieldFactor moved to AliTracker class and
	25	// LoadTPCParam moved to AliTPC class
	26	// TString replaced by Char_t*
	27	//
73f6cbc2	28	// 20.11.2002 ... Changes due to a changed interface of AliTPCtracker.
	29	// Use Riostream.h instead of iostream.h
	30	//
171e422d	31	// 22.08.2002 ... first version
	32	//
	33	////////////////////////////////////////////////////////////////////////
	34
	35	#if !defined(__CINT__) \|\| defined(__MAKECINT__)
73f6cbc2	36	#include "Riostream.h"
171e422d	37	#include "TTree.h"
	38	#include "TSystem.h"
	39	#include "TArrayF.h"
	40	#include "TPC/alles.h"
	41	#include "TPC/AliTPCtracker.h"
	42	#include "TPC/AliTPCclusterer.h"
2d813ef1	43	#include "TPC/AliTPC.h"
171e422d	44	#include "STEER/AliRun.h"
	45	#include "STEER/AliHeader.h"
	46	#include "STEER/AliGenEventHeader.h"
	47	#include "STEER/AliMagF.h"
2d813ef1	48	#include "STEER/AliTracker.h"
171e422d	49
	50	#endif
	51
	52	Int_t gDEBUG = 2;
	53
	54	Int_t TPCFindClusters(Int_t nEvents=1, Int_t firstEvent=0,
2d813ef1	55	Char_t* fileNameDigits="rfio:galiceSDR.root",
2d813ef1	56	Char_t* fileNameClusters="tpc.clusters.root");
171e422d	57	Int_t TPCFindClusters(Int_t nEvents, Int_t firstEvent,
	58	TFile* fileDigits, TFile* fileClusters,
	59	AliTPCParam* paramTPC=0);
	60	Int_t TPCFindTracks(Int_t nEvents=1, Int_t firstEvent=0,
2d813ef1	61	Char_t* fileNameClusters="tpc.clusters.root",
2d813ef1	62	Char_t* fileNameTracks="tpc.tracks.root");
171e422d	63	Int_t TPCFindTracks(Int_t nEvents, Int_t firstEvent,
	64	TFile* fileClusters, TFile* fileTracks,
	65	AliTPCParam* paramTPC=0);
	66
171e422d	67	void FindVertex(Int_t iEvent, Double_t *vertex);
171e422d	68	void PrintVertex(TArrayF &primaryVertex);
171e422d	69
171e422d	70	Int_t AliTPCTracking(Int_t nEvents=1, Int_t firstEvent=0,
2d813ef1	71	Char_t* fileNameHits="rfio:galice.root",
	72	Char_t* fileNameDigits="rfio:galiceSDR.root",
	73	Char_t* fileNameClusters="tpc.clusters.root",
	74	Char_t* fileNameTracks="tpc.tracks.root");
171e422d	75
	76	////////////////////////////////////////////////////////////////////////
	77	Int_t AliTPCTracking( Int_t nEvents, Int_t firstEvent,
2d813ef1	78	Char_t* fileNameHits,
	79	Char_t* fileNameDigits,
	80	Char_t* fileNameClusters,
	81	Char_t* fileNameTracks) {
171e422d	82
2d813ef1	83	AliTracker::SetFieldFactor(fileNameHits,kFALSE);
171e422d	84
	85	// ******** Find TPC clusters ********* //
	86	if (TPCFindClusters(nEvents,firstEvent,fileNameDigits,fileNameClusters)) {
	87	cerr<<"Failed to get TPC clusters: !\n";
	88	return 1;
	89	}
	90
	91	// ******** Find TPC tracks ********* //
	92	if (TPCFindTracks(nEvents,firstEvent,fileNameClusters,fileNameTracks)) {
	93	cerr<<"Failed to get TPC tracks !\n";
	94	return 2;
	95	}
	96
	97	return 0;
	98	}
	99
	100	////////////////////////////////////////////////////////////////////////
	101	Int_t TPCFindClusters(Int_t nEvents, Int_t firstEvent,
2d813ef1	102	Char_t* fileNameDigits, Char_t* fileNameClusters) {
171e422d	103
	104	Int_t rc;
	105	const Char_t *name="TPCFindClusters";
	106	if (gDEBUG>1) cout<<name<<" starts...\n";
	107	if (gDEBUG>1) gBenchmark->Start(name);
	108	TFile *fileClusters = TFile::Open(fileNameClusters,"recreate");
	109	TFile *fileDigits = TFile::Open(fileNameDigits);
	110	if (!fileDigits->IsOpen()) {
	111	cerr<<"Cannnot open "<<fileNameDigits<<" !\n";
	112	return 1;
	113	}
	114	if (!fileClusters->IsOpen()) {
	115	cerr<<"Cannnot open "<<fileNameClusters<<" !\n";
	116	return 1;
	117	}
	118
	119	rc = TPCFindClusters(nEvents,firstEvent,fileDigits,fileClusters);
	120
	121	fileDigits->Close();
	122	fileClusters->Close();
	123	delete fileDigits;
	124	delete fileClusters;
	125	if (gDEBUG>1) gBenchmark->Stop(name);
	126	if (gDEBUG>1) gBenchmark->Show(name);
	127
	128	return rc;
	129	}
	130	////////////////////////////////////////////////////////////////////////
	131	Int_t TPCFindClusters(Int_t nEvents, Int_t firstEvent,
	132	TFile* fileDigits, TFile* fileClusters,
	133	AliTPCParam* paramTPC) {
	134
	135	fileDigits->cd();
2d813ef1	136	if (!paramTPC) paramTPC = AliTPC::LoadTPCParam(fileDigits);
171e422d	137	if (!paramTPC) return 1;
	138
	139	for (Int_t iEvent = firstEvent; iEvent < firstEvent+nEvents; iEvent++){
	140	if (gDEBUG > 2) cout<<"TPCFindClusters: event "<<iEvent<<endl;
	141	AliTPCclusterer::Digits2Clusters(paramTPC, fileClusters, iEvent);
	142	}
	143	return 0;
	144	}
	145	////////////////////////////////////////////////////////////////////////
	146	Int_t TPCFindTracks(Int_t nEvents, Int_t firstEvent,
2d813ef1	147	Char_t* fileNameClusters, Char_t* fileNameTracks) {
171e422d	148
	149	Int_t rc = 0;
	150	const Char_t *name="TPCFindTracks";
	151	if (gDEBUG>1) cout<<name<<" starts"<<endl;
	152	if (gDEBUG>1) gBenchmark->Start(name);
	153	TFile *fileTracks = TFile::Open(fileNameTracks,"recreate");
	154	TFile *fileClusters =TFile::Open(fileNameClusters);
	155
	156	rc = TPCFindTracks(nEvents, firstEvent, fileClusters, fileTracks);
	157
	158	fileClusters->Close();
	159	fileTracks->Close();
	160	delete fileClusters;
	161	delete fileTracks;
	162	if (gDEBUG>1) gBenchmark->Stop(name);
	163	if (gDEBUG>1) gBenchmark->Show(name);
	164	return rc;
	165
	166	}
	167	////////////////////////////////////////////////////////////////////////
	168	Int_t TPCFindTracks(Int_t nEvents, Int_t firstEvent,
	169	TFile fileClusters, TFile fileTracks,
	170	AliTPCParam* paramTPC) {
	171
	172	Int_t rc = 0;
2d813ef1	173	if (!paramTPC) paramTPC = AliTPC::LoadTPCParam(fileClusters);
171e422d	174	if (!paramTPC) return 1;
	175
	176	for (Int_t iEvent = firstEvent; iEvent < firstEvent+nEvents; iEvent++){
	177	if (gDEBUG > 2) cout<<"TPCFindTracks: event "<<iEvent<<endl;
73f6cbc2	178	AliTPCtracker *tracker = new AliTPCtracker(paramTPC);
73f6cbc2	179	tracker->SetEventNumber(iEvent);
171e422d	180	Double_t vertex[3];
	181	FindVertex(iEvent,vertex);
	182	tracker->SetVertex(vertex);
	183	fileClusters->cd();
	184	rc = tracker->Clusters2Tracks(0,fileTracks);
	185	delete tracker;
	186	}
	187	return rc;
	188	}
	189
171e422d	190	////////////////////////////////////////////////////////////////////////
	191	void FindVertex(Int_t eventNr, Double_t *vertex) {
	192
	193	vertex[0] = vertex[1] = vertex[2] = 0.;
	194	if (!gAlice) {
	195	cerr<<"gAlice was not found! Using vertex position (0,0,0).\n";
	196	return;
	197	}
	198
	199	gAlice->GetEvent(eventNr);
	200	AliHeader *header = gAlice->GetHeader();
	201	if (!header) {
	202	cerr<<"header was not found!\n";
	203	return;
	204	}
	205	AliGenEventHeader* genEventHeader = header->GenEventHeader();
	206	if (!genEventHeader) {
	207	cerr<<"AliGenEventHeader was not found!\n";
	208	return;
	209	}
	210
	211	TArrayF primaryVertex(3);
	212	genEventHeader->PrimaryVertex(primaryVertex);
	213	PrintVertex(primaryVertex);
	214	vertex[0] = static_cast<Double_t>(primaryVertex[0]);
	215	vertex[1] = static_cast<Double_t>(primaryVertex[1]);
	216	vertex[2] = static_cast<Double_t>(primaryVertex[2]);
	217	// delete header;
	218	delete genEventHeader;
	219	return;
	220
	221	}
	222	////////////////////////////////////////////////////////////////////////
	223	void PrintVertex(TArrayF &primaryVertex)
	224	{
	225	cout <<"Vertex: "
	226	<<primaryVertex[0]<<" "
	227	<<primaryVertex[1]<<" "
	228	<<primaryVertex[2]<<" "<<endl;
	229	exit;
	230	}
	231	////////////////////////////////////////////////////////////////////////