int* pMaxSlice)
{
// see header file for class documentation
+
int iResult=0;
int iAddedDataBlocks=0;
if (pESD && blocks) {
- pESD->Reset();
-
+ pESD->Reset();
pESD->SetMagneticField(fSolenoidBz);
const AliHLTComponentBlockData* iter = NULL;
AliHLTTPCTrackletData* inPtr=NULL;
if (pMinSlice && (*pMinSlice==-1 || *pMinSlice>minslice)) *pMinSlice=minslice;
if (pMaxSlice && (*pMaxSlice==-1 || *pMaxSlice<maxslice)) *pMaxSlice=maxslice;
}
- //HLTDebug("dataspec %#x minslice %d", iter->fSpecification, minslice);
+ HLTDebug("AliHLTTPCEsdWriterComponent::ProcessBlocks: dataspec %#x minslice %d", iter->fSpecification, minslice);
if (minslice >=-1 && minslice<AliHLTTPCTransform::GetNSlice()) {
if (minslice!=maxslice) {
HLTWarning("data from multiple sectors in one block: "
}
}
}
- if (iAddedDataBlocks>0 && pTree) {
+ if ( iAddedDataBlocks>0 && pTree) {
pTree->Fill();
}
{
// see header file for class documentation
int iResult=0;
- if (pTracks && pESD) {
- HLTDebug("converting %d tracks from track array", pTracks->GetNTracks());
+ if (pTracks && pESD) {
for (int i=0; i<pTracks->GetNTracks() && iResult>=0; i++) {
AliHLTTPCTrack* pTrack=(*pTracks)[i];
if (pTrack) {
//pTrack->Print();
int iLocal=pTrack->Convert2AliKalmanTrack();
if (iLocal>=0) {
+
AliESDtrack iotrack;
iotrack.UpdateTrackParams(pTrack,AliESDtrack::kTPCin);
- iotrack.SetTPCPoints(pTrack->GetPoints());
+
+ Float_t points[4] = {pTrack->GetFirstPointX(), pTrack->GetFirstPointY(), pTrack->GetLastPointX(), pTrack->GetLastPointY() };
+
+ if(pTrack->GetSector() == -1){ // Set first and last points for global tracks
+ Double_t s = TMath::Sin( pTrack->GetAlpha() );
+ Double_t c = TMath::Cos( pTrack->GetAlpha() );
+ points[0] = pTrack->GetFirstPointX()*c + pTrack->GetFirstPointY()*s;
+ points[1] = -pTrack->GetFirstPointX()*s + pTrack->GetFirstPointY()*c;
+ points[2] = pTrack->GetLastPointX() *c + pTrack->GetLastPointY() *s;
+ points[3] = -pTrack->GetLastPointX() *s + pTrack->GetLastPointY() *c;
+ }
+ iotrack.SetTPCPoints(points);
+ //iotrack.SetTPCPoints(pTrack->GetPoints());
pESD->AddTrack(&iotrack);
} else {
HLTError("conversion to AliKalmanTrack failed for track %d of %d", i, pTracks->GetNTracks());
// or
// refer to README to build package
// or
- // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
-}
+ // visit http://web.ift.uib.no/~kjeks
+ }
AliHLTTPCEsdWriterComponent::AliConverter::~AliConverter()
{
} else if (argument.CompareTo("-tree")==0) {
fWriteTree=1;
- } else {
+ }else if (argument.CompareTo("-solenoidBz")==0) {
+ TString tmp="-solenoidBz ";
+ tmp+=argv[++i];
+ fBase->Configure(tmp.Data());
+ }
+ else {
HLTError("unknown argument %s", argument.Data());
break;
}
if (iResult>=0) {
iResult=fBase->Reconfigure(NULL, NULL);
- }
+ }
return iResult;
}
AliHLTUInt8_t* /*outputPtr*/,
AliHLTUInt32_t& size,
AliHLTComponentBlockDataList& /*outputBlocks*/ )
-{
+{
// see header file for class documentation
int iResult=0;
// no direct writing to the output buffer
}
AliESDEvent* pESD = fESD;
+
if (pESD && fBase) {
+
TTree* pTree = NULL;
// TODO: Matthias 06.12.2007
// Tried to write the ESD directly instead to a tree, but this did not work
// out. Information in the ESD is different, needs investigation.
+
if (fWriteTree)
pTree = new TTree("esdTree", "Tree with HLT ESD objects");
+
if (pTree) {
pTree->SetDirectory(0);
pESD->WriteToTree(pTree);
}
if ((iResult=fBase->ProcessBlocks(pTree, pESD, blocks, (int)evtData.fBlockCnt))>0) {
- // TODO: set the specification correctly
+ // TODO: set the specification correctly
if (pTree) {
// the esd structure is written to the user info and is
// needed in te ReadFromTree method to read all objects correctly
Double_t dr = .5*(slice->Param().RMax()-slice->Param().RMin());
Double_t cx = 0;
Double_t cy = r0;
- fYX->Range(cx-dr, cy-dr*1.05, cx+dr, cy+dr);
Double_t cz = .5*(slice->Param().ZMax()+slice->Param().ZMin());
Double_t dz = .5*(slice->Param().ZMax()-slice->Param().ZMin())*1.2;
- fZX->Range(cz-dz, cy-dr*1.05, cz+dz, cy+dr);//+dr);
+ fYX->Range(cx-dr, cy-dr*1.05, cx+dr, cy+dr);
+ fZX->Range(cz-dz, cy-dr*1.05, cz+dz, cy+dr);
+
+ //fYX->Range(cx-dr/6, cy-dr, cx+dr/8, cy-dr + dr/8);
+ //fZX->Range(cz+dz/2+dz/6, cy-dr , cz+dz-dz/8, cy-dr + dr/8);
//fYX->Range(cx-dr/3, cy-dr/3, cx+dr, cy+dr);
//fZX->Range(cz-dz/3, cy-dr/3, cz+dz, cy+dr);//+dr);
Double_t x12 = row1.X();
Double_t y11 = h11.Y() - h11.ErrY()*3;
Double_t y12 = h12.Y() + h12.ErrY()*3;
- Double_t z11 = h11.Z();
- Double_t z12 = h12.Z();
+ Double_t z11 = h11.Z() - h11.ErrZ()*3;
+ Double_t z12 = h12.Z() + h12.ErrZ()*3;
Double_t x21 = row2.X();
Double_t x22 = row2.X();
Double_t y21 = h21.Y() - h21.ErrY()*3;
Double_t y22 = h22.Y() + h22.ErrY()*3;
- Double_t z21 = h21.Z();
- Double_t z22 = h22.Z();
+ Double_t z21 = h21.Z() - h21.ErrZ()*3;
+ Double_t z22 = h22.Z() + h22.ErrZ()*3;
Double_t vx11, vx12, vy11, vy12, vx21, vx22, vy21, vy22;
-void AliHLTTPCCADisplay::DrawTrack1( AliHLTTPCCATrack &track, Int_t color, Bool_t DrawCells )
+void AliHLTTPCCADisplay::DrawTrack( AliHLTTPCCATrack &track, Int_t color, Bool_t DrawCells )
{
// draw track
if( track.NCells()<2 ) return;
- int width = 3;
+ int width = 1;
AliHLTTPCCADisplayTmpCell *vCells = new AliHLTTPCCADisplayTmpCell[track.NCells()];
AliHLTTPCCATrackParam &t = fSlice->ID2Point(track.PointID()[0]).Param();
delete[] vCells;
}
+
void AliHLTTPCCADisplay::DrawTrackletPoint( AliHLTTPCCATrackParam &t, Int_t color )
{
// draw tracklet point
void ConnectCells( Int_t iRow1, AliHLTTPCCACell &cell1, Int_t iRow2, AliHLTTPCCACell &cell2, Int_t color=-1 );
- void DrawTrack1( AliHLTTPCCATrack &track, Int_t color=-1, Bool_t DrawCells=1 );
+ void DrawTrack( AliHLTTPCCATrack &track, Int_t color=-1, Bool_t DrawCells=1 );
void DrawTrackletPoint( AliHLTTPCCATrackParam &t, Int_t color=-1 );
void SetSliceTransform( Double_t alpha );
{
public:
AliHLTTPCCAGBHit()
- :fX(0),fY(0),fZ(0),fErrX(0),fErrY(0),fErrZ(0),
+ :fX(0),fY(0),fZ(0),fErrX(0),fErrY(0),fErrZ(0),fAmp(0),
fISlice(0), fIRow(0), fID(0), fIsUsed(0),fSliceHit(){}
virtual ~AliHLTTPCCAGBHit(){}
Float_t &ErrX(){ return fErrX; }
Float_t &ErrY(){ return fErrY; }
Float_t &ErrZ(){ return fErrZ; }
+ Float_t &Amp() { return fAmp; }
Int_t &ISlice(){ return fISlice; }
Int_t &IRow(){ return fIRow; }
Float_t fErrY; //* Y position error
Float_t fErrZ; //* Z position error
+ Float_t fAmp; //* Maximal amplitude
Int_t fISlice; //* slice number
Int_t fIRow; //* row number
Int_t fID; //* external ID (id of AliTPCcluster)
{
public:
- AliHLTTPCCAGBTrack():fFirstHitRef(0),fNHits(0),fParam(),fAlpha(0){ ; }
+ AliHLTTPCCAGBTrack():fFirstHitRef(0),fNHits(0),fParam(),fAlpha(0),fDeDx(0){ ; }
virtual ~AliHLTTPCCAGBTrack(){ ; }
Int_t &NHits() { return fNHits; }
Int_t &FirstHitRef() { return fFirstHitRef; }
AliHLTTPCCATrackParam &Param() { return fParam; }
- Double_t &Alpha() { return fAlpha; }
-
+ Float_t &Alpha() { return fAlpha; }
+ Float_t &DeDx() { return fDeDx; }
static Bool_t ComparePNClusters( const AliHLTTPCCAGBTrack *a, const AliHLTTPCCAGBTrack *b){
return (a->fNHits > b->fNHits);
}
Int_t fFirstHitRef; // index of the first hit reference in track->hit reference array
Int_t fNHits; // number of track hits
AliHLTTPCCATrackParam fParam;// fitted track parameters
- Double_t fAlpha; //* Alpha angle of the parametrerisation
+ Float_t fAlpha; //* Alpha angle of the parametrerisation
+ Float_t fDeDx; //* DE/DX
private:
#include "TMath.h"
#include "TStopwatch.h"
#include "Riostream.h"
-#include "TROOT.h"
//#define DRAW
fTracks(0),
fNTracks(0),
fSliceTrackInfos(0),
+ fTime(0),
fStatNEvents(0)
{
//* constructor
- fStatTime[0] = 0;
- fStatTime[1] = 0;
- fStatTime[2] = 0;
- fStatTime[3] = 0;
- fStatTime[4] = 0;
- fStatTime[5] = 0;
- fStatTime[6] = 0;
- fStatTime[7] = 0;
- fStatTime[8] = 0;
- fStatTime[9] = 0;
+ for( Int_t i=0; i<20; i++ ) fStatTime[i] = 0;
}
AliHLTTPCCAGBTracker::AliHLTTPCCAGBTracker(const AliHLTTPCCAGBTracker&)
fTracks(0),
fNTracks(0),
fSliceTrackInfos(0),
+ fTime(0),
fStatNEvents(0)
{
//* dummy
fNHits = 0;
}
-void AliHLTTPCCAGBTracker::ReadHit( Double_t x, Double_t y, Double_t z,
- Double_t errY, Double_t errZ,
+void AliHLTTPCCAGBTracker::ReadHit( Float_t x, Float_t y, Float_t z,
+ Float_t errY, Float_t errZ, Float_t amp,
Int_t ID, Int_t iSlice, Int_t iRow )
{
//* read the hit to local array
hit.ErrX() = 1.e-4;//fSlices[iSlice].Param().ErrX();
hit.ErrY() = errY;
hit.ErrZ() = errZ;
+ hit.Amp() = amp;
hit.ID() = ID;
hit.ISlice()=iSlice;
hit.IRow() = iRow;
void AliHLTTPCCAGBTracker::FindTracks()
{
//* main tracking routine
-
+ fTime = 0;
fStatNEvents++;
#ifdef DRAW
AliHLTTPCCATracker &slice = fSlices[iSlice];
slice.Reconstruct();
timer.Stop();
+ fTime = timer.CpuTime();
fStatTime[0] += timer.CpuTime();
fStatTime[1]+=slice.Timers()[0];
fStatTime[2]+=slice.Timers()[1];
fStatTime[5]+=slice.Timers()[4];
fStatTime[6]+=slice.Timers()[5];
fStatTime[7]+=slice.Timers()[6];
-#ifdef DRAW
- //SlicePerformance( iSlice );
- //if( slice.NTracks()>0 ) AliHLTTPCCADisplay::Instance().Ask();
-#endif
+ fStatTime[8]+=slice.Timers()[7];
}
//cout<<"End CA reconstruction"<<endl;
TStopwatch timerMerge;
Merging();
timerMerge.Stop();
- fStatTime[8]+=timerMerge.CpuTime();
+ fStatTime[9]+=timerMerge.CpuTime();
+ fTime+=timerMerge.CpuTime();
//cout<<"End CA merging"<<endl;
#ifdef DRAW
{
//* track merging between slices
- Double_t dalpha = fSlices[1].Param().Alpha() - fSlices[0].Param().Alpha();
+ Float_t dalpha = fSlices[1].Param().Alpha() - fSlices[0].Param().Alpha();
+ Int_t nextSlice[100], prevSlice[100];
+ for( Int_t iSlice=0; iSlice<fNSlices; iSlice++ ){
+ nextSlice[iSlice] = iSlice + 1;
+ prevSlice[iSlice] = iSlice - 1;
+ }
+ nextSlice[ fNSlices/2 - 1 ] = 0;
+ prevSlice[ 0 ] = fNSlices/2 - 1;
+ nextSlice[ fNSlices - 1 ] = fNSlices/2;
+ prevSlice[ fNSlices/2 ] = fNSlices - 1;
+
+ TStopwatch timerMerge1;
+
Int_t maxNSliceTracks = 0;
for( Int_t iSlice=0; iSlice<fNSlices; iSlice++ ){
AliHLTTPCCATracker &iS = fSlices[iSlice];
}
for( Int_t iSlice=0; iSlice<fNSlices; iSlice++ ){
+ //cout<<"\nMerge slice "<<iSlice<<endl<<endl;
AliHLTTPCCATracker &iS = fSlices[iSlice];
- Int_t jSlice = iSlice+1;
- if( iSlice==fNSlices/2-1 ) jSlice = 0;
- else if( iSlice==fNSlices-1 ) jSlice = fNSlices/2;
+ Int_t jSlice = nextSlice[iSlice];
AliHLTTPCCATracker &jS = fSlices[jSlice];
int iNTracks = iS.NOutTracks();
int jNTracks = jS.NOutTracks();
for (Int_t itr=0; itr<iNTracks; itr++) {
iOK[0][itr] = 0;
iOK[1][itr] = 0;
- if( iS.OutTracks()[itr].NHits()<30 ) continue;
+ if( iS.OutTracks()[itr].NHits()<10 ) continue;
AliHLTTPCCATrackParam &iT1 = iTrParams[0][itr];
AliHLTTPCCATrackParam &iT2 = iTrParams[1][itr];
iT1 = iS.OutTracks()[itr].StartPoint();
for (Int_t jtr=0; jtr<jNTracks; jtr++) {
jOK[0][jtr] = 0;
jOK[1][jtr] = 0;
- if( jS.OutTracks()[jtr].NHits()<30 ) continue;
+ if( jS.OutTracks()[jtr].NHits()<10 ) continue;
AliHLTTPCCATrackParam &jT1 = jTrParams[0][jtr];
AliHLTTPCCATrackParam &jT2 = jTrParams[1][jtr];
jT1 = jS.OutTracks()[jtr].StartPoint();
}
//* start merging
-
+ //cout<<"Start slice merging.."<<endl;
for (Int_t itr=0; itr<iNTracks; itr++) {
if( !iOK[0][itr] && !iOK[1][itr] ) continue;
Int_t jBest = -1;
fSliceTrackInfos[iSlice][oldi].fNextNeighbour = -1;
} else continue;
}
-
+ //SG!!!
fSliceTrackInfos[iSlice][itr].fNextNeighbour = jBest;
fSliceTrackInfos[jSlice][jBest].fPrevNeighbour = itr;
}
if( iOK[i] ) delete[] iOK[i];
if( jOK[i] ) delete[] jOK[i];
}
+ timerMerge1.Stop();
+ fStatTime[10]+=timerMerge1.CpuTime();
+
+ TStopwatch timerMerge2;
Int_t nTracksTot = 0;
for( Int_t iSlice = 0; iSlice<fNSlices; iSlice++ ){
Int_t nTrackHits = 0;
+ //cout<<"\nStart global track creation...\n"<<endl;
+
+ static int nRejected = 0;
+
+ Int_t maxNRows = fSlices[0].Param().NRows();
+
for( Int_t iSlice = 0; iSlice<fNSlices; iSlice++ ){
AliHLTTPCCATracker &slice = fSlices[iSlice];
for( Int_t itr=0; itr<slice.NOutTracks(); itr++ ){
if( fSliceTrackInfos[iSlice][itr].fUsed ) continue;
- fSliceTrackInfos[iSlice][itr].fUsed = 1;
+ //cout<<"\n slice "<<iSlice<<", track "<<itr<<"\n"<<endl;
AliHLTTPCCAOutTrack &tCA = slice.OutTracks()[itr];
AliHLTTPCCAGBTrack &t = fTracks[fNTracks];
- t.Param() = tCA.StartPoint();
+ //t.Param() = tCA.StartPoint();
+ //t.Alpha() = slice.Param().Alpha();
t.NHits() = 0;
t.FirstHitRef() = nTrackHits;
- t.Alpha() = slice.Param().Alpha();
-
- Int_t ihit = 0;
+
+ struct FitPoint{
+ Int_t fISlice;
+ Int_t fHitID;
+ Float_t fX, fY, fZ, fErr2Y, fErr2Z, fAmp;
+ } fitPoints[300];
+ for( Int_t i=0; i<maxNRows; i++ ) fitPoints[i].fISlice = -1;
+
+ Int_t nHits = 0;
Int_t jSlice = iSlice;
Int_t jtr = itr;
- for( int jhit=0; jhit<tCA.NHits(); jhit++){
- int index = slice.OutTrackHits()[tCA.FirstHitRef()+jhit];
- fTrackHits[nTrackHits+ihit] = index;
- ihit++;
- }
- while( fSliceTrackInfos[jSlice][jtr].fNextNeighbour >=0 ){
- jtr = fSliceTrackInfos[jSlice][jtr].fNextNeighbour;
- if( jSlice==fNSlices/2-1 ) jSlice = 0;
- else if( jSlice==fNSlices-1 ) jSlice = fNSlices/2;
- else jSlice = jSlice + 1;
+ do{
if( fSliceTrackInfos[jSlice][jtr].fUsed ) break;
fSliceTrackInfos[jSlice][jtr].fUsed = 1;
- AliHLTTPCCAOutTrack &jTr = fSlices[jSlice].OutTracks()[jtr];
+ AliHLTTPCCATracker &jslice = fSlices[jSlice];
+ AliHLTTPCCAOutTrack &jTr = jslice.OutTracks()[jtr];
for( int jhit=0; jhit<jTr.NHits(); jhit++){
- int index = fSlices[jSlice].OutTrackHits()[jTr.FirstHitRef()+jhit];
- fTrackHits[nTrackHits+ihit] = index;
- ihit++;
+ int id = jslice.OutTrackHits()[jTr.FirstHitRef()+jhit];
+ AliHLTTPCCAGBHit &h = fHits[id];
+ FitPoint &p = fitPoints[h.IRow()];
+ if( p.fISlice >=0 ) continue;
+ p.fISlice = h.ISlice();
+ p.fHitID = id;
+ p.fX = jslice.Rows()[h.IRow()].X();
+ p.fY = h.Y();
+ p.fZ = h.Z();
+ p.fErr2Y = h.ErrY()*h.ErrY();
+ p.fErr2Z = h.ErrZ()*h.ErrZ();
+ p.fAmp = h.Amp();
+ nHits++;
}
+ jtr = fSliceTrackInfos[jSlice][jtr].fNextNeighbour;
+ jSlice = nextSlice[jSlice];
+ } while( jtr >=0 );
+
+ if( nHits < 30 ) continue; //SG!!!
+
+
+ Int_t firstRow = 0, lastRow = maxNRows-1;
+ for( firstRow=0; firstRow<maxNRows; firstRow++ ){
+ if( fitPoints[firstRow].fISlice>=0 ) break;
+ }
+ for( lastRow=maxNRows-1; lastRow>=0; lastRow-- ){
+ if( fitPoints[lastRow].fISlice>=0 ) break;
}
- t.NHits() = ihit;
- if( t.NHits()<50 ) continue;
- Int_t nHitsOld = t.NHits();
+ Int_t mmidRow = (firstRow + lastRow )/2;
+ Int_t midRow = firstRow;
+ for( Int_t i=firstRow+1; i<=lastRow; i++ ){
+ if( fitPoints[i].fISlice<0 ) continue;
+ if( TMath::Abs(i-mmidRow)>=TMath::Abs(midRow-mmidRow) ) continue;
+ midRow = i;
+ }
+ if( midRow==firstRow || midRow==lastRow ) continue;
+
+ Int_t searchRows[300];
+ Int_t nSearchRows = 0;
+ for( Int_t i=firstRow; i<=lastRow; i++ ) searchRows[nSearchRows++] = i;
+ for( Int_t i=lastRow+1; i<maxNRows; i++ ) searchRows[nSearchRows++] = i;
+ for( Int_t i=firstRow-1; i>=0; i-- ) searchRows[nSearchRows++] = i;
+
// refit
- {
- if( t.NHits()<10 ) continue;
+ AliHLTTPCCATrackParam t0;
- AliHLTTPCCAGBHit* *vhits = new AliHLTTPCCAGBHit*[nHitsOld];
-
- for( Int_t ih=0; ih<nHitsOld; ih++ ){
- vhits[ih] = &(fHits[fTrackHits[t.FirstHitRef() + ih]]);
- }
-
- sort(vhits, vhits+nHitsOld, AliHLTTPCCAGBHit::ComparePRowDown );
+ {
- AliHLTTPCCATrackParam t0;
-
{
- AliHLTTPCCAGBHit &h0 = *(vhits[0]);
- AliHLTTPCCAGBHit &h1 = *(vhits[nHitsOld/2]);
- AliHLTTPCCAGBHit &h2 = *(vhits[nHitsOld-1]);
- if( h0.IRow()==h1.IRow() || h0.IRow()==h2.IRow() || h1.IRow()==h2.IRow() ) continue;
- Double_t x1,y1,z1, x2, y2, z2, x3, y3, z3;
- fSlices[h0.ISlice()].Param().Slice2Global(h0.X(), h0.Y(), h0.Z(), &x1,&y1,&z1);
- fSlices[h1.ISlice()].Param().Slice2Global(h1.X(), h1.Y(), h1.Z(), &x2,&y2,&z2);
- fSlices[h2.ISlice()].Param().Slice2Global(h2.X(), h2.Y(), h2.Z(), &x3,&y3,&z3);
- fSlices[h0.ISlice()].Param().Global2Slice(x1, y1, z1, &x1,&y1,&z1);
- fSlices[h0.ISlice()].Param().Global2Slice(x2, y2, z2, &x2,&y2,&z2);
- fSlices[h0.ISlice()].Param().Global2Slice(x3, y3, z3, &x3,&y3,&z3);
- Float_t sp0[5] = {x1, y1, z1, .5, .5 };
- Float_t sp1[5] = {x2, y2, z2, .5, .5 };
- Float_t sp2[5] = {x3, y3, z3, .5, .5 };
- t0.ConstructXYZ3(sp0,sp1,sp2,1., 0);
- }
-
- Int_t currslice = vhits[0]->ISlice();
- Int_t currrow = fSlices[currslice].Param().NRows()-1;
- Double_t currd2 = 1.e10;
- AliHLTTPCCAGBHit *currhit = 0;
-
- for( Int_t ih=0; ih<nHitsOld; ih++ ){
- Double_t y0 = t0.GetY();
- Double_t z0 = t0.GetZ();
- AliHLTTPCCAGBHit &h = *(vhits[ih]);
- //cout<<"hit,slice,row="<<ih<<" "<<h.slice<<" "<<h.row<<endl;
- if( h.ISlice() == currslice && h.IRow() == currrow ){
- //* select the best hit in the row
- Double_t dy = h.Y() - y0;
- Double_t dz = h.Z() - z0;
- Double_t d2 = dy*dy + dz*dz;
- if( d2<currd2 ){
- currhit = &h;
- currd2 = d2;
- }
- continue;
- }else{
- if( currhit != 0 ){
- //* update track
- t0.Filter(currhit->Y(), currhit->Z(), currhit->ErrY(), currhit->ErrZ());
- currhit = 0;
- }
- if( h.ISlice() != currslice ){
- //* Rotate to the new slice
- currhit = 0;
- if( !t0.Rotate( -fSlices[currslice].Param().Alpha() +fSlices[h.ISlice()].Param().Alpha() ) ) break;
- currslice = h.ISlice();
- //currrow = 300;
- currd2 = 1.e10;
- }
- //* search for missed hits in rows
- {
- Double_t factor2 = 3.5*3.5;
- AliHLTTPCCATracker &cslice = fSlices[currslice];
- for( Int_t srow=currrow-1; srow>h.IRow(); srow--){
- AliHLTTPCCARow &row = cslice.Rows()[srow];
- if( !t0.TransportToX( row.X() ) ) continue;
- Int_t bestsh = -1;
- Double_t ds = 1.e10;
- for( Int_t ish=0; ish<row.NHits(); ish++ ){
- AliHLTTPCCAHit &sh = row.Hits()[ish];
- Double_t dy = sh.Y() - t0.GetY();
- Double_t dz = sh.Z() - t0.GetZ();
- Double_t dds = dy*dy+dz*dz;
- if( dds<ds ){
- ds = dds;
- bestsh = ish;
- }
- }
- if( bestsh<0 ) continue;
- AliHLTTPCCAHit &sh = row.Hits()[bestsh];
- Double_t dy = sh.Y() - t0.GetY();
- Double_t dz = sh.Z() - t0.GetZ();
- Double_t s2z = /*t0.GetErr2Z() + */ sh.ErrZ()*sh.ErrZ();
- if( dz*dz>factor2*s2z ) continue;
- Double_t s2y = /*t0.GetErr2Y() + */ sh.ErrY()*sh.ErrY();
- if( dy*dy>factor2*s2y ) continue;
- //* update track
- t0.Filter(sh.Y(), sh.Z(), sh.ErrY()/.33, sh.ErrZ()/.33);
- fTrackHits[nTrackHits+t.NHits()] = sh.ID();
- t.NHits()++;
- }
- }
- //* transport to the new row
- currrow = h.IRow();
- Bool_t ok = t0.TransportToX( h.X() );
- if( ok ){
- currrow = h.IRow();
- Double_t dy = h.Y() - t0.GetY();
- Double_t dz = h.Z() - t0.GetZ();
- currd2 = dy*dy + dz*dz;
- currhit = &h;
- }
+ FitPoint &p0 = fitPoints[firstRow];
+ FitPoint &p1 = fitPoints[midRow];
+ FitPoint &p2 = fitPoints[lastRow];
+ Float_t x0=p0.fX, y0=p0.fY, z0=p0.fZ;
+ Float_t x1=p1.fX, y1=p1.fY, z1=p1.fZ;
+ Float_t x2=p2.fX, y2=p2.fY, z2=p2.fZ;
+ if( p1.fISlice!=p0.fISlice ){
+ Float_t dAlpha = fSlices[p0.fISlice].Param().Alpha() - fSlices[p1.fISlice].Param().Alpha();
+ Float_t c = TMath::Cos(dAlpha);
+ Float_t s = TMath::Sin(dAlpha);
+ x1 = p1.fX*c + p1.fY*s;
+ y1 = p1.fY*c - p1.fX*s;
}
- }
- if( currhit != 0 ){ // update track
-
- t0.Filter(currhit->Y(), currhit->Z(), currhit->ErrY(), currhit->ErrZ());
+ if( p2.fISlice!=p0.fISlice ){
+ Float_t dAlpha = fSlices[p0.fISlice].Param().Alpha() - fSlices[p2.fISlice].Param().Alpha();
+ Float_t c = TMath::Cos(dAlpha);
+ Float_t s = TMath::Sin(dAlpha);
+ x2 = p2.fX*c + p2.fY*s;
+ y2 = p2.fY*c - p2.fX*s;
+ }
+
+ Float_t sp0[5] = {x0, y0, z0, .5, .5 };
+ Float_t sp1[5] = {x1, y1, z1, .5, .5 };
+ Float_t sp2[5] = {x2, y2, z2, .5, .5 };
+ t0.ConstructXYZ3(sp0,sp1,sp2,1., 0);
}
- //* search for missed hits in rows
- {
- Double_t factor2 = 3.5*3.5;
- for( Int_t srow=currrow-1; srow>=0; srow--){
+ Int_t currslice = fitPoints[firstRow].fISlice;
+
+ for( Int_t rowID=0; rowID<nSearchRows; rowID++ ){
+ Int_t iRow = searchRows[rowID];
+ FitPoint &p = fitPoints[iRow];
+
+ if( p.fISlice>=0 ){
+
+ //* Existing hit
+
+ //* Rotate to the new slice
+
+ if( p.fISlice!=currslice ){
+ if( !t0.Rotate( fSlices[p.fISlice].Param().Alpha() - fSlices[currslice].Param().Alpha() ) ) continue;
+ currslice = p.fISlice;
+ }
+ //* Transport to the new row
+
+ if( !t0.TransportToX( p.fX ) ) continue;
+
+ //* Calculate hit errors
+
+ GetErrors2( p.fISlice, iRow, t0, p.fErr2Y, p.fErr2Z );
+
+ } else {
+ //continue; //SG!!
+ //* Search for the missed hit
+
+ Float_t factor2 = 3.5*3.5;
+
AliHLTTPCCATracker *cslice = &(fSlices[currslice]);
- AliHLTTPCCARow *row = &(cslice->Rows()[srow]);
+ AliHLTTPCCARow *row = &(cslice->Rows()[iRow]);
if( !t0.TransportToX( row->X() ) ) continue;
+
if( t0.GetY() > row->MaxY() ){ //next slice
- Int_t j = currslice+1;
- if( currslice==fNSlices/2-1 ) j = 0;
- else if( currslice==fNSlices-1 ) j = fNSlices/2;
+
+ Int_t j = nextSlice[currslice];
+
//* Rotate to the new slice
- if( !t0.Rotate( -fSlices[currslice].Param().Alpha() +fSlices[j].Param().Alpha() ) ) break;
+
+ if( !t0.Rotate( -fSlices[currslice].Param().Alpha() +fSlices[j].Param().Alpha() ) ) continue;
currslice = j;
cslice = &(fSlices[currslice]);
- row = &(cslice->Rows()[srow]);
- if( !t0.TransportToX( row->X() ) ) continue;
+ row = &(cslice->Rows()[iRow]);
+ if( !t0.TransportToX( row->X() ) ) continue;
+ if( TMath::Abs(t0.GetY()) > row->MaxY() ) continue;
+
}else if( t0.GetY() < -row->MaxY() ){ //prev slice
- Int_t j = currslice-1;
- if( currslice==0 ) j = fNSlices/2-1;
- else if( currslice==fNSlices/2 ) j = fNSlices-1;
+ Int_t j = prevSlice[currslice];
//* Rotate to the new slice
if( !t0.Rotate( -fSlices[currslice].Param().Alpha() +fSlices[j].Param().Alpha() ) ) break;
currslice = j;
cslice = &(fSlices[currslice]);
- row = &(cslice->Rows()[srow]);
+ row = &(cslice->Rows()[iRow]);
if( !t0.TransportToX( row->X() ) ) continue;
+ if( TMath::Abs(t0.GetY()) > row->MaxY() ) continue;
}
+
Int_t bestsh = -1;
- Double_t ds = 1.e10;
+ Float_t ds = 1.e10;
for( Int_t ish=0; ish<row->NHits(); ish++ ){
AliHLTTPCCAHit &sh = row->Hits()[ish];
- Double_t dy = sh.Y() - t0.GetY();
- Double_t dz = sh.Z() - t0.GetZ();
- Double_t dds = dy*dy+dz*dz;
+ Float_t dy = sh.Y() - t0.GetY();
+ Float_t dz = sh.Z() - t0.GetZ();
+ Float_t dds = dy*dy+dz*dz;
if( dds<ds ){
ds = dds;
bestsh = ish;
}
}
if( bestsh<0 ) continue;
+
+ //* Calculate hit errors
+
+ GetErrors2( currslice, iRow, t0, p.fErr2Y, p.fErr2Z );
+
AliHLTTPCCAHit &sh = row->Hits()[bestsh];
- Double_t dy = sh.Y() - t0.GetY();
- Double_t dz = sh.Z() - t0.GetZ();
- Double_t s2z = /*t0.GetErr2Z() + */ sh.ErrZ()*sh.ErrZ();
+ Float_t dy = sh.Y() - t0.GetY();
+ Float_t dz = sh.Z() - t0.GetZ();
+ Float_t s2z = /*t0.GetErr2Z() + */ p.fErr2Z;
if( dz*dz>factor2*s2z ) continue;
- Double_t s2y = /*t0.GetErr2Y() + */ sh.ErrY()*sh.ErrY();
+ Float_t s2y = /*t0.GetErr2Y() + */ p.fErr2Y;
if( dy*dy>factor2*s2y ) continue;
- //* update track
- t0.Filter(sh.Y(), sh.Z(), sh.ErrY()/.33, sh.ErrZ()/.33);
- fTrackHits[nTrackHits+t.NHits()] = sh.ID();
- t.NHits()++;
- }
+
+ p.fISlice = currslice;
+ p.fHitID = sh.ID();
+ p.fX = row->X();
+ p.fY = sh.Y();
+ p.fZ = sh.Z();
+ p.fAmp = fHits[p.fHitID].Amp();
+ }
+
+ //* Update the track
+
+ t0.Filter2( p.fY, p.fZ, p.fErr2Y, p.fErr2Z );
}
- if( vhits ) delete[] vhits;
- //* search for missed hits in rows
+ //* final refit, dE/dx calculation
+ //cout<<"\n\nstart refit..\n"<<endl;
{
- Double_t factor2 = 3.5*3.5;
- AliHLTTPCCAGBHit &h0 = fHits[fTrackHits[t.FirstHitRef()]];
-
- Bool_t ok = t0.Rotate( -fSlices[currslice].Param().Alpha() +fSlices[h0.ISlice()].Param().Alpha() );
-
- currslice = h0.ISlice();
- currrow = h0.IRow();
-
- for( Int_t srow=currrow+1; srow<fSlices[0].Param().NRows()&&ok; srow++){
- AliHLTTPCCATracker *cslice = &(fSlices[currslice]);
- AliHLTTPCCARow *row = &(cslice->Rows()[srow]);
- if( !t0.TransportToX( row->X() ) ) continue;
- if( t0.GetY() > row->MaxY() ){ //next slice
- Int_t j = currslice+1;
- if( currslice==fNSlices/2-1 ) j = 0;
- else if( currslice==fNSlices-1 ) j = fNSlices/2;
- //* Rotate to the new slice
- if( !t0.Rotate( -fSlices[currslice].Param().Alpha() +fSlices[j].Param().Alpha() ) ) break;
- currslice = j;
- cslice = &(fSlices[currslice]);
- row = &(cslice->Rows()[srow]);
- if( !t0.TransportToX( row->X() ) ) continue;
- }else if( t0.GetY() < -row->MaxY() ){ //prev slice
- Int_t j = currslice-1;
- if( currslice==0 ) j = fNSlices/2-1;
- else if( currslice==fNSlices/2 ) j = fNSlices-1;
- //* Rotate to the new slice
- if( !t0.Rotate( -fSlices[currslice].Param().Alpha() +fSlices[j].Param().Alpha() ) ) break;
- currslice = j;
- cslice = &(fSlices[currslice]);
- row = &(cslice->Rows()[srow]);
- if( !t0.TransportToX( row->X() ) ) continue;
+ Double_t sumDeDx = 0;
+ Int_t nDeDx = 0;
+ t.NHits() = 0;
+
+ AliHLTTPCCATrackFitParam fitPar;
+ t0.CalculateFitParameters( fitPar, fSlices[0].Param().Bz() );
+
+ t0.Cov()[ 0] = .1;
+ t0.Cov()[ 1] = 0;
+ t0.Cov()[ 2] = .1;
+ t0.Cov()[ 3] = 0;
+ t0.Cov()[ 4] = 0;
+ t0.Cov()[ 5] = .1;
+ t0.Cov()[ 6] = 0;
+ t0.Cov()[ 7] = 0;
+ t0.Cov()[ 8] = 0;
+ t0.Cov()[ 9] = .1;
+ t0.Cov()[10] = 0;
+ t0.Cov()[11] = 0;
+ t0.Cov()[12] = 0;
+ t0.Cov()[13] = 0;
+ t0.Cov()[14] = .1;
+ t0.Chi2() = 0;
+ t0.NDF() = -5;
+ bool first = 1;
+ for( Int_t iRow = maxNRows-1; iRow>=0; iRow-- ){
+ FitPoint &p = fitPoints[iRow];
+ if( p.fISlice<0 ) continue;
+ fTrackHits[nTrackHits+t.NHits()] = p.fHitID;
+ t.NHits()++;
+
+ //* Rotate to the new slice
+
+ if( p.fISlice!=currslice ){
+ //cout<<"rotate..."<<endl;
+ //cout<<" before rotation:"<<endl;
+ //t0.Print();
+ if( !t0.Rotate( fSlices[p.fISlice].Param().Alpha() - fSlices[currslice].Param().Alpha() ) ) continue;
+ //cout<<" after rotation:"<<endl;
+ //t0.Print();
+ currslice = p.fISlice;
}
- Int_t bestsh = -1;
- Double_t ds = 1.e10;
- for( Int_t ish=0; ish<row->NHits(); ish++ ){
- AliHLTTPCCAHit &sh = row->Hits()[ish];
- Double_t dy = sh.Y() - t0.GetY();
- Double_t dz = sh.Z() - t0.GetZ();
- Double_t dds = dy*dy+dz*dz;
- if( dds<ds ){
- ds = dds;
- bestsh = ish;
- }
+ //* Transport to the new row
+
+ //cout<<" before transport:"<<endl;
+ //t0.Print();
+
+ //if( !t0.TransportToX( p.fX ) ) continue;
+ if( !t0.TransportToXWithMaterial( p.fX, fitPar ) ) continue;
+ //if( !t0.TransportToX( p.fX ) ) continue;
+ //cout<<" after transport:"<<endl;
+ //t0.Print();
+
+ //* Update the track
+
+ if( first ){
+ t0.Cov()[ 0] = .5*.5;
+ t0.Cov()[ 1] = 0;
+ t0.Cov()[ 2] = .5*.5;
+ t0.Cov()[ 3] = 0;
+ t0.Cov()[ 4] = 0;
+ t0.Cov()[ 5] = .2*.2;
+ t0.Cov()[ 6] = 0;
+ t0.Cov()[ 7] = 0;
+ t0.Cov()[ 8] = 0;
+ t0.Cov()[ 9] = .2*.2;
+ t0.Cov()[10] = 0;
+ t0.Cov()[11] = 0;
+ t0.Cov()[12] = 0;
+ t0.Cov()[13] = 0;
+ t0.Cov()[14] = .5*.5;
+ t0.Chi2() = 0;
+ t0.NDF() = -5;
}
- if( bestsh<0 ) continue;
- AliHLTTPCCAHit &sh = row->Hits()[bestsh];
- Double_t dy = sh.Y() - t0.GetY();
- Double_t dz = sh.Z() - t0.GetZ();
- Double_t s2z = /*t0.GetErr2Z() + */ sh.ErrZ()*sh.ErrZ();
- if( dz*dz>factor2*s2z ) continue;
- Double_t s2y = /*t0.GetErr2Y() + */ sh.ErrY()*sh.ErrY();
- if( dy*dy>factor2*s2y ) continue;
- //* update track
- t0.Filter(sh.Y(), sh.Z(), sh.ErrY()/.33, sh.ErrZ()/.33);
- fTrackHits[nTrackHits+t.NHits()] = sh.ID();
- t.NHits()++;
- }
+
+ //cout<<" before filtration:"<<endl;
+ //t0.Print();
+
+ if( !t0.Filter2( p.fY, p.fZ, p.fErr2Y, p.fErr2Z ) ) continue;
+ //cout<<" after filtration:"<<endl;
+ //t0.Print();
+ first = 0;
+
+ if( TMath::Abs( t0.CosPhi() )>1.e-4 ){
+ Float_t dLdX = TMath::Sqrt(1.+t0.DzDs()*t0.DzDs())/TMath::Abs(t0.CosPhi());
+ sumDeDx+=p.fAmp/dLdX;
+ nDeDx++;
+ }
+ }
+ t.DeDx() = 0;
+ if( nDeDx >0 ) t.DeDx() = sumDeDx/nDeDx;
+ if( t0.GetErr2Y()<=0 ){
+ cout<<"nhits = "<<t.NHits()<<", t0.GetErr2Y() = "<<t0.GetErr2Y()<<endl;
+ t0.Print();
+ //exit(1);
+ }
}
-
- Bool_t ok=1;
+
+ if( t.NHits()<30 ) continue;//SG!!
+
{
- for( int i=0; i<15; i++ ) ok = ok && finite(t0.Cov()[i]);
+ Bool_t ok=1;
+
+ Float_t *c = t0.Cov();
+ for( int i=0; i<15; i++ ) ok = ok && finite(c[i]);
for( int i=0; i<5; i++ ) ok = ok && finite(t0.Par()[i]);
ok = ok && (t0.GetX()>50);
- }
- if(!ok) continue;
- if( TMath::Abs(t0.Kappa())<1.e-8 ) t0.Kappa() = 1.e-8;
- if( nHitsOld != t.NHits() ){
- //cout<<"N matched hits = "<<(t.NHits() - nHitsOld )<<" / "<<nHitsOld<<endl;
+
+ if( c[0]<=0 || c[2]<=0 || c[5]<=0 || c[9]<=0 || c[14]<=0 ) ok = 0;
+ //if( c[0]>5. || c[2]>5. || c[5]>2. || c[9]>2 || c[14]>2 ) ok = 0;
+ if(!ok){
+ nRejected++;
+ //cout<<"\n\nRejected: "<<nRejected<<"\n"<<endl;
+ continue;
+ }
}
+ if( TMath::Abs(t0.Kappa())<1.e-8 ) t0.Kappa() = 1.e-8;
t.Param() = t0;
t.Alpha() = fSlices[currslice].Param().Alpha();
- if( t.NHits()<50 ) continue;
nTrackHits+= t.NHits();
- fNTracks++;
+ fNTracks++;
}
}
}
-
- //* selection
+ cout<<"\n\nRejected: "<<nRejected<<"\n"<<endl;
+ timerMerge2.Stop();
+ fStatTime[11]+=timerMerge2.CpuTime();
+ TStopwatch timerMerge3;
+
+ //* selection
+ //cout<<"Selection..."<<endl;
{
AliHLTTPCCAGBTrack *vtracks = new AliHLTTPCCAGBTrack [fNTracks];
Int_t *vhits = new Int_t [fNHits];
to.NHits()++;
h.IsUsed() = 1;
}
- if( to.NHits()<50 ) continue;
+ if( to.NHits()<30 ) continue;//SG!!!
nHits+=to.NHits();
nTracks++;
+ //cout<<to.Param().GetErr2Y()<<" "<<to.Param().GetErr2Z()<<endl;
}
fNTracks = nTracks;
if( fTrackHits ) delete[] fTrackHits;
fTracks = vtracks;
delete[] vptracks;
}
+ timerMerge3.Stop();
+ fStatTime[12]+=timerMerge3.CpuTime();
+}
+
+void AliHLTTPCCAGBTracker::GetErrors2( Int_t iSlice, Int_t iRow, AliHLTTPCCATrackParam &t, Float_t &Err2Y, Float_t &Err2Z )
+{
+ //
+ // Use calibrated cluster error from OCDB
+ //
+
+ Float_t z = TMath::Abs((250.-0.275)-TMath::Abs(t.GetZ()));
+ Int_t type = (iRow<63) ? 0: (iRow>126) ? 1:2;
+ Float_t cosPhiInv = TMath::Abs(t.GetCosPhi())>1.e-2 ?1./t.GetCosPhi() :0;
+ Float_t angleY = t.GetSinPhi()*cosPhiInv ;
+ Float_t angleZ = t.GetDzDs()*cosPhiInv ;
+
+ AliHLTTPCCATracker &slice = fSlices[iSlice];
+
+ Err2Y = slice.Param().GetClusterError2(0,type, z,angleY);
+ Err2Z = slice.Param().GetClusterError2(1,type, z,angleZ);
+}
+
+void AliHLTTPCCAGBTracker::GetErrors2( AliHLTTPCCAGBHit &h, AliHLTTPCCATrackParam &t, Float_t &Err2Y, Float_t &Err2Z )
+{
+ //
+ // Use calibrated cluster error from OCDB
+ //
+
+ GetErrors2( h.ISlice(), h.IRow(), t, Err2Y, Err2Z );
+}
+
+void AliHLTTPCCAGBTracker::WriteSettings( ostream &out )
+{
+ out<< NSlices()<<endl;
+ for( int iSlice=0; iSlice<NSlices(); iSlice++ ){
+ fSlices[iSlice].Param().WriteSettings( out );
+ }
+}
+
+void AliHLTTPCCAGBTracker::ReadSettings( istream &in )
+{
+ Int_t nSlices=0;
+ in >> nSlices;
+ SetNSlices( nSlices );
+ for( int iSlice=0; iSlice<NSlices(); iSlice++ ){
+ AliHLTTPCCAParam param;
+ param.ReadSettings ( in );
+ fSlices[iSlice].Initialize( param );
+ }
+}
+
+void AliHLTTPCCAGBTracker::WriteEvent( ostream &out )
+{
+ out<<NHits()<<endl;
+ for (Int_t ih=0; ih<NHits(); ih++) {
+ AliHLTTPCCAGBHit &h = fHits[ih];
+ out<<h.X()<<" ";
+ out<<h.Y()<<" ";
+ out<<h.Z()<<" ";
+ out<<h.ErrY()<<" ";
+ out<<h.ErrZ()<<" ";
+ out<<h.Amp()<<" ";
+ out<<h.ID()<<" ";
+ out<<h.ISlice()<<" ";
+ out<<h.IRow()<<endl;
+ }
+}
+
+void AliHLTTPCCAGBTracker::ReadEvent( istream &in )
+{
+ StartEvent();
+ Int_t nHits;
+ in >> nHits;
+ SetNHits(nHits);
+ for (Int_t i=0; i<nHits; i++) {
+ Float_t x, y, z, errY, errZ;
+ Float_t amp;
+ Int_t id, iSlice, iRow;
+ in>>x>>y>>z>>errY>>errZ>>amp>>id>>iSlice>>iRow;
+ ReadHit( x, y, z, errY, errZ, amp, id, iSlice, iRow );
+ }
+}
+
+void AliHLTTPCCAGBTracker::WriteTracks( ostream &out )
+{
+ out<<fTime<<endl;
+ Int_t nTrackHits = 0;
+ for( Int_t itr=0; itr<NTracks(); itr++ ){
+ AliHLTTPCCAGBTrack &t = Tracks()[itr];
+ nTrackHits+=t.NHits();
+ }
+ out<<nTrackHits<<endl;
+ for( Int_t ih=0; ih<nTrackHits; ih++ ){
+ out<< TrackHits()[ih]<<" ";
+ }
+ out<<endl;
+
+ out<<NTracks()<<endl;
+ for( Int_t itr=0; itr<NTracks(); itr++ ){
+ AliHLTTPCCAGBTrack &t = Tracks()[itr];
+ AliHLTTPCCATrackParam &p = t.Param();
+ out<< t.NHits()<<" ";
+ out<< t.FirstHitRef()<<" ";
+ out<< t.Alpha()<<" ";
+ out<< t.DeDx()<<endl;
+ out<< p.X()<<" ";
+ out<< p.CosPhi()<<" ";
+ out<< p.Chi2()<<" ";
+ out<< p.NDF()<<endl;
+ for( Int_t i=0; i<5; i++ ) out<<p.Par()[i]<<" ";
+ out<<endl;
+ for( Int_t i=0; i<15; i++ ) out<<p.Cov()[i]<<" ";
+ out<<endl;
+ }
+}
+
+void AliHLTTPCCAGBTracker::ReadTracks( istream &in )
+{
+ in>>fTime;
+ fStatTime[0]+=fTime;
+ fStatNEvents++;
+ if( fTrackHits ) delete[] fTrackHits;
+ fTrackHits = 0;
+ Int_t nTrackHits = 0;
+ in >> nTrackHits;
+ fTrackHits = new Int_t [nTrackHits];
+ for( Int_t ih=0; ih<nTrackHits; ih++ ){
+ in >> TrackHits()[ih];
+ }
+ if( fTracks ) delete[] fTracks;
+ fTracks = 0;
+ in >> fNTracks;
+ fTracks = new AliHLTTPCCAGBTrack[fNTracks];
+ for( Int_t itr=0; itr<NTracks(); itr++ ){
+ AliHLTTPCCAGBTrack &t = Tracks()[itr];
+ AliHLTTPCCATrackParam &p = t.Param();
+ in>> t.NHits();
+ in>> t.FirstHitRef();
+ in>> t.Alpha();
+ in>> t.DeDx();
+ in>> p.X();
+ in>> p.CosPhi();
+ in>> p.Chi2();
+ in>> p.NDF();
+ for( Int_t i=0; i<5; i++ ) in>>p.Par()[i];
+ for( Int_t i=0; i<15; i++ ) in>>p.Cov()[i];
+ }
}
class AliHLTTPCCAGBHit;
class TParticle;
class TProfile;
-
+class AliHLTTPCCATrackParam;
/**
* @class AliHLTTPCCAGBTracker
void SetNSlices( Int_t N );
void SetNHits( Int_t nHits );
- void ReadHit( Double_t x, Double_t y, Double_t z,
- Double_t ErrY, Double_t ErrZ, Int_t ID,
- Int_t iSlice, Int_t iRow );
+ void ReadHit( Float_t x, Float_t y, Float_t z,
+ Float_t ErrY, Float_t ErrZ, Float_t amp,
+ Int_t ID, Int_t iSlice, Int_t iRow );
void FindTracks();
void Merging();
AliHLTTPCCAGBHit *Hits(){ return fHits; }
Int_t NHits() const { return fNHits; }
Int_t NSlices() const { return fNSlices; }
+ Double_t Time() const { return fTime; }
Double_t StatTime( Int_t iTimer ) const { return fStatTime[iTimer]; }
Int_t StatNEvents() const { return fStatNEvents; }
Int_t NTracks() const { return fNTracks; }
AliHLTTPCCAGBTrack *Tracks(){ return fTracks; }
Int_t *TrackHits() {return fTrackHits; }
+ void GetErrors2( AliHLTTPCCAGBHit &h, AliHLTTPCCATrackParam &t, Float_t &Err2Y, Float_t &Err2Z );
+ void GetErrors2( Int_t iSlice, Int_t iRow, AliHLTTPCCATrackParam &t, Float_t &Err2Y, Float_t &Err2Z );
+
+ void WriteSettings( ostream &out );
+ void ReadSettings( istream &in );
+ void WriteEvent( ostream &out );
+ void ReadEvent( istream &in );
+ void WriteTracks( ostream &out );
+ void ReadTracks( istream &in );
protected:
};
AliHLTTPCCAGBSliceTrackInfo **fSliceTrackInfos; //* additional information for slice tracks
-
- Double_t fStatTime[10]; //* time measured
+ Double_t fTime;
+ Double_t fStatTime[20]; //* timers
Int_t fStatNEvents; //* n events proceed
ClassDef(AliHLTTPCCAGBTracker,1) //Base class for conformal mapping tracking
--- /dev/null
+// $Id: AliHLTTPCCAMCPoint.cxx 27042 2008-07-02 12:06:02Z richterm $
+//***************************************************************************
+// This file is property of and copyright by the ALICE HLT Project *
+// ALICE Experiment at CERN, All rights reserved. *
+// *
+// Primary Authors: Sergey Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de> *
+// Ivan Kisel <kisel@kip.uni-heidelberg.de> *
+// for The ALICE HLT Project. *
+// *
+// 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. *
+//***************************************************************************
+
+#include "AliHLTTPCCAMCPoint.h"
+
+
+AliHLTTPCCAMCPoint::AliHLTTPCCAMCPoint()
+ : fX(0), fY(0), fZ(0), fSx(0), fSy(0), fSz(0), fTime(0), fISlice(0), fTrackID(0)
+{
+ //* Default constructor
+}
--- /dev/null
+//-*- Mode: C++ -*-
+
+//* This file is property of and copyright by the ALICE HLT Project *
+//* ALICE Experiment at CERN, All rights reserved. *
+//* See cxx source for full Copyright notice *
+
+#ifndef ALIHLTTPCCAMCPOINT_H
+#define ALIHLTTPCCAMCPOINT_H
+
+#include "Rtypes.h"
+
+
+/**
+ * @class AliHLTTPCCAMCPoint
+ * store MC point information for AliHLTTPCCAPerformance
+ */
+class AliHLTTPCCAMCPoint
+{
+ public:
+
+ AliHLTTPCCAMCPoint();
+
+ Float_t &X() { return fX; }
+ Float_t &Y() { return fY; }
+ Float_t &Z() { return fZ; }
+ Float_t &Sx() { return fSx; }
+ Float_t &Sy() { return fSy; }
+ Float_t &Sz() { return fSz; }
+ Float_t &Time() { return fTime; }
+ Int_t &ISlice() { return fISlice; }
+ Int_t &TrackID() { return fTrackID; }
+
+ static Bool_t Compare( const AliHLTTPCCAMCPoint &p1, const AliHLTTPCCAMCPoint &p2 )
+ {
+ return (p1.fTrackID < p2.fTrackID);
+ }
+
+ protected:
+
+ Float_t fX; //* global X position
+ Float_t fY; //* global Y position
+ Float_t fZ; //* global Z position
+ Float_t fSx; //* slice X position
+ Float_t fSy; //* slice Y position
+ Float_t fSz; //* slice Z position
+ Float_t fTime; //* time
+ Int_t fISlice; //* slice number
+ Int_t fTrackID; //* mc track number
+};
+
+#endif
AliHLTTPCCAMCTrack::AliHLTTPCCAMCTrack()
- : fP(0), fPt(0), fNHits(0), fNReconstructed(0), fSet(0), fNTurns(0)
+ : fPDG(0), fP(0), fPt(0), fNHits(0), fNMCPoints(0), fFirstMCPointID(0), fNReconstructed(0), fSet(0), fNTurns(0)
{
//* Default constructor
}
AliHLTTPCCAMCTrack::AliHLTTPCCAMCTrack( const TParticle *part )
- : fP(0), fPt(0), fNHits(0), fNReconstructed(0), fSet(0), fNTurns(0)
+ : fPDG(0), fP(0), fPt(0), fNHits(0), fNMCPoints(0), fFirstMCPointID(0), fNReconstructed(0), fSet(0), fNTurns(0)
{
//* Constructor from TParticle
for( Int_t i=0; i<7; i++ ) fPar[i] = 0;
+ for( Int_t i=0; i<7; i++ ) fTPCPar[i] = 0;
fP = 0;
fPt = 0;
fPar[4] = part->Py()*pi;
fPar[5] = part->Pz()*pi;
fPar[6] = 0;
- Int_t pdg = part->GetPdgCode();
- if ( TMath::Abs(pdg) < 100000 ){
- TParticlePDG *pPDG = TDatabasePDG::Instance()->GetParticle(pdg);
+ fPDG = part->GetPdgCode();
+ if ( TMath::Abs(fPDG) < 100000 ){
+ TParticlePDG *pPDG = TDatabasePDG::Instance()->GetParticle(fPDG);
if( pPDG ) fPar[6] = pPDG->Charge()/3.0*pi;
}
}
+
+void AliHLTTPCCAMCTrack::SetTPCPar( Float_t X, Float_t Y, Float_t Z,
+ Float_t Px, Float_t Py, Float_t Pz )
+{
+ //* Set parameters at TPC entrance
+
+ for( Int_t i=0; i<7; i++ ) fTPCPar[i] = 0;
+
+ fTPCPar[0] = X;
+ fTPCPar[1] = Y;
+ fTPCPar[2] = Z;
+ Double_t p = TMath::Sqrt(Px*Px + Py*Py + Pz*Pz);
+ Double_t pi = ( p >1.e-4 ) ?1./p :0;
+ fTPCPar[3] = Px*pi;
+ fTPCPar[4] = Py*pi;
+ fTPCPar[5] = Pz*pi;
+ fTPCPar[6] = 0;
+ if ( TMath::Abs(fPDG) < 100000 ){
+ TParticlePDG *pPDG = TDatabasePDG::Instance()->GetParticle(fPDG);
+ if( pPDG ) fTPCPar[6] = pPDG->Charge()/3.0*pi;
+ }
+}
AliHLTTPCCAMCTrack();
AliHLTTPCCAMCTrack( const TParticle *part );
+ void SetTPCPar( Float_t X, Float_t Y, Float_t Z, Float_t Px, Float_t Py, Float_t Pz );
+ Int_t &PDG() { return fPDG;}
Double_t *Par() { return fPar; }
+ Double_t *TPCPar() { return fTPCPar; }
Double_t &P() { return fP; }
Double_t &Pt() { return fPt; }
+
Int_t &NHits() { return fNHits;}
+ Int_t &NMCPoints() { return fNMCPoints;}
+ Int_t &FirstMCPointID(){ return fFirstMCPointID;}
Int_t &NReconstructed(){ return fNReconstructed; }
Int_t &Set() { return fSet; }
Int_t &NTurns() { return fNTurns; }
protected:
+ Int_t fPDG; //* particle pdg code
Double_t fPar[7]; //* x,y,z,ex,ey,ez,q/p
+ Double_t fTPCPar[7]; //* x,y,z,ex,ey,ez,q/p at TPC entrance (x=y=0 means no information)
Double_t fP, fPt; //* momentum and transverse momentum
Int_t fNHits; //* N TPC clusters
+ Int_t fNMCPoints; //* N MC points
+ Int_t fFirstMCPointID; //* id of the first MC point in the points array
Int_t fNReconstructed; //* how many times is reconstructed
Int_t fSet; //* set of tracks 0-OutSet, 1-ExtraSet, 2-RefSet
Int_t fNTurns; //* N of turns in the current sector
-
+
};
#endif
fCosAlpha(0), fSinAlpha(0), fAngleMin(0), fAngleMax(0), fRMin(83.65), fRMax(133.3),
fZMin(0.0529937), fZMax(249.778), fErrX(0), fErrY(0), fErrZ(0.228808),fPadPitch(0.4),fBz(-5.),
fYErrorCorrection(0.33), fZErrorCorrection(0.45),
- fCellConnectionAngleXY(35./180.*TMath::Pi()),
- fCellConnectionAngleXZ(35./180.*TMath::Pi()),
+ fCellConnectionAngleXY(45./180.*TMath::Pi()),
+ fCellConnectionAngleXZ(45./180.*TMath::Pi()),
fMaxTrackMatchDRow(4), fTrackConnectionFactor(3.5), fTrackChiCut(3.5), fTrackChi2Cut(10)
{
+ fParamS0Par[0][0][0] = 0.00047013;
+ fParamS0Par[0][0][1] = 2.00135e-05;
+ fParamS0Par[0][0][2] = 0.0106533;
+ fParamS0Par[0][0][3] = 5.27104e-08;
+ fParamS0Par[0][0][4] = 0.012829;
+ fParamS0Par[0][0][5] = 0.000147125;
+ fParamS0Par[0][0][6] = 4.99432;
+ fParamS0Par[0][1][0] = 0.000883342;
+ fParamS0Par[0][1][1] = 1.07011e-05;
+ fParamS0Par[0][1][2] = 0.0103187;
+ fParamS0Par[0][1][3] = 4.25141e-08;
+ fParamS0Par[0][1][4] = 0.0224292;
+ fParamS0Par[0][1][5] = 8.27274e-05;
+ fParamS0Par[0][1][6] = 4.17233;
+ fParamS0Par[0][2][0] = 0.000745399;
+ fParamS0Par[0][2][1] = 5.62408e-06;
+ fParamS0Par[0][2][2] = 0.0151562;
+ fParamS0Par[0][2][3] = 5.08757e-08;
+ fParamS0Par[0][2][4] = 0.0601004;
+ fParamS0Par[0][2][5] = 7.97129e-05;
+ fParamS0Par[0][2][6] = 4.84913;
+ fParamS0Par[1][0][0] = 0.00215126;
+ fParamS0Par[1][0][1] = 6.82233e-05;
+ fParamS0Par[1][0][2] = 0.0221867;
+ fParamS0Par[1][0][3] = -6.27825e-09;
+ fParamS0Par[1][0][4] = -0.00745378;
+ fParamS0Par[1][0][5] = 0.000172629;
+ fParamS0Par[1][0][6] = 6.24987;
+ fParamS0Par[1][1][0] = 0.00181667;
+ fParamS0Par[1][1][1] = -4.17772e-06;
+ fParamS0Par[1][1][2] = 0.0253429;
+ fParamS0Par[1][1][3] = 1.3011e-07;
+ fParamS0Par[1][1][4] = -0.00362827;
+ fParamS0Par[1][1][5] = 0.00030406;
+ fParamS0Par[1][1][6] = 17.7775;
+ fParamS0Par[1][2][0] = 0.00158251;
+ fParamS0Par[1][2][1] = -3.55911e-06;
+ fParamS0Par[1][2][2] = 0.0247899;
+ fParamS0Par[1][2][3] = 7.20604e-08;
+ fParamS0Par[1][2][4] = 0.0179946;
+ fParamS0Par[1][2][5] = 0.000425504;
+ fParamS0Par[1][2][6] = 20.9294;
+
Update();
}
void AliHLTTPCCAParam::Initialize( Int_t iSlice,
- Int_t nRows, Double_t rowX[],
- Double_t alpha, Double_t dAlpha,
- Double_t rMin, Double_t rMax,
- Double_t zMin, Double_t zMax,
- Double_t padPitch, Double_t zSigma,
- Double_t bz
+ Int_t nRows, Float_t rowX[],
+ Float_t alpha, Float_t dAlpha,
+ Float_t rMin, Float_t rMax,
+ Float_t zMin, Float_t zMax,
+ Float_t padPitch, Float_t zSigma,
+ Float_t bz
)
{
// initialization
fTrackChi2Cut = fTrackChiCut * fTrackChiCut;
}
-void AliHLTTPCCAParam::Slice2Global( Double_t x, Double_t y, Double_t z,
- Double_t *X, Double_t *Y, Double_t *Z ) const
+void AliHLTTPCCAParam::Slice2Global( Float_t x, Float_t y, Float_t z,
+ Float_t *X, Float_t *Y, Float_t *Z ) const
{
// conversion of coorinates sector->global
*X = x*fCosAlpha - y*fSinAlpha;
*Z = z;
}
-void AliHLTTPCCAParam::Global2Slice( Double_t X, Double_t Y, Double_t Z,
- Double_t *x, Double_t *y, Double_t *z ) const
+void AliHLTTPCCAParam::Global2Slice( Float_t X, Float_t Y, Float_t Z,
+ Float_t *x, Float_t *y, Float_t *z ) const
{
// conversion of coorinates global->sector
*x = X*fCosAlpha + Y*fSinAlpha;
*y = Y*fCosAlpha - X*fSinAlpha;
*z = Z;
}
+
+Float_t AliHLTTPCCAParam::GetClusterError2( Int_t yz, Int_t type, Float_t z, Float_t angle )
+{
+ //* recalculate the cluster error wih respect to the track slope
+ Float_t angle2 = angle*angle;
+ Float_t *c = fParamS0Par[yz][type];
+ Float_t v = c[0] + z*(c[1] + c[3]*z) + angle2*(c[2] + angle2*c[4] + c[5]*z );
+ return TMath::Abs(v);
+}
+
+void AliHLTTPCCAParam::WriteSettings( ostream &out )
+{
+ out << fISlice<<endl;
+ out << fNRows<<endl;
+ out << fAlpha<<endl;
+ out << fDAlpha<<endl;
+ out << fCosAlpha<<endl;
+ out << fSinAlpha<<endl;
+ out << fAngleMin<<endl;
+ out << fAngleMax<<endl;
+ out << fRMin<<endl;
+ out << fRMax<<endl;
+ out << fZMin<<endl;
+ out << fZMax<<endl;
+ out << fErrX<<endl;
+ out << fErrY<<endl;
+ out << fErrZ<<endl;
+ out << fPadPitch<<endl;
+ out << fBz<<endl;
+ out << fYErrorCorrection<<endl;
+ out << fZErrorCorrection<<endl;
+ out << fCellConnectionAngleXY<<endl;
+ out << fCellConnectionAngleXZ<<endl;
+ out << fMaxTrackMatchDRow<<endl;
+ out << fTrackConnectionFactor<<endl;
+ out << fTrackChiCut<<endl;
+ out << fTrackChi2Cut<<endl;
+ for( Int_t iRow = 0; iRow<fNRows; iRow++ ){
+ out << fRowX[iRow]<<endl;
+ }
+ out<<endl;
+ for( Int_t i=0; i<2; i++ )
+ for( Int_t j=0; j<3; j++ )
+ for( Int_t k=0; k<7; k++ )
+ out << fParamS0Par[i][j][k]<<endl;
+ out<<endl;
+}
+
+void AliHLTTPCCAParam::ReadSettings( istream &in )
+{
+ in >> fISlice;
+ in >> fNRows;
+ in >> fAlpha;
+ in >> fDAlpha;
+ in >> fCosAlpha;
+ in >> fSinAlpha;
+ in >> fAngleMin;
+ in >> fAngleMax;
+ in >> fRMin;
+ in >> fRMax;
+ in >> fZMin;
+ in >> fZMax;
+ in >> fErrX;
+ in >> fErrY;
+ in >> fErrZ;
+ in >> fPadPitch;
+ in >> fBz;
+ in >> fYErrorCorrection;
+ in >> fZErrorCorrection;
+ in >> fCellConnectionAngleXY;
+ in >> fCellConnectionAngleXZ;
+ in >> fMaxTrackMatchDRow;
+ in >> fTrackConnectionFactor;
+ in >> fTrackChiCut;
+ in >> fTrackChi2Cut;
+ for( Int_t iRow = 0; iRow<fNRows; iRow++ ){
+ in >> fRowX[iRow];
+ }
+ for( Int_t i=0; i<2; i++ )
+ for( Int_t j=0; j<3; j++ )
+ for( Int_t k=0; k<7; k++ )
+ in >> fParamS0Par[i][j][k];
+}
#define ALIHLTTPCCAPARAM_H
#include "Rtypes.h"
+#include "Riostream.h"
/**
* @class ALIHLTTPCCAParam
AliHLTTPCCAParam();
virtual ~AliHLTTPCCAParam(){;}
- void Initialize( Int_t iSlice, Int_t nRows, Double_t rowX[],
- Double_t alpha, Double_t dAlpha,
- Double_t rMin, Double_t rMax, Double_t zMin, Double_t zMax,
- Double_t padPitch, Double_t zSigma, Double_t bz );
+ void Initialize( Int_t iSlice, Int_t nRows, Float_t rowX[],
+ Float_t alpha, Float_t dAlpha,
+ Float_t rMin, Float_t rMax, Float_t zMin, Float_t zMax,
+ Float_t padPitch, Float_t zSigma, Float_t bz );
void Update();
- void Slice2Global( Double_t x, Double_t y, Double_t z,
- Double_t *X, Double_t *Y, Double_t *Z ) const;
- void Global2Slice( Double_t x, Double_t y, Double_t z,
- Double_t *X, Double_t *Y, Double_t *Z ) const;
+ void Slice2Global( Float_t x, Float_t y, Float_t z,
+ Float_t *X, Float_t *Y, Float_t *Z ) const;
+ void Global2Slice( Float_t x, Float_t y, Float_t z,
+ Float_t *X, Float_t *Y, Float_t *Z ) const;
Int_t &ISlice(){ return fISlice;}
Int_t &NRows(){ return fNRows;}
- Double_t &RowX( Int_t iRow ){ return fRowX[iRow]; }
+ Float_t &RowX( Int_t iRow ){ return fRowX[iRow]; }
- Double_t &Alpha(){ return fAlpha;}
- Double_t &DAlpha(){ return fDAlpha;}
- Double_t &CosAlpha(){ return fCosAlpha;}
- Double_t &SinAlpha(){ return fSinAlpha;}
- Double_t &AngleMin(){ return fAngleMin;}
- Double_t &AngleMax(){ return fAngleMax;}
- Double_t &RMin(){ return fRMin;}
- Double_t &RMax(){ return fRMax;}
- Double_t &ZMin(){ return fZMin;}
- Double_t &ZMax(){ return fZMax;}
- Double_t &ErrZ(){ return fErrZ;}
- Double_t &ErrX(){ return fErrX;}
- Double_t &ErrY(){ return fErrY;}
- Double_t &Bz(){ return fBz;}
-
- Double_t &TrackConnectionFactor(){ return fTrackConnectionFactor; }
- Double_t &TrackChiCut() { return fTrackChiCut; }
- Double_t &TrackChi2Cut(){ return fTrackChi2Cut; }
- Int_t &MaxTrackMatchDRow(){ return fMaxTrackMatchDRow; }
- Double_t &YErrorCorrection(){ return fYErrorCorrection; }
- Double_t &ZErrorCorrection(){ return fZErrorCorrection; }
- Double_t &CellConnectionAngleXY(){ return fCellConnectionAngleXY; }
- Double_t &CellConnectionAngleXZ(){ return fCellConnectionAngleXZ; }
+ Float_t &Alpha(){ return fAlpha;}
+ Float_t &DAlpha(){ return fDAlpha;}
+ Float_t &CosAlpha(){ return fCosAlpha;}
+ Float_t &SinAlpha(){ return fSinAlpha;}
+ Float_t &AngleMin(){ return fAngleMin;}
+ Float_t &AngleMax(){ return fAngleMax;}
+ Float_t &RMin(){ return fRMin;}
+ Float_t &RMax(){ return fRMax;}
+ Float_t &ZMin(){ return fZMin;}
+ Float_t &ZMax(){ return fZMax;}
+ Float_t &ErrZ(){ return fErrZ;}
+ Float_t &ErrX(){ return fErrX;}
+ Float_t &ErrY(){ return fErrY;}
+ Float_t &Bz(){ return fBz;}
+
+ Float_t &TrackConnectionFactor(){ return fTrackConnectionFactor; }
+ Float_t &TrackChiCut() { return fTrackChiCut; }
+ Float_t &TrackChi2Cut(){ return fTrackChi2Cut; }
+ Int_t &MaxTrackMatchDRow(){ return fMaxTrackMatchDRow; }
+ Float_t &YErrorCorrection(){ return fYErrorCorrection; }
+ Float_t &ZErrorCorrection(){ return fZErrorCorrection; }
+ Float_t &CellConnectionAngleXY(){ return fCellConnectionAngleXY; }
+ Float_t &CellConnectionAngleXZ(){ return fCellConnectionAngleXZ; }
+
+ Float_t GetClusterError2(Int_t yz, Int_t type, Float_t z, Float_t angle );
+
+ void WriteSettings( ostream &out );
+ void ReadSettings( istream &in );
protected:
Int_t fISlice; // slice number
Int_t fNRows; // number of rows
- Double_t fAlpha, fDAlpha; // slice angle and angular size
- Double_t fCosAlpha, fSinAlpha;// sign and cosine of the slice angle
- Double_t fAngleMin, fAngleMax; // minimal and maximal angle
- Double_t fRMin, fRMax;// slice R range
- Double_t fZMin, fZMax;// slice Z range
- Double_t fErrX, fErrY, fErrZ;// default cluster errors
- Double_t fPadPitch; // pad pitch
- Double_t fBz; // magnetic field value (only constant field can be used)
-
- Double_t fYErrorCorrection;// correction factor for Y error of input clusters
- Double_t fZErrorCorrection;// correction factor for Z error of input clusters
-
- Double_t fCellConnectionAngleXY; // max phi angle between connected cells
- Double_t fCellConnectionAngleXZ; // max psi angle between connected cells
- Int_t fMaxTrackMatchDRow;// maximal jump in TPC row for connecting track segments
- Double_t fTrackConnectionFactor; // allowed distance in Chi^2/3.5 for neighbouring tracks
- Double_t fTrackChiCut; // cut for track Sqrt(Chi2/NDF);
- Double_t fTrackChi2Cut;// cut for track Chi^2/NDF
-
- Double_t fRowX[200];// X-coordinate of rows
+ Float_t fAlpha, fDAlpha; // slice angle and angular size
+ Float_t fCosAlpha, fSinAlpha;// sign and cosine of the slice angle
+ Float_t fAngleMin, fAngleMax; // minimal and maximal angle
+ Float_t fRMin, fRMax;// slice R range
+ Float_t fZMin, fZMax;// slice Z range
+ Float_t fErrX, fErrY, fErrZ;// default cluster errors
+ Float_t fPadPitch; // pad pitch
+ Float_t fBz; // magnetic field value (only constant field can be used)
+
+ Float_t fYErrorCorrection;// correction factor for Y error of input clusters
+ Float_t fZErrorCorrection;// correction factor for Z error of input clusters
+
+ Float_t fCellConnectionAngleXY; // max phi angle between connected cells
+ Float_t fCellConnectionAngleXZ; // max psi angle between connected cells
+ Int_t fMaxTrackMatchDRow;// maximal jump in TPC row for connecting track segments
+ Float_t fTrackConnectionFactor; // allowed distance in Chi^2/3.5 for neighbouring tracks
+ Float_t fTrackChiCut; // cut for track Sqrt(Chi2/NDF);
+ Float_t fTrackChi2Cut;// cut for track Chi^2/NDF
+
+ Float_t fRowX[200];// X-coordinate of rows
+ Float_t fParamS0Par[2][3][7]; // cluster error parameterization coeficients
ClassDef(AliHLTTPCCAParam,1);
};
#include "AliHLTTPCCAPerformance.h"
#include "AliHLTTPCCAGBHit.h"
#include "AliHLTTPCCAMCTrack.h"
+#include "AliHLTTPCCAMCPoint.h"
#include "AliHLTTPCCAOutTrack.h"
+#include "AliHLTTPCCAGBTrack.h"
#include "AliHLTTPCCAGBTracker.h"
#include "AliHLTTPCCATracker.h"
+
#include "TMath.h"
#include "TROOT.h"
#include "Riostream.h"
#include "TFile.h"
#include "TH1.h"
+#include "TH2.h"
#include "TProfile.h"
fNHits(0),
fMCTracks(0),
fNMCTracks(0),
+ fMCPoints(0),
+ fNMCPoints(0),
+ fDoClusterPulls(1),
fStatNEvents(0),
fStatNRecTot(0),
fStatNRecOut(0),
fStatNMCRef(0),
fStatNRecRef(0),
fStatNClonesRef(0),
+ fStatGBNRecTot(0),
+ fStatGBNRecOut(0),
+ fStatGBNGhost(0),
+ fStatGBNMCAll(0),
+ fStatGBNRecAll(0),
+ fStatGBNClonesAll(0),
+ fStatGBNMCRef(0),
+ fStatGBNRecRef(0),
+ fStatGBNClonesRef(0),
fHistoDir(0),
+ fhResY(0),
+ fhResZ(0),
+ fhResSinPhi(0),
+ fhResDzDs(0),
+ fhResPt(0),
+ fhPullY(0),
+ fhPullZ(0),
+ fhPullSinPhi(0),
+ fhPullDzDs(0),
+ fhPullQPt(0),
fhHitErrY(0),
fhHitErrZ(0),
- fhHitResX(0),
fhHitResY(0),
fhHitResZ(0),
- fhHitPullX(0),
fhHitPullY(0),
fhHitPullZ(0),
+ fhHitResY1(0),
+ fhHitResZ1(0),
+ fhHitPullY1(0),
+ fhHitPullZ1(0),
fhCellPurity(0),
fhCellNHits(0),
fhCellPurityVsN(0),
fhCellPurityVsPt(0),
- fhEffVsP(0)
+ fhEffVsP(0),
+ fhGBEffVsP(0)
{
//* constructor
}
fNHits(0),
fMCTracks(0),
fNMCTracks(0),
+ fMCPoints(0),
+ fNMCPoints(0),
+ fDoClusterPulls(1),
fStatNEvents(0),
fStatNRecTot(0),
fStatNRecOut(0),
fStatNMCRef(0),
fStatNRecRef(0),
fStatNClonesRef(0),
+ fStatGBNRecTot(0),
+ fStatGBNRecOut(0),
+ fStatGBNGhost(0),
+ fStatGBNMCAll(0),
+ fStatGBNRecAll(0),
+ fStatGBNClonesAll(0),
+ fStatGBNMCRef(0),
+ fStatGBNRecRef(0),
+ fStatGBNClonesRef(0),
fHistoDir(0),
+ fhResY(0),
+ fhResZ(0),
+ fhResSinPhi(0),
+ fhResDzDs(0),
+ fhResPt(0),
+ fhPullY(0),
+ fhPullZ(0),
+ fhPullSinPhi(0),
+ fhPullDzDs(0),
+ fhPullQPt(0),
fhHitErrY(0),
- fhHitErrZ(0),
- fhHitResX(0),
+ fhHitErrZ(0),
fhHitResY(0),
fhHitResZ(0),
- fhHitPullX(0),
fhHitPullY(0),
fhHitPullZ(0),
+ fhHitResY1(0),
+ fhHitResZ1(0),
+ fhHitPullY1(0),
+ fhHitPullZ1(0),
fhCellPurity(0),
fhCellNHits(0),
fhCellPurityVsN(0),
fhCellPurityVsPt(0),
- fhEffVsP(0)
+ fhEffVsP(0),
+ fhGBEffVsP(0)
{
//* dummy
}
if( fMCTracks ) delete[] fMCTracks;
fMCTracks = 0;
fNMCTracks = 0;
+ if( fMCPoints ) delete[] fMCPoints;
+ fMCPoints = 0;
+ fNMCPoints = 0;
}
void AliHLTTPCCAPerformance::SetNHits( Int_t NHits )
fNMCTracks = NMCTracks;
}
+void AliHLTTPCCAPerformance::SetNMCPoints( Int_t NMCPoints )
+{
+ //* set number of MC points
+ if( fMCPoints ) delete[] fMCPoints;
+ fMCPoints = 0;
+ fMCPoints = new AliHLTTPCCAMCPoint[ NMCPoints ];
+ fNMCPoints = 0;
+}
+
void AliHLTTPCCAPerformance::ReadHitLabel( Int_t HitID,
Int_t lab0, Int_t lab1, Int_t lab2 )
{
void AliHLTTPCCAPerformance::ReadMCTrack( Int_t index, const TParticle *part )
{
- //* read mc track to local array
+ //* read mc track to the local array
fMCTracks[index] = AliHLTTPCCAMCTrack(part);
}
+void AliHLTTPCCAPerformance::ReadMCTPCTrack( Int_t index, Float_t X, Float_t Y, Float_t Z,
+ Float_t Px, Float_t Py, Float_t Pz )
+{
+ //* read mc track parameters at TPC
+ fMCTracks[index].SetTPCPar(X,Y,Z,Px,Py,Pz);
+}
+
+void AliHLTTPCCAPerformance::ReadMCPoint( Int_t TrackID, Float_t X, Float_t Y, Float_t Z, Float_t Time, Int_t iSlice )
+{
+ //* read mc point to the local array
+ AliHLTTPCCAMCPoint &p = fMCPoints[fNMCPoints];
+ p.TrackID() = TrackID;
+ p.X() = X;
+ p.Y() = Y;
+ p.Z() = Z;
+ p.Time() = Time;
+ p.ISlice() = iSlice;
+ fTracker->Slices()[iSlice].Param().Global2Slice( X, Y, Z,
+ &p.Sx(), &p.Sy(), &p.Sz() );
+ if( X*X + Y*Y>10.) fNMCPoints++;
+}
+
void AliHLTTPCCAPerformance::CreateHistos()
{
//* create performance histogramms
TDirectory *curdir = gDirectory;
fHistoDir = gROOT->mkdir("HLTTPCCATrackerPerformance");
fHistoDir->cd();
- gDirectory->mkdir("Fit");
- gDirectory->cd("Fit");
- /*
- fhResY = new TH1D("resY", "Y resoltion [cm]", 100, -5., 5.);
- fhResZ = new TH1D("resZ", "Z resoltion [cm]", 100, -5., 5.);
- fhResTy = new TH1D("resTy", "Ty resoltion ", 100, -1., 1.);
- fhResTz = new TH1D("resTz", "Tz resoltion ", 100, -1., 1.);
- fhResP = new TH1D("resP", "P resoltion ", 100, -10., 10.);
- fhPullY = new TH1D("pullY", "Y pull", 100, -10., 10.);
- fhPullZ = new TH1D("pullZ", "Z pull", 100, -10., 10.);
- fhPullTy = new TH1D("pullTy", "Ty pull", 100, -10., 10.);
- fhPullTz = new TH1D("pullTz", "Tz pull", 100, -10., 10.);
- fhPullQp = new TH1D("pullQp", "Qp pull", 100, -10., 10.);
- */
+ gDirectory->mkdir("TrackFit");
+ gDirectory->cd("TrackFit");
+
+ fhResY = new TH1D("resY", "track Y resoltion [cm]", 30, -.5, .5);
+ fhResZ = new TH1D("resZ", "track Z resoltion [cm]", 30, -.5, .5);
+ fhResSinPhi = new TH1D("resSinPhi", "track SinPhi resoltion ", 30, -.03, .03);
+ fhResDzDs = new TH1D("resDzDs", "track DzDs resoltion ", 30, -.01, .01);
+ fhResPt = new TH1D("resPt", "track telative Pt resoltion", 30, -.2, .2);
+ fhPullY = new TH1D("pullY", "track Y pull", 30, -10., 10.);
+ fhPullZ = new TH1D("pullZ", "track Z pull", 30, -10., 10.);
+ fhPullSinPhi = new TH1D("pullSinPhi", "track SinPhi pull", 30, -10., 10.);
+ fhPullDzDs = new TH1D("pullDzDs", "track DzDs pull", 30, -10., 10.);
+ fhPullQPt = new TH1D("pullQPt", "track Q/Pt pull", 30, -10., 10.);
+
gDirectory->cd("..");
fhEffVsP = new TProfile("EffVsP", "Eff vs P", 100, 0., 5.);
+ fhGBEffVsP = new TProfile("GBEffVsP", "Global tracker: Eff vs P", 100, 0., 5.);
- fhHitResX = new TH1D("resHitX", "X cluster resoltion [cm]", 100, -2., 2.);
+ gDirectory->mkdir("Clusters");
+ gDirectory->cd("Clusters");
fhHitResY = new TH1D("resHitY", "Y cluster resoltion [cm]", 100, -2., 2.);
fhHitResZ = new TH1D("resHitZ", "Z cluster resoltion [cm]", 100, -2., 2.);
- fhHitPullX = new TH1D("pullHitX", "X cluster pull", 100, -10., 10.);
- fhHitPullY = new TH1D("pullHitY", "Y cluster pull", 100, -10., 10.);
- fhHitPullZ = new TH1D("pullHitZ", "Z cluster pull", 100, -10., 10.);
+ fhHitPullY = new TH1D("pullHitY", "Y cluster pull", 50, -10., 10.);
+ fhHitPullZ = new TH1D("pullHitZ", "Z cluster pull", 50, -10., 10.);
+ fhHitResY1 = new TH1D("resHitY1", "Y cluster resoltion [cm]", 100, -2., 2.);
+ fhHitResZ1 = new TH1D("resHitZ1", "Z cluster resoltion [cm]", 100, -2., 2.);
+ fhHitPullY1 = new TH1D("pullHitY1", "Y cluster pull", 50, -10., 10.);
+ fhHitPullZ1 = new TH1D("pullHitZ1", "Z cluster pull", 50, -10., 10.);
fhHitErrY = new TH1D("HitErrY", "Y cluster error [cm]", 100, 0., 1.);
fhHitErrZ = new TH1D("HitErrZ", "Z cluster error [cm]", 100, 0., 1.);
+ gDirectory->cd("..");
+
+ gDirectory->mkdir("Cells");
+ gDirectory->cd("Cells");
fhCellPurity = new TH1D("CellPurity", "Cell Purity", 100, -0.1, 1.1);
fhCellNHits = new TH1D("CellNHits", "Cell NHits", 40, 0., 40.);
fhCellPurityVsN = new TProfile("CellPurityVsN", "Cell purity Vs N hits", 40, 2., 42.);
fhCellPurityVsPt = new TProfile("CellPurityVsPt", "Cell purity Vs Pt", 100, 0., 5.);
+ gDirectory->cd("..");
curdir->cd();
}
{
//* calculate slice tracker performance
if( !fTracker ) return;
-
- int nRecTot = 0, nGhost=0, nRecOut=0;
- int nMCAll = 0, nRecAll=0, nClonesAll=0;
- int nMCRef = 0, nRecRef=0, nClonesRef=0;
+
+ Int_t nRecTot = 0, nGhost=0, nRecOut=0;
+ Int_t nMCAll = 0, nRecAll=0, nClonesAll=0;
+ Int_t nMCRef = 0, nRecRef=0, nClonesRef=0;
AliHLTTPCCATracker &slice = fTracker->Slices()[iSlice];
- int firstSliceHit = 0;
+ Int_t firstSliceHit = 0;
for( ; firstSliceHit<fTracker->NHits(); firstSliceHit++){
if( fTracker->Hits()[firstSliceHit].ISlice()==iSlice ) break;
}
- int endSliceHit = firstSliceHit;
+ Int_t endSliceHit = firstSliceHit;
for( ; endSliceHit<fTracker->NHits(); endSliceHit++){
if( fTracker->Hits()[endSliceHit].ISlice()!=iSlice ) break;
for (Int_t ic = 0; ic<row.NCells(); ic++){
AliHLTTPCCACell &c = row.Cells()[ic];
Int_t *lb = new Int_t[c.NHits()*3];
- Int_t nla = 0;
- //cout<<11<<" "<<c.NHits()<<endl;
+ Int_t nla = 0;
for( Int_t j=0; j<c.NHits(); j++){
AliHLTTPCCAHit &h = row.GetCellHit(c,j);
//cout<<"hit ID="<<h.ID()<<" of"<<fTracker->NHits()<<endl;
if( l.fLab[1]>=0 ) lb[nla++]= l.fLab[1];
if( l.fLab[2]>=0 ) lb[nla++]= l.fLab[2];
}
- //cout<<12<<endl;
sort( lb, lb+nla );
- int labmax = -1, labcur=-1, lmax = 0, lcurr=0;
- for( int i=0; i<nla; i++ ){
+ Int_t labmax = -1, labcur=-1, lmax = 0, lcurr=0;
+ for( Int_t i=0; i<nla; i++ ){
if( lb[i]!=labcur ){
if( labcur>=0 && lmax<lcurr ){
lmax = lcurr;
labmax = labcur;
}
- int label = labmax;
+ Int_t label = labmax;
lmax = 0;
for( Int_t j=0; j<c.NHits(); j++){
AliHLTTPCCAHit &h = row.GetCellHit(c,j);
{
for (Int_t imc=0; imc<fNMCTracks; imc++) fMCTracks[imc].NHits() = 0;
- for( int ih=firstSliceHit; ih<endSliceHit; ih++){
+ for( Int_t ih=firstSliceHit; ih<endSliceHit; ih++){
AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[ih].ID()];
if( l.fLab[0]>=0 ) fMCTracks[l.fLab[0]].NHits()++;
if( l.fLab[1]>=0 ) fMCTracks[l.fLab[1]].NHits()++;
}
}
- int traN = slice.NOutTracks();
- Int_t *traLabels = new Int_t[traN];
- Double_t *traPurity = new Double_t[traN];
+ Int_t traN = slice.NOutTracks();
+ Int_t *traLabels = 0;
+ Double_t *traPurity = 0;
+ traLabels = new Int_t[traN];
+ traPurity = new Double_t[traN];
{
for (Int_t itr=0; itr<traN; itr++) {
traLabels[itr]=-1;
traPurity[itr]= 0;
AliHLTTPCCAOutTrack &tCA = slice.OutTracks()[itr];
- int nhits = tCA.NHits();
- int *lb = new Int_t[nhits*3];
- int nla=0;
- for( int ihit=0; ihit<nhits; ihit++){
- int index = slice.OutTrackHits()[tCA.FirstHitRef()+ihit];
+ Int_t nhits = tCA.NHits();
+ Int_t *lb = new Int_t[nhits*3];
+ Int_t nla=0;
+ //cout<<"\nHit labels:"<<endl;
+ for( Int_t ihit=0; ihit<nhits; ihit++){
+ Int_t index = slice.OutTrackHits()[tCA.FirstHitRef()+ihit];
AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[index].ID()];
+ //cout<<l.fLab[0]<<" "<<l.fLab[1]<<" "<<l.fLab[2]<<endl;
if(l.fLab[0]>=0 ) lb[nla++]= l.fLab[0];
if(l.fLab[1]>=0 ) lb[nla++]= l.fLab[1];
if(l.fLab[2]>=0 ) lb[nla++]= l.fLab[2];
}
sort( lb, lb+nla );
- int labmax = -1, labcur=-1, lmax = 0, lcurr=0;
- for( int i=0; i<nla; i++ ){
+ Int_t labmax = -1, labcur=-1, lmax = 0, lcurr=0;
+ for( Int_t i=0; i<nla; i++ ){
if( lb[i]!=labcur ){
if( labcur>=0 && lmax<lcurr ){
lmax = lcurr;
labmax = labcur;
}
lmax = 0;
- for( int ihit=0; ihit<nhits; ihit++){
- int index = slice.OutTrackHits()[tCA.FirstHitRef()+ihit];
+ for( Int_t ihit=0; ihit<nhits; ihit++){
+ Int_t index = slice.OutTrackHits()[tCA.FirstHitRef()+ihit];
AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[index].ID()];
if( l.fLab[0] == labmax || l.fLab[1] == labmax || l.fLab[2] == labmax
) lmax++;
}
traLabels[itr] = labmax;
traPurity[itr] = ( (nhits>0) ?double(lmax)/double(nhits) :0 );
+ //cout<<"perf track "<<itr<<": "<<nhits<<" "<<labmax<<" "<<traPurity[itr]<<endl;
if( lb ) delete[] lb;
}
}
nRecTot+= traN;
- for(int itr=0; itr<traN; itr++){
- if( traPurity[itr]<.7 || traLabels[itr]<0 || traLabels[itr]>=fNMCTracks){
+
+ for(Int_t itr=0; itr<traN; itr++){
+ if( traPurity[itr]<.9 || traLabels[itr]<0 || traLabels[itr]>=fNMCTracks){
nGhost++;
continue;
}
if( mc.Set()>0 ) fhEffVsP->Fill(mc.P(), ( mc.NReconstructed()>0 ?1 :0));
}
-
if( traLabels ) delete[] traLabels;
if( traPurity ) delete[] traPurity;
void AliHLTTPCCAPerformance::Performance()
{
// main routine for performance calculation
+ //SG!!!
+ /*
+ fStatNEvents=0;
+ fStatNRecTot=0;
+ fStatNRecOut=0;
+ fStatNGhost=0;
+ fStatNMCAll=0;
+ fStatNRecAll=0;
+ fStatNClonesAll=0;
+ fStatNMCRef=0;
+ fStatNRecRef=0;
+ fStatNClonesRef=0;
+ */
fStatNEvents++;
- for( int islice=0; islice<fTracker->NSlices(); islice++){
+ for( Int_t islice=0; islice<fTracker->NSlices(); islice++){
SlicePerformance(islice,0);
}
+
+ // global tracker performance
+ {
+ if( !fTracker ) return;
+
+ Int_t nRecTot = 0, nGhost=0, nRecOut=0;
+ Int_t nMCAll = 0, nRecAll=0, nClonesAll=0;
+ Int_t nMCRef = 0, nRecRef=0, nClonesRef=0;
+
+ // Select reconstructable MC tracks
+
+ {
+ for (Int_t imc=0; imc<fNMCTracks; imc++) fMCTracks[imc].NHits() = 0;
+
+ for( Int_t ih=0; ih<fNHits; ih++){
+ AliHLTTPCCAHitLabel &l = fHitLabels[ih];
+ if( l.fLab[0]>=0 ) fMCTracks[l.fLab[0]].NHits()++;
+ if( l.fLab[1]>=0 ) fMCTracks[l.fLab[1]].NHits()++;
+ if( l.fLab[2]>=0 ) fMCTracks[l.fLab[2]].NHits()++;
+ }
+
+ for (Int_t imc=0; imc<fNMCTracks; imc++) {
+ AliHLTTPCCAMCTrack &mc = fMCTracks[imc];
+ mc.Set() = 0;
+ mc.NReconstructed() = 0;
+ mc.NTurns() = 1;
+ if( mc.NHits() >= 50 && mc.P()>=.05 ){
+ mc.Set() = 1;
+ nMCAll++;
+ if( mc.P()>=1. ){
+ mc.Set() = 2;
+ nMCRef++;
+ }
+ }
+ }
+ }
+
+ Int_t traN = fTracker->NTracks();
+ Int_t *traLabels = 0;
+ Double_t *traPurity = 0;
+ traLabels = new Int_t[traN];
+ traPurity = new Double_t[traN];
+ {
+ for (Int_t itr=0; itr<traN; itr++) {
+ traLabels[itr]=-1;
+ traPurity[itr]= 0;
+ AliHLTTPCCAGBTrack &tCA = fTracker->Tracks()[itr];
+ Int_t nhits = tCA.NHits();
+ Int_t *lb = new Int_t[nhits*3];
+ Int_t nla=0;
+ for( Int_t ihit=0; ihit<nhits; ihit++){
+ Int_t index = fTracker->TrackHits()[tCA.FirstHitRef()+ihit];
+ AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[index].ID()];
+ if(l.fLab[0]>=0 ) lb[nla++]= l.fLab[0];
+ if(l.fLab[1]>=0 ) lb[nla++]= l.fLab[1];
+ if(l.fLab[2]>=0 ) lb[nla++]= l.fLab[2];
+ }
+ sort( lb, lb+nla );
+ Int_t labmax = -1, labcur=-1, lmax = 0, lcurr=0;
+ for( Int_t i=0; i<nla; i++ ){
+ if( lb[i]!=labcur ){
+ if( labcur>=0 && lmax<lcurr ){
+ lmax = lcurr;
+ labmax = labcur;
+ }
+ labcur = lb[i];
+ lcurr = 0;
+ }
+ lcurr++;
+ }
+ if( labcur>=0 && lmax<lcurr ){
+ lmax = lcurr;
+ labmax = labcur;
+ }
+ lmax = 0;
+ for( Int_t ihit=0; ihit<nhits; ihit++){
+ Int_t index = fTracker->TrackHits()[tCA.FirstHitRef()+ihit];
+ AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[index].ID()];
+ if( l.fLab[0] == labmax || l.fLab[1] == labmax || l.fLab[2] == labmax
+ ) lmax++;
+ }
+ traLabels[itr] = labmax;
+ traPurity[itr] = ( (nhits>0) ?double(lmax)/double(nhits) :0 );
+ if( lb ) delete[] lb;
+ }
+ }
+
+ nRecTot+= traN;
+ for(Int_t itr=0; itr<traN; itr++){
+ if( traPurity[itr]<.9 || traLabels[itr]<0 || traLabels[itr]>=fNMCTracks){
+ nGhost++;
+ continue;
+ }
+
+ AliHLTTPCCAMCTrack &mc = fMCTracks[traLabels[itr]];
+
+ mc.NReconstructed()++;
+ if( mc.Set()== 0 ) nRecOut++;
+ else{
+ if( mc.NReconstructed()==1 ) nRecAll++;
+ else if(mc.NReconstructed() > mc.NTurns() ) nClonesAll++;
+ if( mc.Set()==2 ){
+ if( mc.NReconstructed()==1 ) nRecRef++;
+ else if(mc.NReconstructed() > mc.NTurns() ) nClonesRef++;
+ }
+ }
+
+ // track resolutions
+ while( TMath::Abs(mc.TPCPar()[0]) + TMath::Abs(mc.TPCPar()[1])>1 ){
+ if( traPurity[itr]<.90 ) break;
+ AliHLTTPCCAGBTrack &t = fTracker->Tracks()[itr];
+ AliHLTTPCCATrackParam p = t.Param();
+ Double_t cosA = TMath::Cos( t.Alpha() );
+ Double_t sinA = TMath::Sin( t.Alpha() );
+ Double_t mcX = mc.TPCPar()[0]*cosA + mc.TPCPar()[1]*sinA;
+ Double_t mcY = -mc.TPCPar()[0]*sinA + mc.TPCPar()[1]*cosA;
+ Double_t mcZ = mc.TPCPar()[2];
+ Double_t mcEx = mc.TPCPar()[3]*cosA + mc.TPCPar()[4]*sinA;
+ Double_t mcEy = -mc.TPCPar()[3]*sinA + mc.TPCPar()[4]*cosA;
+ Double_t mcEz = mc.TPCPar()[5];
+ Double_t mcEt = TMath::Sqrt(mcEx*mcEx + mcEy*mcEy);
+ if( TMath::Abs(mcEt)<1.e-4 ) break;
+ Double_t mcSinPhi = mcEy / mcEt;
+ Double_t mcDzDs = mcEz / mcEt;
+ Double_t mcQPt = mc.TPCPar()[6]/ mcEt;
+ if( TMath::Abs(mcQPt)<1.e-4 ) break;
+ Double_t mcPt = 1./TMath::Abs(mcQPt);
+ if( mcPt<1. ) break;
+ if( t.NHits() < 50 ) break;
+ Double_t bz = fTracker->Slices()[0].Param().Bz();
+ if( !p.TransportToXWithMaterial( mcX, bz ) ) break;
+ if( p.GetCosPhi()*mcEx < 0 ){ // change direction
+ mcSinPhi = -mcSinPhi;
+ mcDzDs = -mcDzDs;
+ mcQPt = -mcQPt;
+ }
+ const Double_t kCLight = 0.000299792458;
+ Double_t k2QPt = 100;
+ if( TMath::Abs(bz)>1.e-4 ) k2QPt= 1./(bz*kCLight);
+ Double_t qPt = p.GetKappa()*k2QPt;
+ Double_t pt = 100;
+ if( TMath::Abs(qPt) >1.e-4 ) pt = 1./TMath::Abs(qPt);
+
+ fhResY->Fill( p.GetY() - mcY );
+ fhResZ->Fill( p.GetZ() - mcZ );
+ fhResSinPhi->Fill( p.GetSinPhi() - mcSinPhi );
+ fhResDzDs->Fill( p.GetDzDs() - mcDzDs );
+ fhResPt->Fill( ( pt - mcPt )/mcPt );
+
+ if( p.GetErr2Y()>0 ) fhPullY->Fill( (p.GetY() - mcY)/TMath::Sqrt(p.GetErr2Y()) );
+ if( p.GetErr2Z()>0 ) fhPullZ->Fill( (p.GetZ() - mcZ)/TMath::Sqrt(p.GetErr2Z()) );
+ if( p.GetErr2SinPhi()>0 ) fhPullSinPhi->Fill( (p.GetSinPhi() - mcSinPhi)/TMath::Sqrt(p.GetErr2SinPhi()) );
+ if( p.GetErr2DzDs()>0 ) fhPullDzDs->Fill( (p.DzDs() - mcDzDs)/TMath::Sqrt(p.GetErr2DzDs()) );
+ if( p.GetErr2Kappa()>0 ) fhPullQPt->Fill( (qPt - mcQPt)/TMath::Sqrt(p.GetErr2Kappa()*k2QPt*k2QPt) );
+ break;
+ }
+ }
+
+ for (Int_t ipart=0; ipart<fNMCTracks; ipart++) {
+ AliHLTTPCCAMCTrack &mc = fMCTracks[ipart];
+ if( mc.Set()>0 ) fhGBEffVsP->Fill(mc.P(), ( mc.NReconstructed()>0 ?1 :0));
+ }
+
+ if( traLabels ) delete[] traLabels;
+ if( traPurity ) delete[] traPurity;
+
+ fStatGBNRecTot += nRecTot;
+ fStatGBNRecOut += nRecOut;
+ fStatGBNGhost += nGhost;
+ fStatGBNMCAll += nMCAll;
+ fStatGBNRecAll += nRecAll;
+ fStatGBNClonesAll += nClonesAll;
+ fStatGBNMCRef += nMCRef;
+ fStatGBNRecRef += nRecRef;
+ fStatGBNClonesRef += nClonesRef;
+ }
+
// distribution of cluster errors
fhHitErrZ->Fill(hit.ErrZ());
}
}
+
+ // cluster pulls
+
+ if( fDoClusterPulls && fNMCPoints>0 ) {
+
+ {
+ for (Int_t ipart=0; ipart<fNMCTracks; ipart++) {
+ AliHLTTPCCAMCTrack &mc = fMCTracks[ipart];
+ mc.NMCPoints() = 0;
+ }
+ sort(fMCPoints, fMCPoints+fNMCPoints, AliHLTTPCCAMCPoint::Compare );
+
+ for( Int_t ip=0; ip<fNMCPoints; ip++ ){
+ AliHLTTPCCAMCPoint &p = fMCPoints[ip];
+ AliHLTTPCCAMCTrack &t = fMCTracks[p.TrackID()];
+ if( t.NMCPoints()==0 ) t.FirstMCPointID() = ip;
+ t.NMCPoints()++;
+ }
+ }
+
+ for( Int_t ih=0; ih<fNHits; ih++ ){
+
+ AliHLTTPCCAGBHit &hit = fTracker->Hits()[ih];
+ AliHLTTPCCAHitLabel &l = fHitLabels[ih];
+
+ if( l.fLab[0]<0 || l.fLab[0]>=fNMCTracks
+ || l.fLab[1]>=0 || l.fLab[2]>=0 ) continue;
+
+ Int_t lab = l.fLab[0];
+
+ AliHLTTPCCAMCTrack &track = fMCTracks[lab];
+ //if( track.Pt()<1. ) continue;
+ Int_t ip1=-1, ip2=-1;
+ Double_t d1 = 1.e20, d2=1.e20;
+ for( Int_t ip=0; ip<track.NMCPoints(); ip++ ){
+ AliHLTTPCCAMCPoint &p = fMCPoints[track.FirstMCPointID() + ip];
+ if( p.ISlice() != hit.ISlice() ) continue;
+ Double_t dx = p.Sx()-hit.X();
+ Double_t dy = p.Sy()-hit.Y();
+ Double_t dz = p.Sz()-hit.Z();
+ Double_t d = dx*dx + dy*dy + dz*dz;
+ if( p.Sx()< hit.X() ){
+ if( d<d1 ){
+ d1 = d;
+ ip1 = ip;
+ }
+ }else{
+ if( d<d2 ){
+ d2 = d;
+ ip2 = ip;
+ }
+ }
+ }
+
+ if( ip1<0 || ip2<0 ) continue;
+
+ AliHLTTPCCAMCPoint &p1 = fMCPoints[track.FirstMCPointID() + ip1];
+ AliHLTTPCCAMCPoint &p2 = fMCPoints[track.FirstMCPointID() + ip2];
+ Double_t dx = p2.Sx() - p1.Sx();
+ Double_t dy = p2.Sy() - p1.Sy();
+ Double_t dz = p2.Sz() - p1.Sz();
+ if( TMath::Abs(dx)>1.e-8 && TMath::Abs(p1.Sx()-hit.X())<2. && TMath::Abs(p2.Sx()-hit.X())<2. ){
+ Double_t sx = hit.X();
+ Double_t sy = p1.Sy() + dy/dx*(sx-p1.Sx());
+ Double_t sz = p1.Sz() + dz/dx*(sx-p1.Sx());
+
+ Float_t errY, errZ;
+ {
+ AliHLTTPCCATrackParam t;
+ t.Z() = sz;
+ t.SinPhi() = dy/TMath::Sqrt(dx*dx+dy*dy);
+ t.CosPhi() = dx/TMath::Sqrt(dx*dx+dy*dy);
+ t.DzDs() = dz/TMath::Sqrt(dx*dx+dy*dy);
+ fTracker->GetErrors2(hit,t,errY, errZ );
+ errY = TMath::Sqrt(errY);
+ errZ = TMath::Sqrt(errZ);
+ }
+
+ fhHitResY->Fill((hit.Y()-sy));
+ fhHitResZ->Fill((hit.Z()-sz));
+ fhHitPullY->Fill((hit.Y()-sy)/errY);
+ fhHitPullZ->Fill((hit.Z()-sz)/errZ);
+ if( track.Pt()>=1. ){
+ fhHitResY1->Fill((hit.Y()-sy));
+ fhHitResZ1->Fill((hit.Z()-sz));
+ fhHitPullY1->Fill((hit.Y()-sy)/errY);
+ fhHitPullZ1->Fill((hit.Z()-sz)/errZ);
+ }
+ }
+ }
+ }
+
+ {
+ cout<<"\nSlice tracker performance: \n"<<endl;
+ cout<<" N tracks : "
+ <<fStatNMCAll/fStatNEvents<<" mc all, "
+ <<fStatNMCRef/fStatNEvents<<" mc ref, "
+ <<fStatNRecTot/fStatNEvents<<" rec total, "
+ <<fStatNRecAll/fStatNEvents<<" rec all, "
+ <<fStatNClonesAll/fStatNEvents<<" clones all, "
+ <<fStatNRecRef/fStatNEvents<<" rec ref, "
+ <<fStatNClonesRef/fStatNEvents<<" clones ref, "
+ <<fStatNRecOut/fStatNEvents<<" out, "
+ <<fStatNGhost/fStatNEvents<<" ghost"<<endl;
- cout<<" N tracks : "
- <<fStatNMCAll/fStatNEvents<<" mc all, "
- <<fStatNMCRef/fStatNEvents<<" mc ref, "
- <<fStatNRecTot/fStatNEvents<<" rec total, "
- <<fStatNRecAll/fStatNEvents<<" rec all, "
- <<fStatNClonesAll/fStatNEvents<<" clones all, "
- <<fStatNRecRef/fStatNEvents<<" rec ref, "
- <<fStatNClonesRef/fStatNEvents<<" clones ref, "
- <<fStatNRecOut/fStatNEvents<<" out, "
- <<fStatNGhost/fStatNEvents<<" ghost"<<endl;
-
- Int_t nRecExtr = fStatNRecAll - fStatNRecRef;
- Int_t nMCExtr = fStatNMCAll - fStatNMCRef;
- Int_t nClonesExtr = fStatNClonesAll - fStatNClonesRef;
+ Int_t nRecExtr = fStatNRecAll - fStatNRecRef;
+ Int_t nMCExtr = fStatNMCAll - fStatNMCRef;
+ Int_t nClonesExtr = fStatNClonesAll - fStatNClonesRef;
+
+ Double_t dRecTot = (fStatNRecTot>0 ) ? fStatNRecTot :1;
+ Double_t dMCAll = (fStatNMCAll>0 ) ? fStatNMCAll :1;
+ Double_t dMCRef = (fStatNMCRef>0 ) ? fStatNMCRef :1;
+ Double_t dMCExtr = (nMCExtr>0 ) ? nMCExtr :1;
+ Double_t dRecAll = (fStatNRecAll+fStatNClonesAll>0 ) ? fStatNRecAll+fStatNClonesAll :1;
+ Double_t dRecRef = (fStatNRecRef+fStatNClonesRef>0 ) ? fStatNRecRef+fStatNClonesRef :1;
+ Double_t dRecExtr = (nRecExtr+nClonesExtr>0 ) ? nRecExtr+nClonesExtr :1;
+
+ cout<<" EffRef = "<< fStatNRecRef/dMCRef
+ <<", CloneRef = " << fStatNClonesRef/dRecRef <<endl;
+ cout<<" EffExtra = "<<nRecExtr/dMCExtr
+ <<", CloneExtra = " << nClonesExtr/dRecExtr<<endl;
+ cout<<" EffAll = "<<fStatNRecAll/dMCAll
+ <<", CloneAll = " << fStatNClonesAll/dRecAll<<endl;
+ cout<<" Out = "<<fStatNRecOut/dRecTot
+ <<", Ghost = "<<fStatNGhost/dRecTot<<endl;
+ cout<<" Time = "<<fTracker->StatTime(0)/fTracker->StatNEvents()*1.e3<<" msec/event "<<endl;
+ cout<<" Local timers = "
+ <<fTracker->StatTime(1)/fTracker->StatNEvents()*1.e3<<" "
+ <<fTracker->StatTime(2)/fTracker->StatNEvents()*1.e3<<" "
+ <<fTracker->StatTime(3)/fTracker->StatNEvents()*1.e3<<" "
+ <<fTracker->StatTime(4)/fTracker->StatNEvents()*1.e3<<" "
+ <<fTracker->StatTime(5)/fTracker->StatNEvents()*1.e3<<"["
+ <<fTracker->StatTime(6)/fTracker->StatNEvents()*1.e3<<"/"
+ <<fTracker->StatTime(7)/fTracker->StatNEvents()*1.e3<<"] "
+ <<fTracker->StatTime(8)/fTracker->StatNEvents()*1.e3<<" "
+ <<" msec/event "<<endl;
+ }
+
+ {
+ cout<<"\nGlobal tracker performance: \n"<<endl;
+ cout<<" N tracks : "
+ <<fStatGBNMCAll/fStatNEvents<<" mc all, "
+ <<fStatGBNMCRef/fStatNEvents<<" mc ref, "
+ <<fStatGBNRecTot/fStatNEvents<<" rec total, "
+ <<fStatGBNRecAll/fStatNEvents<<" rec all, "
+ <<fStatGBNClonesAll/fStatNEvents<<" clones all, "
+ <<fStatGBNRecRef/fStatNEvents<<" rec ref, "
+ <<fStatGBNClonesRef/fStatNEvents<<" clones ref, "
+ <<fStatGBNRecOut/fStatNEvents<<" out, "
+ <<fStatGBNGhost/fStatNEvents<<" ghost"<<endl;
- Double_t dRecTot = (fStatNRecTot>0 ) ? fStatNRecTot :1;
- Double_t dMCAll = (fStatNMCAll>0 ) ? fStatNMCAll :1;
- Double_t dMCRef = (fStatNMCRef>0 ) ? fStatNMCRef :1;
- Double_t dMCExtr = (nMCExtr>0 ) ? nMCExtr :1;
- Double_t dRecAll = (fStatNRecAll+fStatNClonesAll>0 ) ? fStatNRecAll+fStatNClonesAll :1;
- Double_t dRecRef = (fStatNRecRef+fStatNClonesRef>0 ) ? fStatNRecRef+fStatNClonesRef :1;
- Double_t dRecExtr = (nRecExtr+nClonesExtr>0 ) ? nRecExtr+nClonesExtr :1;
-
- cout<<" EffRef = "<< fStatNRecRef/dMCRef
- <<", CloneRef = " << fStatNClonesRef/dRecRef <<endl;
- cout<<" EffExtra = "<<nRecExtr/dMCExtr
- <<", CloneExtra = " << nClonesExtr/dRecExtr<<endl;
- cout<<" EffAll = "<<fStatNRecAll/dMCAll
- <<", CloneAll = " << fStatNClonesAll/dRecAll<<endl;
- cout<<" Out = "<<fStatNRecOut/dRecTot
- <<", Ghost = "<<fStatNGhost/dRecTot<<endl;
- cout<<" Time = "<<fTracker->StatTime(0)/fTracker->StatNEvents()*1.e3<<" msec/event "<<endl;
- cout<<" Local timers = "
- <<fTracker->StatTime(1)/fTracker->StatNEvents()*1.e3<<" "
- <<fTracker->StatTime(2)/fTracker->StatNEvents()*1.e3<<" "
- <<fTracker->StatTime(3)/fTracker->StatNEvents()*1.e3<<" "
- <<fTracker->StatTime(4)/fTracker->StatNEvents()*1.e3<<" "
- <<fTracker->StatTime(5)/fTracker->StatNEvents()*1.e3<<"["
- <<fTracker->StatTime(6)/fTracker->StatNEvents()*1.e3<<"/"
- <<fTracker->StatTime(7)/fTracker->StatNEvents()*1.e3<<"], merge:"
- <<fTracker->StatTime(8)/fTracker->StatNEvents()*1.e3
- <<" msec/event "<<endl;
+ Int_t nRecExtr = fStatGBNRecAll - fStatGBNRecRef;
+ Int_t nMCExtr = fStatGBNMCAll - fStatGBNMCRef;
+ Int_t nClonesExtr = fStatGBNClonesAll - fStatGBNClonesRef;
+
+ Double_t dRecTot = (fStatGBNRecTot>0 ) ? fStatGBNRecTot :1;
+ Double_t dMCAll = (fStatGBNMCAll>0 ) ? fStatGBNMCAll :1;
+ Double_t dMCRef = (fStatGBNMCRef>0 ) ? fStatGBNMCRef :1;
+ Double_t dMCExtr = (nMCExtr>0 ) ? nMCExtr :1;
+ Double_t dRecAll = (fStatGBNRecAll+fStatGBNClonesAll>0 ) ? fStatGBNRecAll+fStatGBNClonesAll :1;
+ Double_t dRecRef = (fStatGBNRecRef+fStatGBNClonesRef>0 ) ? fStatGBNRecRef+fStatGBNClonesRef :1;
+ Double_t dRecExtr = (nRecExtr+nClonesExtr>0 ) ? nRecExtr+nClonesExtr :1;
+
+ cout<<" EffRef = "<< fStatGBNRecRef/dMCRef
+ <<", CloneRef = " << fStatGBNClonesRef/dRecRef <<endl;
+ cout<<" EffExtra = "<<nRecExtr/dMCExtr
+ <<", CloneExtra = " << nClonesExtr/dRecExtr<<endl;
+ cout<<" EffAll = "<<fStatGBNRecAll/dMCAll
+ <<", CloneAll = " << fStatGBNClonesAll/dRecAll<<endl;
+ cout<<" Out = "<<fStatGBNRecOut/dRecTot
+ <<", Ghost = "<<fStatGBNGhost/dRecTot<<endl;
+ cout<<" Time = "<<( fTracker->StatTime(0)+fTracker->StatTime(9) )/fTracker->StatNEvents()*1.e3<<" msec/event "<<endl;
+ cout<<" Local timers: "<<endl;
+ cout<<" slice tracker "<<fTracker->StatTime(0)/fTracker->StatNEvents()*1.e3<<": "
+ <<fTracker->StatTime(1)/fTracker->StatNEvents()*1.e3<<" "
+ <<fTracker->StatTime(2)/fTracker->StatNEvents()*1.e3<<" "
+ <<fTracker->StatTime(3)/fTracker->StatNEvents()*1.e3<<" "
+ <<fTracker->StatTime(4)/fTracker->StatNEvents()*1.e3<<" "
+ <<fTracker->StatTime(5)/fTracker->StatNEvents()*1.e3<<"["
+ <<fTracker->StatTime(6)/fTracker->StatNEvents()*1.e3<<"/"
+ <<fTracker->StatTime(7)/fTracker->StatNEvents()*1.e3<<"] "
+ <<fTracker->StatTime(8)/fTracker->StatNEvents()*1.e3
+ <<" msec/event "<<endl;
+ cout<<" GB merger "<<fTracker->StatTime(9)/fTracker->StatNEvents()*1.e3<<": "
+ <<fTracker->StatTime(10)/fTracker->StatNEvents()*1.e3<<", "
+ <<fTracker->StatTime(11)/fTracker->StatNEvents()*1.e3<<", "
+ <<fTracker->StatTime(12)/fTracker->StatNEvents()*1.e3<<" "
+ <<" msec/event "<<endl;
+ }
+
WriteHistos();
}
+
+void AliHLTTPCCAPerformance::WriteMCEvent( ostream &out )
+{
+ out<<fNMCTracks<<endl;
+ for( Int_t it=0; it<fNMCTracks; it++ ){
+ AliHLTTPCCAMCTrack &t = fMCTracks[it];
+ out<<it<<" ";
+ out<<t.PDG()<<endl;
+ for( Int_t i=0; i<7; i++ ) out<<t.Par()[i]<<" ";
+ out<<endl<<" ";
+ for( Int_t i=0; i<7; i++ ) out<<t.TPCPar()[i]<<" ";
+ out<<endl<<" ";
+ out<< t.P()<<" ";
+ out<< t.Pt()<<" ";
+ out<< t.NMCPoints()<<" ";
+ out<< t.FirstMCPointID()<<" ";
+ out<< t.NHits()<<" ";
+ out<< t.NReconstructed()<<" ";
+ out<< t.Set()<<" ";
+ out<< t.NTurns()<<endl;
+ }
+
+ out<<fNHits<<endl;
+ for( Int_t ih=0; ih<fNHits; ih++ ){
+ AliHLTTPCCAHitLabel &l = fHitLabels[ih];
+ out<<l.fLab[0]<<" "<<l.fLab[1]<<" "<<l.fLab[2]<<endl;
+ }
+}
+
+void AliHLTTPCCAPerformance::WriteMCPoints( ostream &out )
+{
+ out<<fNMCPoints<<endl;
+ for( Int_t ip=0; ip<fNMCPoints; ip++ ){
+ AliHLTTPCCAMCPoint &p = fMCPoints[ip];
+ out<< p.X()<<" ";
+ out<< p.Y()<<" ";
+ out<< p.Z()<<" ";
+ out<< p.Sx()<<" ";
+ out<< p.Sy()<<" ";
+ out<< p.Sz()<<" ";
+ out<< p.Time()<<" ";
+ out<< p.ISlice()<<" ";
+ out<< p.TrackID()<<endl;
+ }
+}
+
+void AliHLTTPCCAPerformance::ReadMCEvent( istream &in )
+{
+ StartEvent();
+ if( fMCTracks ) delete[] fMCTracks;
+ fMCTracks = 0;
+ fNMCTracks = 0;
+ if( fHitLabels ) delete[] fHitLabels;
+ fHitLabels = 0;
+ fNHits = 0;
+ if( fMCPoints ) delete[] fMCPoints;
+ fMCPoints = 0;
+ fNMCPoints = 0;
+
+ in>>fNMCTracks;
+ fMCTracks = new AliHLTTPCCAMCTrack[fNMCTracks];
+ for( Int_t it=0; it<fNMCTracks; it++ ){
+ AliHLTTPCCAMCTrack &t = fMCTracks[it];
+ Int_t j;
+ in>>j;
+ in>> t.PDG();
+ for( Int_t i=0; i<7; i++ ) in>>t.Par()[i];
+ for( Int_t i=0; i<7; i++ ) in>>t.TPCPar()[i];
+ in>> t.P();
+ in>> t.Pt();
+ in>> t.NHits();
+ in>> t.NMCPoints();
+ in>> t.FirstMCPointID();
+ in>> t.NReconstructed();
+ in>> t.Set();
+ in>> t.NTurns();
+ }
+
+ in>>fNHits;
+ fHitLabels = new AliHLTTPCCAHitLabel[fNHits];
+ for( Int_t ih=0; ih<fNHits; ih++ ){
+ AliHLTTPCCAHitLabel &l = fHitLabels[ih];
+ in>>l.fLab[0]>>l.fLab[1]>>l.fLab[2];
+ }
+}
+
+void AliHLTTPCCAPerformance::ReadMCPoints( istream &in )
+{
+ if( fMCPoints ) delete[] fMCPoints;
+ fMCPoints = 0;
+ fNMCPoints = 0;
+
+ in>>fNMCPoints;
+ fMCPoints = new AliHLTTPCCAMCPoint[fNMCPoints];
+ for( Int_t ip=0; ip<fNMCPoints; ip++ ){
+ AliHLTTPCCAMCPoint &p = fMCPoints[ip];
+ in>> p.X();
+ in>> p.Y();
+ in>> p.Z();
+ in>> p.Sx();
+ in>> p.Sy();
+ in>> p.Sz();
+ in>> p.Time();
+ in>> p.ISlice();
+ in>> p.TrackID();
+ }
+}
class TParticle;
class AliHLTTPCCAMCTrack;
+class AliHLTTPCCAMCPoint;
class AliHLTTPCCAGBTracker;
class TDirectory;
class TH1D;
+class TH2D;
class TProfile;
/**
void StartEvent();
void SetNHits( Int_t NHits );
void SetNMCTracks( Int_t NMCTracks );
+ void SetNMCPoints( Int_t NMCPoints );
+
void ReadHitLabel( Int_t HitID,
Int_t lab0, Int_t lab1, Int_t lab2 );
void ReadMCTrack( Int_t index, const TParticle *part );
+ void ReadMCTPCTrack( Int_t index, Float_t X, Float_t Y, Float_t Z,
+ Float_t Px, Float_t Py, Float_t Pz );
+
+ void ReadMCPoint( Int_t TrackID, Float_t X, Float_t Y, Float_t Z, Float_t Time, Int_t iSlice );
void CreateHistos();
void WriteHistos();
void SlicePerformance( Int_t iSlice, Bool_t PrintFlag );
void Performance();
+ void WriteMCEvent( ostream &out );
+ void ReadMCEvent( istream &in );
+ void WriteMCPoints( ostream &out );
+ void ReadMCPoints( istream &in );
+ Bool_t& DoClusterPulls(){ return fDoClusterPulls; }
+
protected:
AliHLTTPCCAGBTracker *fTracker; //* pointer to the tracker
struct AliHLTTPCCAHitLabel{
Int_t fLab[3]; //* array of 3 MC labels
};
+
AliHLTTPCCAHitLabel *fHitLabels; //* array of hit MC labels
- Int_t fNHits; //* number of hits
- AliHLTTPCCAMCTrack *fMCTracks; //* array of Mc tracks
- Int_t fNMCTracks; //* number of MC tracks
+ Int_t fNHits; //* number of hits
+ AliHLTTPCCAMCTrack *fMCTracks; //* array of MC tracks
+ Int_t fNMCTracks; //* number of MC tracks
+ AliHLTTPCCAMCPoint *fMCPoints; //* array of MC points
+ Int_t fNMCPoints; //* number of MC points
+ Bool_t fDoClusterPulls; //* do cluster pulls (very slow)
Int_t fStatNEvents; //* n of events proceed
Int_t fStatNRecTot; //* total n of reconstructed tracks
Int_t fStatNRecOut; //* n of reconstructed tracks in Out set
Int_t fStatNRecRef; //* n of reconstructed tracks in Ref set
Int_t fStatNClonesRef; //* n of reconstructed clones in Ref set
+ Int_t fStatGBNRecTot; //* global tracker: total n of reconstructed tracks
+ Int_t fStatGBNRecOut; //* global tracker: n of reconstructed tracks in Out set
+ Int_t fStatGBNGhost;//* global tracker: n of reconstructed tracks in Ghost set
+ Int_t fStatGBNMCAll;//* global tracker: n of MC tracks
+ Int_t fStatGBNRecAll; //* global tracker: n of reconstructed tracks in All set
+ Int_t fStatGBNClonesAll;//* global tracker: total n of reconstructed tracks in Clone set
+ Int_t fStatGBNMCRef; //* global tracker: n of MC reference tracks
+ Int_t fStatGBNRecRef; //* global tracker: n of reconstructed tracks in Ref set
+ Int_t fStatGBNClonesRef; //* global tracker: n of reconstructed clones in Ref set
+
TDirectory *fHistoDir; //* ROOT directory with histogramms
+ TH1D
+ *fhResY, //* track Y resolution at the TPC entrance
+ *fhResZ, //* track Z resolution at the TPC entrance
+ *fhResSinPhi, //* track SinPhi resolution at the TPC entrance
+ *fhResDzDs, //* track DzDs resolution at the TPC entrance
+ *fhResPt, //* track Pt relative resolution at the TPC entrance
+ *fhPullY, //* track Y pull at the TPC entrance
+ *fhPullZ, //* track Z pull at the TPC entrance
+ *fhPullSinPhi, //* track SinPhi pull at the TPC entrance
+ *fhPullDzDs, //* track DzDs pull at the TPC entrance
+ *fhPullQPt; //* track Q/Pt pull at the TPC entrance
+
TH1D
*fhHitErrY, //* hit error in Y
- *fhHitErrZ,//* hit error in Z
- *fhHitResX,//* hit resolution X
+ *fhHitErrZ,//* hit error in Z
*fhHitResY,//* hit resolution Y
*fhHitResZ,//* hit resolution Z
- *fhHitPullX,//* hit pull X
*fhHitPullY,//* hit pull Y
- *fhHitPullZ,//* hit pull Z
+ *fhHitPullZ;//* hit pull Z
+
+ TH1D
+ *fhHitResY1,//* hit resolution Y, pt>1GeV
+ *fhHitResZ1,//* hit resolution Z, pt>1GeV
+ *fhHitPullY1,//* hit pull Y, pt>1GeV
+ *fhHitPullZ1;//* hit pull Z, pt>1GeV
+
+ TH1D
*fhCellPurity,//* cell purity
*fhCellNHits//* cell n hits
;
TProfile
*fhCellPurityVsN, //* cell purity vs N hits
*fhCellPurityVsPt,//* cell purity vs MC Pt
- *fhEffVsP; //* reconstruction efficiency vs P plot
+ *fhEffVsP, //* reconstruction efficiency vs P plot
+ *fhGBEffVsP; //* global reconstruction efficiency vs P plot
static void WriteDir2Current( TObject *obj );
--- /dev/null
+// $Id: AliHLTTPCCATrackConvertor.cxx 27042 2008-07-02 12:06:02Z richterm $
+//***************************************************************************
+// This file is property of and copyright by the ALICE HLT Project *
+// ALICE Experiment at CERN, All rights reserved. *
+// *
+// Primary Authors: Sergey Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de> *
+// Ivan Kisel <kisel@kip.uni-heidelberg.de> *
+// for The ALICE HLT Project. *
+// *
+// 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. *
+//***************************************************************************
+
+#include "AliHLTTPCCATrackConvertor.h"
+#include "AliExternalTrackParam.h"
+#include "AliHLTTPCCATrackParam.h"
+#include "TMath.h"
+
+
+void AliHLTTPCCATrackConvertor::GetExtParam( const AliHLTTPCCATrackParam &T1, AliExternalTrackParam &T2, Double_t alpha, Double_t Bz )
+{
+ //* Convert from AliHLTTPCCATrackParam to AliExternalTrackParam parameterisation,
+ //* the angle alpha is the global angle of the local X axis
+
+ Double_t par[5], cov[15];
+ for( Int_t i=0; i<5; i++ ) par[i] = T1.GetPar()[i];
+ for( Int_t i=0; i<15; i++ ) cov[i] = T1.GetCov()[i];
+
+ if(par[2]>.99 ) par[2]=.99;
+ if(par[2]<-.99 ) par[2]=-.99;
+
+ { // kappa => 1/pt
+ const Double_t kCLight = 0.000299792458;
+ Double_t c = 1.e4;
+ if( TMath::Abs(Bz)>1.e-4 ) c = 1./(Bz*kCLight);
+ par[4] *= c;
+ cov[10]*= c;
+ cov[11]*= c;
+ cov[12]*= c;
+ cov[13]*= c;
+ cov[14]*= c*c;
+ }
+ if( T1.GetCosPhi()<0 ){ // change direction
+ par[2] = -par[2]; // sin phi
+ par[3] = -par[3]; // DzDs
+ par[4] = -par[4]; // kappa
+ cov[3] = -cov[3];
+ cov[4] = -cov[4];
+ cov[6] = -cov[6];
+ cov[7] = -cov[7];
+ cov[10] = -cov[10];
+ cov[11] = -cov[11];
+ }
+ T2.Set(T1.GetX(),alpha,par,cov);
+}
+
+void AliHLTTPCCATrackConvertor::SetExtParam( AliHLTTPCCATrackParam &T1, const AliExternalTrackParam &T2, Double_t Bz )
+{
+ //* Convert from AliExternalTrackParam parameterisation
+
+ for( Int_t i=0; i<5; i++ ) T1.Par()[i] = T2.GetParameter()[i];
+ for( Int_t i=0; i<15; i++ ) T1.Cov()[i] = T2.GetCovariance()[i];
+ T1.X() = T2.GetX();
+ if(T1.SinPhi()>.99 ) T1.SinPhi()=.99;
+ if(T1.SinPhi()<-.99 ) T1.SinPhi()=-.99;
+ T1.CosPhi() = TMath::Sqrt(1.-T1.SinPhi()*T1.SinPhi());
+ const Double_t kCLight = 0.000299792458;
+ Double_t c = Bz*kCLight;
+ { // 1/pt -> kappa
+ T1.Par()[4] *= c;
+ T1.Cov()[10]*= c;
+ T1.Cov()[11]*= c;
+ T1.Cov()[12]*= c;
+ T1.Cov()[13]*= c;
+ T1.Cov()[14]*= c*c;
+ }
+}
+
--- /dev/null
+//-*- Mode: C++ -*-
+
+//* This file is property of and copyright by the ALICE HLT Project *
+//* ALICE Experiment at CERN, All rights reserved. *
+//* See cxx source for full Copyright notice *
+
+
+#ifndef ALIHLTTPCCATRACKCONVERTOR_H
+#define ALIHLTTPCCATRACKCONVERTOR_H
+#include "Rtypes.h"
+
+class AliExternalTrackParam;
+class AliHLTTPCCATrackParam;
+
+/**
+ * @class AliHLTTPCCATrackConvertor
+ *
+ * AliHLTTPCCATrackConvertor class converts tracks to different parameterisations
+ * it is used by the AliHLTTPCCATrackerComponent
+ *
+ */
+class AliHLTTPCCATrackConvertor
+{
+ public:
+
+ AliHLTTPCCATrackConvertor(){}
+
+ static void GetExtParam( const AliHLTTPCCATrackParam &T1, AliExternalTrackParam &T2, Double_t alpha, Double_t Bz );
+ static void SetExtParam( AliHLTTPCCATrackParam &T1, const AliExternalTrackParam &T2, Double_t Bz );
+
+};
+
+
+#endif
#include "AliHLTTPCCATrackParam.h"
#include "TMath.h"
-#include "AliExternalTrackParam.h"
+#include "Riostream.h"
//ClassImp(AliHLTTPCCATrackParam)
Float_t s0 = sigmaY2[0];
Float_t s1 = sigmaY2[1];
Float_t s2 = sigmaY2[2];
-
+
fC[0] = s0;
fC[1] = 0;
- fC[2] = 0;
+ fC[2] = 100.;
fC[3] = d0[0]*s0;
fC[4] = 0;
fC[6] = 0;
fC[7] = 0;
fC[8] = 0;
- fC[9] = 0;
+ fC[9] = 100.;
fC[10] = d0[1]*s0;
fC[11] = 0;
Float_t dSin = dl*k/2;
if( dSin > 1 ) dSin = 1;
if( dSin <-1 ) dSin = -1;
- Float_t dS = ( TMath::Abs(k)>1.e-4) ? (2*TMath::ASin(dSin)/k) :dl;
+ Float_t dS = ( TMath::Abs(k)>1.e-4) ? (2*TMath::ASin(dSin)/k) :dl;
Float_t dz = dS*DzDs();
Float_t cci = 0, exi = 0, ex1i = 0;
else ret = 0;
if( !ret ) return ret;
+
X() += dx;
fP[0]+= dy;
fP[1]+= dz;
return ret;
}
+Bool_t AliHLTTPCCATrackParam::TransportToXWithMaterial( Float_t X, Float_t Bz )
+{
+ AliHLTTPCCATrackFitParam par;
+ CalculateFitParameters( par, Bz );
+ return TransportToXWithMaterial(X, par );
+}
-Bool_t AliHLTTPCCATrackParam::Rotate( Float_t alpha )
+Bool_t AliHLTTPCCATrackParam::TransportToXWithMaterial( Float_t x, AliHLTTPCCATrackFitParam &par )
{
- //* Rotate the coordinate system in XY on the angle alpha
+ //* Transport the track parameters to X=x
Bool_t ret = 1;
+ Float_t oldX=GetX();
+
+ Float_t x0 = X();
+ //Float_t y0 = Y();
+ Float_t k = Kappa();
+ Float_t ex = CosPhi();
+ Float_t ey = SinPhi();
+ Float_t dx = x - x0;
+
+ Float_t ey1 = k*dx + ey;
+ Float_t ex1;
+ if( TMath::Abs(ey1)>.99 ){ // no intersection -> check the border
+ ey1 = ( ey1>0 ) ?1 :-1;
+ ex1 = 0;
+ dx = ( TMath::Abs(k)>1.e-4) ? ( (ey1-ey)/k ) :0;
+
+ Float_t ddx = TMath::Abs(x0+dx - x)*k*k;
+ Float_t hx[] = {0, -k, 1+ey };
+ Float_t sx2 = hx[1]*hx[1]*fC[ 3] + hx[2]*hx[2]*fC[ 5];
+ if( ddx*ddx>3.5*3.5*sx2 ) ret = 0; // x not withing the error
+ ret = 0; // any case
+ return ret;
+ }else{
+ ex1 = TMath::Sqrt(1 - ey1*ey1);
+ if( ex<0 ) ex1 = -ex1;
+ }
+
+ Float_t dx2 = dx*dx;
+ CosPhi() = ex1;
+ Float_t ss = ey+ey1;
+ Float_t cc = ex+ex1;
+ Float_t tg = 0;
+ if( TMath::Abs(cc)>1.e-4 ) tg = ss/cc; // tan((phi1+phi)/2)
+ else ret = 0;
+ Float_t dy = dx*tg;
+ Float_t dl = dx*TMath::Sqrt(1+tg*tg);
+
+ if( cc<0 ) dl = -dl;
+ Float_t dSin = dl*k/2;
+ if( dSin > 1 ) dSin = 1;
+ if( dSin <-1 ) dSin = -1;
+ Float_t dS = ( TMath::Abs(k)>1.e-4) ? (2*TMath::ASin(dSin)/k) :dl;
+ Float_t dz = dS*DzDs();
+
+ Float_t cci = 0, exi = 0, ex1i = 0;
+ if( TMath::Abs(cc)>1.e-4 ) cci = 1./cc;
+ else ret = 0;
+ if( TMath::Abs(ex)>1.e-4 ) exi = 1./ex;
+ else ret = 0;
+ if( TMath::Abs(ex1)>1.e-4 ) ex1i = 1./ex1;
+ else ret = 0;
+
+ if( !ret ) return ret;
+
+ X() += dx;
+ fP[0]+= dy;
+ fP[1]+= dz;
+ fP[2] = ey1;
+ fP[3] = fP[3];
+ fP[4] = fP[4];
+
+ Float_t h2 = dx*(1+ ex*ex1 + ey*ey1 )*cci*exi*ex1i;
+ Float_t h4 = dx2*(cc + ss*ey1*ex1i )*cci*cci;
+
+ Float_t c00 = fC[0];
+ Float_t c10 = fC[1];
+ Float_t c11 = fC[2];
+ Float_t c20 = fC[3];
+ Float_t c21 = fC[4];
+ Float_t c22 = fC[5];
+ Float_t c30 = fC[6];
+ Float_t c31 = fC[7];
+ Float_t c32 = fC[8];
+ Float_t c33 = fC[9];
+ Float_t c40 = fC[10];
+ Float_t c41 = fC[11];
+ Float_t c42 = fC[12];
+ Float_t c43 = fC[13];
+ Float_t c44 = fC[14];
+
+ //Float_t H0[5] = { 1,0, h2, 0, h4 };
+ //Float_t H1[5] = { 0, 1, 0, dS, 0 };
+ //Float_t H2[5] = { 0, 0, 1, 0, dx };
+ //Float_t H3[5] = { 0, 0, 0, 1, 0 };
+ //Float_t H4[5] = { 0, 0, 0, 0, 1 };
+
+
+ fC[0]=( c00 + h2*h2*c22 + h4*h4*c44
+ + 2*( h2*c20 + h4*c40 + h2*h4*c42 ) );
+
+ fC[1]= c10 + h2*c21 + h4*c41 + dS*(c30 + h2*c32 + h4*c43);
+ fC[2]= c11 + 2*dS*c31 + dS*dS*c33;
+
+ fC[3]= c20 + h2*c22 + h4*c42 + dx*( c40 + h2*c42 + h4*c44);
+ fC[4]= c21 + dS*c32 + dx*(c41 + dS*c43);
+ fC[5]= c22 +2*dx*c42 + dx2*c44;
+
+ fC[6]= c30 + h2*c32 + h4*c43;
+ fC[7]= c31 + dS*c33;
+ fC[8]= c32 + dx*c43;
+ fC[9]= c33;
+
+ fC[10]= c40 + h2*c42 + h4*c44;
+ fC[11]= c41 + dS*c43;
+ fC[12]= c42 + dx*c44;
+ fC[13]= c43;
+ fC[14]= c44;
+
+ Float_t d = TMath::Sqrt(dS*dS + dz*dz );
+
+ if (oldX > GetX() ) d = -d;
+ {
+ Float_t rho=0.9e-3;
+ Float_t radLen=28.94;
+ CorrectForMeanMaterial(d*rho/radLen,d*rho,par);
+ }
+
+ return ret;
+}
+
+
+
+Float_t AliHLTTPCCATrackParam::ApproximateBetheBloch( Float_t beta2 )
+{
+ //------------------------------------------------------------------
+ // This is an approximation of the Bethe-Bloch formula with
+ // the density effect taken into account at beta*gamma > 3.5
+ // (the approximation is reasonable only for solid materials)
+ //------------------------------------------------------------------
+ if (beta2 >= 1) return 0;
+
+ if (beta2/(1-beta2)>3.5*3.5)
+ return 0.153e-3/beta2*( log(3.5*5940)+0.5*log(beta2/(1-beta2)) - beta2);
+ return 0.153e-3/beta2*(log(5940*beta2/(1-beta2)) - beta2);
+}
+
+
+void AliHLTTPCCATrackParam::CalculateFitParameters( AliHLTTPCCATrackFitParam &par, Float_t Bz, Float_t mass )
+{
+ //*!
+
+ const Float_t kCLight = 0.000299792458;
+ Float_t c = Bz*kCLight;
+ Float_t p2 = (1.+ fP[3]*fP[3])*c*c;
+ Float_t k2 = fP[4]*fP[4];
+ Float_t beta2= p2 / (p2 + mass*mass*k2);
+ Float_t Bethe = ApproximateBetheBloch(beta2);
+
+ Float_t P2 = (k2>1.e-8) ?p2/k2 :10000; // impuls 2
+ par.fBethe = Bethe;
+ par.fE = TMath::Sqrt( P2 + mass*mass);
+ par.fTheta2 = 14.1*14.1/(beta2*p2*1e6)*k2;
+ par.fEP2 = par.fE/p2*k2;
+
+ // Approximate energy loss fluctuation (M.Ivanov)
+
+ const Float_t knst=0.07; // To be tuned.
+ par.fSigmadE2 = knst*par.fEP2*fP[4];
+ par.fSigmadE2 = par.fSigmadE2 * par.fSigmadE2;
+
+ par.fK22 = (1. + fP[3]*fP[3]);
+ par.fK33 = par.fK22*par.fK22;
+ par.fK43 = fP[3]*fP[4]*par.fK22;
+ par.fK44 = fP[3]*fP[3]*fP[4]*fP[4];
+}
+
+
+Bool_t AliHLTTPCCATrackParam::CorrectForMeanMaterial( Float_t xOverX0, Float_t xTimesRho, AliHLTTPCCATrackFitParam &par )
+{
+ //------------------------------------------------------------------
+ // This function corrects the track parameters for the crossed material.
+ // "xOverX0" - X/X0, the thickness in units of the radiation length.
+ // "xTimesRho" - is the product length*density (g/cm^2).
+ //------------------------------------------------------------------
+
+ Float_t &fC22=fC[5];
+ Float_t &fC33=fC[9];
+ Float_t &fC40=fC[10];
+ Float_t &fC41=fC[11];
+ Float_t &fC42=fC[12];
+ Float_t &fC43=fC[13];
+ Float_t &fC44=fC[14];
+
+ //Energy losses************************
+
+ Float_t dE = par.fBethe*xTimesRho;
+ if ( TMath::Abs(dE) > 0.3*par.fE ) return 0; //30% energy loss is too much!
+ Float_t corr = (1.- par.fEP2*dE);
+ if( corr<0.3 ) return 0;
+ fP[4]*= corr;
+ fC40*= corr;
+ fC41*= corr;
+ fC42*= corr;
+ fC43*= corr;
+ fC44*= corr*corr;
+ fC44+= par.fSigmadE2*TMath::Abs(dE);
+
+
+ //Multiple scattering******************
+
+ Float_t theta2 = par.fTheta2*TMath::Abs(xOverX0);
+ fC22 += theta2*par.fK22*(1.- fP[2]*fP[2]);
+ fC33 += theta2*par.fK33;
+ fC43 += theta2*par.fK43;
+ fC44 += theta2*par.fK44;
+
+ return 1;
+}
+
+
+
+#include "Riostream.h"
+
+Bool_t AliHLTTPCCATrackParam::Rotate( Float_t alpha )
+{
+ //* Rotate the coordinate system in XY on the angle alpha
+
Float_t cA = TMath::Cos( alpha );
Float_t sA = TMath::Sin( alpha );
Float_t x = X(), y= Y(), sP= SinPhi(), cP= CosPhi();
-
+ Float_t cosPhi = cP*cA + sP*sA;
+ Float_t sinPhi =-cP*sA + sP*cA;
+
+ if( TMath::Abs(sinPhi)>.99 || TMath::Abs(cosPhi)<1.e-2 || TMath::Abs(cP)<1.e-2 ) return 0;
+
+ Float_t j0 = cP/cosPhi;
+ Float_t j2 = cosPhi/cP;
+
X() = x*cA + y*sA;
Y() = -x*sA + y*cA;
- CosPhi() = cP*cA + sP*sA;
- SinPhi() = -cP*sA + sP*cA;
+ CosPhi() = cosPhi;
+ SinPhi() = sinPhi;
- Float_t j0 = 0, j2 = 0;
-
- if( TMath::Abs(CosPhi())>1.e-4 ) j0 = cP/CosPhi(); else ret = 0;
- if( TMath::Abs(cP)>1.e-4 ) j2 = CosPhi()/cP; else ret = 0;
//Float_t J[5][5] = { { j0, 0, 0, 0, 0 }, // Y
// { 0, 1, 0, 0, 0 }, // Z
// { 0, 0, j2, 0, 0 }, // SinPhi
// { 0, 0, 0, 1, 0 }, // DzDs
// { 0, 0, 0, 0, 1 } }; // Kappa
-
+ //cout<<"alpha="<<alpha<<" "<<x<<" "<<y<<" "<<sP<<" "<<cP<<" "<<j0<<" "<<j2<<endl;
+ //cout<<" "<<fC[0]<<" "<<fC[1]<<" "<<fC[6]<<" "<<fC[10]<<" "<<fC[4]<<" "<<fC[5]<<" "<<fC[8]<<" "<<fC[12]<<endl;
fC[0]*= j0*j0;
fC[1]*= j0;
//fC[3]*= j0;
fC[5]*= j2*j2;
fC[8]*= j2;
fC[12]*= j2;
- return ret;
+ //cout<<" "<<fC[0]<<" "<<fC[1]<<" "<<fC[6]<<" "<<fC[10]<<" "<<fC[4]<<" "<<fC[5]<<" "<<fC[8]<<" "<<fC[12]<<endl;
+ return 1;
}
-void AliHLTTPCCATrackParam::GetExtParam( AliExternalTrackParam &T, Double_t alpha, Double_t Bz ) const
-{
- //* Convert to AliExternalTrackParam parameterisation,
- //* the angle alpha is the global angle of the local X axis
-
- Double_t par[5], cov[15];
- for( Int_t i=0; i<5; i++ ) par[i] = fP[i];
- for( Int_t i=0; i<15; i++ ) cov[i] = fC[i];
-
- if(par[2]>.999 ) par[2]=.999;
- if(par[2]<-.999 ) par[2]=-.999;
-
- const Double_t kCLight = 0.000299792458;
- Double_t c = (TMath::Abs(Bz)>1.e-4) ?1./(Bz*kCLight) :1./(5.*kCLight);
- { // kappa => 1/pt
- par[4] *= c;
- cov[10]*= c;
- cov[11]*= c;
- cov[12]*= c;
- cov[13]*= c;
- cov[14]*= c*c;
- }
- if( GetCosPhi()<0 ){ // change direction
- par[2] = -par[2]; // sin phi
- par[3] = -par[3]; // DzDs
- par[4] = -par[4]; // kappa
- cov[3] = -cov[3];
- cov[4] = -cov[4];
- cov[6] = -cov[6];
- cov[7] = -cov[7];
- cov[10] = -cov[10];
- cov[11] = -cov[11];
- }
- T.Set(GetX(),alpha,par,cov);
-}
-
-void AliHLTTPCCATrackParam::SetExtParam( const AliExternalTrackParam &T, Double_t Bz )
-{
- //* Convert from AliExternalTrackParam parameterisation
-
- for( Int_t i=0; i<5; i++ ) fP[i] = T.GetParameter()[i];
- for( Int_t i=0; i<15; i++ ) fC[i] = T.GetCovariance()[i];
- X() = T.GetX();
- if(SinPhi()>.999 ) SinPhi()=.999;
- if(SinPhi()<-.999 ) SinPhi()=-.999;
- CosPhi() = TMath::Sqrt(1.-SinPhi()*SinPhi());
- const Double_t kCLight = 0.000299792458;
- Double_t c = Bz*kCLight;
- { // 1/pt -> kappa
- fP[4] *= c;
- fC[10]*= c;
- fC[11]*= c;
- fC[12]*= c;
- fC[13]*= c;
- fC[14]*= c*c;
- }
-}
-
-
-void AliHLTTPCCATrackParam::Filter( Float_t y, Float_t z, Float_t erry, Float_t errz )
+Bool_t AliHLTTPCCATrackParam::Filter2( Float_t y, Float_t z, Float_t err2Y, Float_t err2Z )
{
//* Add the y,z measurement with the Kalman filter
z0 = y-fP[0],
z1 = z-fP[1];
- Float_t v[3] = {erry*erry, 0, errz*errz};
+ Float_t v[3] = {err2Y, 0, err2Z};
Float_t mS[3] = { c00+v[0], c10+v[1], c11+v[2] };
Float_t mSi[3];
Float_t det = (mS[0]*mS[2] - mS[1]*mS[1]);
- if( TMath::Abs(det)<1.e-8 ) return;
+ if( det < 1.e-8 ) return 0;
det = 1./det;
mSi[0] = mS[2]*det;
mSi[1] = -mS[1]*det;
mSi[2] = mS[0]*det;
- fNDF += 2;
- fChi2 += ( +(mSi[0]*z0 + mSi[1]*z1 )*z0
- +(mSi[1]*z0 + mSi[2]*z1 )*z1 );
-
// K = CHtS
Float_t k00, k01 , k10, k11, k20, k21, k30, k31, k40, k41;
k40 = c40*mSi[0] + c41*mSi[1]; k41 = c40*mSi[1] + c41*mSi[2] ;
Float_t sinPhi = fP[2] + k20*z0 + k21*z1 ;
- if( TMath::Abs(sinPhi)>=0.99 ) return;
+ if( TMath::Abs(sinPhi)>=0.99 ) return 0;
+
+ fNDF += 2;
+ fChi2 += ( +(mSi[0]*z0 + mSi[1]*z1 )*z0
+ +(mSi[1]*z0 + mSi[2]*z1 )*z1 );
fP[ 0]+= k00*z0 + k01*z1 ;
fP[ 1]+= k10*z0 + k11*z1 ;
- fP[ 2]+= k20*z0 + k21*z1 ;
+ fP[ 2] = sinPhi;
fP[ 3]+= k30*z0 + k31*z1 ;
fP[ 4]+= k40*z0 + k41*z1 ;
}else{
CosPhi() = -TMath::Sqrt(1-SinPhi()*SinPhi());
}
+ return 1;
+}
+
+void AliHLTTPCCATrackParam::FilterY( Float_t y, Float_t erry )
+{
+ //* Add the y measurement with the Kalman filter
+
+ Float_t
+ c00 = fC[ 0],
+ c10 = fC[ 1],
+ c20 = fC[ 3],
+ c30 = fC[ 6],
+ c40 = fC[10];
+
+ Float_t
+ z0 = y-fP[0];
+
+ Float_t s = { c00+erry*erry };
+ if( TMath::Abs(s)<1.e-4 ) return;
+
+ Float_t si = 1/s;
+
+ fNDF += 1;
+ fChi2 += si*z0*z0;
+
+ // K = CHtS
+
+ Float_t k0, k1 , k2, k3, k4;
+ k0 = c00*si;
+ k1 = c10*si;
+ k2 = c20*si;
+ k3 = c30*si;
+ k4 = c40*si;
+
+ Float_t sinPhi = fP[2] + k2*z0 ;
+ if( TMath::Abs(sinPhi)>=0.99 ) return;
+
+ fP[ 0]+= k0*z0 ;
+ fP[ 1]+= k1*z0 ;
+ fP[ 2] = sinPhi;
+ fP[ 3]+= k3*z0 ;
+ fP[ 4]+= k4*z0 ;
+
+ fC[ 0]-= k0*c00;
+
+ fC[ 1]-= k1*c00;
+ fC[ 2]-= k1*c10;
+
+ fC[ 3]-= k2*c00;
+ fC[ 4]-= k2*c10;
+ fC[ 5]-= k2*c20;
+
+ fC[ 6]-= k3*c00;
+ fC[ 7]-= k3*c10;
+ fC[ 8]-= k3*c20;
+ fC[ 9]-= k3*c30;
+
+ fC[10]-= k4*c00;
+ fC[11]-= k4*c10;
+ fC[12]-= k4*c20;
+ fC[13]-= k4*c30;
+ fC[14]-= k4*c40;
+
+ if( CosPhi()>=0 ){
+ CosPhi() = TMath::Sqrt(1-SinPhi()*SinPhi());
+ }else{
+ CosPhi() = -TMath::Sqrt(1-SinPhi()*SinPhi());
+ }
+
+}
+
+void AliHLTTPCCATrackParam::FilterZ( Float_t z, Float_t errz )
+{
+ //* Add the z measurement with the Kalman filter
+
+ Float_t
+ c01 = fC[ 1],
+ c11 = fC[ 2],
+ c21 = fC[ 4],
+ c31 = fC[ 7],
+ c41 = fC[11];
+
+ Float_t
+ z1 = z-fP[1];
+
+ Float_t s = c11 + errz*errz;
+ if( TMath::Abs(s)<1.e-4 ) return;
+
+ Float_t si = 1./s;
+
+ fNDF += 1;
+ fChi2 += si*z1*z1;
+
+ // K = CHtS
+
+ Float_t k0, k1 , k2, k3, k4;
+
+ k0 = 0;//c01*si;
+ k1 = c11*si;
+ k2 = 0;//c21*si;
+ k3 = c31*si;
+ k4 = 0;//c41*si;
+
+ Float_t sinPhi = fP[2] + k2*z1 ;
+ if( TMath::Abs(sinPhi)>=0.99 ) return;
+
+ fP[ 0]+= k0*z1 ;
+ fP[ 1]+= k1*z1 ;
+ fP[ 2] = sinPhi;
+ fP[ 3]+= k3*z1 ;
+ fP[ 4]+= k4*z1 ;
+
+
+ fC[ 0]-= k0*c01 ;
+
+ fC[ 1]-= k1*c01 ;
+ fC[ 2]-= k1*c11 ;
+
+ fC[ 3]-= k2*c01 ;
+ fC[ 4]-= k2*c11 ;
+ fC[ 5]-= k2*c21 ;
+
+ fC[ 6]-= k3*c01 ;
+ fC[ 7]-= k3*c11 ;
+ fC[ 8]-= k3*c21 ;
+ fC[ 9]-= k3*c31 ;
+
+ fC[10]-= k4*c01 ;
+ fC[11]-= k4*c11 ;
+ fC[12]-= k4*c21 ;
+ fC[13]-= k4*c31 ;
+ fC[14]-= k4*c41 ;
+
+ if( CosPhi()>=0 ){
+ CosPhi() = TMath::Sqrt(1-SinPhi()*SinPhi());
+ }else{
+ CosPhi() = -TMath::Sqrt(1-SinPhi()*SinPhi());
+ }
+
+}
+
+void AliHLTTPCCATrackParam::Print() const
+{
+ cout<<"track: "<<GetX()<<" "<<GetY()<<" "<<GetZ()<<" "<<GetSinPhi()<<" "<<GetDzDs()<<" "<<GetKappa()<<endl;
+ cout<<"errs2: "<<GetErr2Y()<<" "<<GetErr2Z()<<" "<<GetErr2SinPhi()<<" "<<GetErr2DzDs()<<" "<<GetErr2Kappa()<<endl;
}
#include "Rtypes.h"
+class AliHLTTPCCATrackFitParam
+{
+public:
+ Float_t fBethe, fE,fTheta2, fEP2, fSigmadE2, fK22,fK33,fK43,fK44;
+};
+
/**
* @class AliHLTTPCCATrackParam
*
* which is used by the AliHLTTPCCATracker slice tracker.
*
*/
-class AliExternalTrackParam;
-
class AliHLTTPCCATrackParam
{
public:
Float_t *Par(){ return fP; }
Float_t *Cov(){ return fC; }
+ const Float_t *GetPar() const { return fP; }
+ const Float_t *GetCov() const { return fC; }
+
void ConstructXY3( const Float_t x[3], const Float_t y[3], const Float_t sigmaY2[3], Float_t CosPhi0 );
void ConstructXYZ3( const Float_t p0[5], const Float_t p1[5], const Float_t p2[5],
Float_t &px, Float_t &py, Float_t &pz ) const;
Bool_t TransportToX( Float_t X );
+ Bool_t TransportToXWithMaterial( Float_t X, AliHLTTPCCATrackFitParam &par );
+ Bool_t TransportToXWithMaterial( Float_t X, Float_t Bz );
Bool_t Rotate( Float_t alpha );
- void GetExtParam( AliExternalTrackParam &T, Double_t alpha, Double_t Bz ) const;
- void SetExtParam( const AliExternalTrackParam &T, Double_t Bz );
- void Filter( Float_t y, Float_t z, Float_t erry, Float_t errz );
+ Bool_t Filter2( Float_t y, Float_t z, Float_t err2Y, Float_t err2Z );
+ void FilterY( Float_t y, Float_t erry );
+ void FilterZ( Float_t z, Float_t errz );
+
+ static Float_t ApproximateBetheBloch( Float_t beta2 );
+ void CalculateFitParameters( AliHLTTPCCATrackFitParam &par, Float_t Bz, Float_t mass = 0.13957 );
+ Bool_t CorrectForMeanMaterial( Float_t xOverX0, Float_t xTimesRho, AliHLTTPCCATrackFitParam &par );
+ void Print() const;
private:
}
//AliHLTTPCCADisplay::Instance().Init();
- AliHLTTPCCADisplay::Instance().SetCurrentSector( this );
- AliHLTTPCCADisplay::Instance().SetSectorView();
- AliHLTTPCCADisplay::Instance().DrawSector( this );
- for( Int_t iRow=0; iRow<fParam.NRows(); iRow++ )
- for (Int_t i = 0; i<fRows[iRow].NHits(); i++)
- AliHLTTPCCADisplay::Instance().DrawHit( iRow, i );
- AliHLTTPCCADisplay::Instance().Ask();
+ AliHLTTPCCADisplay::Instance().SetCurrentSlice( this );
+ AliHLTTPCCADisplay::Instance().SetSliceView();
+ AliHLTTPCCADisplay::Instance().DrawSlice( this );
+ //for( Int_t iRow=0; iRow<fParam.NRows(); iRow++ )
+ //for (Int_t i = 0; i<fRows[iRow].NHits(); i++)
+ //AliHLTTPCCADisplay::Instance().DrawHit( iRow, i );
+ //AliHLTTPCCADisplay::Instance().Ask();
#endif
- //if( fParam.ISec()!=8 ) return;
- //if( fParam.ISec()!=8 ) return;
- //if( fParam.ISec()!=33 ) return;
- //if( fParam.ISec()!=2 ) return;
fTimers[0] = 0;
fTimers[1] = 0;
fTimers[2] = 0;
fTimers[4] = 0;
fTimers[5] = 0;
fTimers[6] = 0;
+ fTimers[7] = 0;
if( fNHitsTotal < 1 ) return;
//cout<<"Find Cells..."<<endl;
// cell.Track = -1
//
- TStopwatch timer1;
+ TStopwatch timer;
Int_t nStartCells = 0;
for( Int_t iRow1=0; iRow1<fParam.NRows(); iRow1++ ){
Float_t deltaY = row1.DeltaY();
Float_t deltaZ = row1.DeltaZ();
-
+ Float_t xStep = 1;
+ if( iRow1 < fParam.NRows()-1 ) xStep = fParam.RowX(iRow1+1) - fParam.RowX(iRow1);
+ Float_t tx = xStep/row1.X();
+
Int_t lastRow2 = iRow1+3;
if( lastRow2>=fParam.NRows() ) lastRow2 = fParam.NRows()-1;
for (Int_t i1 = 0; i1<row1.NCells(); i1++){
AliHLTTPCCACell &c1 = row1.Cells()[i1];
//cout<<"row, cell= "<<iRow1<<" "<<i1<<" "<<c1.Y()<<" "<<c1.ErrY()<<" "<<c1.Z()<<" "<<c1.ErrZ()<<endl;
- //Float_t sy1 = c1.ErrY()*c1.ErrY();
- Float_t yMin = c1.Y() - 3.5*c1.ErrY() - deltaY;
- Float_t yMax = c1.Y() + 3.5*c1.ErrY() + deltaY;
- Float_t zMin = c1.Z() - 3.5*c1.ErrZ() - deltaZ;
- Float_t zMax = c1.Z() + 3.5*c1.ErrZ() + deltaZ;
+ //Float_t sy1 = c1.ErrY()*c1.ErrY();
+ Float_t yy = c1.Y() +tx*c1.Y();
+ Float_t zz = c1.Z() +tx*c1.Z();
+
+ Float_t yMin = yy - 3.5*c1.ErrY() - deltaY;
+ Float_t yMax = yy + 3.5*c1.ErrY() + deltaY;
+ Float_t zMin = zz - 3.5*c1.ErrZ() - deltaZ;
+ Float_t zMax = zz + 3.5*c1.ErrZ() + deltaZ;
//Float_t sz1 = c1.ErrZ()*c1.ErrZ();
if( c1.Status()<=0 ) nStartCells++;
}
}//row1
- timer1.Stop();
- fTimers[1] = timer1.CpuTime();
+ timer.Stop();
+ fTimers[1] = timer.CpuTime();
// Second step: create tracks
// for each sequence of neighbouring Cells create Track object
AliHLTTPCCATrack &iTrack = fTracks[itr];
//if( iTrack.NCells()<3 ) continue;
//cout<<" fit track "<<itr<<", NCells="<<iTrack.NCells()<<endl;
- ID2Point(iTrack.PointID()[0]).Param().CosPhi() = -1;//[2] = -TMath::Pi();
- ID2Point(iTrack.PointID()[1]).Param().CosPhi() = 1;//[2] = 0;
+ ID2Point(iTrack.PointID()[0]).Param().CosPhi() = -1;
+ ID2Point(iTrack.PointID()[1]).Param().CosPhi() = 1;
FitTrack( iTrack );
//if( iTrack.Param().Chi2() > fParam.TrackChi2Cut()*iTrack.Param().NDF() ){
//iTrack.Alive() = 0;
}
//AliHLTTPCCADisplay::Instance().Init();
- AliHLTTPCCADisplay::Instance().SetCurrentSector( this );
- AliHLTTPCCADisplay::Instance().DrawSector( this );
+ AliHLTTPCCADisplay::Instance().SetCurrentSlice( this );
+ AliHLTTPCCADisplay::Instance().DrawSlice( this );
cout<<"draw hits..."<<endl;
for( Int_t iRow=0; iRow<fParam.NRows(); iRow++ )
for (Int_t i = 0; i<fRows[iRow].NHits(); i++)
AliHLTTPCCADisplay::Instance().DrawHit( iRow, i );
AliHLTTPCCADisplay::Instance().Ask();
- //AliHLTTPCCADisplay::Instance().Clear();
- //AliHLTTPCCADisplay::Instance().DrawSector( this );
+ AliHLTTPCCADisplay::Instance().ClearView();
+ AliHLTTPCCADisplay::Instance().DrawSlice( this );
cout<<"draw cells..."<<endl;
for( Int_t iRow=0; iRow<fParam.NRows(); iRow++ )
for (Int_t i = 0; i<fRows[iRow].NCells(); i++)
}
- AliHLTTPCCADisplay::Instance().Clear();
- AliHLTTPCCADisplay::Instance().DrawSector( this );
+ AliHLTTPCCADisplay::Instance().ClearView();
+ AliHLTTPCCADisplay::Instance().DrawSlice( this );
for( Int_t iRow=0; iRow<fParam.NRows(); iRow++ )
for (Int_t i = 0; i<fRows[iRow].NCells(); i++)
AliHLTTPCCADisplay::Instance().DrawCell( iRow, i );
-
+
cout<<"draw initial tracks"<<endl;
for( Int_t itr=0; itr<fNTracks; itr++ ){
- //if( itr!=58 ) continue;
AliHLTTPCCATrack &iTrack = fTracks[itr];
if( iTrack.NCells()<3 ) continue;
- if( !iTrack.Alive() ) continue;
- AliHLTTPCCADisplay::Instance().DrawTrack1( iTrack );
+ if( !iTrack.Alive() ) continue;
+ AliHLTTPCCADisplay::Instance().DrawTrack( iTrack );
}
//if( fNTracks>0 )
AliHLTTPCCADisplay::Instance().Ask();
TStopwatch timer4;
- Bool_t doMerging=1;
+ Bool_t doMerging=1;//SG!!!
Float_t factor2 = fParam.TrackConnectionFactor()*fParam.TrackConnectionFactor();
Int_t *refEndPoints = new Int_t[fNHitsTotal];
Float_t dy2, dz2;
AliHLTTPCCATrackParam &jPar = jp.Param();
-
+
+ // check direction
+ {
+ if( jPar.GetCosPhi()*iPar.GetCosPhi()>=0 ) continue;
+ }
// check for neighbouring
{
float d = jPar.GetY() - iPar.GetY();
if( d*d>factor2*s2 ){
continue;
}
+ //cout<<"\ndy="<<TMath::Sqrt(d*d)<<", err="<<TMath::Sqrt(factor2*s2)<<endl;
d = jPar.GetZ() - iPar.GetZ();
s2 = jPar.GetErr2Z() + iPar.GetErr2Z();
if( d*d>factor2*s2 ){
if( d>0 ) break;
continue;
}
+ //cout<<"dz="<<TMath::Sqrt(d*d)<<", err="<<TMath::Sqrt(factor2*s2)<<endl;
if( iTrack.NCells()>=3 && jTrack.NCells()>=3 ){
d = jPar.GetSinPhi() + iPar.GetSinPhi(); //! phi sign is different
s2 = jPar.GetErr2SinPhi() + iPar.GetErr2SinPhi();
if( d*d>factor2*s2 ) continue;
+ //cout<<"dphi="<<TMath::Sqrt(d*d)<<", err="<<TMath::Sqrt(factor2*s2)<<endl;
d = jPar.GetKappa() + iPar.GetKappa(); // ! kappa sign iz different
s2 = jPar.GetErr2Kappa() + iPar.GetErr2Kappa();
if( d*d>factor2*s2 ) continue;
+ //cout<<"dk="<<TMath::Sqrt(d*d)<<", err="<<TMath::Sqrt(factor2*s2)<<endl;
d = jPar.GetDzDs() + iPar.GetDzDs(); // ! DzDs sign iz different
s2 = jPar.GetErr2DzDs() + iPar.GetErr2DzDs();
if( d*d>factor2*s2 ) continue;
+ //cout<<"dlam="<<TMath::Sqrt(d*d)<<", err="<<TMath::Sqrt(factor2*s2)<<endl;
}
}
cout<<"merged points..."<<ipID<<"/"<<jpID<<endl;
//AliHLTTPCCADisplay::Instance().ConnectCells( iRow,ic,ID2IRow(jcID),jc,kRed );
AliHLTTPCCADisplay::Instance().ConnectEndPoints( ipID,jpID,1.,2,kRed );
- //AliHLTTPCCADisplay::Instance().DrawEndPoint( ipID,1.,2,kRed );
- //AliHLTTPCCADisplay::Instance().DrawEndPoint( jpID,1.,2,kRed );
+ AliHLTTPCCADisplay::Instance().DrawEndPoint( ipID,1.,2,kRed );
+ AliHLTTPCCADisplay::Instance().DrawEndPoint( jpID,1.,2,kRed );
AliHLTTPCCADisplay::Instance().Ask();
cout<<"merged track"<<endl;
- AliHLTTPCCADisplay::Instance().DrawTrack1(iTrack);
+ AliHLTTPCCADisplay::Instance().DrawTrack(iTrack);
AliHLTTPCCADisplay::Instance().Ask();
#endif
/*
fTimers[4] = timer4.CpuTime();
#ifdef DRAW
- if( fNTracks>0 ) AliHLTTPCCADisplay::Instance().Ask();
+ //if( fNTracks>0 ) AliHLTTPCCADisplay::Instance().Ask();
#endif
-#ifdef DRAW
- AliHLTTPCCADisplay::Instance().Clear();
- AliHLTTPCCADisplay::Instance().DrawSector( this );
+#ifdef DRAWXX
+ AliHLTTPCCADisplay::Instance().ClearView();
+ AliHLTTPCCADisplay::Instance().DrawSlice( this );
for( Int_t iRow=0; iRow<fParam.NRows(); iRow++ )
for (Int_t i = 0; i<fRows[iRow].NCells(); i++)
AliHLTTPCCADisplay::Instance().DrawCell( iRow, i );
for( Int_t itr=0; itr<fNTracks; itr++ ){
AliHLTTPCCATrack &iTrack = fTracks[itr];
- //if( itr!=3 ) continue;
if( iTrack.NCells()<3 ) continue;
if( !iTrack.Alive() ) continue;
- AliHLTTPCCADisplay::Instance().DrawTrack1( iTrack );
+ AliHLTTPCCADisplay::Instance().DrawTrack( iTrack );
+ cout<<"final track "<<itr<<", ncells="<<iTrack.NCells()<<endl;
+ AliHLTTPCCADisplay::Instance().Ask();
}
AliHLTTPCCADisplay::Instance().Ask();
#endif
// write output
+ TStopwatch timer7;
+
//cout<<"write output"<<endl;
+#ifdef DRAW
+ AliHLTTPCCADisplay::Instance().ClearView();
+ AliHLTTPCCADisplay::Instance().DrawSlice( this );
+ for( Int_t iRow=0; iRow<fParam.NRows(); iRow++ )
+ for (Int_t i = 0; i<fRows[iRow].NCells(); i++)
+ AliHLTTPCCADisplay::Instance().DrawCell( iRow, i );
+
+ cout<<"draw out tracks"<<endl;
+#endif
+
fOutTrackHits = new Int_t[fNHitsTotal];
fOutTracks = new AliHLTTPCCAOutTrack[fNTracks];
fNOutTrackHits = 0;
out.StartPoint() = p0.Param();
out.EndPoint() = p2.Param();
}
- AliHLTTPCCATrackParam t = out.StartPoint();
+ AliHLTTPCCATrackParam &t = out.StartPoint();//SG!!!
+ AliHLTTPCCATrackParam t0 = t;
+
+ t.Chi2() = 0;
+ t.NDF() = -5;
+ Bool_t first = 1;
+
Int_t iID = iTrack.FirstCellID();
Int_t fNOutTrackHitsOld = fNOutTrackHits;
for( AliHLTTPCCACell *ic = &ID2Cell(iID); ic->Link()>=0; iID = ic->Link(), ic = &ID2Cell(iID) ){
//cout<<"itr="<<iTr<<", cell ="<<ID2IRow(iID)<<" "<<ID2ICell(iID)<<endl;
AliHLTTPCCARow &row = ID2Row(iID);
+ if( !t0.TransportToX( row.X() ) ) continue;
+ Int_t jHit = -1;
+ Float_t dy, dz, d = 1.e10;
for( Int_t iHit=0; iHit<ic->NHits(); iHit++ ){
AliHLTTPCCAHit &h = row.GetCellHit(*ic,iHit);
- // check for wrong hits
-
+ // check for wrong hits
{
- if( !t.TransportToX( row.X() ) ) continue;
- Float_t dy = t.GetY() - h.Y();
- Float_t dz = t.GetZ() - h.Z();
- if( dy*dy > 3.5*3.5*(t.GetErr2Y() + h.ErrY()*h.ErrY() ) ) continue;
- if( dz*dz > 3.5*3.5*(t.GetErr2Z() + h.ErrZ()*h.ErrZ() ) ) continue;
- /*
- Float_t m[3] = {row.X(), h.Y(), h.Z() };
- Float_t mV[6] = {fParam.ErrX()*fParam.ErrX(), 0, h.ErrY()*h.ErrY(), 0, 0, h.ErrZ()*h.ErrZ() };
- Float_t mG[6];
- t.TransportBz(fParam.Bz(),m);
- t.GetConnectionMatrix(fParam.Bz(),m, mG);
- Float_t chi2 = t.GetChi2( m, mV, mG )/2.;
- if( chi2>10. ) continue;
- */
- }
-
- fOutTrackHits[fNOutTrackHits] = h.ID();
- fNOutTrackHits++;
- if( fNOutTrackHits>fNHitsTotal ){
- //cout<<"fNOutTrackHits>fNHitsTotal"<<endl;
- //exit(0);
- return;
+ Float_t ddy = t0.GetY() - h.Y();
+ Float_t ddz = t0.GetZ() - h.Z();
+ Float_t dd = ddy*ddy+ddz*ddz;
+ if( dd<d ){
+ d = dd;
+ dy = ddy;
+ dz = ddz;
+ jHit = iHit;
+ }
}
- out.NHits()++;
}
+ if( jHit<0 ) continue;
+ AliHLTTPCCAHit &h = row.GetCellHit(*ic,jHit);
+ //if( dy*dy > 3.5*3.5*(/*t0.GetErr2Y() + */h.ErrY()*h.ErrY() ) ) continue;//SG!!!
+ //if( dz*dz > 3.5*3.5*(/*t0.GetErr2Z() + */h.ErrZ()*h.ErrZ() ) ) continue;
+ //if( !t0.Filter2( h.Y(), h.Z(), h.ErrY()*h.ErrY(), h.ErrZ()*h.ErrZ() ) ) continue;
+
+
+ if( !t.TransportToX( row.X() ) ) continue;
+
+ //* Update the track
+
+ if( first ){
+ t.Cov()[ 0] = .5*.5;
+ t.Cov()[ 1] = 0;
+ t.Cov()[ 2] = .5*.5;
+ t.Cov()[ 3] = 0;
+ t.Cov()[ 4] = 0;
+ t.Cov()[ 5] = .2*.2;
+ t.Cov()[ 6] = 0;
+ t.Cov()[ 7] = 0;
+ t.Cov()[ 8] = 0;
+ t.Cov()[ 9] = .2*.2;
+ t.Cov()[10] = 0;
+ t.Cov()[11] = 0;
+ t.Cov()[12] = 0;
+ t.Cov()[13] = 0;
+ t.Cov()[14] = .2*.2;
+ t.Chi2() = 0;
+ t.NDF() = -5;
+ }
+
+ if( t.Filter2( h.Y(), h.Z(), h.ErrY()*h.ErrY(), h.ErrZ()*h.ErrZ() ) ) first = 0;
+ else continue;
+
+ fOutTrackHits[fNOutTrackHits] = h.ID();
+ fNOutTrackHits++;
+ if( fNOutTrackHits>fNHitsTotal ){
+ cout<<"fNOutTrackHits>fNHitsTotal"<<endl;
+ exit(0);//SG!!!
+ return;
+ }
+ out.NHits()++;
}
+ //cout<<fNOutTracks<<": itr = "<<iTr<<", n outhits = "<<out.NHits()<<endl;
if( out.NHits() > 3 ){
fNOutTracks++;
+#ifdef DRAW
+ AliHLTTPCCADisplay::Instance().DrawTrack( iTrack );
+ cout<<"out track "<<(fNOutTracks-1)<<", orig = "<<iTr<<", nhits="<<out.NHits()<<endl;
+ AliHLTTPCCADisplay::Instance().Ask();
+#endif
+
}else {
fNOutTrackHits = fNOutTrackHitsOld;
}
}
//cout<<"end writing"<<endl;
+ timer7.Stop();
+ fTimers[7] = timer7.CpuTime();
+
#ifdef DRAW
AliHLTTPCCADisplay::Instance().Ask();
//AliHLTTPCCADisplay::Instance().DrawMCTracks(fParam.fISec);
#include "AliHLTTPCCAHit.h"
#include "AliHLTTPCCAOutTrack.h"
#include "AliHLTTPCCAParam.h"
+#include "AliHLTTPCCATrackConvertor.h"
#include "AliHLTTPCSpacePointData.h"
#include "AliHLTTPCClusterDataFormat.h"
AliHLTTPCCATrackerComponent::AliHLTTPCCATrackerComponent()
:
fTracker(NULL),
- fBField(0),
- fMinNTrackClusters(30),
+ fSolenoidBz(0),
+ fMinNTrackClusters(0),
+ fCellConnectionAngleXY(45),
+ fCellConnectionAngleXZ(45),
+ fClusterZCut(500.),
fFullTime(0),
fRecoTime(0),
fNEvents(0)
:
AliHLTProcessor(),
fTracker(NULL),
- fBField(0),
+ fSolenoidBz(0),
fMinNTrackClusters(30),
- fFullTime(0),
+ fCellConnectionAngleXY(35),
+ fCellConnectionAngleXZ(35),
+ fClusterZCut(500.),
+ fFullTime(0),
fRecoTime(0),
fNEvents(0)
{
// Initialize the CA tracker component
//
// arguments could be:
- // bfield - the magnetic field value
+ // solenoidBz - the magnetic field value
+ // minNTrackClusters - required minimum of clusters on the track
//
if ( fTracker ) return EINPROGRESS;
int i = 0;
char* cpErr;
while ( i < argc ){
- if ( !strcmp( argv[i], "bfield" ) ){
+ if ( !strcmp( argv[i], "solenoidBz" ) || !strcmp( argv[i], "-solenoidBz" ) ){
if ( i+1 >= argc )
{
- Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing B-field", "Missing B-field specifier." );
+ Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing solenoidBz", "Missing solenoidBz specifier." );
return ENOTSUP;
}
- fBField = strtod( argv[i+1], &cpErr );
+ fSolenoidBz = strtod( argv[i+1], &cpErr );
if ( *cpErr )
{
- Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing multiplicity", "Cannot convert B-field specifier '%s'.", argv[i+1] );
+ Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing multiplicity", "Cannot convert solenoidBz specifier '%s'.", argv[i+1] );
return EINVAL;
}
Logging( kHLTLogInfo, "HLT::TPCCATracker::DoInit", "Reading command line",
- "Magnetic field value is set to %f kG", fBField );
+ "Magnetic field value is set to %f kG", fSolenoidBz );
i += 2;
continue;
}
- if ( !strcmp( argv[i], "MinNTrackClusters" ) ){
+ if ( !strcmp( argv[i], "minNTrackClusters" ) || !strcmp( argv[i], "-minNTrackClusters" ) ){
if ( i+1 >= argc )
{
- Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing MinNTrackClusters", "Missing MinNTrackClusters specifier." );
+ Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing minNTrackClusters", "Missing minNTrackClusters specifier." );
return ENOTSUP;
}
fMinNTrackClusters = (Int_t ) strtod( argv[i+1], &cpErr );
if ( *cpErr )
{
- Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing multiplicity", "Cannot convert MinNTrackClusters '%s'.", argv[i+1] );
+ Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing multiplicity", "Cannot convert minNTrackClusters '%s'.", argv[i+1] );
return EINVAL;
}
Logging( kHLTLogInfo, "HLT::TPCCATracker::DoInit", "Reading command line",
- "MinNTrackClusters is set to %i ", fMinNTrackClusters );
+ "minNTrackClusters is set to %i ", fMinNTrackClusters );
i += 2;
continue;
}
-
+
+ if ( !strcmp( argv[i], "cellConnectionAngleXY" ) || !strcmp( argv[i], "-cellConnectionAngleXY" ) ){
+ if ( i+1 >= argc )
+ {
+ Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing cellConnectionAngleXY", "Missing cellConnectionAngleXY specifier." );
+ return ENOTSUP;
+ }
+ fCellConnectionAngleXY = strtod( argv[i+1], &cpErr );
+ if ( *cpErr )
+ {
+ Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing multiplicity", "Cannot convert cellConnectionAngleXY '%s'.", argv[i+1] );
+ return EINVAL;
+ }
+
+ Logging( kHLTLogInfo, "HLT::TPCCATracker::DoInit", "Reading command line",
+ "cellConnectionAngleXY is set to %f ", fCellConnectionAngleXY );
+
+ i += 2;
+ continue;
+ }
+ if ( !strcmp( argv[i], "cellConnectionAngleXZ" ) || !strcmp( argv[i], "-cellConnectionAngleXZ" ) ){
+ if ( i+1 >= argc )
+ {
+ Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing cellConnectionAngleXZ", "Missing cellConnectionAngleXZ specifier." );
+ return ENOTSUP;
+ }
+ fCellConnectionAngleXZ = strtod( argv[i+1], &cpErr );
+ if ( *cpErr )
+ {
+ Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing multiplicity", "Cannot convert cellConnectionAngleXZ '%s'.", argv[i+1] );
+ return EINVAL;
+ }
+
+ Logging( kHLTLogInfo, "HLT::TPCCATracker::DoInit", "Reading command line",
+ "cellConnectionAngleXZ is set to %f ", fCellConnectionAngleXZ );
+
+ i += 2;
+ continue;
+ }
+ if ( !strcmp( argv[i], "clusterZCut" ) || !strcmp( argv[i], "-clusterZCut" ) ){
+ if ( i+1 >= argc )
+ {
+ Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing clusterZCut", "Missing clusterZCut specifier." );
+ return ENOTSUP;
+ }
+ fClusterZCut = TMath::Abs(strtod( argv[i+1], &cpErr ));
+ if ( *cpErr )
+ {
+ Logging( kHLTLogError, "HLT::TPCCATracker::DoInit", "Missing multiplicity", "Cannot convert clusterZCut '%s'.", argv[i+1] );
+ return EINVAL;
+ }
+
+ Logging( kHLTLogInfo, "HLT::TPCCATracker::DoInit", "Reading command line",
+ "clusterZCut is set to %f ", fClusterZCut );
+
+ i += 2;
+ continue;
+ }
+
Logging(kHLTLogError, "HLT::TPCCATracker::DoInit", "Unknown Option", "Unknown option '%s'", argv[i] );
return EINVAL;
}
return 0;
}
+Bool_t AliHLTTPCCATrackerComponent::CompareClusters(AliHLTTPCSpacePointData *a, AliHLTTPCSpacePointData *b)
+{
+ //* Comparison function for sorting clusters
+ if( a->fPadRow<b->fPadRow ) return 1;
+ if( a->fPadRow>b->fPadRow ) return 0;
+ return (a->fZ < b->fZ);
+}
+
int AliHLTTPCCATrackerComponent::DoEvent
(
const AliHLTComponentEventData& evtData,
AliHLTUInt32_t MaxBufferSize = size;
size = 0; // output size
+ if(GetFirstInputBlock( kAliHLTDataTypeSOR ) || GetFirstInputBlock( kAliHLTDataTypeEOR )){
+ return 0;
+ }
+
TStopwatch timer;
// Event reconstruction in one TPC slice with CA Tracker
-
+
//Logging( kHLTLogWarning, "HLT::TPCCATracker::DoEvent", "DoEvent", "CA::DoEvent()" );
if ( evtData.fBlockCnt<=0 )
{
// Initialize the tracker
- Double_t Bz = fBField;
+ Float_t Bz = fSolenoidBz;
{
if( !fTracker ) fTracker = new AliHLTTPCCATracker;
Int_t iSec = slice;
- Double_t inRmin = 83.65;
- // Double_t inRmax = 133.3;
- // Double_t outRmin = 133.5;
- Double_t outRmax = 247.7;
- Double_t plusZmin = 0.0529937;
- Double_t plusZmax = 249.778;
- Double_t minusZmin = -249.645;
- Double_t minusZmax = -0.0799937;
- Double_t dalpha = 0.349066;
- Double_t alpha = 0.174533 + dalpha*iSec;
+ Float_t inRmin = 83.65;
+ // Float_t inRmax = 133.3;
+ // Float_t outRmin = 133.5;
+ Float_t outRmax = 247.7;
+ Float_t plusZmin = 0.0529937;
+ Float_t plusZmax = 249.778;
+ Float_t minusZmin = -249.645;
+ Float_t minusZmax = -0.0799937;
+ Float_t dalpha = 0.349066;
+ Float_t alpha = 0.174533 + dalpha*iSec;
Bool_t zPlus = (iSec<18 );
- Double_t zMin = zPlus ?plusZmin :minusZmin;
- Double_t zMax = zPlus ?plusZmax :minusZmax;
+ Float_t zMin = zPlus ?plusZmin :minusZmin;
+ Float_t zMax = zPlus ?plusZmax :minusZmax;
//TPCZmin = -249.645, ZMax = 249.778
- // Double_t rMin = inRmin;
- // Double_t rMax = outRmax;
+ // Float_t rMin = inRmin;
+ // Float_t rMax = outRmax;
Int_t NRows = AliHLTTPCTransform::GetNRows();
- Double_t padPitch = 0.4;
- Double_t sigmaZ = 0.228808;
+ Float_t padPitch = 0.4;
+ Float_t sigmaZ = 0.228808;
- Double_t *rowX = new Double_t [NRows];
+ Float_t *rowX = new Float_t [NRows];
for( Int_t irow=0; irow<NRows; irow++){
rowX[irow] = AliHLTTPCTransform::Row2X( irow );
}
inRmin, outRmax, zMin, zMax, padPitch, sigmaZ, Bz );
param.YErrorCorrection() = 1;
param.ZErrorCorrection() = 2;
-
+ param.CellConnectionAngleXY() = fCellConnectionAngleXY/180.*TMath::Pi();
+ param.CellConnectionAngleXZ() = fCellConnectionAngleXZ/180.*TMath::Pi();
+ param.Update();
fTracker->Initialize( param );
delete[] rowX;
}
fTracker->StartEvent();
- AliHLTTPCCAHit *vHits = new AliHLTTPCCAHit [nHitsTotal]; // CA hit array
- Double_t *vHitStoreX = new Double_t [nHitsTotal]; // hit X coordinates
- Int_t *vHitStoreID = new Int_t [nHitsTotal]; // hit ID's
- Int_t *vHitRowID = new Int_t [nHitsTotal]; // hit ID's
-
- Int_t nHits = 0;
-
Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Reading hits",
"Total %d hits to read for slice %d", nHitsTotal, slice );
+
+ AliHLTTPCSpacePointData** vOrigClusters = new AliHLTTPCSpacePointData* [ nHitsTotal];
+
Int_t nClusters=0;
for( std::vector<unsigned long>::iterator pIter = patchIndices.begin(); pIter!=patchIndices.end(); pIter++ ){
Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Reading hits",
"Reading %d hits for slice %d - patch %d", inPtrSP->fSpacePointCnt, slice, patch );
- // Read patch hits, row by row
+ for (UInt_t i=0; i<inPtrSP->fSpacePointCnt; i++ ){
+ vOrigClusters[nClusters++] = &(inPtrSP->fSpacePoints[i]);
+ }
+ }
+ // sort clusters since they are not sorted fore some reason
+
+ sort( vOrigClusters, vOrigClusters+nClusters, CompareClusters );
+
+ AliHLTTPCCAHit *vHits = new AliHLTTPCCAHit [nHitsTotal]; // CA hit array
+ Float_t *vHitStoreX = new Float_t [nHitsTotal]; // hit X coordinates
+ Float_t *vHitStoreY = new Float_t [nHitsTotal]; // hit Y coordinates
+ Float_t *vHitStoreZ = new Float_t [nHitsTotal]; // hit Z coordinates
+ Int_t *vHitStoreID = new Int_t [nHitsTotal]; // hit ID's
+ Int_t *vHitRowID = new Int_t [nHitsTotal]; // hit ID's
+
+ Int_t nHits = 0;
+
+ {
Int_t oldRow = -1;
Int_t nRowHits = 0;
Int_t firstRowHit = 0;
- for (UInt_t i=0; i<inPtrSP->fSpacePointCnt; i++ ){
- AliHLTTPCSpacePointData* pSP = &(inPtrSP->fSpacePoints[i]);
+ for (Int_t i=0; i<nClusters; i++ ){
+ AliHLTTPCSpacePointData* pSP = vOrigClusters[i];
if( pSP->fPadRow != oldRow ){
- if( oldRow>=0 ) fTracker->ReadHitRow( oldRow, vHits+firstRowHit, nRowHits );
+ if( oldRow>=0 ){
+ if( fTracker->Rows()[oldRow].NHits()!=0 ) HLTError("CA: clusters from row %d are readed twice",oldRow);
+ fTracker->ReadHitRow( oldRow, vHits+firstRowHit, nRowHits );
+ }
oldRow = pSP->fPadRow;
firstRowHit = nHits;
nRowHits = 0;
}
AliHLTTPCCAHit &h = vHits[nHits];
- if( TMath::Abs(pSP->fX- fTracker->Rows()[pSP->fPadRow].X() )>1.e-4 ) cout<<"row "<<(Int_t)pSP->fPadRow<<" "<<fTracker->Rows()[pSP->fPadRow].X()-pSP->fX <<endl;
+ //if( TMath::Abs(pSP->fX- fTracker->Rows()[pSP->fPadRow].X() )>1.e-4 ) cout<<"row "<<(Int_t)pSP->fPadRow<<" "<<fTracker->Rows()[pSP->fPadRow].X()-pSP->fX <<endl;
+ if( TMath::Abs(pSP->fX- fTracker->Rows()[pSP->fPadRow].X() )>1.e-4 ) HLTError( "row %d, %f",(Int_t)pSP->fPadRow, fTracker->Rows()[pSP->fPadRow].X()-pSP->fX );
h.Y() = pSP->fY;
h.Z() = pSP->fZ;
- if( TMath::Abs(h.Z())>230.) continue;
+ if( TMath::Abs(h.Z())>fClusterZCut) continue;
h.ErrY() = TMath::Sqrt(TMath::Abs(pSP->fSigmaY2));
h.ErrZ() = TMath::Sqrt(TMath::Abs(pSP->fSigmaZ2));
if( h.ErrY()<.1 ) h.ErrY() = .1;
if( h.ErrY()>1. ) h.ErrY() = 1.;
if( h.ErrZ()>1. ) h.ErrZ() = 1.;
h.ID() = nHits;
- vHitStoreX[nHits] = pSP->fX;
+ vHitStoreX[nHits] = pSP->fX;
+ vHitStoreY[nHits] = h.Y();
+ vHitStoreZ[nHits] = h.Z();
vHitStoreID[nHits] = pSP->fID;
vHitRowID[nHits] = pSP->fPadRow;
nHits++;
nRowHits++;
- nClusters++;
}
- if( oldRow>=0 ) fTracker->ReadHitRow( oldRow, vHits+firstRowHit, nRowHits );
+ if( oldRow>=0 ){
+ if( fTracker->Rows()[oldRow].NHits()!=0 ) HLTError("CA: clusters from row %d are readed twice",oldRow);
+ fTracker->ReadHitRow( oldRow, vHits+firstRowHit, nRowHits );
+ }
}
+ if( vOrigClusters ) delete[] vOrigClusters;
// reconstruct the event
}
// convert CA track parameters to HLT Track Segment
-
+
+ Int_t iFirstRow = 1000;
+ Int_t iLastRow = -1;
Int_t iFirstHit = fTracker->OutTrackHits()[t.FirstHitRef()];
- Int_t iLastHit = fTracker->OutTrackHits()[t.FirstHitRef()+t.NHits()-1];
-
+ Int_t iLastHit = iFirstHit;
+ for( Int_t ih=0; ih<t.NHits(); ih++ ){
+ Int_t hitID = fTracker->OutTrackHits()[t.FirstHitRef() + ih ];
+ Int_t iRow = vHitRowID[hitID];
+ if( iRow<iFirstRow ){ iFirstRow = iRow; iFirstHit = hitID; }
+ if( iRow>iLastRow ){ iLastRow = iRow; iLastHit = hitID; }
+ }
+
AliHLTTPCCATrackParam par = t.StartPoint();
par.TransportToX( vHitStoreX[iFirstHit] );
AliExternalTrackParam tp;
- par.GetExtParam( tp, 0, fBField );
+ AliHLTTPCCATrackConvertor::GetExtParam( par, tp, 0, fSolenoidBz );
currOutTracklet->fX = tp.GetX();
currOutTracklet->fY = tp.GetY();
currOutTracklet->fZ = tp.GetZ();
currOutTracklet->fCharge = (Int_t ) tp.GetSign();
currOutTracklet->fPt = TMath::Abs(tp.GetSignedPt());
- Double_t snp = tp.GetSnp() ;
+ Float_t snp = tp.GetSnp() ;
if( snp>.999 ) snp=.999;
if( snp<-.999 ) snp=-.999;
currOutTracklet->fPsi = TMath::ASin( snp );
currOutTracklet->fTgl = tp.GetTgl();
- Double_t h = -currOutTracklet->fPt*currOutTracklet->fPt;
+
+ currOutTracklet->fY0err = tp.GetSigmaY2();
+ currOutTracklet->fZ0err = tp.GetSigmaZ2();
+ Float_t h = -currOutTracklet->fPt*currOutTracklet->fPt;
currOutTracklet->fPterr = h*h*tp.GetSigma1Pt2();
h = 1./TMath::Sqrt(1-snp*snp);
currOutTracklet->fPsierr = h*h*tp.GetSigmaSnp2();
currOutTracklet->fTglerr = tp.GetSigmaTgl2();
-
- par.TransportToX( vHitStoreX[iLastHit] );
- currOutTracklet->fLastX = par.GetX();
- currOutTracklet->fLastY = par.GetY();
- currOutTracklet->fLastZ = par.GetZ();
-
+
+ if( par.TransportToX( vHitStoreX[iLastHit] ) ){
+ currOutTracklet->fLastX = par.GetX();
+ currOutTracklet->fLastY = par.GetY();
+ currOutTracklet->fLastZ = par.GetZ();
+ } else {
+ currOutTracklet->fLastX = vHitStoreX[iLastHit];
+ currOutTracklet->fLastY = vHitStoreY[iLastHit];
+ currOutTracklet->fLastZ = vHitStoreZ[iLastHit];
+ }
+ //if( currOutTracklet->fLastX<10. ) {
+ //HLTError("CA last point: hitxyz=%f,%f,%f, track=%f,%f,%f, tracklet=%f,%f,%f, nhits=%d",vHitStoreX[iLastHit],vHitStoreY[iLastHit],vHitStoreZ[iLastHit],
+ //par.GetX(), par.GetY(),par.GetZ(),currOutTracklet->fLastX,currOutTracklet->fLastY ,currOutTracklet->fLastZ, t.NHits());
+ //}
#ifdef INCLUDE_TPC_HOUGH
#ifdef ROWHOUGHPARAMS
currOutTracklet->fTrackID = 0;
outPtr->fTrackletCnt++;
}
- delete[] vHits;
- delete[] vHitStoreX;
- delete[] vHitStoreID;
- delete[] vHitRowID;
+ if( vHits ) delete[] vHits;
+ if( vHitStoreX ) delete[] vHitStoreX;
+ if( vHitStoreY ) delete[] vHitStoreY;
+ if( vHitStoreZ ) delete[] vHitStoreZ;
+ if( vHitStoreID ) delete[] vHitStoreID;
+ if( vHitRowID ) delete[] vHitRowID;
AliHLTComponentBlockData bd;
FillBlockData( bd );
timer.Stop();
- fFullTime+= timer.CpuTime();
- fRecoTime+= timerReco.CpuTime();
+ fFullTime+= timer.RealTime();
+ fRecoTime+= timerReco.RealTime();
fNEvents++;
// Set log level to "Warning" for on-line system monitoring
-
+ Int_t hz = (Int_t) (fFullTime>1.e-10 ?fNEvents/fFullTime :100000);
+ Int_t hz1 = (Int_t) (fRecoTime>1.e-10 ?fNEvents/fRecoTime :100000);
Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Tracks",
- "CATracker slice %d: output %d tracks; input %d clusters, patches %d..%d, rows %d..%d; reco time %d/%d us",
- slice, ntracks, nClusters, minPatch, maxPatch, row[0], row[1], (Int_t) (fFullTime/fNEvents*1.e6), (Int_t) (fRecoTime/fNEvents*1.e6) );
+ "CATracker slice %d: output %d tracks; input %d clusters, patches %d..%d, rows %d..%d; reco time %d/%d Hz",
+ slice, ntracks, nClusters, minPatch, maxPatch, row[0], row[1], hz, hz1 );
return ret;
#include "AliHLTProcessor.h"
class AliHLTTPCCATracker;
+class AliHLTTPCSpacePointData;
/**
* @class AliHLTTPCCATrackerComponent
AliHLTTPCCATracker* fTracker; //! transient
/** magnetic field */
- Double_t fBField; // see above
+ Double_t fSolenoidBz; // see above
Int_t fMinNTrackClusters; //* required min number of clusters on the track
+ Double_t fCellConnectionAngleXY; //* max phi angle between connected cells (deg)
+ Double_t fCellConnectionAngleXZ; //* max psi angle between connected cells (deg)
+ Double_t fClusterZCut; //* cut on cluster Z position (for noise rejection at the age of TPC)
+ Double_t fFullTime; //* total time for DoEvent() [s]
+ Double_t fRecoTime; //* total reconstruction time [s]
+ Long_t fNEvents; //* number of reconstructed events
- Double_t fFullTime; //! total time for DoEvent() [s]
- Double_t fRecoTime; //! total reconstruction time [s]
- Long_t fNEvents; //! number of reconstructed events
+ static Bool_t CompareClusters(AliHLTTPCSpacePointData *a, AliHLTTPCSpacePointData *b);
ClassDef(AliHLTTPCCATrackerComponent, 0);
#include "AliHLTTPCCAOutTrack.h"
#include "AliHLTTPCCAPerformance.h"
#include "AliHLTTPCCAParam.h"
+#include "AliHLTTPCCATrackConvertor.h"
#include "TMath.h"
#include "AliTPCLoader.h"
#include "AliTPCtrack.h"
#include "AliESDtrack.h"
#include "AliESDEvent.h"
+#include "AliTrackReference.h"
+#include <fstream.h>
ClassImp(AliTPCtrackerCA)
AliTPCtrackerCA::AliTPCtrackerCA()
- :AliTracker(),fParam(0), fClusters(0), fNClusters(0), fHLTTracker(0),fHLTPerformance(0),fDoHLTPerformance(0)
+ :AliTracker(),fParam(0), fClusters(0), fNClusters(0), fHLTTracker(0),fHLTPerformance(0),fDoHLTPerformance(0),fDoHLTPerformanceClusters(0),fStatNEvents(0)
{
//* default constructor
}
AliTPCtrackerCA::AliTPCtrackerCA(const AliTPCtrackerCA &):
- AliTracker(),fParam(0), fClusters(0), fNClusters(0), fHLTTracker(0),fHLTPerformance(0),fDoHLTPerformance(0)
+ AliTracker(),fParam(0), fClusters(0), fNClusters(0), fHLTTracker(0),fHLTPerformance(0),fDoHLTPerformance(0),fDoHLTPerformanceClusters(0),fStatNEvents(0)
{
//* dummy
}
}
AliTPCtrackerCA::AliTPCtrackerCA(const AliTPCParam *par):
- AliTracker(),fParam(par), fClusters(0), fNClusters(0), fHLTTracker(0), fHLTPerformance(0),fDoHLTPerformance(0)
+ AliTracker(),fParam(par), fClusters(0), fNClusters(0), fHLTTracker(0), fHLTPerformance(0),fDoHLTPerformance(0),fDoHLTPerformanceClusters(0),fStatNEvents(0)
{
//* constructor
- DoHLTPerformance() = 0;
+ DoHLTPerformance() = 1;
+ DoHLTPerformanceClusters() = 1;
fHLTTracker = new AliHLTTPCCAGBTracker;
fHLTTracker->SetNSlices( fParam->GetNSector()/2 );
for( int iSlice=0; iSlice<fHLTTracker->NSlices(); iSlice++ ){
- Double_t bz = AliTracker::GetBz();
-
- Double_t inRmin = fParam->GetInnerRadiusLow();
- //Double_t inRmax = fParam->GetInnerRadiusUp();
- //Double_t outRmin = fParam->GetOuterRadiusLow();
- Double_t outRmax = fParam->GetOuterRadiusUp();
- Double_t plusZmin = 0.0529937;
- Double_t plusZmax = 249.778;
- Double_t minusZmin = -249.645;
- Double_t minusZmax = -0.0799937;
- Double_t dalpha = 0.349066;
- Double_t alpha = 0.174533 + dalpha*iSlice;
+ Float_t bz = AliTracker::GetBz();
+
+ Float_t inRmin = fParam->GetInnerRadiusLow();
+ //Float_t inRmax = fParam->GetInnerRadiusUp();
+ //Float_t outRmin = fParam->GetOuterRadiusLow();
+ Float_t outRmax = fParam->GetOuterRadiusUp();
+ Float_t plusZmin = 0.0529937;
+ Float_t plusZmax = 249.778;
+ Float_t minusZmin = -249.645;
+ Float_t minusZmax = -0.0799937;
+ Float_t dalpha = 0.349066;
+ Float_t alpha = 0.174533 + dalpha*iSlice;
Bool_t zPlus = (iSlice<18 );
- Double_t zMin = zPlus ?plusZmin :minusZmin;
- Double_t zMax = zPlus ?plusZmax :minusZmax;
+ Float_t zMin = zPlus ?plusZmin :minusZmin;
+ Float_t zMax = zPlus ?plusZmax :minusZmax;
//TPCZmin = -249.645, ZMax = 249.778
- //Double_t rMin = inRmin;
- //Double_t rMax = outRmax;
+ //Float_t rMin = inRmin;
+ //Float_t rMax = outRmax;
- Double_t padPitch = 0.4;
- Double_t sigmaZ = 0.228808;
+ Float_t padPitch = 0.4;
+ Float_t sigmaZ = 0.228808;
Int_t NRows = fParam->GetNRowLow()+fParam->GetNRowUp();
- Double_t rowX[200];
+ Float_t rowX[200];
for( Int_t irow=0; irow<fParam->GetNRowLow(); irow++){
rowX[irow] = fParam->GetPadRowRadiiLow(irow);
}
AliHLTTPCCAParam param;
param.Initialize( iSlice, NRows, rowX, alpha, dalpha,
inRmin, outRmax, zMin, zMax, padPitch, sigmaZ, bz );
- param.YErrorCorrection() = .33;//1;
- param.ZErrorCorrection() = .33;//2;
+ param.YErrorCorrection() = 1.;//.33;
+ param.ZErrorCorrection() = 1.;//.33;
param.MaxTrackMatchDRow() = 5;
- param.TrackConnectionFactor() = 5.;
+ param.TrackConnectionFactor() = 3.5;
fHLTTracker->Slices()[iSlice].Initialize( param );
}
}
if( !fParam ) return 1;
// load mc tracks
- if( fDoHLTPerformance ){
- if( !gAlice ) return 0;
+ while( fDoHLTPerformance ){
+ if( !gAlice ) break;
AliRunLoader *rl = gAlice->GetRunLoader();
- if( !rl ) return 0;
+ if( !rl ) break;
rl->LoadKinematics();
AliStack *stack = rl->Stack();
- if( !stack ) return 0 ;
+ if( !stack ) break;
fHLTPerformance->SetNMCTracks( stack->GetNtrack() );
TParticle *part = stack->Particle(itr);
fHLTPerformance->ReadMCTrack( itr, part );
}
- }
+
+ { // check for MC tracks at the TPC entrance
+ Bool_t *isTPC = 0;
+ isTPC = new Bool_t [stack->GetNtrack()];
+ for( Int_t i=0; i<stack->GetNtrack(); i++ ) isTPC[i] = 0;
+ rl->LoadTrackRefs();
+ TTree *TR = rl->TreeTR();
+ if( !TR ) break;
+ TBranch *branch=TR->GetBranch("TrackReferences");
+ if (!branch ) break;
+ TClonesArray tpcdummy("AliTrackReference",1000), *tpcRefs=&tpcdummy;
+ branch->SetAddress(&tpcRefs);
+ Int_t nr=(Int_t)TR->GetEntries();
+ for (Int_t r=0; r<nr; r++) {
+ TR->GetEvent(r);
+ Int_t nref = tpcRefs->GetEntriesFast();
+ if (!nref) continue;
+ AliTrackReference *tpcRef= 0x0;
+ for (Int_t iref=0; iref<nref; ++iref) {
+ tpcRef = (AliTrackReference*)tpcRefs->UncheckedAt(iref);
+ if (tpcRef->DetectorId() == AliTrackReference::kTPC) break;
+ tpcRef = 0x0;
+ }
+ if (!tpcRef) continue;
+
+ if( isTPC[tpcRef->Label()] ) continue;
+
+ fHLTPerformance->ReadMCTPCTrack(tpcRef->Label(),
+ tpcRef->X(),tpcRef->Y(),tpcRef->Z(),
+ tpcRef->Px(),tpcRef->Py(),tpcRef->Pz() );
+ isTPC[tpcRef->Label()] = 1;
+ tpcRefs->Clear();
+ }
+ if( isTPC ) delete[] isTPC;
+ }
+
+ while( fDoHLTPerformanceClusters ){
+ AliTPCLoader *tpcl = (AliTPCLoader*) rl->GetDetectorLoader("TPC");
+ if( !tpcl ) break;
+ if( tpcl->TreeH() == 0x0 ){
+ if( tpcl->LoadHits() ) break;
+ }
+ if( tpcl->TreeH() == 0x0 ) break;
+
+ AliTPC *tpc = (AliTPC*) gAlice->GetDetector("TPC");
+ Int_t nEnt=(Int_t)tpcl->TreeH()->GetEntries();
+ Int_t nPoints = 0;
+ for (Int_t iEnt=0; iEnt<nEnt; iEnt++) {
+ tpc->ResetHits();
+ tpcl->TreeH()->GetEvent(iEnt);
+ AliTPChit *phit = (AliTPChit*)tpc->FirstHit(-1);
+ for ( ; phit; phit=(AliTPChit*)tpc->NextHit() ) nPoints++;
+ }
+ fHLTPerformance->SetNMCPoints( nPoints );
+
+ for (Int_t iEnt=0; iEnt<nEnt; iEnt++) {
+ tpc->ResetHits();
+ tpcl->TreeH()->GetEvent(iEnt);
+ AliTPChit *phit = (AliTPChit*)tpc->FirstHit(-1);
+ for ( ; phit; phit=(AliTPChit*)tpc->NextHit() ){
+ fHLTPerformance->ReadMCPoint( phit->GetTrack(),phit->X(), phit->Y(),phit->Z(),phit->Time(), phit->fSector%36);
+ }
+ }
+ break;
+ }
+ break;
+ }
TBranch * br = tree->GetBranch("Segment");
if( !br ) return 1;
Int_t sec,row;
fParam->AdjustSectorRow(clrow->GetID(),sec,row);
int NClu = clrow->GetArray()->GetEntriesFast();
- Double_t x = fParam->GetPadRowRadii(sec,row);
+ Float_t x = fParam->GetPadRowRadii(sec,row);
for (Int_t icl=0; icl<NClu; icl++){
Int_t lab0 = -1;
Int_t lab1 = -1;
cluster->SetY(posC[1]);
cluster->SetZ(posC[2]);
- Double_t y = cluster->GetY();
- Double_t z = cluster->GetZ();
+ Float_t y = cluster->GetY();
+ Float_t z = cluster->GetZ();
if( sec>=36 ){
sec = sec - 36;
Int_t index = ind++;
fClusters[index] = *cluster;
fHLTTracker->ReadHit( x, y, z,
- TMath::Sqrt(cluster->GetSigmaY2()), TMath::Sqrt(cluster->GetSigmaZ2()),
- index, sec, row );
+ TMath::Sqrt(cluster->GetSigmaY2()), TMath::Sqrt(cluster->GetSigmaZ2()),
+ cluster->GetMax(), index, sec, row );
if( fDoHLTPerformance ) fHLTPerformance->ReadHitLabel(index, lab0, lab1, lab2 );
}
}
fHLTTracker->FindTracks();
if( fDoHLTPerformance ) fHLTPerformance->Performance();
+ if( 0 ) {// Write Event
+ if( fStatNEvents == 0 ){
+ fstream geo;
+ geo.open("CAEvents/settings.dat", ios::out);
+ if( geo.is_open() ){
+ fHLTTracker->WriteSettings(geo);
+ }
+ geo.close();
+ }
+
+ fstream event;
+ char name[255];
+ sprintf( name,"CAEvents/%i.event.dat",fStatNEvents );
+ event.open(name, ios::out);
+ if( event.is_open() ){
+ fHLTTracker->WriteEvent(event);
+ fstream tracks;
+ sprintf( name,"CAEvents/%i.tracks.dat",fStatNEvents );
+ tracks.open(name, ios::out);
+ fHLTTracker->WriteTracks(tracks);
+ }
+ event.close();
+ if( fDoHLTPerformance ){
+ fstream mcevent, mcpoints;
+ char mcname[255];
+ sprintf( mcname,"CAEvents/%i.mcevent.dat",fStatNEvents );
+ mcevent.open(mcname, ios::out);
+ if( mcevent.is_open() ){
+ fHLTPerformance->WriteMCEvent(mcevent);
+ }
+ if(fDoHLTPerformanceClusters ){
+ sprintf( mcname,"CAEvents/%i.mcpoints.dat",fStatNEvents );
+ mcpoints.open(mcname, ios::out);
+ if( mcpoints.is_open() ){
+ fHLTPerformance->WriteMCPoints(mcpoints);
+ }
+ mcpoints.close();
+ }
+ mcevent.close();
+ }
+ }
+ fStatNEvents++;
+
if( event ){
+ Float_t bz = fHLTTracker->Slices()[0].Param().Bz();
+
for( Int_t itr=0; itr<fHLTTracker->NTracks(); itr++ ){
AliTPCtrack tTPC;
AliHLTTPCCAGBTrack &tCA = fHLTTracker->Tracks()[itr];
AliHLTTPCCATrackParam &par = tCA.Param();
-
- par.GetExtParam( tTPC, tCA.Alpha(), fHLTTracker->Slices()[0].Param().Bz() );
-
+ AliHLTTPCCATrackConvertor::GetExtParam( par, tTPC, tCA.Alpha(), bz );
tTPC.SetMass(0.13957);
+ tTPC.SetdEdx( tCA.DeDx() );
+ if( TMath::Abs(tTPC.GetSigned1Pt())>1./0.1 ) continue;
int nhits = tCA.NHits();
if( nhits>kMaxRow ) nhits = kMaxRow;
- tTPC.SetNumberOfClusters(nhits);
+ tTPC.SetNumberOfClusters(nhits);
for( Int_t ih=0; ih<nhits; ih++ ){
Int_t index = fHLTTracker->TrackHits()[tCA.FirstHitRef()+ih];
Int_t ext_index = fHLTTracker->Hits()[index].ID();
tTPC.SetClusterIndex(ih, ext_index);
}
- CookLabel(&tTPC,0.1);
+ CookLabel(&tTPC,0.1);
{
- Double_t xTPC=83.65;
- if (tTPC.AliExternalTrackParam::PropagateTo(xTPC,5)) {
+ Double_t xTPC=fParam->GetInnerRadiusLow();
+ Double_t dAlpha = fParam->GetInnerAngle()/180.*TMath::Pi();
+ if (tTPC.AliExternalTrackParam::PropagateTo(xTPC,bz)) {
Double_t y=tTPC.GetY();
- Double_t ymax=xTPC*TMath::Tan(1.74532920122146606e-01);
+ Double_t ymax=xTPC*TMath::Tan(dAlpha/2.);
if (y > ymax) {
- if (tTPC.Rotate(2*1.74532920122146606e-01)) tTPC.AliExternalTrackParam::PropagateTo(xTPC,5);
+ if (tTPC.Rotate(dAlpha)) tTPC.AliExternalTrackParam::PropagateTo(xTPC,bz);
} else if (y <-ymax) {
- if (tTPC.Rotate(-2*1.74532920122146606e-01)) tTPC.AliExternalTrackParam::PropagateTo(xTPC,5);
+ if (tTPC.Rotate(-dAlpha)) tTPC.AliExternalTrackParam::PropagateTo(xTPC,bz);
}
}
}
tESD.UpdateTrackParams( &(tTPC),AliESDtrack::kTPCin);
//tESD.SetStatus( AliESDtrack::kTPCrefit );
//tESD.SetTPCPoints(tTPC.GetPoints());
- //tESD.myTPC = tTPC;
+ //tESD.myTPC = tTPC;
event->AddTrack(&tESD);
}
}
}
-Int_t AliTPCtrackerCA::RefitInward (AliESDEvent *)
+Int_t AliTPCtrackerCA::RefitInward (AliESDEvent *event)
{
- //* not implemented yet
- return 0;
+ //* back propagation of ESD tracks (not fully functional)
+
+ Float_t bz = fHLTTracker->Slices()[0].Param().Bz();
+ Float_t xTPC = fParam->GetInnerRadiusLow();
+ Float_t dAlpha = fParam->GetInnerAngle()/180.*TMath::Pi();
+ Float_t yMax = xTPC*TMath::Tan(dAlpha/2.);
+
+ Int_t nentr=event->GetNumberOfTracks();
+
+ for (Int_t i=0; i<nentr; i++) {
+ AliESDtrack *esd=event->GetTrack(i);
+ ULong_t status=esd->GetStatus();
+ if (!(status&AliESDtrack::kTPCin)) continue;
+ AliHLTTPCCATrackParam t0;
+ AliHLTTPCCATrackConvertor::SetExtParam(t0,*esd, bz );
+ Float_t alpha = esd->GetAlpha();
+ if( t0.TransportToXWithMaterial( xTPC, bz) ){
+ if (t0.GetY() > yMax) {
+ if (t0.Rotate(dAlpha)){
+ alpha+=dAlpha;
+ t0.TransportToXWithMaterial( xTPC, bz);
+ }
+ } else if (t0.GetY() <-yMax) {
+ if (t0.Rotate(-dAlpha)){
+ alpha+=-dAlpha;
+ t0.TransportToXWithMaterial( xTPC, bz);
+ }
+ }
+ }
+ AliTPCtrack tt(*esd);
+ AliHLTTPCCATrackConvertor::GetExtParam(t0,tt,alpha,bz);
+ esd->UpdateTrackParams( &tt,AliESDtrack::kTPCrefit);
+ }
+ return 0;
}
Int_t AliTPCtrackerCA::PropagateBack(AliESDEvent *)
void UnloadClusters(){ return ; }
AliCluster * GetCluster(Int_t index) const;
Bool_t &DoHLTPerformance(){ return fDoHLTPerformance; }
+ Bool_t &DoHLTPerformanceClusters(){ return fDoHLTPerformanceClusters; }
//
protected:
AliHLTTPCCAGBTracker *fHLTTracker; //* pointer to the HLT tracker
AliHLTTPCCAPerformance *fHLTPerformance; //* performance calculations
Bool_t fDoHLTPerformance; //* flag for call AliHLTTPCCAPerformance
+ Bool_t fDoHLTPerformanceClusters; //* flag for call AliHLTTPCCAPerformance with cluster pulls (takes some time to load TPC MC points)
+ Int_t fStatNEvents; //* N of reconstructed events
ClassDef(AliTPCtrackerCA,1)
};
tracking-ca/AliHLTTPCCAHit.h \
tracking-ca/AliHLTTPCCAOutTrack.h \
tracking-ca/AliHLTTPCCATrackParam.h \
+ tracking-ca/AliHLTTPCCATrackConvertor.h \
tracking-ca/AliHLTTPCCAParam.h \
tracking-ca/AliHLTTPCCARow.h \
tracking-ca/AliHLTTPCCAGrid.h \
tracking-ca/AliHLTTPCCAGBHit.h \
tracking-ca/AliHLTTPCCAMCTrack.h \
+ tracking-ca/AliHLTTPCCAMCPoint.h \
tracking-ca/AliHLTTPCCAGBTrack.h \
tracking-ca/AliHLTTPCCAGBTracker.h \
tracking-ca/AliHLTTPCCAPerformance.h \
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff