[u/mrichter/AliRoot.git] / PWG2 / FEMTOSCOPY / AliFemto / Reader / AliFemtoEventReaderESDChain.cxx

#include "AliFemtoEventReaderESDChain.h"

#include "TFile.h"
#include "TTree.h"
#include "AliESD.h"
#include "AliESDtrack.h"

#include "Infrastructure/AliFmPhysicalHelixD.h"
#include "Infrastructure/AliFmThreeVectorF.h"

#include "Base/SystemOfUnits.h"

#include "Infrastructure/AliFemtoEvent.h"

ClassImp(AliFemtoEventReaderESDChain)

#if !(ST_NO_NAMESPACES)
  using namespace units;
#endif

using namespace std;
//____________________________
//constructor with 0 parameters , look at default settings 
AliFemtoEventReaderESDChain::AliFemtoEventReaderESDChain():
  fFileName(" "),
  fConstrained(true),
  fNumberofEvent(0),
  fCurEvent(0),
  fCurFile(0),
  fEvent(0x0),
  fEventFriend(0),
  fSharedList(0x0),
  fClusterPerPadrow(0x0)
{
  fClusterPerPadrow = (list<Int_t> **) malloc(sizeof(list<Int_t> *) * AliESDfriendTrack::kMaxTPCcluster);
  for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
    fClusterPerPadrow[tPad] = new list<Int_t>();
  }
  fSharedList = (list<Int_t> **) malloc(sizeof(list<Int_t> *) * AliESDfriendTrack::kMaxTPCcluster);
  for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
    fSharedList[tPad] = new list<Int_t>();
  }
}

//__________________
// Copy constructor
AliFemtoEventReaderESDChain::AliFemtoEventReaderESDChain(const AliFemtoEventReaderESDChain& aReader):
  fFileName(" "),
  fConstrained(true),
  fNumberofEvent(0),
  fCurEvent(0),
  fCurFile(0),
  fEvent(0x0),
  fEventFriend(0),
  fSharedList(0x0),
  fClusterPerPadrow(0x0)
{
  fConstrained = aReader.fConstrained;
  fNumberofEvent = aReader.fNumberofEvent;
  fCurEvent = aReader.fCurEvent;
  fCurFile = aReader.fCurFile;
  fEvent = new AliESD(*aReader.fEvent);
  fEventFriend = aReader.fEventFriend;
  fClusterPerPadrow = (list<Int_t> **) malloc(sizeof(list<Int_t> *) * AliESDfriendTrack::kMaxTPCcluster);
  for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
    fClusterPerPadrow[tPad] = new list<Int_t>();
    list<Int_t>::iterator iter;
    for (iter=aReader.fClusterPerPadrow[tPad]->begin(); iter!=aReader.fClusterPerPadrow[tPad]->end(); iter++) {
      fClusterPerPadrow[tPad]->push_back(*iter);
    }
  }
  fSharedList = (list<Int_t> **) malloc(sizeof(list<Int_t> *) * AliESDfriendTrack::kMaxTPCcluster);
  for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
    fSharedList[tPad] = new list<Int_t>();
    list<Int_t>::iterator iter;
    for (iter=aReader.fSharedList[tPad]->begin(); iter!=aReader.fSharedList[tPad]->end(); iter++) {
      fSharedList[tPad]->push_back(*iter);
    }
  }
}
//__________________
//Destructor
AliFemtoEventReaderESDChain::~AliFemtoEventReaderESDChain()
{
  delete fEvent;

  for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
    fClusterPerPadrow[tPad]->clear();
    delete fClusterPerPadrow[tPad];
  }
  delete [] fClusterPerPadrow;
  for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
    fSharedList[tPad]->clear();
    delete fSharedList[tPad];
  }
  delete [] fSharedList;
}

//__________________
// Assignment operator
AliFemtoEventReaderESDChain& AliFemtoEventReaderESDChain::operator=(const AliFemtoEventReaderESDChain& aReader)
{
  if (this == &aReader)
    return *this;

  fConstrained = aReader.fConstrained;
  fNumberofEvent = aReader.fNumberofEvent;
  fCurEvent = aReader.fCurEvent;
  fCurFile = aReader.fCurFile;
  if (fEvent) delete fEvent;
  fEvent = new AliESD(*aReader.fEvent);

  fEventFriend = aReader.fEventFriend;
  
  if (fClusterPerPadrow) {
    for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
      fClusterPerPadrow[tPad]->clear();
      delete fClusterPerPadrow[tPad];
    }
    delete [] fClusterPerPadrow;
  }
  
  if (fSharedList) {
    for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
      fSharedList[tPad]->clear();
      delete fSharedList[tPad];
    }
    delete [] fSharedList;
  }

  fClusterPerPadrow = (list<Int_t> **) malloc(sizeof(list<Int_t> *) * AliESDfriendTrack::kMaxTPCcluster);
  for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
    fClusterPerPadrow[tPad] = new list<Int_t>();
    list<Int_t>::iterator iter;
    for (iter=aReader.fClusterPerPadrow[tPad]->begin(); iter!=aReader.fClusterPerPadrow[tPad]->end(); iter++) {
      fClusterPerPadrow[tPad]->push_back(*iter);
    }
  }
  fSharedList = (list<Int_t> **) malloc(sizeof(list<Int_t> *) * AliESDfriendTrack::kMaxTPCcluster);
  for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
    fSharedList[tPad] = new list<Int_t>();
    list<Int_t>::iterator iter;
    for (iter=aReader.fSharedList[tPad]->begin(); iter!=aReader.fSharedList[tPad]->end(); iter++) {
      fSharedList[tPad]->push_back(*iter);
    }
  }
  
  return *this;
}
//__________________
// Simple report
AliFemtoString AliFemtoEventReaderESDChain::Report()
{
  AliFemtoString temp = "\n This is the AliFemtoEventReaderESDChain\n";
  return temp;
}

//__________________
// Select whether to read constrained or not constrained momentum
void AliFemtoEventReaderESDChain::SetConstrained(const bool constrained)
{
  fConstrained=constrained;
}
//__________________
// Check whether we read constrained or not constrained momentum
bool AliFemtoEventReaderESDChain::GetConstrained() const
{
  return fConstrained;
}
//__________________
AliFemtoEvent* AliFemtoEventReaderESDChain::ReturnHbtEvent()
{
  // Get the event, read all the relevant information
  // and fill the AliFemtoEvent class
  // Returns a valid AliFemtoEvent
  AliFemtoEvent *hbtEvent = 0;
  string tFriendFileName;

  // Get the friend information
  cout<<"starting to read event "<<fCurEvent<<endl;
  fEvent->SetESDfriend(fEventFriend);
  vector<int> label_table;//to check labels
	
  hbtEvent = new AliFemtoEvent;
  //setting basic things
  //  hbtEvent->SetEventNumber(fEvent->GetEventNumber());
  hbtEvent->SetRunNumber(fEvent->GetRunNumber());
  //hbtEvent->SetNumberOfTracks(fEvent->GetNumberOfTracks());
  hbtEvent->SetMagneticField(fEvent->GetMagneticField()*kilogauss);//to check if here is ok
  hbtEvent->SetZDCN1Energy(fEvent->GetZDCN1Energy());
  hbtEvent->SetZDCP1Energy(fEvent->GetZDCP1Energy());
  hbtEvent->SetZDCN2Energy(fEvent->GetZDCN2Energy());
  hbtEvent->SetZDCP2Energy(fEvent->GetZDCP2Energy());
  hbtEvent->SetZDCEMEnergy(fEvent->GetZDCEMEnergy());
  hbtEvent->SetZDCParticipants(fEvent->GetZDCParticipants());
  hbtEvent->SetTriggerMask(fEvent->GetTriggerMask());
  hbtEvent->SetTriggerCluster(fEvent->GetTriggerCluster());
	
  //Vertex
  double fV1[3];
  fEvent->GetVertex()->GetXYZ(fV1);

  AliFmThreeVectorF vertex(fV1[0],fV1[1],fV1[2]);
  hbtEvent->SetPrimVertPos(vertex);
	
  //starting to reading tracks
  int nofTracks=0;  //number of reconstructed tracks in event
  nofTracks=fEvent->GetNumberOfTracks();
  int realnofTracks=0;//number of track which we use ina analysis

  // Clear the shared cluster list
  for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
    fClusterPerPadrow[tPad]->clear();
  }
  for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
    fSharedList[tPad]->clear();
  }


  for (int i=0;i<nofTracks;i++) {
    const AliESDtrack *esdtrack=fEvent->GetTrack(i);//getting next track

    list<Int_t>::iterator tClustIter;

    Int_t tTrackIndices[AliESDfriendTrack::kMaxTPCcluster];
    Int_t tNClusters = esdtrack->GetTPCclusters(tTrackIndices);
    for (int tNcl=0; tNcl<AliESDfriendTrack::kMaxTPCcluster; tNcl++) {
      if (tTrackIndices[tNcl] >= 0) {
	tClustIter = find(fClusterPerPadrow[tNcl]->begin(), fClusterPerPadrow[tNcl]->end(), tTrackIndices[tNcl]);
	if (tClustIter == fClusterPerPadrow[tNcl]->end()) {
	  fClusterPerPadrow[tNcl]->push_back(tTrackIndices[tNcl]);
	}
	else {
	  fSharedList[tNcl]->push_back(tTrackIndices[tNcl]);
	}
      }
    }
      
  }

  for (int i=0;i<nofTracks;i++)
    {
      bool  good_momentum=true; //flaga to chcek if we can read momentum of this track
		
      AliFemtoTrack* trackCopy = new AliFemtoTrack();	
      const AliESDtrack *esdtrack=fEvent->GetTrack(i);//getting next track
      //      const AliESDfriendTrack *tESDfriendTrack = esdtrack->GetFriendTrack();

      trackCopy->SetCharge((short)esdtrack->GetSign());

      //in aliroot we have AliPID 
      //0-electron 1-muon 2-pion 3-kaon 4-proton 5-photon 6-pi0 7-neutron 8-kaon0 9-eleCon   
      //we use only 5 first
      double esdpid[5];
      esdtrack->GetESDpid(esdpid);
      trackCopy->SetPidProbElectron(esdpid[0]);
      trackCopy->SetPidProbMuon(esdpid[1]);
      trackCopy->SetPidProbPion(esdpid[2]);
      trackCopy->SetPidProbKaon(esdpid[3]);
      trackCopy->SetPidProbProton(esdpid[4]);
						
      double pxyz[3];
      if (fConstrained==true)		    
	good_momentum=esdtrack->GetConstrainedPxPyPz(pxyz); //reading constrained momentum
      else
	good_momentum=esdtrack->GetPxPyPz(pxyz);//reading noconstarined momentum

      AliFemtoThreeVector v(pxyz[0],pxyz[1],pxyz[2]);
      if (v.mag() == 0) {
	//	cout << "Found 0 momentum ???? " <<endl;
	delete trackCopy;
	continue;
      }
      trackCopy->SetP(v);//setting momentum
      trackCopy->SetPt(sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]));
      const AliFmThreeVectorD p(pxyz[0],pxyz[1],pxyz[2]);
      const AliFmThreeVectorD origin(fV1[0],fV1[1],fV1[2]);
      //setting helix I do not if it is ok
      AliFmPhysicalHelixD helix(p,origin,(double)(fEvent->GetMagneticField())*kilogauss,(double)(trackCopy->Charge())); 
      trackCopy->SetHelix(helix);
	    	
      trackCopy->SetTrackId(esdtrack->GetID());
      trackCopy->SetFlags(esdtrack->GetStatus());
      //trackCopy->SetLabel(esdtrack->GetLabel());
		
      //some stuff which could be useful 
      float impact[2];
      float covimpact[3];
      esdtrack->GetImpactParameters(impact,covimpact);
      trackCopy->SetImpactD(impact[0]);
      trackCopy->SetImpactZ(impact[1]);
      trackCopy->SetCdd(covimpact[0]);
      trackCopy->SetCdz(covimpact[1]);
      trackCopy->SetCzz(covimpact[2]);
      trackCopy->SetITSchi2(esdtrack->GetITSchi2());    
      trackCopy->SetITSncls(esdtrack->GetNcls(0));     
      trackCopy->SetTPCchi2(esdtrack->GetTPCchi2());       
      trackCopy->SetTPCncls(esdtrack->GetTPCNcls());       
      trackCopy->SetTPCnclsF(esdtrack->GetTPCNclsF());      
      trackCopy->SetTPCsignalN((short)esdtrack->GetTPCsignalN()); //due to bug in aliesdtrack class   
      trackCopy->SetTPCsignalS(esdtrack->GetTPCsignalSigma()); 

      // Fill cluster per padrow information
      Int_t tTrackIndices[AliESDfriendTrack::kMaxTPCcluster];
      Int_t tNClusters = esdtrack->GetTPCclusters(tTrackIndices);
      for (int tNcl=0; tNcl<AliESDfriendTrack::kMaxTPCcluster; tNcl++) {
	if (tTrackIndices[tNcl] > 0)
	  trackCopy->SetTPCcluster(tNcl, 1);
	else
	  trackCopy->SetTPCcluster(tNcl, 0);
      }
      
      // Fill shared cluster information
      list<Int_t>::iterator tClustIter;

      for (int tNcl=0; tNcl<AliESDfriendTrack::kMaxTPCcluster; tNcl++) {
	if (tTrackIndices[tNcl] > 0) {
	  tClustIter = find(fSharedList[tNcl]->begin(), fSharedList[tNcl]->end(), tTrackIndices[tNcl]);
	  if (tClustIter != fSharedList[tNcl]->end()) {
	    trackCopy->SetTPCshared(tNcl, 1);
	    cout << "Event next" <<  endl;
	    cout << "Track: " << i << endl;
	    cout << "Shared cluster: " << tNcl << " " << tTrackIndices[tNcl] << endl;
	  }
	  else {
	    trackCopy->SetTPCshared(tNcl, 0);
	  }
	}
      }

      if (good_momentum==true)
	{
	  hbtEvent->TrackCollection()->push_back(trackCopy);//adding track to analysis
	  realnofTracks++;//real number of tracks
	  //	  delete trackCopy;
	}
      else
	{
	  delete  trackCopy;
	}
		
    }

  hbtEvent->SetNumberOfTracks(realnofTracks);//setting number of track which we read in event	
  fCurEvent++;	
  cout<<"end of reading nt "<<nofTracks<<" real number "<<realnofTracks<<endl;
  return hbtEvent; 
}
//___________________
// The chain loads the ESD for us
// You must provide the address where it can be found
void AliFemtoEventReaderESDChain::SetESDSource(AliESD *aESD)
{
  fEvent = aESD;
}
//___________________
// We need the ESD tree to obtain
// information about the friend file location
void AliFemtoEventReaderESDChain::SetESDfriendSource(AliESDfriend *aFriend)
{
  fEventFriend = aFriend;
}
Commit	Line	Data
67427ff7	1	#include "AliFemtoEventReaderESDChain.h"
	2
	3	#include "TFile.h"
	4	#include "TTree.h"
	5	#include "AliESD.h"
	6	#include "AliESDtrack.h"
	7
b2f60a91	8	#include "Infrastructure/AliFmPhysicalHelixD.h"
b2f60a91	9	#include "Infrastructure/AliFmThreeVectorF.h"
67427ff7	10
b2f60a91	11	#include "Base/SystemOfUnits.h"
67427ff7	12
b2f60a91	13	#include "Infrastructure/AliFemtoEvent.h"
67427ff7	14
	15	ClassImp(AliFemtoEventReaderESDChain)
	16
	17	#if !(ST_NO_NAMESPACES)
	18	using namespace units;
	19	#endif
	20
	21	using namespace std;
	22	//____________________________
	23	//constructor with 0 parameters , look at default settings
	24	AliFemtoEventReaderESDChain::AliFemtoEventReaderESDChain():
	25	fFileName(" "),
	26	fConstrained(true),
	27	fNumberofEvent(0),
	28	fCurEvent(0),
	29	fCurFile(0),
	30	fEvent(0x0),
0215f606	31	fEventFriend(0),
	32	fSharedList(0x0),
	33	fClusterPerPadrow(0x0)
67427ff7	34	{
	35	fClusterPerPadrow = (list<Int_t> *) malloc(sizeof(list<Int_t> ) * AliESDfriendTrack::kMaxTPCcluster);
	36	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	37	fClusterPerPadrow[tPad] = new list<Int_t>();
	38	}
	39	fSharedList = (list<Int_t> *) malloc(sizeof(list<Int_t> ) * AliESDfriendTrack::kMaxTPCcluster);
	40	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	41	fSharedList[tPad] = new list<Int_t>();
	42	}
	43	}
	44
0215f606	45	//__________________
	46	// Copy constructor
	47	AliFemtoEventReaderESDChain::AliFemtoEventReaderESDChain(const AliFemtoEventReaderESDChain& aReader):
	48	fFileName(" "),
	49	fConstrained(true),
	50	fNumberofEvent(0),
	51	fCurEvent(0),
	52	fCurFile(0),
	53	fEvent(0x0),
	54	fEventFriend(0),
	55	fSharedList(0x0),
	56	fClusterPerPadrow(0x0)
	57	{
	58	fConstrained = aReader.fConstrained;
	59	fNumberofEvent = aReader.fNumberofEvent;
	60	fCurEvent = aReader.fCurEvent;
	61	fCurFile = aReader.fCurFile;
	62	fEvent = new AliESD(*aReader.fEvent);
	63	fEventFriend = aReader.fEventFriend;
	64	fClusterPerPadrow = (list<Int_t> *) malloc(sizeof(list<Int_t> ) * AliESDfriendTrack::kMaxTPCcluster);
	65	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	66	fClusterPerPadrow[tPad] = new list<Int_t>();
	67	list<Int_t>::iterator iter;
	68	for (iter=aReader.fClusterPerPadrow[tPad]->begin(); iter!=aReader.fClusterPerPadrow[tPad]->end(); iter++) {
	69	fClusterPerPadrow[tPad]->push_back(*iter);
	70	}
	71	}
	72	fSharedList = (list<Int_t> *) malloc(sizeof(list<Int_t> ) * AliESDfriendTrack::kMaxTPCcluster);
	73	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	74	fSharedList[tPad] = new list<Int_t>();
	75	list<Int_t>::iterator iter;
	76	for (iter=aReader.fSharedList[tPad]->begin(); iter!=aReader.fSharedList[tPad]->end(); iter++) {
	77	fSharedList[tPad]->push_back(*iter);
	78	}
	79	}
	80	}
67427ff7	81	//__________________
	82	//Destructor
	83	AliFemtoEventReaderESDChain::~AliFemtoEventReaderESDChain()
	84	{
	85	delete fEvent;
	86
	87	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	88	fClusterPerPadrow[tPad]->clear();
	89	delete fClusterPerPadrow[tPad];
	90	}
	91	delete [] fClusterPerPadrow;
	92	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	93	fSharedList[tPad]->clear();
	94	delete fSharedList[tPad];
	95	}
	96	delete [] fSharedList;
	97	}
	98
0215f606	99	//__________________
	100	// Assignment operator
	101	AliFemtoEventReaderESDChain& AliFemtoEventReaderESDChain::operator=(const AliFemtoEventReaderESDChain& aReader)
	102	{
	103	if (this == &aReader)
	104	return *this;
	105
	106	fConstrained = aReader.fConstrained;
	107	fNumberofEvent = aReader.fNumberofEvent;
	108	fCurEvent = aReader.fCurEvent;
	109	fCurFile = aReader.fCurFile;
	110	if (fEvent) delete fEvent;
	111	fEvent = new AliESD(*aReader.fEvent);
	112
	113	fEventFriend = aReader.fEventFriend;
	114
	115	if (fClusterPerPadrow) {
	116	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	117	fClusterPerPadrow[tPad]->clear();
	118	delete fClusterPerPadrow[tPad];
	119	}
	120	delete [] fClusterPerPadrow;
	121	}
	122
	123	if (fSharedList) {
	124	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	125	fSharedList[tPad]->clear();
	126	delete fSharedList[tPad];
	127	}
	128	delete [] fSharedList;
	129	}
	130
	131	fClusterPerPadrow = (list<Int_t> *) malloc(sizeof(list<Int_t> ) * AliESDfriendTrack::kMaxTPCcluster);
	132	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	133	fClusterPerPadrow[tPad] = new list<Int_t>();
	134	list<Int_t>::iterator iter;
	135	for (iter=aReader.fClusterPerPadrow[tPad]->begin(); iter!=aReader.fClusterPerPadrow[tPad]->end(); iter++) {
	136	fClusterPerPadrow[tPad]->push_back(*iter);
	137	}
	138	}
	139	fSharedList = (list<Int_t> *) malloc(sizeof(list<Int_t> ) * AliESDfriendTrack::kMaxTPCcluster);
	140	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	141	fSharedList[tPad] = new list<Int_t>();
	142	list<Int_t>::iterator iter;
	143	for (iter=aReader.fSharedList[tPad]->begin(); iter!=aReader.fSharedList[tPad]->end(); iter++) {
	144	fSharedList[tPad]->push_back(*iter);
	145	}
	146	}
	147
	148	return *this;
	149	}
67427ff7	150	//__________________
	151	// Simple report
	152	AliFemtoString AliFemtoEventReaderESDChain::Report()
	153	{
	154	AliFemtoString temp = "\n This is the AliFemtoEventReaderESDChain\n";
	155	return temp;
	156	}
	157
	158	//__________________
	159	// Select whether to read constrained or not constrained momentum
	160	void AliFemtoEventReaderESDChain::SetConstrained(const bool constrained)
	161	{
	162	fConstrained=constrained;
	163	}
	164	//__________________
	165	// Check whether we read constrained or not constrained momentum
	166	bool AliFemtoEventReaderESDChain::GetConstrained() const
	167	{
	168	return fConstrained;
	169	}
	170	//__________________
	171	AliFemtoEvent* AliFemtoEventReaderESDChain::ReturnHbtEvent()
	172	{
	173	// Get the event, read all the relevant information
	174	// and fill the AliFemtoEvent class
	175	// Returns a valid AliFemtoEvent
	176	AliFemtoEvent *hbtEvent = 0;
	177	string tFriendFileName;
	178
	179	// Get the friend information
	180	cout<<"starting to read event "<<fCurEvent<<endl;
	181	fEvent->SetESDfriend(fEventFriend);
	182	vector<int> label_table;//to check labels
	183
	184	hbtEvent = new AliFemtoEvent;
	185	//setting basic things
a531c6fc	186	// hbtEvent->SetEventNumber(fEvent->GetEventNumber());
67427ff7	187	hbtEvent->SetRunNumber(fEvent->GetRunNumber());
	188	//hbtEvent->SetNumberOfTracks(fEvent->GetNumberOfTracks());
	189	hbtEvent->SetMagneticField(fEvent->GetMagneticField()*kilogauss);//to check if here is ok
	190	hbtEvent->SetZDCN1Energy(fEvent->GetZDCN1Energy());
	191	hbtEvent->SetZDCP1Energy(fEvent->GetZDCP1Energy());
	192	hbtEvent->SetZDCN2Energy(fEvent->GetZDCN2Energy());
	193	hbtEvent->SetZDCP2Energy(fEvent->GetZDCP2Energy());
	194	hbtEvent->SetZDCEMEnergy(fEvent->GetZDCEMEnergy());
	195	hbtEvent->SetZDCParticipants(fEvent->GetZDCParticipants());
	196	hbtEvent->SetTriggerMask(fEvent->GetTriggerMask());
	197	hbtEvent->SetTriggerCluster(fEvent->GetTriggerCluster());
	198
	199	//Vertex
	200	double fV1[3];
	201	fEvent->GetVertex()->GetXYZ(fV1);
	202
	203	AliFmThreeVectorF vertex(fV1[0],fV1[1],fV1[2]);
	204	hbtEvent->SetPrimVertPos(vertex);
	205
	206	//starting to reading tracks
	207	int nofTracks=0; //number of reconstructed tracks in event
	208	nofTracks=fEvent->GetNumberOfTracks();
	209	int realnofTracks=0;//number of track which we use ina analysis
	210
	211	// Clear the shared cluster list
	212	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	213	fClusterPerPadrow[tPad]->clear();
	214	}
	215	for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
	216	fSharedList[tPad]->clear();
	217	}
	218
	219
	220	for (int i=0;i<nofTracks;i++) {
	221	const AliESDtrack *esdtrack=fEvent->GetTrack(i);//getting next track
	222
	223	list<Int_t>::iterator tClustIter;
	224
	225	Int_t tTrackIndices[AliESDfriendTrack::kMaxTPCcluster];
	226	Int_t tNClusters = esdtrack->GetTPCclusters(tTrackIndices);
	227	for (int tNcl=0; tNcl<AliESDfriendTrack::kMaxTPCcluster; tNcl++) {
	228	if (tTrackIndices[tNcl] >= 0) {
	229	tClustIter = find(fClusterPerPadrow[tNcl]->begin(), fClusterPerPadrow[tNcl]->end(), tTrackIndices[tNcl]);
	230	if (tClustIter == fClusterPerPadrow[tNcl]->end()) {
	231	fClusterPerPadrow[tNcl]->push_back(tTrackIndices[tNcl]);
	232	}
	233	else {
	234	fSharedList[tNcl]->push_back(tTrackIndices[tNcl]);
	235	}
	236	}
	237	}
	238
	239	}
	240
	241	for (int i=0;i<nofTracks;i++)
	242	{
	243	bool good_momentum=true; //flaga to chcek if we can read momentum of this track
	244
	245	AliFemtoTrack* trackCopy = new AliFemtoTrack();
	246	const AliESDtrack *esdtrack=fEvent->GetTrack(i);//getting next track
0215f606	247	// const AliESDfriendTrack *tESDfriendTrack = esdtrack->GetFriendTrack();
67427ff7	248
	249	trackCopy->SetCharge((short)esdtrack->GetSign());
	250
	251	//in aliroot we have AliPID
	252	//0-electron 1-muon 2-pion 3-kaon 4-proton 5-photon 6-pi0 7-neutron 8-kaon0 9-eleCon
	253	//we use only 5 first
	254	double esdpid[5];
	255	esdtrack->GetESDpid(esdpid);
	256	trackCopy->SetPidProbElectron(esdpid[0]);
	257	trackCopy->SetPidProbMuon(esdpid[1]);
	258	trackCopy->SetPidProbPion(esdpid[2]);
	259	trackCopy->SetPidProbKaon(esdpid[3]);
	260	trackCopy->SetPidProbProton(esdpid[4]);
	261
	262	double pxyz[3];
	263	if (fConstrained==true)
	264	good_momentum=esdtrack->GetConstrainedPxPyPz(pxyz); //reading constrained momentum
	265	else
	266	good_momentum=esdtrack->GetPxPyPz(pxyz);//reading noconstarined momentum
	267
	268	AliFemtoThreeVector v(pxyz[0],pxyz[1],pxyz[2]);
	269	if (v.mag() == 0) {
	270	// cout << "Found 0 momentum ???? " <<endl;
	271	delete trackCopy;
	272	continue;
	273	}
	274	trackCopy->SetP(v);//setting momentum
	275	trackCopy->SetPt(sqrt(pxyz[0]pxyz[0]+pxyz[1]pxyz[1]));
	276	const AliFmThreeVectorD p(pxyz[0],pxyz[1],pxyz[2]);
	277	const AliFmThreeVectorD origin(fV1[0],fV1[1],fV1[2]);
	278	//setting helix I do not if it is ok
	279	AliFmPhysicalHelixD helix(p,origin,(double)(fEvent->GetMagneticField())*kilogauss,(double)(trackCopy->Charge()));
	280	trackCopy->SetHelix(helix);
	281
	282	trackCopy->SetTrackId(esdtrack->GetID());
	283	trackCopy->SetFlags(esdtrack->GetStatus());
	284	//trackCopy->SetLabel(esdtrack->GetLabel());
	285
	286	//some stuff which could be useful
	287	float impact[2];
	288	float covimpact[3];
	289	esdtrack->GetImpactParameters(impact,covimpact);
	290	trackCopy->SetImpactD(impact[0]);
	291	trackCopy->SetImpactZ(impact[1]);
	292	trackCopy->SetCdd(covimpact[0]);
	293	trackCopy->SetCdz(covimpact[1]);
	294	trackCopy->SetCzz(covimpact[2]);
	295	trackCopy->SetITSchi2(esdtrack->GetITSchi2());
	296	trackCopy->SetITSncls(esdtrack->GetNcls(0));
	297	trackCopy->SetTPCchi2(esdtrack->GetTPCchi2());
	298	trackCopy->SetTPCncls(esdtrack->GetTPCNcls());
	299	trackCopy->SetTPCnclsF(esdtrack->GetTPCNclsF());
	300	trackCopy->SetTPCsignalN((short)esdtrack->GetTPCsignalN()); //due to bug in aliesdtrack class
	301	trackCopy->SetTPCsignalS(esdtrack->GetTPCsignalSigma());
	302
	303	// Fill cluster per padrow information
	304	Int_t tTrackIndices[AliESDfriendTrack::kMaxTPCcluster];
	305	Int_t tNClusters = esdtrack->GetTPCclusters(tTrackIndices);
	306	for (int tNcl=0; tNcl<AliESDfriendTrack::kMaxTPCcluster; tNcl++) {
	307	if (tTrackIndices[tNcl] > 0)
	308	trackCopy->SetTPCcluster(tNcl, 1);
	309	else
	310	trackCopy->SetTPCcluster(tNcl, 0);
	311	}
312
313	// Fill shared cluster information
314	list<Int_t>::iterator tClustIter;
315
316	for (int tNcl=0; tNcl<AliESDfriendTrack::kMaxTPCcluster; tNcl++) {
317	if (tTrackIndices[tNcl] > 0) {
318	tClustIter = find(fSharedList[tNcl]->begin(), fSharedList[tNcl]->end(), tTrackIndices[tNcl]);
319	if (tClustIter != fSharedList[tNcl]->end()) {
320	trackCopy->SetTPCshared(tNcl, 1);
321	cout << "Event next" << endl;
322	cout << "Track: " << i << endl;
323	cout << "Shared cluster: " << tNcl << " " << tTrackIndices[tNcl] << endl;
324	}
325	else {
326	trackCopy->SetTPCshared(tNcl, 0);
327	}
328	}
329	}
330
331	if (good_momentum==true)
332	{
333	hbtEvent->TrackCollection()->push_back(trackCopy);//adding track to analysis
334	realnofTracks++;//real number of tracks
335	// delete trackCopy;
336	}
337	else
338	{
339	delete trackCopy;
340	}
341
342	}
343
344	hbtEvent->SetNumberOfTracks(realnofTracks);//setting number of track which we read in event
345	fCurEvent++;
346	cout<<"end of reading nt "<<nofTracks<<" real number "<<realnofTracks<<endl;
347	return hbtEvent;
348	}
349	//___________________
350	// The chain loads the ESD for us
351	// You must provide the address where it can be found
352	void AliFemtoEventReaderESDChain::SetESDSource(AliESD *aESD)
353	{
354	fEvent = aESD;
355	}
356	//___________________
357	// We need the ESD tree to obtain
358	// information about the friend file location
359	void AliFemtoEventReaderESDChain::SetESDfriendSource(AliESDfriend *aFriend)
360	{
361	fEventFriend = aFriend;
362	}
363
364
365
366
367
368
369
370
371