// Origin: Marian Ivanov Marian.Ivanov@cern.ch
//
// AliTPC parallel tracker
+//
+// The track fitting is based on Kalaman filtering approach
+
+// The track finding steps:
+// 1. Seeding - with and without vertex constraint
+// - seeding with vertex constain done at first n^2 proble
+// - seeding without vertex constraint n^3 problem
+// 2. Tracking - follow prolongation road - find cluster - update kalman track
+
+// The seeding and tracking is repeated several times, in different seeding region.
+// This approach enables to find the track which cannot be seeded in some region of TPC
+// This can happen because of low momenta (track do not reach outer radius), or track is currently in the ded region between sectors, or the track is for the moment overlapped with other track (seed quality is poor) ...
+
+// With this approach we reach almost 100 % efficiency also for high occupancy events.
+// (If the seeding efficiency in a region is about 90 % than with logical or of several
+// regions we will reach 100% (in theory - supposing independence)
+
+// Repeating several seeding - tracking procedures some of the tracks can be find
+// several times.
+
+// The procedures to remove multi find tacks are impremented:
+// RemoveUsed2 - fast procedure n problem -
+// Algorithm - Sorting tracks according quality
+// remove tracks with some shared fraction
+// Sharing in respect to all tacks
+// Signing clusters in gold region
+// FindSplitted - slower algorithm n^2
+// Sort the tracks according quality
+// Loop over pair of tracks
+// If overlap with other track bigger than threshold - remove track
+//
+// FindCurling - Finds the pair of tracks which are curling
+// - About 10% of tracks can be find with this procedure
+// The combinatorial background is too big to be used in High
+// multiplicity environment
+// - n^2 problem - Slow procedure - currently it is disabled because of
+// low efficiency
+//
+// The number of splitted tracks can be reduced disabling the sharing of the cluster.
+// tpcRecoParam-> SetClusterSharing(kFALSE);
+// IT IS HIGHLY non recomended to use it in high flux enviroonment
+// Even using this switch some tracks can be found more than once
+// (because of multiple seeding and low quality tracks which will not cross full chamber)
+//
+//
+// The tracker itself can be debugged - the information about tracks can be stored in several // phases of the reconstruction
+// To enable storage of the TPC tracks in the ESD friend track
+// use AliTPCReconstructor::SetStreamLevel(n); where nis bigger 0
+//
+// The debug level - different procedure produce tree for numerical debugging
+// To enable them set AliTPCReconstructor::SetStreamLevel(n); where nis bigger 1
+//
+
+//
+// Adding systematic errors to the covariance:
+//
+// The systematic errors due to the misalignment and miscalibration are added to the covariance matrix
+// of the tracks (not to the clusters as they are dependent):
+// The parameters form AliTPCRecoParam are used AliTPCRecoParam::GetSystematicError
+// The systematic errors are expressed there in RMS - position (cm), angle (rad), curvature (1/cm)
+// The default values are 0.
+//
+// The sytematic errors are added to the covariance matrix in following places:
+//
+// 1. During fisrt itteration - AliTPCtrackerMI::FillESD
+// 2. Second iteration -
+// 2.a ITS->TPC - AliTPCtrackerMI::ReadSeeds
+// 2.b TPC->TRD - AliTPCtrackerMI::PropagateBack
+// 3. Third iteration -
+// 3.a TRD->TPC - AliTPCtrackerMI::ReadSeeds
+// 3.b TPC->ITS - AliTPCtrackerMI::RefitInward
+//
+// There are several places in the code which can be numerically debuged
+// This code is keeped in order to enable code development and to check the calibration implementtion
+//
+// 1. ErrParam stream (Log level 9) - dump information about
+// 1.a) cluster
+// 2.a) cluster error estimate
+// 3.a) cluster shape estimate
+//
+//
//-------------------------------------------------------
#include <TObjArray.h>
#include <TTree.h>
#include "AliLog.h"
-
#include "AliComplexCluster.h"
-#include "AliESD.h"
+#include "AliESDEvent.h"
+#include "AliESDtrack.h"
+#include "AliESDVertex.h"
#include "AliKink.h"
#include "AliV0.h"
#include "AliHelix.h"
#include "AliTPCReconstructor.h"
#include "AliTPCpolyTrack.h"
#include "AliTPCreco.h"
-#include "AliTPCseed.h"
+#include "AliTPCseed.h"
+
+#include "AliTPCtrackerSector.h"
#include "AliTPCtrackerMI.h"
#include "TStopwatch.h"
#include "AliTPCReconstructor.h"
#include "TTreeStream.h"
#include "AliAlignObj.h"
#include "AliTrackPointArray.h"
+#include "TRandom.h"
+#include "AliTPCcalibDB.h"
+#include "AliTPCTransform.h"
+#include "AliTPCClusterParam.h"
//
ClassImp(AliTPCtrackerMI)
+
class AliTPCFastMath {
public:
AliTPCFastMath();
fNtracks(0),
fSeeds(0),
fIteration(0),
- fParam(0),
+ fkParam(0),
fDebugStreamer(0)
{
//
track->SetRow((i&0x00ff0000)>>16);
track->SetSector(sec);
// Int_t index = i&0xFFFF;
- if (sec>=fParam->GetNInnerSector()) track->SetRow(track->GetRow()+fParam->GetNRowLow());
+ if (sec>=fkParam->GetNInnerSector()) track->SetRow(track->GetRow()+fkParam->GetNRowLow());
track->SetClusterIndex2(track->GetRow(), i);
//track->fFirstPoint = row;
//if ( track->fLastPoint<row) track->fLastPoint =row;
Double_t chi2 = track->GetPredictedChi2(track->GetCurrentCluster());
//
- track->SetErrorY2(track->GetErrorY2()*1.3);
- track->SetErrorY2(track->GetErrorY2()+0.01);
- track->SetErrorZ2(track->GetErrorZ2()*1.3);
- track->SetErrorZ2(track->GetErrorZ2()+0.005);
+// track->SetErrorY2(track->GetErrorY2()*1.3);
+// track->SetErrorY2(track->GetErrorY2()+0.01);
+// track->SetErrorZ2(track->GetErrorZ2()*1.3);
+// track->SetErrorZ2(track->GetErrorZ2()+0.005);
//}
if (accept>0) return 0;
if (track->GetNumberOfClusters()%20==0){
-Int_t AliTPCtrackerMI::AcceptCluster(AliTPCseed * seed, AliTPCclusterMI * cluster, Float_t factor,
- Float_t cory, Float_t corz)
+Int_t AliTPCtrackerMI::AcceptCluster(AliTPCseed * seed, AliTPCclusterMI * cluster)
{
//
// decide according desired precision to accept given
// cluster for tracking
- Double_t sy2=ErrY2(seed,cluster)*cory;
- Double_t sz2=ErrZ2(seed,cluster)*corz;
- //sy2=ErrY2(seed,cluster)*cory;
- //sz2=ErrZ2(seed,cluster)*cory;
+ Double_t sy2=ErrY2(seed,cluster);
+ Double_t sz2=ErrZ2(seed,cluster);
Double_t sdistancey2 = sy2+seed->GetSigmaY2();
Double_t sdistancez2 = sz2+seed->GetSigmaZ2();
Double_t rdistance2 = rdistancey2+rdistancez2;
//Int_t accept =0;
+ if (AliTPCReconstructor::StreamLevel()>5 && seed->GetNumberOfClusters()>20) {
+ Float_t rmsy2 = seed->GetCurrentSigmaY2();
+ Float_t rmsz2 = seed->GetCurrentSigmaZ2();
+ Float_t rmsy2p30 = seed->GetCMeanSigmaY2p30();
+ Float_t rmsz2p30 = seed->GetCMeanSigmaZ2p30();
+ Float_t rmsy2p30R = seed->GetCMeanSigmaY2p30R();
+ Float_t rmsz2p30R = seed->GetCMeanSigmaZ2p30R();
+
+ AliExternalTrackParam param(*seed);
+ static TVectorD gcl(3),gtr(3);
+ Float_t gclf[3];
+ param.GetXYZ(gcl.GetMatrixArray());
+ cluster->GetGlobalXYZ(gclf);
+ gcl[0]=gclf[0]; gcl[1]=gclf[1]; gcl[2]=gclf[2];
+
+ if (AliTPCReconstructor::StreamLevel()>0) {
+ (*fDebugStreamer)<<"ErrParam"<<
+ "Cl.="<<cluster<<
+ "T.="<<¶m<<
+ "gcl.="<<&gcl<<
+ "gtr.="<<>r<<
+ "erry2="<<sy2<<
+ "errz2="<<sz2<<
+ "rmsy2="<<rmsy2<<
+ "rmsz2="<<rmsz2<<
+ "rmsy2p30="<<rmsy2p30<<
+ "rmsz2p30="<<rmsz2p30<<
+ "rmsy2p30R="<<rmsy2p30R<<
+ "rmsz2p30R="<<rmsz2p30R<<
+ "\n";
+ }
+ }
+
if (rdistance2>16) return 3;
- if ((rdistancey2>9.*factor || rdistancez2>9.*factor) && cluster->GetType()==0)
+ if ((rdistancey2>9. || rdistancez2>9.) && cluster->GetType()==0)
return 2; //suspisiouce - will be changed
- if ((rdistancey2>6.25*factor || rdistancez2>6.25*factor) && cluster->GetType()>0)
+ if ((rdistancey2>6.25 || rdistancez2>6.25) && cluster->GetType()>0)
// strict cut on overlaped cluster
return 2; //suspisiouce - will be changed
- if ( (rdistancey2>1.*factor || rdistancez2>6.25*factor )
+ if ( (rdistancey2>1. || rdistancez2>6.25 )
&& cluster->GetType()<0){
seed->SetNFoundable(seed->GetNFoundable()-1);
return 2;
fNtracks(0),
fSeeds(0),
fIteration(0),
- fParam(0),
+ fkParam(0),
fDebugStreamer(0)
{
//---------------------------------------------------------------------
// The main TPC tracker constructor
//---------------------------------------------------------------------
- fInnerSec=new AliTPCSector[fkNIS];
- fOuterSec=new AliTPCSector[fkNOS];
+ fInnerSec=new AliTPCtrackerSector[fkNIS];
+ fOuterSec=new AliTPCtrackerSector[fkNOS];
Int_t i;
for (i=0; i<fkNIS; i++) fInnerSec[i].Setup(par,0);
for (i=0; i<fkNOS; i++) fOuterSec[i].Setup(par,1);
- fParam = par;
+ fkParam = par;
Int_t nrowlow = par->GetNRowLow();
Int_t nrowup = par->GetNRowUp();
- for (Int_t i=0;i<nrowlow;i++){
+ for (i=0;i<nrowlow;i++){
fXRow[i] = par->GetPadRowRadiiLow(i);
fPadLength[i]= par->GetPadPitchLength(0,i);
fYMax[i] = fXRow[i]*TMath::Tan(0.5*par->GetInnerAngle());
}
- for (Int_t i=0;i<nrowup;i++){
+ for (i=0;i<nrowup;i++){
fXRow[i+nrowlow] = par->GetPadRowRadiiUp(i);
fPadLength[i+nrowlow] = par->GetPadPitchLength(60,i);
fYMax[i+nrowlow] = fXRow[i+nrowlow]*TMath::Tan(0.5*par->GetOuterAngle());
}
- fDebugStreamer = new TTreeSRedirector("TPCdebug.root");
+ if (AliTPCReconstructor::StreamLevel()>0) {
+ fDebugStreamer = new TTreeSRedirector("TPCdebug.root");
+ }
}
//________________________________________________________________________
AliTPCtrackerMI::AliTPCtrackerMI(const AliTPCtrackerMI &t):
fNtracks(0),
fSeeds(0),
fIteration(0),
- fParam(0),
+ fkParam(0),
fDebugStreamer(0)
{
//------------------------------------
//------------------------------------------------------------------
fOutput=t.fOutput;
}
-AliTPCtrackerMI & AliTPCtrackerMI::operator=(const AliTPCtrackerMI& /*r*/){
+AliTPCtrackerMI & AliTPCtrackerMI::operator=(const AliTPCtrackerMI& /*r*/)
+{
//------------------------------
// dummy
//--------------------------------------------------------------
if (fDebugStreamer) delete fDebugStreamer;
}
-void AliTPCtrackerMI::SetIO()
-{
- //
- fNewIO = kTRUE;
- fInput = AliRunLoader::GetTreeR("TPC", kFALSE,AliConfig::GetDefaultEventFolderName());
-
- fOutput = AliRunLoader::GetTreeT("TPC", kTRUE,AliConfig::GetDefaultEventFolderName());
- if (fOutput){
- AliTPCtrack *iotrack= new AliTPCtrack;
- fOutput->Branch("tracks","AliTPCtrack",&iotrack,32000,100);
- delete iotrack;
- }
-}
-
-
-void AliTPCtrackerMI::SetIO(TTree * input, TTree * output, AliESD * event)
-{
-
- // set input
- fNewIO = kFALSE;
- fInput = 0;
- fOutput = 0;
- fSeedTree = 0;
- fTreeDebug =0;
- fInput = input;
- if (input==0){
- return;
- }
- //set output
- fOutput = output;
- if (output){
- AliTPCtrack *iotrack= new AliTPCtrack;
- // iotrack->fHelixIn = new TClonesArray("AliHelix");
- //iotrack->fHelixOut = new TClonesArray("AliHelix");
- fOutput->Branch("tracks","AliTPCtrack",&iotrack,32000,100);
- delete iotrack;
- }
- if (output && (fDebug&2)){
- //write the full seed information if specified in debug mode
- //
- fSeedTree = new TTree("Seeds","Seeds");
- AliTPCseed * vseed = new AliTPCseed;
- //
- TClonesArray * arrtr = new TClonesArray("AliTPCTrackPoint",160);
- arrtr->ExpandCreateFast(160);
- TClonesArray * arre = new TClonesArray("AliTPCExactPoint",160);
- //
- vseed->SetPoints(arrtr);
- vseed->SetEPoints(arre);
- // vseed->fClusterPoints = arrcl;
- fSeedTree->Branch("seeds","AliTPCseed",&vseed,32000,99);
- delete arrtr;
- delete arre;
- fTreeDebug = new TTree("trackDebug","trackDebug");
- TClonesArray * arrd = new TClonesArray("AliTPCTrackPoint2",0);
- fTreeDebug->Branch("debug",&arrd,32000,99);
- }
-
- //set ESD event
- fEvent = event;
-}
-
-void AliTPCtrackerMI::FillESD(TObjArray* arr)
+void AliTPCtrackerMI::FillESD(const TObjArray* arr)
{
//
//
AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
if (!pt) continue;
pt->UpdatePoints();
- // pt->PropagateTo(fParam->GetInnerRadiusLow());
+ AddCovariance(pt);
+ if (AliTPCReconstructor::StreamLevel()>1) {
+ (*fDebugStreamer)<<"Track0"<<
+ "Tr0.="<<pt<<
+ "\n";
+ }
+ // pt->PropagateTo(fkParam->GetInnerRadiusLow());
if (pt->GetKinkIndex(0)<=0){ //don't propagate daughter tracks
- pt->PropagateTo(fParam->GetInnerRadiusLow());
+ pt->PropagateTo(fkParam->GetInnerRadiusLow());
}
if (( pt->GetPoints()[2]- pt->GetPoints()[0])>5 && pt->GetPoints()[3]>0.8){
printf("Number of filled ESDs-\t%d\n",fEvent->GetNumberOfTracks());
}
-void AliTPCtrackerMI::WriteTracks(TTree * tree)
-{
- //
- // write tracks from seed array to selected tree
- //
- fOutput = tree;
- if (fOutput){
- AliTPCtrack *iotrack= new AliTPCtrack;
- fOutput->Branch("tracks","AliTPCtrack",&iotrack,32000,100);
- }
- WriteTracks();
-}
-void AliTPCtrackerMI::WriteTracks()
-{
- //
- // write tracks to the given output tree -
- // output specified with SetIO routine
- if (!fSeeds) return;
- if (!fOutput){
- SetIO();
- }
-
- if (fOutput){
- AliTPCtrack *iotrack= 0;
- Int_t nseed=fSeeds->GetEntriesFast();
- //for (Int_t i=0; i<nseed; i++) {
- // iotrack= (AliTPCtrack*)fSeeds->UncheckedAt(i);
- // if (iotrack) break;
- //}
- //TBranch * br = fOutput->Branch("tracks","AliTPCtrack",&iotrack,32000,100);
- TBranch * br = fOutput->GetBranch("tracks");
- br->SetAddress(&iotrack);
- //
- for (Int_t i=0; i<nseed; i++) {
- AliTPCseed *pt=(AliTPCseed*)fSeeds->UncheckedAt(i);
- if (!pt) continue;
- AliTPCtrack * track = new AliTPCtrack(*pt);
- iotrack = track;
- pt->SetLab2(i);
- // br->SetAddress(&iotrack);
- fOutput->Fill();
- delete track;
- iotrack =0;
- }
- //fOutput->GetDirectory()->cd();
- //fOutput->Write();
- }
- // delete iotrack;
- //
- if (fSeedTree){
- //write the full seed information if specified in debug mode
-
- AliTPCseed * vseed = new AliTPCseed;
- //
- TClonesArray * arrtr = new TClonesArray("AliTPCTrackPoint",160);
- arrtr->ExpandCreateFast(160);
- //TClonesArray * arrcl = new TClonesArray("AliTPCclusterMI",160);
- //arrcl->ExpandCreateFast(160);
- TClonesArray * arre = new TClonesArray("AliTPCExactPoint",160);
- //
- vseed->SetPoints(arrtr);
- vseed->SetEPoints(arre);
- // vseed->fClusterPoints = arrcl;
- //TBranch * brseed = seedtree->Branch("seeds","AliTPCseed",&vseed,32000,99);
- TBranch * brseed = fSeedTree->GetBranch("seeds");
-
- Int_t nseed=fSeeds->GetEntriesFast();
-
- for (Int_t i=0; i<nseed; i++) {
- AliTPCseed *pt=(AliTPCseed*)fSeeds->UncheckedAt(i);
- if (!pt) continue;
- pt->SetPoints(arrtr);
- // pt->fClusterPoints = arrcl;
- pt->SetEPoints(arre);
- pt->RebuildSeed();
- vseed = pt;
- brseed->SetAddress(&vseed);
- fSeedTree->Fill();
- pt->SetPoints(0);
- pt->SetEPoints(0);
- // pt->fClusterPoints = 0;
- }
- fSeedTree->Write();
- if (fTreeDebug) fTreeDebug->Write();
- }
-
-}
-
-Double_t AliTPCtrackerMI::ErrY2(AliTPCseed* seed, AliTPCclusterMI * cl){
+Double_t AliTPCtrackerMI::ErrY2(AliTPCseed* seed, const AliTPCclusterMI * cl){
//
//
- //seed->SetErrorY2(0.1);
- //return 0.1;
- //calculate look-up table at the beginning
- static Bool_t ginit = kFALSE;
- static Float_t gnoise1,gnoise2,gnoise3;
- static Float_t ggg1[10000];
- static Float_t ggg2[10000];
- static Float_t ggg3[10000];
- static Float_t glandau1[10000];
- static Float_t glandau2[10000];
- static Float_t glandau3[10000];
+ // Use calibrated cluster error from OCDB
//
- static Float_t gcor01[500];
- static Float_t gcor02[500];
- static Float_t gcorp[500];
+ AliTPCClusterParam * clparam = AliTPCcalibDB::Instance()->GetClusterParam();
//
+ Float_t z = TMath::Abs(fkParam->GetZLength(0)-TMath::Abs(seed->GetZ()));
+ Int_t ctype = cl->GetType();
+ Int_t type = (cl->GetRow()<63) ? 0: (cl->GetRow()>126) ? 1:2;
+ Double_t angle = seed->GetSnp()*seed->GetSnp();
+ angle = TMath::Sqrt(TMath::Abs(angle/(1.-angle)));
+ Double_t erry2 = clparam->GetError0Par(0,type, z,angle);
+ if (ctype<0) {
+ erry2+=0.5; // edge cluster
+ }
+ erry2*=erry2;
+ seed->SetErrorY2(erry2);
+ //
+ return erry2;
+
+//calculate look-up table at the beginning
+// static Bool_t ginit = kFALSE;
+// static Float_t gnoise1,gnoise2,gnoise3;
+// static Float_t ggg1[10000];
+// static Float_t ggg2[10000];
+// static Float_t ggg3[10000];
+// static Float_t glandau1[10000];
+// static Float_t glandau2[10000];
+// static Float_t glandau3[10000];
+// //
+// static Float_t gcor01[500];
+// static Float_t gcor02[500];
+// static Float_t gcorp[500];
+// //
- //
- if (ginit==kFALSE){
- for (Int_t i=1;i<500;i++){
- Float_t rsigma = float(i)/100.;
- gcor02[i] = TMath::Max(0.78 +TMath::Exp(7.4*(rsigma-1.2)),0.6);
- gcor01[i] = TMath::Max(0.72 +TMath::Exp(3.36*(rsigma-1.2)),0.6);
- gcorp[i] = TMath::Max(TMath::Power((rsigma+0.5),1.5),1.2);
- }
+// //
+// if (ginit==kFALSE){
+// for (Int_t i=1;i<500;i++){
+// Float_t rsigma = float(i)/100.;
+// gcor02[i] = TMath::Max(0.78 +TMath::Exp(7.4*(rsigma-1.2)),0.6);
+// gcor01[i] = TMath::Max(0.72 +TMath::Exp(3.36*(rsigma-1.2)),0.6);
+// gcorp[i] = TMath::Max(TMath::Power((rsigma+0.5),1.5),1.2);
+// }
- //
- for (Int_t i=3;i<10000;i++){
- //
- //
- // inner sector
- Float_t amp = float(i);
- Float_t padlength =0.75;
- gnoise1 = 0.0004/padlength;
- Float_t nel = 0.268*amp;
- Float_t nprim = 0.155*amp;
- ggg1[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.001*nel/(padlength*padlength))/nel;
- glandau1[i] = (2.+0.12*nprim)*0.5* (2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
- if (glandau1[i]>1) glandau1[i]=1;
- glandau1[i]*=padlength*padlength/12.;
- //
- // outer short
- padlength =1.;
- gnoise2 = 0.0004/padlength;
- nel = 0.3*amp;
- nprim = 0.133*amp;
- ggg2[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
- glandau2[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
- if (glandau2[i]>1) glandau2[i]=1;
- glandau2[i]*=padlength*padlength/12.;
- //
- //
- // outer long
- padlength =1.5;
- gnoise3 = 0.0004/padlength;
- nel = 0.3*amp;
- nprim = 0.133*amp;
- ggg3[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
- glandau3[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
- if (glandau3[i]>1) glandau3[i]=1;
- glandau3[i]*=padlength*padlength/12.;
- //
- }
- ginit = kTRUE;
- }
- //
- //
- //
- Int_t amp = int(TMath::Abs(cl->GetQ()));
- if (amp>9999) {
- seed->SetErrorY2(1.);
- return 1.;
- }
- Float_t snoise2;
- Float_t z = TMath::Abs(fParam->GetZLength()-TMath::Abs(seed->GetZ()));
- Int_t ctype = cl->GetType();
- Float_t padlength= GetPadPitchLength(seed->GetRow());
- Double_t angle2 = seed->GetSnp()*seed->GetSnp();
- angle2 = angle2/(1-angle2);
- //
- //cluster "quality"
- Int_t rsigmay = int(100.*cl->GetSigmaY2()/(seed->GetCurrentSigmaY2()));
- Float_t res;
- //
- if (fSectors==fInnerSec){
- snoise2 = gnoise1;
- res = ggg1[amp]*z+glandau1[amp]*angle2;
- if (ctype==0) res *= gcor01[rsigmay];
- if ((ctype>0)){
- res+=0.002;
- res*= gcorp[rsigmay];
- }
- }
- else {
- if (padlength<1.1){
- snoise2 = gnoise2;
- res = ggg2[amp]*z+glandau2[amp]*angle2;
- if (ctype==0) res *= gcor02[rsigmay];
- if ((ctype>0)){
- res+=0.002;
- res*= gcorp[rsigmay];
- }
- }
- else{
- snoise2 = gnoise3;
- res = ggg3[amp]*z+glandau3[amp]*angle2;
- if (ctype==0) res *= gcor02[rsigmay];
- if ((ctype>0)){
- res+=0.002;
- res*= gcorp[rsigmay];
- }
- }
- }
+// //
+// for (Int_t i=3;i<10000;i++){
+// //
+// //
+// // inner sector
+// Float_t amp = float(i);
+// Float_t padlength =0.75;
+// gnoise1 = 0.0004/padlength;
+// Float_t nel = 0.268*amp;
+// Float_t nprim = 0.155*amp;
+// ggg1[i] = fkParam->GetDiffT()*fkParam->GetDiffT()*(2+0.001*nel/(padlength*padlength))/nel;
+// glandau1[i] = (2.+0.12*nprim)*0.5* (2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
+// if (glandau1[i]>1) glandau1[i]=1;
+// glandau1[i]*=padlength*padlength/12.;
+// //
+// // outer short
+// padlength =1.;
+// gnoise2 = 0.0004/padlength;
+// nel = 0.3*amp;
+// nprim = 0.133*amp;
+// ggg2[i] = fkParam->GetDiffT()*fkParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
+// glandau2[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
+// if (glandau2[i]>1) glandau2[i]=1;
+// glandau2[i]*=padlength*padlength/12.;
+// //
+// //
+// // outer long
+// padlength =1.5;
+// gnoise3 = 0.0004/padlength;
+// nel = 0.3*amp;
+// nprim = 0.133*amp;
+// ggg3[i] = fkParam->GetDiffT()*fkParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
+// glandau3[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
+// if (glandau3[i]>1) glandau3[i]=1;
+// glandau3[i]*=padlength*padlength/12.;
+// //
+// }
+// ginit = kTRUE;
+// }
+// //
+// //
+// //
+// Int_t amp = int(TMath::Abs(cl->GetQ()));
+// if (amp>9999) {
+// seed->SetErrorY2(1.);
+// return 1.;
+// }
+// Float_t snoise2;
+// Float_t z = TMath::Abs(fkParam->GetZLength(0)-TMath::Abs(seed->GetZ()));
+// Int_t ctype = cl->GetType();
+// Float_t padlength= GetPadPitchLength(seed->GetRow());
+// Double_t angle2 = seed->GetSnp()*seed->GetSnp();
+// angle2 = angle2/(1-angle2);
+// //
+// //cluster "quality"
+// Int_t rsigmay = int(100.*cl->GetSigmaY2()/(seed->GetCurrentSigmaY2()));
+// Float_t res;
+// //
+// if (fSectors==fInnerSec){
+// snoise2 = gnoise1;
+// res = ggg1[amp]*z+glandau1[amp]*angle2;
+// if (ctype==0) res *= gcor01[rsigmay];
+// if ((ctype>0)){
+// res+=0.002;
+// res*= gcorp[rsigmay];
+// }
+// }
+// else {
+// if (padlength<1.1){
+// snoise2 = gnoise2;
+// res = ggg2[amp]*z+glandau2[amp]*angle2;
+// if (ctype==0) res *= gcor02[rsigmay];
+// if ((ctype>0)){
+// res+=0.002;
+// res*= gcorp[rsigmay];
+// }
+// }
+// else{
+// snoise2 = gnoise3;
+// res = ggg3[amp]*z+glandau3[amp]*angle2;
+// if (ctype==0) res *= gcor02[rsigmay];
+// if ((ctype>0)){
+// res+=0.002;
+// res*= gcorp[rsigmay];
+// }
+// }
+// }
- if (ctype<0){
- res+=0.005;
- res*=2.4; // overestimate error 2 times
- }
- res+= snoise2;
+// if (ctype<0){
+// res+=0.005;
+// res*=2.4; // overestimate error 2 times
+// }
+// res+= snoise2;
- if (res<2*snoise2)
- res = 2*snoise2;
+// if (res<2*snoise2)
+// res = 2*snoise2;
- seed->SetErrorY2(res);
- return res;
+// seed->SetErrorY2(res);
+// return res;
}
-Double_t AliTPCtrackerMI::ErrZ2(AliTPCseed* seed, AliTPCclusterMI * cl){
- //
+Double_t AliTPCtrackerMI::ErrZ2(AliTPCseed* seed, const AliTPCclusterMI * cl){
//
- //seed->SetErrorY2(0.1);
- //return 0.1;
- //calculate look-up table at the beginning
- static Bool_t ginit = kFALSE;
- static Float_t gnoise1,gnoise2,gnoise3;
- static Float_t ggg1[10000];
- static Float_t ggg2[10000];
- static Float_t ggg3[10000];
- static Float_t glandau1[10000];
- static Float_t glandau2[10000];
- static Float_t glandau3[10000];
//
- static Float_t gcor01[1000];
- static Float_t gcor02[1000];
- static Float_t gcorp[1000];
+ // Use calibrated cluster error from OCDB
//
-
+ AliTPCClusterParam * clparam = AliTPCcalibDB::Instance()->GetClusterParam();
//
- if (ginit==kFALSE){
- for (Int_t i=1;i<1000;i++){
- Float_t rsigma = float(i)/100.;
- gcor02[i] = TMath::Max(0.81 +TMath::Exp(6.8*(rsigma-1.2)),0.6);
- gcor01[i] = TMath::Max(0.72 +TMath::Exp(2.04*(rsigma-1.2)),0.6);
- gcorp[i] = TMath::Max(TMath::Power((rsigma+0.5),1.5),1.2);
- }
-
- //
- for (Int_t i=3;i<10000;i++){
- //
- //
- // inner sector
- Float_t amp = float(i);
- Float_t padlength =0.75;
- gnoise1 = 0.0004/padlength;
- Float_t nel = 0.268*amp;
- Float_t nprim = 0.155*amp;
- ggg1[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.001*nel/(padlength*padlength))/nel;
- glandau1[i] = (2.+0.12*nprim)*0.5* (2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
- if (glandau1[i]>1) glandau1[i]=1;
- glandau1[i]*=padlength*padlength/12.;
- //
- // outer short
- padlength =1.;
- gnoise2 = 0.0004/padlength;
- nel = 0.3*amp;
- nprim = 0.133*amp;
- ggg2[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
- glandau2[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
- if (glandau2[i]>1) glandau2[i]=1;
- glandau2[i]*=padlength*padlength/12.;
- //
- //
- // outer long
- padlength =1.5;
- gnoise3 = 0.0004/padlength;
- nel = 0.3*amp;
- nprim = 0.133*amp;
- ggg3[i] = fParam->GetDiffT()*fParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
- glandau3[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
- if (glandau3[i]>1) glandau3[i]=1;
- glandau3[i]*=padlength*padlength/12.;
- //
- }
- ginit = kTRUE;
- }
- //
- //
- //
- Int_t amp = int(TMath::Abs(cl->GetQ()));
- if (amp>9999) {
- seed->SetErrorY2(1.);
- return 1.;
- }
- Float_t snoise2;
- Float_t z = TMath::Abs(fParam->GetZLength()-TMath::Abs(seed->GetZ()));
+ Float_t z = TMath::Abs(fkParam->GetZLength(0)-TMath::Abs(seed->GetZ()));
Int_t ctype = cl->GetType();
- Float_t padlength= GetPadPitchLength(seed->GetRow());
+ Int_t type = (cl->GetRow()<63) ? 0: (cl->GetRow()>126) ? 1:2;
//
Double_t angle2 = seed->GetSnp()*seed->GetSnp();
- // if (angle2<0.6) angle2 = 0.6;
angle2 = seed->GetTgl()*seed->GetTgl()*(1+angle2/(1-angle2));
- //
- //cluster "quality"
- Int_t rsigmaz = int(100.*cl->GetSigmaZ2()/(seed->GetCurrentSigmaZ2()));
- Float_t res;
- //
- if (fSectors==fInnerSec){
- snoise2 = gnoise1;
- res = ggg1[amp]*z+glandau1[amp]*angle2;
- if (ctype==0) res *= gcor01[rsigmaz];
- if ((ctype>0)){
- res+=0.002;
- res*= gcorp[rsigmaz];
- }
+ Double_t angle = TMath::Sqrt(TMath::Abs(angle2));
+ Double_t errz2 = clparam->GetError0Par(1,type, z,angle);
+ if (ctype<0) {
+ errz2+=0.5; // edge cluster
}
- else {
- if (padlength<1.1){
- snoise2 = gnoise2;
- res = ggg2[amp]*z+glandau2[amp]*angle2;
- if (ctype==0) res *= gcor02[rsigmaz];
- if ((ctype>0)){
- res+=0.002;
- res*= gcorp[rsigmaz];
- }
- }
- else{
- snoise2 = gnoise3;
- res = ggg3[amp]*z+glandau3[amp]*angle2;
- if (ctype==0) res *= gcor02[rsigmaz];
- if ((ctype>0)){
- res+=0.002;
- res*= gcorp[rsigmaz];
- }
- }
- }
-
- if (ctype<0){
- res+=0.002;
- res*=1.3;
- }
- if ((ctype<0) &&<70){
- res+=0.002;
- res*=1.3;
- }
- res += snoise2;
- if (res<2*snoise2)
- res = 2*snoise2;
- if (res>3) res =3;
- seed->SetErrorZ2(res);
- return res;
-}
-
-
-
-/*
-Double_t AliTPCtrackerMI::ErrZ2(AliTPCseed* seed, AliTPCclusterMI * cl){
- //
+ errz2*=errz2;
+ seed->SetErrorZ2(errz2);
//
- //seed->SetErrorZ2(0.1);
- //return 0.1;
+ return errz2;
- Float_t snoise2;
- Float_t z = TMath::Abs(fParam->GetZLength()-TMath::Abs(seed->GetZ()));
- //
- Float_t rsigmaz = cl->GetSigmaZ2()/(seed->fCurrentSigmaZ2);
- Int_t ctype = cl->GetType();
- Float_t amp = TMath::Abs(cl->GetQ());
-
- Float_t nel;
- Float_t nprim;
- //
- Float_t landau=2 ; //landau fluctuation part
- Float_t gg=2; // gg fluctuation part
- Float_t padlength= GetPadPitchLength(seed->GetX());
-
- if (fSectors==fInnerSec){
- snoise2 = 0.0004/padlength;
- nel = 0.268*amp;
- nprim = 0.155*amp;
- gg = (2+0.001*nel/(padlength*padlength))/nel;
- landau = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
- if (landau>1) landau=1;
- }
- else {
- snoise2 = 0.0004/padlength;
- nel = 0.3*amp;
- nprim = 0.133*amp;
- gg = (2+0.0008*nel/(padlength*padlength))/nel;
- landau = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
- if (landau>1) landau=1;
- }
- Float_t sdiff = gg*fParam->GetDiffT()*fParam->GetDiffT()*z;
- //
- Float_t angle2 = seed->GetSnp()*seed->GetSnp();
- angle2 = TMath::Sqrt((1-angle2));
- if (angle2<0.6) angle2 = 0.6;
- //angle2 = 1;
- Float_t angle = seed->GetTgl()/angle2;
- Float_t angular = landau*angle*angle*padlength*padlength/12.;
- Float_t res = sdiff + angular;
+// //seed->SetErrorY2(0.1);
+// //return 0.1;
+// //calculate look-up table at the beginning
+// static Bool_t ginit = kFALSE;
+// static Float_t gnoise1,gnoise2,gnoise3;
+// static Float_t ggg1[10000];
+// static Float_t ggg2[10000];
+// static Float_t ggg3[10000];
+// static Float_t glandau1[10000];
+// static Float_t glandau2[10000];
+// static Float_t glandau3[10000];
+// //
+// static Float_t gcor01[1000];
+// static Float_t gcor02[1000];
+// static Float_t gcorp[1000];
+// //
-
- if ((ctype==0) && (fSectors ==fOuterSec))
- res *= 0.81 +TMath::Exp(6.8*(rsigmaz-1.2));
+// //
+// if (ginit==kFALSE){
+// for (Int_t i=1;i<1000;i++){
+// Float_t rsigma = float(i)/100.;
+// gcor02[i] = TMath::Max(0.81 +TMath::Exp(6.8*(rsigma-1.2)),0.6);
+// gcor01[i] = TMath::Max(0.72 +TMath::Exp(2.04*(rsigma-1.2)),0.6);
+// gcorp[i] = TMath::Max(TMath::Power((rsigma+0.5),1.5),1.2);
+// }
- if ((ctype==0) && (fSectors ==fInnerSec))
- res *= 0.72 +TMath::Exp(2.04*(rsigmaz-1.2));
-
- if ((ctype>0)){
- res+=0.005;
- res*= TMath::Power(rsigmaz+0.5,1.5); //0.31+0.147*ctype;
- }
- if (ctype<0){
- res+=0.002;
- res*=1.3;
- }
- if ((ctype<0) &&<70){
- res+=0.002;
- res*=1.3;
- }
- res += snoise2;
- if (res<2*snoise2)
- res = 2*snoise2;
+// //
+// for (Int_t i=3;i<10000;i++){
+// //
+// //
+// // inner sector
+// Float_t amp = float(i);
+// Float_t padlength =0.75;
+// gnoise1 = 0.0004/padlength;
+// Float_t nel = 0.268*amp;
+// Float_t nprim = 0.155*amp;
+// ggg1[i] = fkParam->GetDiffT()*fkParam->GetDiffT()*(2+0.001*nel/(padlength*padlength))/nel;
+// glandau1[i] = (2.+0.12*nprim)*0.5* (2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
+// if (glandau1[i]>1) glandau1[i]=1;
+// glandau1[i]*=padlength*padlength/12.;
+// //
+// // outer short
+// padlength =1.;
+// gnoise2 = 0.0004/padlength;
+// nel = 0.3*amp;
+// nprim = 0.133*amp;
+// ggg2[i] = fkParam->GetDiffT()*fkParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
+// glandau2[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
+// if (glandau2[i]>1) glandau2[i]=1;
+// glandau2[i]*=padlength*padlength/12.;
+// //
+// //
+// // outer long
+// padlength =1.5;
+// gnoise3 = 0.0004/padlength;
+// nel = 0.3*amp;
+// nprim = 0.133*amp;
+// ggg3[i] = fkParam->GetDiffT()*fkParam->GetDiffT()*(2+0.0008*nel/(padlength*padlength))/nel;
+// glandau3[i] = (2.+0.12*nprim)*0.5*(2.+nprim*nprim*0.001/(padlength*padlength))/nprim;
+// if (glandau3[i]>1) glandau3[i]=1;
+// glandau3[i]*=padlength*padlength/12.;
+// //
+// }
+// ginit = kTRUE;
+// }
+// //
+// //
+// //
+// Int_t amp = int(TMath::Abs(cl->GetQ()));
+// if (amp>9999) {
+// seed->SetErrorY2(1.);
+// return 1.;
+// }
+// Float_t snoise2;
+// Float_t z = TMath::Abs(fkParam->GetZLength(0)-TMath::Abs(seed->GetZ()));
+// Int_t ctype = cl->GetType();
+// Float_t padlength= GetPadPitchLength(seed->GetRow());
+// //
+// Double_t angle2 = seed->GetSnp()*seed->GetSnp();
+// // if (angle2<0.6) angle2 = 0.6;
+// angle2 = seed->GetTgl()*seed->GetTgl()*(1+angle2/(1-angle2));
+// //
+// //cluster "quality"
+// Int_t rsigmaz = int(100.*cl->GetSigmaZ2()/(seed->GetCurrentSigmaZ2()));
+// Float_t res;
+// //
+// if (fSectors==fInnerSec){
+// snoise2 = gnoise1;
+// res = ggg1[amp]*z+glandau1[amp]*angle2;
+// if (ctype==0) res *= gcor01[rsigmaz];
+// if ((ctype>0)){
+// res+=0.002;
+// res*= gcorp[rsigmaz];
+// }
+// }
+// else {
+// if (padlength<1.1){
+// snoise2 = gnoise2;
+// res = ggg2[amp]*z+glandau2[amp]*angle2;
+// if (ctype==0) res *= gcor02[rsigmaz];
+// if ((ctype>0)){
+// res+=0.002;
+// res*= gcorp[rsigmaz];
+// }
+// }
+// else{
+// snoise2 = gnoise3;
+// res = ggg3[amp]*z+glandau3[amp]*angle2;
+// if (ctype==0) res *= gcor02[rsigmaz];
+// if ((ctype>0)){
+// res+=0.002;
+// res*= gcorp[rsigmaz];
+// }
+// }
+// }
- seed->SetErrorZ2(res);
- return res;
+// if (ctype<0){
+// res+=0.002;
+// res*=1.3;
+// }
+// if ((ctype<0) &&<70){
+// res+=0.002;
+// res*=1.3;
+// }
+// res += snoise2;
+// if (res<2*snoise2)
+// res = 2*snoise2;
+// if (res>3) res =3;
+// seed->SetErrorZ2(res);
+// return res;
}
-*/
+
//_____________________________________________________________________________
Double_t AliTPCtrackerMI::F1old(Double_t x1,Double_t y1,
Double_t x2,Double_t y2,
- Double_t x3,Double_t y3)
+ Double_t x3,Double_t y3) const
{
//-----------------------------------------------------------------
// Initial approximation of the track curvature
//_____________________________________________________________________________
Double_t AliTPCtrackerMI::F1(Double_t x1,Double_t y1,
Double_t x2,Double_t y2,
- Double_t x3,Double_t y3)
+ Double_t x3,Double_t y3) const
{
//-----------------------------------------------------------------
// Initial approximation of the track curvature
Double_t AliTPCtrackerMI::F2(Double_t x1,Double_t y1,
Double_t x2,Double_t y2,
- Double_t x3,Double_t y3)
+ Double_t x3,Double_t y3) const
{
//-----------------------------------------------------------------
// Initial approximation of the track curvature
//_____________________________________________________________________________
Double_t AliTPCtrackerMI::F2old(Double_t x1,Double_t y1,
Double_t x2,Double_t y2,
- Double_t x3,Double_t y3)
+ Double_t x3,Double_t y3) const
{
//-----------------------------------------------------------------
// Initial approximation of the track curvature times center of curvature
//_____________________________________________________________________________
Double_t AliTPCtrackerMI::F3(Double_t x1,Double_t y1,
Double_t x2,Double_t y2,
- Double_t z1,Double_t z2)
+ Double_t z1,Double_t z2) const
{
//-----------------------------------------------------------------
// Initial approximation of the tangent of the track dip angle
Double_t AliTPCtrackerMI::F3n(Double_t x1,Double_t y1,
Double_t x2,Double_t y2,
- Double_t z1,Double_t z2, Double_t c)
+ Double_t z1,Double_t z2, Double_t c) const
{
//-----------------------------------------------------------------
// Initial approximation of the tangent of the track dip angle
return kFALSE;
}
- Double_t c1=x[4]*x1 - x[2], r1=sqrt(1.- c1*c1);
- Double_t c2=x[4]*x2 - x[2], r2=sqrt(1.- c2*c2);
+ Double_t c1=x[4]*x1 - x[2], r1=sqrt((1.-c1)*(1.+c1));
+ Double_t c2=x[4]*x2 - x[2], r2=sqrt((1.-c2)*(1.+c2));
y = x[0];
z = x[1];
return LoadClusters();
}
+
+Int_t AliTPCtrackerMI::LoadClusters(const TObjArray *arr)
+{
+ //
+ // load clusters to the memory
+ AliTPCClustersRow *clrow = 0x0;
+ Int_t lower = arr->LowerBound();
+ Int_t entries = arr->GetEntriesFast();
+
+ for (Int_t i=lower; i<entries; i++) {
+ clrow = (AliTPCClustersRow*) arr->At(i);
+ if(!clrow) continue;
+ if(!clrow->GetArray()) continue;
+
+ //
+ Int_t sec,row;
+ fkParam->AdjustSectorRow(clrow->GetID(),sec,row);
+
+ for (Int_t icl=0; icl<clrow->GetArray()->GetEntriesFast(); icl++){
+ Transform((AliTPCclusterMI*)(clrow->GetArray()->At(icl)));
+ }
+ //
+ if (clrow->GetArray()->GetEntriesFast()<=0) continue;
+ AliTPCtrackerRow * tpcrow=0;
+ Int_t left=0;
+ if (sec<fkNIS*2){
+ tpcrow = &(fInnerSec[sec%fkNIS][row]);
+ left = sec/fkNIS;
+ }
+ else{
+ tpcrow = &(fOuterSec[(sec-fkNIS*2)%fkNOS][row]);
+ left = (sec-fkNIS*2)/fkNOS;
+ }
+ if (left ==0){
+ tpcrow->SetN1(clrow->GetArray()->GetEntriesFast());
+ tpcrow->SetClusters1(new AliTPCclusterMI[tpcrow->GetN1()]);
+ for (Int_t j=0;j<tpcrow->GetN1();++j)
+ tpcrow->SetCluster1(j, *(AliTPCclusterMI*)(clrow->GetArray()->At(j)));
+ }
+ if (left ==1){
+ tpcrow->SetN2(clrow->GetArray()->GetEntriesFast());
+ tpcrow->SetClusters2(new AliTPCclusterMI[tpcrow->GetN2()]);
+ for (Int_t j=0;j<tpcrow->GetN2();++j)
+ tpcrow->SetCluster2(j,*(AliTPCclusterMI*)(clrow->GetArray()->At(j)));
+ }
+ }
+ //
+ delete clrow;
+ LoadOuterSectors();
+ LoadInnerSectors();
+ return 0;
+}
+
+Int_t AliTPCtrackerMI::LoadClusters(const TClonesArray *arr)
+{
+ //
+ // load clusters to the memory from one
+ // TClonesArray
+ //
+ AliTPCclusterMI *clust=0;
+ Int_t count[72][96] = { {0} , {0} };
+
+ // loop over clusters
+ for (Int_t icl=0; icl<arr->GetEntriesFast(); icl++) {
+ clust = (AliTPCclusterMI*)arr->At(icl);
+ if(!clust) continue;
+ //printf("cluster: det %d, row %d \n", clust->GetDetector(),clust->GetRow());
+
+ // transform clusters
+ Transform(clust);
+
+ // count clusters per pad row
+ count[clust->GetDetector()][clust->GetRow()]++;
+ }
+
+ // insert clusters to sectors
+ for (Int_t icl=0; icl<arr->GetEntriesFast(); icl++) {
+ clust = (AliTPCclusterMI*)arr->At(icl);
+ if(!clust) continue;
+
+ Int_t sec = clust->GetDetector();
+ Int_t row = clust->GetRow();
+
+ // filter overlapping pad rows needed by HLT
+ if(sec<fkNIS*2) { //IROCs
+ if(row == 30) continue;
+ }
+ else { // OROCs
+ if(row == 27 || row == 76) continue;
+ }
+
+ Int_t left=0;
+ if (sec<fkNIS*2){
+ left = sec/fkNIS;
+ fInnerSec[sec%fkNIS].InsertCluster(clust, count[sec][row], fkParam);
+ }
+ else{
+ left = (sec-fkNIS*2)/fkNOS;
+ fOuterSec[(sec-fkNIS*2)%fkNOS].InsertCluster(clust, count[sec][row], fkParam);
+ }
+ }
+
+ // Load functions must be called behind LoadCluster(TClonesArray*)
+ // needed by HLT
+ //LoadOuterSectors();
+ //LoadInnerSectors();
+
+ return 0;
+}
+
+
Int_t AliTPCtrackerMI::LoadClusters()
{
//
br->GetEntry(i);
//
Int_t sec,row;
- fParam->AdjustSectorRow(clrow->GetID(),sec,row);
+ fkParam->AdjustSectorRow(clrow->GetID(),sec,row);
for (Int_t icl=0; icl<clrow->GetArray()->GetEntriesFast(); icl++){
Transform((AliTPCclusterMI*)(clrow->GetArray()->At(icl)));
}
//
- AliTPCRow * tpcrow=0;
+ AliTPCtrackerRow * tpcrow=0;
Int_t left=0;
if (sec<fkNIS*2){
tpcrow = &(fInnerSec[sec%fkNIS][row]);
if (left ==0){
tpcrow->SetN1(clrow->GetArray()->GetEntriesFast());
tpcrow->SetClusters1(new AliTPCclusterMI[tpcrow->GetN1()]);
- for (Int_t i=0;i<tpcrow->GetN1();i++)
- tpcrow->SetCluster1(i, *(AliTPCclusterMI*)(clrow->GetArray()->At(i)));
+ for (Int_t k=0;k<tpcrow->GetN1();++k)
+ tpcrow->SetCluster1(k, *(AliTPCclusterMI*)(clrow->GetArray()->At(k)));
}
if (left ==1){
tpcrow->SetN2(clrow->GetArray()->GetEntriesFast());
tpcrow->SetClusters2(new AliTPCclusterMI[tpcrow->GetN2()]);
- for (Int_t i=0;i<tpcrow->GetN2();i++)
- tpcrow->SetCluster2(i,*(AliTPCclusterMI*)(clrow->GetArray()->At(i)));
+ for (Int_t k=0;k<tpcrow->GetN2();++k)
+ tpcrow->SetCluster2(k,*(AliTPCclusterMI*)(clrow->GetArray()->At(k)));
}
}
//
Int_t nrows = fOuterSec->GetNRows();
for (Int_t sec = 0;sec<fkNOS;sec++)
for (Int_t row = 0;row<nrows;row++){
- AliTPCRow* tpcrow = &(fOuterSec[sec%fkNOS][row]);
+ AliTPCtrackerRow* tpcrow = &(fOuterSec[sec%fkNOS][row]);
// if (tpcrow){
// if (tpcrow->fClusters1) delete []tpcrow->fClusters1;
// if (tpcrow->fClusters2) delete []tpcrow->fClusters2;
nrows = fInnerSec->GetNRows();
for (Int_t sec = 0;sec<fkNIS;sec++)
for (Int_t row = 0;row<nrows;row++){
- AliTPCRow* tpcrow = &(fInnerSec[sec%fkNIS][row]);
+ AliTPCtrackerRow* tpcrow = &(fInnerSec[sec%fkNIS][row]);
//if (tpcrow){
// if (tpcrow->fClusters1) delete []tpcrow->fClusters1;
//if (tpcrow->fClusters2) delete []tpcrow->fClusters2;
return ;
}
+void AliTPCtrackerMI::FillClusterArray(TObjArray* array) const{
+ //
+ // Filling cluster to the array - For visualization purposes
+ //
+ Int_t nrows=0;
+ nrows = fOuterSec->GetNRows();
+ for (Int_t sec = 0;sec<fkNOS;sec++)
+ for (Int_t row = 0;row<nrows;row++){
+ AliTPCtrackerRow* tpcrow = &(fOuterSec[sec%fkNOS][row]);
+ if (!tpcrow) continue;
+ for (Int_t icl = 0;icl<tpcrow->GetN();icl++){
+ array->AddLast((TObject*)((*tpcrow)[icl]));
+ }
+ }
+ nrows = fInnerSec->GetNRows();
+ for (Int_t sec = 0;sec<fkNIS;sec++)
+ for (Int_t row = 0;row<nrows;row++){
+ AliTPCtrackerRow* tpcrow = &(fInnerSec[sec%fkNIS][row]);
+ if (!tpcrow) continue;
+ for (Int_t icl = 0;icl<tpcrow->GetN();icl++){
+ array->AddLast((TObject*)(*tpcrow)[icl]);
+ }
+ }
+}
+
+
void AliTPCtrackerMI::Transform(AliTPCclusterMI * cluster){
+ //
+ //
+ //
+ AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
+ if (!transform) {
+ AliFatal("Tranformations not in calibDB");
+ }
+ transform->SetCurrentRecoParam((AliTPCRecoParam*)AliTPCReconstructor::GetRecoParam());
+ Double_t x[3]={cluster->GetRow(),cluster->GetPad(),cluster->GetTimeBin()};
+ Int_t i[1]={cluster->GetDetector()};
+ transform->Transform(x,i,0,1);
+ // if (cluster->GetDetector()%36>17){
+ // x[1]*=-1;
+ //}
+
+ //
+ // in debug mode check the transformation
+ //
+ if (AliTPCReconstructor::StreamLevel()>1) {
+ Float_t gx[3];
+ cluster->GetGlobalXYZ(gx);
+ Int_t event = (fEvent==NULL)? 0: fEvent->GetEventNumberInFile();
+ TTreeSRedirector &cstream = *fDebugStreamer;
+ cstream<<"Transform"<<
+ "event="<<event<<
+ "x0="<<x[0]<<
+ "x1="<<x[1]<<
+ "x2="<<x[2]<<
+ "gx0="<<gx[0]<<
+ "gx1="<<gx[1]<<
+ "gx2="<<gx[2]<<
+ "Cl.="<<cluster<<
+ "\n";
+ }
+ cluster->SetX(x[0]);
+ cluster->SetY(x[1]);
+ cluster->SetZ(x[2]);
+ // The old stuff:
+
//
//
//
- //if (!fParam->IsGeoRead()) fParam->ReadGeoMatrices();
- TGeoHMatrix *mat = fParam->GetClusterMatrix(cluster->GetDetector());
+ //if (!fkParam->IsGeoRead()) fkParam->ReadGeoMatrices();
+ TGeoHMatrix *mat = fkParam->GetClusterMatrix(cluster->GetDetector());
//TGeoHMatrix mat;
Double_t pos[3]= {cluster->GetX(),cluster->GetY(),cluster->GetZ()};
Double_t posC[3]={cluster->GetX(),cluster->GetY(),cluster->GetZ()};
else{
// chack Loading of Geo matrices from GeoManager - TEMPORARY FIX
}
- //cluster->SetX(posC[0]);
- //cluster->SetY(posC[1]);
- //cluster->SetZ(posC[2]);
+ cluster->SetX(posC[0]);
+ cluster->SetY(posC[1]);
+ cluster->SetZ(posC[2]);
}
//_____________________________________________________________________________
UInt_t index=0;
for (Int_t sec = 0;sec<fkNOS;sec++)
for (Int_t row = 0;row<nrows;row++){
- AliTPCRow* tpcrow = &(fOuterSec[sec%fkNOS][row]);
+ AliTPCtrackerRow* tpcrow = &(fOuterSec[sec%fkNOS][row]);
Int_t sec2 = sec+2*fkNIS;
//left
Int_t ncl = tpcrow->GetN1();
UInt_t index=0;
for (Int_t sec = 0;sec<fkNIS;sec++)
for (Int_t row = 0;row<nrows;row++){
- AliTPCRow* tpcrow = &(fInnerSec[sec%fkNIS][row]);
+ AliTPCtrackerRow* tpcrow = &(fInnerSec[sec%fkNIS][row]);
//
//left
Int_t ncl = tpcrow->GetN1();
Int_t row=(index&0x00ff0000)>>16;
Int_t ncl=(index&0x00007fff)>>00;
- const AliTPCRow * tpcrow=0;
+ const AliTPCtrackerRow * tpcrow=0;
AliTPCclusterMI * clrow =0;
if (sec<0 || sec>=fkNIS*4) {
}
if (sec<fkNIS*2){
- tpcrow = &(fInnerSec[sec%fkNIS][row]);
+ AliTPCtrackerSector& tracksec = fInnerSec[sec%fkNIS];
+ if (tracksec.GetNRows()<=row) return 0;
+ tpcrow = &(tracksec[row]);
if (tpcrow==0) return 0;
if (sec<fkNIS) {
}
}
else {
- tpcrow = &(fOuterSec[(sec-fkNIS*2)%fkNOS][row]);
+ AliTPCtrackerSector& tracksec = fOuterSec[(sec-fkNIS*2)%fkNOS];
+ if (tracksec.GetNRows()<=row) return 0;
+ tpcrow = &(tracksec[row]);
if (tpcrow==0) return 0;
if (sec-2*fkNIS<fkNOS) {
//-----------------------------------------------------------------
//
Double_t x= GetXrow(nr), ymax=GetMaxY(nr);
+ //
+ //
AliTPCclusterMI *cl=0;
Int_t tpcindex= t.GetClusterIndex2(nr);
//
//
t.SetCurrentCluster(cl);
t.SetRow(nr);
- Int_t accept = AcceptCluster(&t,t.GetCurrentCluster(),1.);
+ Int_t accept = AcceptCluster(&t,t.GetCurrentCluster());
if ((tpcindex&0x8000)==0) accept =0;
if (accept<3) {
//if founded cluster is acceptible
}
}
if (TMath::Abs(t.GetSnp())>AliTPCReconstructor::GetMaxSnpTracker()) return 0; // cut on angle
- if (fIteration>1){
+ if (fIteration>1 && IsFindable(t)){
// not look for new cluster during refitting
t.SetNFoundable(t.GetNFoundable()+1);
return 0;
if (!isActive || !isActive2) return 0;
- const AliTPCRow &krow=GetRow(t.GetRelativeSector(),nr);
+ const AliTPCtrackerRow &krow=GetRow(t.GetRelativeSector(),nr);
if ( (t.GetSigmaY2()<0) || t.GetSigmaZ2()<0) return 0;
Double_t roady =1.;
Double_t roadz = 1.;
}
else
{
- if (TMath::Abs(z)<(AliTPCReconstructor::GetCtgRange()*x+10) && TMath::Abs(z)<fParam->GetZLength() && (TMath::Abs(t.GetSnp())<AliTPCReconstructor::GetMaxSnpTracker()))
+ if (IsFindable(t))
+ // if (TMath::Abs(z)<(AliTPCReconstructor::GetCtgRange()*x+10) && TMath::Abs(z)<fkParam->GetZLength(0) && (TMath::Abs(t.GetSnp())<AliTPCReconstructor::GetMaxSnpTracker()))
t.SetNFoundable(t.GetNFoundable()+1);
else
return 0;
t.SetCurrentCluster(cl);
t.SetRow(nr);
if (fIteration==2&&cl->IsUsed(10)) return 0;
- Int_t accept = AcceptCluster(&t,t.GetCurrentCluster(),1.);
+ Int_t accept = AcceptCluster(&t,t.GetCurrentCluster());
if (fIteration==2&&cl->IsUsed(11)) {
t.SetErrorY2(t.GetErrorY2()+0.03);
t.SetErrorZ2(t.GetErrorZ2()+0.03);
return 1;
}
-Int_t AliTPCtrackerMI::FollowToNextFast(AliTPCseed& t, Int_t nr) {
- //-----------------------------------------------------------------
- // This function tries to find a track prolongation to next pad row
- //-----------------------------------------------------------------
- //
- Double_t x= GetXrow(nr), ymax=GetMaxY(nr);
- Double_t y,z;
- if (!t.GetProlongation(x,y,z)) {
- t.SetRemoval(10);
- return 0;
- }
- //
- //
- if (TMath::Abs(y)>ymax){
-
- if (y > ymax) {
- t.SetRelativeSector((t.GetRelativeSector()+1) % fN);
- if (!t.Rotate(fSectors->GetAlpha()))
- return 0;
- } else if (y <-ymax) {
- t.SetRelativeSector((t.GetRelativeSector()-1+fN) % fN);
- if (!t.Rotate(-fSectors->GetAlpha()))
- return 0;
- }
- if (!t.PropagateTo(x)) {
- return 0;
- }
- t.GetProlongation(x,y,z);
- }
- //
- // update current shape info every 3 pad-row
- if ( (nr%6==0) || t.GetNumberOfClusters()<2 || (t.GetCurrentSigmaY2()<0.0001) ){
- // t.fCurrentSigmaY = GetSigmaY(&t);
- //t.fCurrentSigmaZ = GetSigmaZ(&t);
- GetShape(&t,nr);
- }
- //
- AliTPCclusterMI *cl=0;
- UInt_t index=0;
-
-
- //Int_t nr2 = nr;
- const AliTPCRow &krow=GetRow(t.GetRelativeSector(),nr);
- if ( (t.GetSigmaY2()<0) || t.GetSigmaZ2()<0) return 0;
- Double_t roady =1.;
- Double_t roadz = 1.;
- //
- Int_t row = nr;
- if (TMath::Abs(TMath::Abs(y)-ymax)<krow.GetDeadZone()){
- t.SetInDead(kTRUE);
- t.SetClusterIndex2(row,-1);
- return 0;
- }
- else
- {
- if (TMath::Abs(z)>(AliTPCReconstructor::GetCtgRange()*x+10)) t.SetClusterIndex2(row,-1);
- }
- //calculate
-
- if ((cl==0)&&(krow)) {
- // cl = krow.FindNearest2(y+10,z,roady,roadz,index);
- cl = krow.FindNearest2(y,z,roady,roadz,index);
-
- if (cl) t.SetCurrentClusterIndex1(krow.GetIndex(index));
- }
-
- if (cl) {
- t.SetCurrentCluster(cl);
- // Int_t accept = AcceptCluster(&t,t.fCurrentCluster,1.);
- //if (accept<3){
- t.SetClusterIndex2(row,index);
- t.SetClusterPointer(row, cl);
- //}
- }
- return 1;
-}
-
//_________________________________________________________________________
Int_t sector = (index&0xff000000)>>24;
// Int_t row = (index&0x00ff0000)>>16;
Float_t xyz[3];
- // xyz[0] = fParam->GetPadRowRadii(sector,row);
+ // xyz[0] = fkParam->GetPadRowRadii(sector,row);
xyz[0] = cl->GetX();
xyz[1] = cl->GetY();
xyz[2] = cl->GetZ();
Float_t sin,cos;
- fParam->AdjustCosSin(sector,cos,sin);
+ fkParam->AdjustCosSin(sector,cos,sin);
Float_t x = cos*xyz[0]-sin*xyz[1];
Float_t y = cos*xyz[1]+sin*xyz[0];
Float_t cov[6];
Float_t sigmaY2 = 0.027*cl->GetSigmaY2();
- if (sector < fParam->GetNInnerSector()) sigmaY2 *= 2.07;
+ if (sector < fkParam->GetNInnerSector()) sigmaY2 *= 2.07;
Float_t sigmaZ2 = 0.066*cl->GetSigmaZ2();
- if (sector < fParam->GetNInnerSector()) sigmaZ2 *= 1.77;
+ if (sector < fkParam->GetNInnerSector()) sigmaZ2 *= 1.77;
cov[0] = sin*sin*sigmaY2;
cov[1] = -sin*cos*sigmaY2;
cov[2] = 0.;
cov[4] = 0.;
cov[5] = sigmaZ2;
p.SetXYZ(x,y,xyz[2],cov);
- AliAlignObj::ELayerID iLayer;
+ AliGeomManager::ELayerID iLayer;
Int_t idet;
- if (sector < fParam->GetNInnerSector()) {
- iLayer = AliAlignObj::kTPC1;
+ if (sector < fkParam->GetNInnerSector()) {
+ iLayer = AliGeomManager::kTPC1;
idet = sector;
}
else {
- iLayer = AliAlignObj::kTPC2;
- idet = sector - fParam->GetNInnerSector();
+ iLayer = AliGeomManager::kTPC2;
+ idet = sector - fkParam->GetNInnerSector();
}
- UShort_t volid = AliAlignObj::LayerToVolUID(iLayer,idet);
+ UShort_t volid = AliGeomManager::LayerToVolUID(iLayer,idet);
p.SetVolumeID(volid);
return kTRUE;
}
}
//AliInfo(Form("A - Sector%d phi %f - alpha %f", t.fRelativeSector,y/x, t.GetAlpha()));
- AliTPCRow &krow=GetRow(t.GetRelativeSector(),nr);
+ AliTPCtrackerRow &krow=GetRow(t.GetRelativeSector(),nr);
if (TMath::Abs(TMath::Abs(y)-ymax)<krow.GetDeadZone()){
t.SetInDead(kTRUE);
}
else
{
- if (TMath::Abs(t.GetZ())<(AliTPCReconstructor::GetCtgRange()*t.GetX()+10) && (TMath::Abs(t.GetSnp())<AliTPCReconstructor::GetMaxSnpTracker())) t.SetNFoundable(t.GetNFoundable()+1);
+
+ // if (TMath::Abs(t.GetZ())<(AliTPCReconstructor::GetCtgRange()*t.GetX()+10) && (TMath::Abs(t.GetSnp())<AliTPCReconstructor::GetMaxSnpTracker()))
+ if (IsFindable(t)) t.SetNFoundable(t.GetNFoundable()+1);
else
return 0;
}
// update current
- if ( (nr%6==0) || t.GetNumberOfClusters()<2){
+ if ( (nr%2==0) || t.GetNumberOfClusters()<2){
// t.fCurrentSigmaY = GetSigmaY(&t);
//t.fCurrentSigmaZ = GetSigmaZ(&t);
GetShape(&t,nr);
if (t.GetCurrentCluster()) {
t.SetRow(nr);
- Int_t accept = AcceptCluster(&t,t.GetCurrentCluster(),1.);
+ Int_t accept = AcceptCluster(&t,t.GetCurrentCluster());
if (t.GetCurrentCluster()->IsUsed(10)){
//
}
-//_____________________________________________________________________________
-Int_t AliTPCtrackerMI::FollowProlongationFast(AliTPCseed& t, Int_t rf, Int_t step) {
- //-----------------------------------------------------------------
- // This function tries to find a track prolongation.
- //-----------------------------------------------------------------
- Double_t xt=t.GetX();
- //
- Double_t alpha=t.GetAlpha() - fSectors->GetAlphaShift();
- if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
- if (alpha < 0. ) alpha += 2.*TMath::Pi();
- t.SetRelativeSector(Int_t(alpha/fSectors->GetAlpha()+0.0001)%fN);
-
- for (Int_t nr=GetRowNumber(xt)-1; nr>=rf; nr-=step) {
-
- if (FollowToNextFast(t,nr)==0)
- if (!t.IsActive()) return 0;
-
- }
- return 1;
-}
-
{
//
//
- if (TMath::Abs(s1->GetC()-s2->GetC())>0.004) return;
- if (TMath::Abs(s1->GetTgl()-s2->GetTgl())>0.6) return;
-
- Float_t dz2 =(s1->GetZ() - s2->GetZ());
- dz2*=dz2;
- Float_t dy2 =(s1->GetY() - s2->GetY());
- dy2*=dy2;
- Float_t distance = dz2+dy2;
- if (distance>325.) return ; // if there are far away - not overlap - to reduce combinatorics
+ Float_t thetaCut = 0.2;//+10.*TMath::Sqrt(s1->GetSigmaTglZ()+ s2->GetSigmaTglZ());
+ if (TMath::Abs(s1->GetTgl()-s2->GetTgl())>thetaCut) return;
+ Float_t minCl = TMath::Min(s1->GetNumberOfClusters(),s2->GetNumberOfClusters());
+ Int_t cutN0 = TMath::Max(5,TMath::Nint(0.1*minCl));
//
Int_t sumshared=0;
//
- Int_t firstpoint = TMath::Max(s1->GetFirstPoint(),s2->GetFirstPoint());
- Int_t lastpoint = TMath::Min(s1->GetLastPoint(),s2->GetLastPoint());
+ //Int_t firstpoint = TMath::Max(s1->GetFirstPoint(),s2->GetFirstPoint());
+ //Int_t lastpoint = TMath::Min(s1->GetLastPoint(),s2->GetLastPoint());
+ Int_t firstpoint = 0;
+ Int_t lastpoint = 160;
//
- if (firstpoint>=lastpoint-5) return;;
+ // if (firstpoint>=lastpoint-5) return;;
for (Int_t i=firstpoint;i<lastpoint;i++){
// if ( (s1->GetClusterIndex2(i)&0xFFFF8FFF)==(s2->GetClusterIndex2(i)&0xFFFF8FFF) && s1->GetClusterIndex2(i)>0) {
if ( (s1->GetClusterIndex2(i))==(s2->GetClusterIndex2(i)) && s1->GetClusterIndex2(i)>0) {
sumshared++;
+ s1->SetSharedMapBit(i, kTRUE);
+ s2->SetSharedMapBit(i, kTRUE);
}
+ if (s1->GetClusterIndex2(i)>0)
+ s1->SetClusterMapBit(i, kTRUE);
}
- if (sumshared>4){
+ if (sumshared>cutN0){
// sign clusters
//
for (Int_t i=firstpoint;i<lastpoint;i++){
}
}
//
- if (sumshared>10){
+ if (sumshared>cutN0){
for (Int_t i=0;i<4;i++){
if (s1->GetOverlapLabel(3*i)==0){
s1->SetOverlapLabel(3*i, s2->GetLabel());
}
}
}
-
}
void AliTPCtrackerMI::SignShared(TObjArray * arr)
pt->SetSort(0);
}
arr->UnSort();
- arr->Sort(); // sorting according z
+ arr->Sort(); // sorting according relative sectors
arr->Expand(arr->GetEntries());
//
//
if (pt->GetRemoval()>10) continue;
for (Int_t j=i+1; j<nseed; j++){
AliTPCseed *pt2=(AliTPCseed*)arr->UncheckedAt(j);
+ if (TMath::Abs(pt->GetRelativeSector()-pt2->GetRelativeSector())>1) continue;
// if (pt2){
if (pt2->GetRemoval()<=10) {
- if ( TMath::Abs(pt->GetRelativeSector()-pt2->GetRelativeSector())>0) break;
+ //if ( TMath::Abs(pt->GetRelativeSector()-pt2->GetRelativeSector())>0) break;
SignShared(pt,pt2);
}
}
}
}
-void AliTPCtrackerMI::RemoveDouble(TObjArray * arr, Float_t factor1, Float_t factor2, Int_t removalindex)
-{
- //
- //sort trackss according sectors
- //
- if (fDebug&1) {
- Info("RemoveDouble","Number of tracks before double removal- %d\n",arr->GetEntries());
- }
- //
- for (Int_t i=0; i<arr->GetEntriesFast(); i++) {
- AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
- if (!pt) continue;
- pt->SetSort(0);
- }
- arr->UnSort();
- arr->Sort(); // sorting according z
- arr->Expand(arr->GetEntries());
- //
- //reset overlap labels
- //
- Int_t nseed=arr->GetEntriesFast();
- for (Int_t i=0; i<nseed; i++) {
- AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
- if (!pt) continue;
- pt->SetUniqueID(i);
- for (Int_t j=0;j<=12;j++){
- pt->SetOverlapLabel(j,0);
- }
- }
- //
- //sign shared tracks
- for (Int_t i=0; i<nseed; i++) {
- AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
- if (!pt) continue;
- if (pt->GetRemoval()>10) continue;
- Float_t deltac = pt->GetC()*0.1;
- for (Int_t j=i+1; j<nseed; j++){
- AliTPCseed *pt2=(AliTPCseed*)arr->UncheckedAt(j);
- // if (pt2){
- if (pt2->GetRemoval()<=10) {
- if ( TMath::Abs(pt->GetRelativeSector()-pt2->GetRelativeSector())>0) break;
- if (TMath::Abs(pt->GetC() -pt2->GetC())>deltac) continue;
- if (TMath::Abs(pt->GetTgl()-pt2->GetTgl())>0.05) continue;
- //
- SignShared(pt,pt2);
- }
- }
- }
- //
- // remove highly shared tracks
- for (Int_t i=0; i<nseed; i++) {
- AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
- if (!pt) continue;
- if (pt->GetRemoval()>10) continue;
- //
- Int_t sumshared =0;
- for (Int_t j=0;j<4;j++){
- sumshared = pt->GetOverlapLabel(j*3+1);
- }
- Float_t factor = factor1;
- if (pt->GetRemoval()>0) factor = factor2;
- if (sumshared/pt->GetNumberOfClusters()>factor){
- for (Int_t j=0;j<4;j++){
- if (pt->GetOverlapLabel(3*j)==0) continue;
- if (pt->GetOverlapLabel(3*j+1)<5) continue;
- if (pt->GetRemoval()==removalindex) continue;
- AliTPCseed * pt2 = (AliTPCseed*)arr->UncheckedAt(pt->GetOverlapLabel(3*j+2));
- if (!pt2) continue;
- if (pt2->GetSigma2C()<pt->GetSigma2C()){
- // pt->fRemoval = removalindex;
- delete arr->RemoveAt(i);
- break;
- }
- }
- }
- }
- arr->Compress();
- if (fDebug&1) {
- Info("RemoveDouble","Number of tracks after double removal- %d\n",arr->GetEntries());
- }
-}
-
-
-
-
-
void AliTPCtrackerMI::SortTracks(TObjArray * arr, Int_t mode) const
{
arr->Sort();
}
-void AliTPCtrackerMI::RemoveUsed(TObjArray * arr, Float_t factor1, Float_t factor2, Int_t removalindex)
-{
-
- //Loop over all tracks and remove "overlaps"
- //
- //
- Int_t nseed = arr->GetEntriesFast();
- Int_t good =0;
-
- for (Int_t i=0; i<nseed; i++) {
- AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
- if (!pt) {
- delete arr->RemoveAt(i);
- }
- else{
- pt->SetSort(1);
- pt->SetBSigned(kFALSE);
- }
- }
- arr->Compress();
- nseed = arr->GetEntriesFast();
- arr->UnSort();
- arr->Sort();
- //
- //unsign used
- UnsignClusters();
- //
- for (Int_t i=0; i<nseed; i++) {
- AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
- if (!pt) {
- continue;
- }
- Int_t found,foundable,shared;
- if (pt->IsActive())
- pt->GetClusterStatistic(0,160,found, foundable,shared,kFALSE);
- else
- pt->GetClusterStatistic(0,160,found, foundable,shared,kTRUE);
- //
- Double_t factor = factor2;
- if (pt->GetBConstrain()) factor = factor1;
-
- if ((Float_t(shared)/Float_t(found))>factor){
- pt->Desactivate(removalindex);
- continue;
- }
-
- good++;
- for (Int_t i=0; i<160; i++) {
- Int_t index=pt->GetClusterIndex2(i);
- if (index<0 || index&0x8000 ) continue;
- AliTPCclusterMI *c= pt->GetClusterPointer(i);
- if (!c) continue;
- // if (!c->IsUsed(10)) c->Use(10);
- //if (pt->IsActive())
- c->Use(10);
- //else
- // c->Use(5);
- }
-
- }
- fNtracks = good;
- if (fDebug>0){
- Info("RemoveUsed","\n*****\nNumber of good tracks after shared removal\t%d\n",fNtracks);
- }
-}
-
void AliTPCtrackerMI::RemoveUsed2(TObjArray * arr, Float_t factor1, Float_t factor2, Int_t minimal)
{
-
- //Loop over all tracks and remove "overlaps"
//
+ // Loop over all tracks and remove overlaped tracks (with lower quality)
+ // Algorithm:
+ // 1. Unsign clusters
+ // 2. Sort tracks according quality
+ // Quality is defined by the number of cluster between first and last points
+ //
+ // 3. Loop over tracks - decreasing quality order
+ // a.) remove - If the fraction of shared cluster less than factor (1- n or 2)
+ // b.) remove - If the minimal number of clusters less than minimal and not ITS
+ // c.) if track accepted - sign clusters
+ //
+ //Called in - AliTPCtrackerMI::Clusters2Tracks()
+ // - AliTPCtrackerMI::PropagateBack()
+ // - AliTPCtrackerMI::RefitInward()
+ //
+ // Arguments:
+ // factor1 - factor for constrained
+ // factor2 - for non constrained tracks
+ // if (Float_t(shared+1)/Float_t(found+1)>factor) - DELETE
//
UnsignClusters();
//
pt->UpdatePoints(); //select first last max dens points
Float_t * points = pt->GetPoints();
if (points[3]<0.8) quality[i] =-1;
- //
quality[i] = (points[2]-points[0])+pt->GetNumberOfClusters();
+ //prefer high momenta tracks if overlaps
+ quality[i] *= TMath::Sqrt(TMath::Abs(pt->Pt())+0.5);
}
TMath::Sort(nseed,quality,indexes);
//
//
for (Int_t itrack=0; itrack<nseed; itrack++) {
Int_t trackindex = indexes[itrack];
- AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(trackindex);
+ AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(trackindex);
+ if (!pt) continue;
+ //
if (quality[trackindex]<0){
if (pt) {
delete arr->RemoveAt(trackindex);
continue;
}
//
+ //
Int_t first = Int_t(pt->GetPoints()[0]);
Int_t last = Int_t(pt->GetPoints()[2]);
Double_t factor = (pt->GetBConstrain()) ? factor1: factor2;
//
Int_t found,foundable,shared;
- pt->GetClusterStatistic(first,last, found, foundable,shared,kFALSE);
- Float_t sharedfactor = Float_t(shared+1)/Float_t(found+1);
+ pt->GetClusterStatistic(first,last, found, foundable,shared,kFALSE); // better to get statistic in "high-dens" region do't use full track as in line bellow
+ // pt->GetClusterStatistic(0,160, found, foundable,shared,kFALSE);
Bool_t itsgold =kFALSE;
if (pt->GetESD()){
Int_t dummy[12];
//
if (Float_t(shared+1)/Float_t(found+1)>factor){
if (pt->GetKinkIndexes()[0]!=0) continue; //don't remove tracks - part of the kinks
+ if( AliTPCReconstructor::StreamLevel()>15){
+ TTreeSRedirector &cstream = *fDebugStreamer;
+ cstream<<"RemoveUsed"<<
+ "iter="<<fIteration<<
+ "pt.="<<pt<<
+ "\n";
+ }
delete arr->RemoveAt(trackindex);
continue;
}
if (pt->GetNumberOfClusters()<50&&(found-0.5*shared)<minimal){ //remove short tracks
if (pt->GetKinkIndexes()[0]!=0) continue; //don't remove tracks - part of the kinks
+ if( AliTPCReconstructor::StreamLevel()>15){
+ TTreeSRedirector &cstream = *fDebugStreamer;
+ cstream<<"RemoveShort"<<
+ "iter="<<fIteration<<
+ "pt.="<<pt<<
+ "\n";
+ }
delete arr->RemoveAt(trackindex);
continue;
}
}
good++;
- if (sharedfactor>0.4) continue;
+ //if (sharedfactor>0.4) continue;
if (pt->GetKinkIndexes()[0]>0) continue;
+ //Remove tracks with undefined properties - seems
+ if (pt->GetSigmaY2()<kAlmost0) continue; // ? what is the origin ?
+ //
for (Int_t i=first; i<last; i++) {
Int_t index=pt->GetClusterIndex2(i);
// if (index<0 || index&0x8000 ) continue;
-void AliTPCtrackerMI::SignClusters(TObjArray * arr, Float_t fnumber, Float_t fdensity)
+void AliTPCtrackerMI::SignClusters(const TObjArray * arr, Float_t fnumber, Float_t fdensity)
{
//
//sign clusters to be "used"
meanchi = sumchi/sum;
//
sdensity = sumdens2/sum-mdensity*mdensity;
- sdensity = TMath::Sqrt(sdensity);
+ if (sdensity >= 0)
+ sdensity = TMath::Sqrt(sdensity);
+ else
+ sdensity = 0.1;
//
smeann = sumn2/sum-meann*meann;
- smeann = TMath::Sqrt(smeann);
+ if (smeann >= 0)
+ smeann = TMath::Sqrt(smeann);
+ else
+ smeann = 10;
//
smeanchi = sumchi2/sum - meanchi*meanchi;
- smeanchi = TMath::Sqrt(smeanchi);
+ if (smeanchi >= 0)
+ smeanchi = TMath::Sqrt(smeanchi);
+ else
+ smeanchi = 0.4;
}
isok =kTRUE;
if (isok)
- for (Int_t i=0; i<160; i++) {
- Int_t index=pt->GetClusterIndex2(i);
+ for (Int_t j=0; j<160; ++j) {
+ Int_t index=pt->GetClusterIndex2(j);
if (index<0) continue;
- AliTPCclusterMI *c= pt->GetClusterPointer(i);
+ AliTPCclusterMI *c= pt->GetClusterPointer(j);
if (!c) continue;
//if (!(c->IsUsed(10))) c->Use();
c->Use(10);
if ((dens>mindens && pt->GetNumberOfClusters()>minn) && chi<maxchi ){
//Int_t noc=pt->GetNumberOfClusters();
pt->SetBSigned(kTRUE);
- for (Int_t i=0; i<160; i++) {
+ for (Int_t j=0; j<160; ++j) {
- Int_t index=pt->GetClusterIndex2(i);
+ Int_t index=pt->GetClusterIndex2(j);
if (index<0) continue;
- AliTPCclusterMI *c= pt->GetClusterPointer(i);
+ AliTPCclusterMI *c= pt->GetClusterPointer(j);
if (!c) continue;
// if (!(c->IsUsed(10))) c->Use();
c->Use(10);
}
-void AliTPCtrackerMI::StopNotActive(TObjArray * arr, Int_t row0, Float_t th0, Float_t th1, Float_t th2) const
+void AliTPCtrackerMI::StopNotActive(const TObjArray * arr, Int_t row0, Float_t th0, Float_t th1, Float_t th2) const
{
// stop not active tracks
// take th1 as threshold for number of founded to number of foundable on last 10 active rows
}
-Int_t AliTPCtrackerMI::RefitInward(AliESD *event)
+Int_t AliTPCtrackerMI::RefitInward(AliESDEvent *event)
{
//
// back propagation of ESD tracks
//
//return 0;
+ const Int_t kMaxFriendTracks=2000;
fEvent = event;
ReadSeeds(event,2);
fIteration=2;
//PrepareForProlongation(fSeeds,1);
PropagateForward2(fSeeds);
+ RemoveUsed2(fSeeds,0.4,0.4,20);
+
+ TObjArray arraySeed(fSeeds->GetEntries());
+ for (Int_t i=0;i<fSeeds->GetEntries();i++) {
+ arraySeed.AddAt(fSeeds->At(i),i);
+ }
+ SignShared(&arraySeed);
+ // FindCurling(fSeeds, event,2); // find multi found tracks
+ FindSplitted(fSeeds, event,2); // find multi found tracks
+ if (AliTPCReconstructor::StreamLevel()>2) FindMultiMC(fSeeds, fEvent,2); // find multi found tracks
Int_t ntracks=0;
Int_t nseed = fSeeds->GetEntriesFast();
AliTPCseed * seed = (AliTPCseed*) fSeeds->UncheckedAt(i);
if (!seed) continue;
if (seed->GetKinkIndex(0)>0) UpdateKinkQualityD(seed); // update quality informations for kinks
+ AliESDtrack *esd=event->GetTrack(i);
+
+ if (seed->GetNumberOfClusters()<60 && seed->GetNumberOfClusters()<(esd->GetTPCclusters(0) -5)*0.8){
+ AliExternalTrackParam paramIn;
+ AliExternalTrackParam paramOut;
+ Int_t ncl = seed->RefitTrack(seed,¶mIn,¶mOut);
+ if (AliTPCReconstructor::StreamLevel()>0) {
+ (*fDebugStreamer)<<"RecoverIn"<<
+ "seed.="<<seed<<
+ "esd.="<<esd<<
+ "pin.="<<¶mIn<<
+ "pout.="<<¶mOut<<
+ "ncl="<<ncl<<
+ "\n";
+ }
+ if (ncl>15) {
+ seed->Set(paramIn.GetX(),paramIn.GetAlpha(),paramIn.GetParameter(),paramIn.GetCovariance());
+ seed->SetNumberOfClusters(ncl);
+ }
+ }
- seed->PropagateTo(fParam->GetInnerRadiusLow());
+ seed->PropagateTo(fkParam->GetInnerRadiusLow());
seed->UpdatePoints();
- AliESDtrack *esd=event->GetTrack(i);
+ AddCovariance(seed);
+ MakeBitmaps(seed);
seed->CookdEdx(0.02,0.6);
CookLabel(seed,0.1); //For comparison only
+ esd->SetTPCClusterMap(seed->GetClusterMap());
+ esd->SetTPCSharedMap(seed->GetSharedMap());
//
- if (AliTPCReconstructor::StreamLevel()>0 && seed!=0&&esd!=0) {
+ if (AliTPCReconstructor::StreamLevel()>1 && seed!=0&&esd!=0) {
TTreeSRedirector &cstream = *fDebugStreamer;
cstream<<"Crefit"<<
"Esd.="<<esd<<
esd->UpdateTrackParams(seed,AliESDtrack::kTPCrefit);
esd->SetTPCPoints(seed->GetPoints());
esd->SetTPCPointsF(seed->GetNFoundable());
- Int_t ndedx = seed->GetNCDEDX(0)+seed->GetNCDEDX(1)+seed->GetNCDEDX(2)+seed->GetNCDEDX(3);
- Float_t sdedx = (seed->GetSDEDX(0)+seed->GetSDEDX(1)+seed->GetSDEDX(2)+seed->GetSDEDX(3))*0.25;
+ Int_t ndedx = seed->GetNCDEDX(0);
+ Float_t sdedx = seed->GetSDEDX(0);
Float_t dedx = seed->GetdEdx();
esd->SetTPCsignal(dedx, sdedx, ndedx);
//
// add seed to the esd track in Calib level
//
- if (AliTPCReconstructor::StreamLevel()>0){
+ Bool_t storeFriend = gRandom->Rndm()<(kMaxFriendTracks)/Float_t(nseed);
+ if (AliTPCReconstructor::StreamLevel()>0 &&storeFriend){
AliTPCseed * seedCopy = new AliTPCseed(*seed, kTRUE);
esd->AddCalibObject(seedCopy);
}
}
//FindKinks(fSeeds,event);
Info("RefitInward","Number of refitted tracks %d",ntracks);
- fEvent =0;
- //WriteTracks();
return 0;
}
-Int_t AliTPCtrackerMI::PropagateBack(AliESD *event)
+Int_t AliTPCtrackerMI::PropagateBack(AliESDEvent *event)
{
//
// back propagation of ESD tracks
ReadSeeds(event,1);
PropagateBack(fSeeds);
RemoveUsed2(fSeeds,0.4,0.4,20);
+ //FindCurling(fSeeds, fEvent,1);
+ FindSplitted(fSeeds, event,1); // find multi found tracks
+ if (AliTPCReconstructor::StreamLevel()>2) FindMultiMC(fSeeds, fEvent,1); // find multi found tracks
+
//
Int_t nseed = fSeeds->GetEntriesFast();
Int_t ntracks=0;
if (!seed) continue;
if (seed->GetKinkIndex(0)<0) UpdateKinkQualityM(seed); // update quality informations for kinks
seed->UpdatePoints();
+ AddCovariance(seed);
AliESDtrack *esd=event->GetTrack(i);
+ if (seed->GetNumberOfClusters()<60 && seed->GetNumberOfClusters()<(esd->GetTPCclusters(0) -5)*0.8){
+ AliExternalTrackParam paramIn;
+ AliExternalTrackParam paramOut;
+ Int_t ncl = seed->RefitTrack(seed,¶mIn,¶mOut);
+ if (AliTPCReconstructor::StreamLevel()>0) {
+ (*fDebugStreamer)<<"RecoverBack"<<
+ "seed.="<<seed<<
+ "esd.="<<esd<<
+ "pin.="<<¶mIn<<
+ "pout.="<<¶mOut<<
+ "ncl="<<ncl<<
+ "\n";
+ }
+ if (ncl>15) {
+ seed->Set(paramOut.GetX(),paramOut.GetAlpha(),paramOut.GetParameter(),paramOut.GetCovariance());
+ seed->SetNumberOfClusters(ncl);
+ }
+ }
seed->CookdEdx(0.02,0.6);
CookLabel(seed,0.1); //For comparison only
if (seed->GetNumberOfClusters()>15){
esd->UpdateTrackParams(seed,AliESDtrack::kTPCout);
esd->SetTPCPoints(seed->GetPoints());
esd->SetTPCPointsF(seed->GetNFoundable());
- Int_t ndedx = seed->GetNCDEDX(0)+seed->GetNCDEDX(1)+seed->GetNCDEDX(2)+seed->GetNCDEDX(3);
- Float_t sdedx = (seed->GetSDEDX(0)+seed->GetSDEDX(1)+seed->GetSDEDX(2)+seed->GetSDEDX(3))*0.25;
+ Int_t ndedx = seed->GetNCDEDX(0);
+ Float_t sdedx = seed->GetSDEDX(0);
Float_t dedx = seed->GetdEdx();
esd->SetTPCsignal(dedx, sdedx, ndedx);
ntracks++;
Int_t eventnumber = event->GetEventNumberInFile();// patch 28 fev 06
// This is most likely NOT the event number you'd like to use. It has nothing to do with the 'real' event number
- (*fDebugStreamer)<<"Cback"<<
- "Tr0.="<<seed<<
- "EventNrInFile="<<eventnumber<<
- "\n"; // patch 28 fev 06
+ if (AliTPCReconstructor::StreamLevel()>1 && esd) {
+ (*fDebugStreamer)<<"Cback"<<
+ "Tr0.="<<seed<<
+ "esd.="<<esd<<
+ "EventNrInFile="<<eventnumber<<
+ "\n";
+ }
}
}
//FindKinks(fSeeds,event);
Info("PropagateBack","Number of back propagated tracks %d",ntracks);
fEvent =0;
- //WriteTracks();
-
+
return 0;
}
fSeeds =0;
}
-void AliTPCtrackerMI::ReadSeeds(AliESD *event, Int_t direction)
+void AliTPCtrackerMI::ReadSeeds(AliESDEvent *event, Int_t direction)
{
//
//read seeds from the event
t.SetNumberOfClusters(0);
// AliTPCseed *seed = new AliTPCseed(t,t.GetAlpha());
AliTPCseed *seed = new AliTPCseed(t/*,t.GetAlpha()*/);
+ seed->SetUniqueID(esd->GetID());
+ AddCovariance(seed); //add systematic ucertainty
for (Int_t ikink=0;ikink<3;ikink++) {
Int_t index = esd->GetKinkIndex(ikink);
seed->GetKinkIndexes()[ikink] = index;
}
if (((status&AliESDtrack::kITSout)==0)&&(direction==1)) seed->ResetCovariance(10.);
if ( direction ==2 &&(status & AliESDtrack::kTRDrefit) == 0 ) seed->ResetCovariance(10.);
- if ( direction ==2 && ((status & AliESDtrack::kTPCout) == 0) ) {
- fSeeds->AddAt(0,i);
- delete seed;
- continue;
- }
+ //if ( direction ==2 && ((status & AliESDtrack::kTPCout) == 0) ) {
+ // fSeeds->AddAt(0,i);
+ // delete seed;
+ // continue;
+ //}
if ( direction ==2 &&(status & AliESDtrack::kTRDrefit) > 0 ) {
Double_t par0[5],par1[5],alpha,x;
esd->GetInnerExternalParameters(alpha,x,par0);
Double_t x[5], c[15];
// Int_t di = i1-i2;
//
- AliTPCseed * seed = new AliTPCseed;
+ AliTPCseed * seed = new AliTPCseed();
Double_t alpha=fSectors->GetAlpha(), shift=fSectors->GetAlphaShift();
Double_t cs=cos(alpha), sn=sin(alpha);
//
Int_t imiddle = (i2+i1)/2; //middle pad row index
Double_t xm = GetXrow(imiddle); // radius of middle pad-row
- const AliTPCRow& krm=GetRow(sec,imiddle); //middle pad -row
+ const AliTPCtrackerRow& krm=GetRow(sec,imiddle); //middle pad -row
//
Int_t ns =sec;
- const AliTPCRow& kr1=GetRow(ns,i1);
+ const AliTPCtrackerRow& kr1=GetRow(ns,i1);
Double_t ymax = GetMaxY(i1)-kr1.GetDeadZone()-1.5;
Double_t ymaxm = GetMaxY(imiddle)-kr1.GetDeadZone()-1.5;
for (Int_t dsec = dsec1; dsec<=dsec2;dsec++){
Int_t sec2 = sec + dsec;
//
- // AliTPCRow& kr2 = fOuterSec[(sec2+fkNOS)%fkNOS][i2];
- //AliTPCRow& kr2m = fOuterSec[(sec2+fkNOS)%fkNOS][imiddle];
- AliTPCRow& kr2 = GetRow((sec2+fkNOS)%fkNOS,i2);
- AliTPCRow& kr2m = GetRow((sec2+fkNOS)%fkNOS,imiddle);
+ // AliTPCtrackerRow& kr2 = fOuterSec[(sec2+fkNOS)%fkNOS][i2];
+ //AliTPCtrackerRow& kr2m = fOuterSec[(sec2+fkNOS)%fkNOS][imiddle];
+ AliTPCtrackerRow& kr2 = GetRow((sec2+fkNOS)%fkNOS,i2);
+ AliTPCtrackerRow& kr2m = GetRow((sec2+fkNOS)%fkNOS,imiddle);
Int_t index1 = TMath::Max(kr2.Find(extraz-0.6-dddz1*TMath::Abs(z1)*0.05)-1,0);
Int_t index2 = TMath::Min(kr2.Find(extraz+0.6+dddz2*TMath::Abs(z1)*0.05)+1,kr2);
// if (!BuildSeed(kr1[is],kcl,0,x1,x2,x3,x,c)) continue;
UInt_t index=kr1.GetIndex(is);
+ seed->~AliTPCseed(); // this does not set the pointer to 0...
AliTPCseed *track=new(seed) AliTPCseed(x1, ns*alpha+shift, x, c, index);
track->SetIsSeeding(kTRUE);
}
}
}
-
+
track->SetSeedType(0);
arr->AddLast(track);
seed = new AliTPCseed;
// first 3 padrows
Double_t x1 = GetXrow(i1-1);
- const AliTPCRow& kr1=GetRow(sec,i1-1);
+ const AliTPCtrackerRow& kr1=GetRow(sec,i1-1);
Double_t y1max = GetMaxY(i1-1)-kr1.GetDeadZone()-1.5;
//
Double_t x1p = GetXrow(i1);
- const AliTPCRow& kr1p=GetRow(sec,i1);
+ const AliTPCtrackerRow& kr1p=GetRow(sec,i1);
//
Double_t x1m = GetXrow(i1-2);
- const AliTPCRow& kr1m=GetRow(sec,i1-2);
+ const AliTPCtrackerRow& kr1m=GetRow(sec,i1-2);
//
//last 3 padrow for seeding
- AliTPCRow& kr3 = GetRow((sec+fkNOS)%fkNOS,i1-7);
+ AliTPCtrackerRow& kr3 = GetRow((sec+fkNOS)%fkNOS,i1-7);
Double_t x3 = GetXrow(i1-7);
// Double_t y3max= GetMaxY(i1-7)-kr3.fDeadZone-1.5;
//
- AliTPCRow& kr3p = GetRow((sec+fkNOS)%fkNOS,i1-6);
+ AliTPCtrackerRow& kr3p = GetRow((sec+fkNOS)%fkNOS,i1-6);
Double_t x3p = GetXrow(i1-6);
//
- AliTPCRow& kr3m = GetRow((sec+fkNOS)%fkNOS,i1-8);
+ AliTPCtrackerRow& kr3m = GetRow((sec+fkNOS)%fkNOS,i1-8);
Double_t x3m = GetXrow(i1-8);
//
// middle padrow
Int_t im = i1-4; //middle pad row index
Double_t xm = GetXrow(im); // radius of middle pad-row
- const AliTPCRow& krm=GetRow(sec,im); //middle pad -row
+ const AliTPCtrackerRow& krm=GetRow(sec,im); //middle pad -row
// Double_t ymmax = GetMaxY(im)-kr1.fDeadZone-1.5;
//
//
// if (!BuildSeed(kr1[is],kcl,0,x1,x2,x3,x,c)) continue;
- UInt_t index=kr1.GetIndex(is);
+ index=kr1.GetIndex(is);
+ seed->~AliTPCseed();
AliTPCseed *track=new(seed) AliTPCseed(x1, sec*alpha+shift, x, c, index);
track->SetIsSeeding(kTRUE);
// Double_t cs=cos(alpha), sn=sin(alpha);
Int_t row0 = (i1+i2)/2;
Int_t drow = (i1-i2)/2;
- const AliTPCRow& kr0=fSectors[sec][row0];
- AliTPCRow * kr=0;
+ const AliTPCtrackerRow& kr0=fSectors[sec][row0];
+ AliTPCtrackerRow * kr=0;
AliTPCpolyTrack polytrack;
Int_t nclusters=fSectors[sec][row0];
Int_t nfound =0;
Int_t nfoundable =0;
for (Int_t iter =1; iter<2; iter++){ //iterations
- const AliTPCRow& krm=fSectors[sec][row0-iter];
- const AliTPCRow& krp=fSectors[sec][row0+iter];
+ const AliTPCtrackerRow& krm=fSectors[sec][row0-iter];
+ const AliTPCtrackerRow& krp=fSectors[sec][row0+iter];
const AliTPCclusterMI * cl= kr0[is];
if (cl->IsUsed(10)) {
Int_t row = row0+ddrow*delta;
kr = &(fSectors[sec][row]);
Double_t xn = kr->GetX();
- Double_t ymax = fSectors->GetMaxY(row)-kr->GetDeadZone()-1.5;
+ Double_t ymax1 = fSectors->GetMaxY(row)-kr->GetDeadZone()-1.5;
polytrack.GetFitPoint(xn,yn,zn);
- if (TMath::Abs(yn)>ymax) continue;
+ if (TMath::Abs(yn)>ymax1) continue;
nfoundable++;
AliTPCclusterMI * cln = kr->FindNearest(yn,zn,roady,roadz);
if (cln) {
UInt_t index=0;
//kr0.GetIndex(is);
- AliTPCseed *track=new (seed) AliTPCseed(x1,sec*alpha+shift,x,c,index);
+ seed->~AliTPCseed();
+ AliTPCseed *track=new(seed) AliTPCseed(x1,sec*alpha+shift,x,c,index);
track->SetIsSeeding(kTRUE);
Int_t rc=FollowProlongation(*track, i2);
if (constrain) track->SetBConstrain(1);
}
-AliTPCseed *AliTPCtrackerMI::ReSeed(AliTPCseed *track, Float_t r0, Float_t r1, Float_t r2)
+AliTPCseed *AliTPCtrackerMI::ReSeed(const AliTPCseed *track, Float_t r0, Float_t r1, Float_t r2)
{
//
//
//
//
//
- // Calculate seed state vector and covariance matrix
+ // Calculate seed state vector and covariance matrix
+
+ Double_t alpha, cs,sn, xx2,yy2;
+ //
+ alpha = (sec[1]-sec[2])*fSectors->GetAlpha();
+ cs = TMath::Cos(alpha);
+ sn = TMath::Sin(alpha);
+ xx2= xyz[1][0]*cs-xyz[1][1]*sn;
+ yy2= xyz[1][0]*sn+xyz[1][1]*cs;
+ xyz[1][0] = xx2;
+ xyz[1][1] = yy2;
+ //
+ alpha = (sec[0]-sec[2])*fSectors->GetAlpha();
+ cs = TMath::Cos(alpha);
+ sn = TMath::Sin(alpha);
+ xx2= xyz[0][0]*cs-xyz[0][1]*sn;
+ yy2= xyz[0][0]*sn+xyz[0][1]*cs;
+ xyz[0][0] = xx2;
+ xyz[0][1] = yy2;
+ //
+ //
+ //
+ Double_t x[5],c[15];
+ //
+ x[0]=xyz[2][1];
+ x[1]=xyz[2][2];
+ x[4]=F1(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1]);
+ x[2]=F2(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1]);
+ x[3]=F3n(xyz[2][0],xyz[2][1],xyz[0][0],xyz[0][1],xyz[2][2],xyz[0][2],x[4]);
+ //
+ Double_t sy =0.1, sz =0.1;
+ //
+ Double_t sy1=0.2, sz1=0.2;
+ Double_t sy2=0.2, sz2=0.2;
+ Double_t sy3=0.2;
+ //
+ Double_t f40=(F1(xyz[2][0],xyz[2][1]+sy,xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1])-x[4])/sy;
+ Double_t f42=(F1(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1]+sy,xyz[0][0],xyz[0][1])-x[4])/sy;
+ Double_t f43=(F1(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1]+sy)-x[4])/sy;
+ Double_t f20=(F2(xyz[2][0],xyz[2][1]+sy,xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1])-x[2])/sy;
+ Double_t f22=(F2(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1]+sy,xyz[0][0],xyz[0][1])-x[2])/sy;
+ Double_t f23=(F2(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1]+sy)-x[2])/sy;
+ //
+ Double_t f30=(F3(xyz[2][0],xyz[2][1]+sy,xyz[0][0],xyz[0][1],xyz[2][2],xyz[0][2])-x[3])/sy;
+ Double_t f31=(F3(xyz[2][0],xyz[2][1],xyz[0][0],xyz[0][1],xyz[2][2]+sz,xyz[0][2])-x[3])/sz;
+ Double_t f32=(F3(xyz[2][0],xyz[2][1],xyz[0][0],xyz[0][1]+sy,xyz[2][2],xyz[0][2])-x[3])/sy;
+ Double_t f34=(F3(xyz[2][0],xyz[2][1],xyz[0][0],xyz[0][1],xyz[2][2],xyz[0][2]+sz)-x[3])/sz;
+
+
+ c[0]=sy1;
+ c[1]=0.; c[2]=sz1;
+ c[3]=f20*sy1; c[4]=0.; c[5]=f20*sy1*f20+f22*sy2*f22+f23*sy3*f23;
+ c[6]=f30*sy1; c[7]=f31*sz1; c[8]=f30*sy1*f20+f32*sy2*f22;
+ c[9]=f30*sy1*f30+f31*sz1*f31+f32*sy2*f32+f34*sz2*f34;
+ c[10]=f40*sy1; c[11]=0.; c[12]=f40*sy1*f20+f42*sy2*f22+f43*sy3*f23;
+ c[13]=f30*sy1*f40+f32*sy2*f42;
+ c[14]=f40*sy1*f40+f42*sy2*f42+f43*sy3*f43;
+
+ // Int_t row1 = fSectors->GetRowNumber(xyz[2][0]);
+ AliTPCseed *seed=new AliTPCseed(xyz[2][0], sec[2]*fSectors->GetAlpha()+fSectors->GetAlphaShift(), x, c, 0);
+ seed->SetLastPoint(row[2]);
+ seed->SetFirstPoint(row[2]);
+ for (Int_t i=row[0];i<row[2];i++){
+ seed->SetClusterIndex(i, track->GetClusterIndex(i));
+ }
+
+ return seed;
+}
+
+
+
+void AliTPCtrackerMI::FindMultiMC(const TObjArray * array, AliESDEvent */*esd*/, Int_t iter)
+{
+ //
+ // find multi tracks - THIS FUNCTION IS ONLY FOR DEBUG PURPOSES
+ // USES MC LABELS
+ // Use AliTPCReconstructor::StreamLevel()>2 if you want to tune parameters - cuts
+ //
+ // Two reasons to have multiple find tracks
+ // 1. Curling tracks can be find more than once
+ // 2. Splitted tracks
+ // a.) Multiple seeding to increase tracking efficiency - (~ 100% reached)
+ // b.) Edge effect on the sector boundaries
+ //
+ //
+ // Algorithm done in 2 phases - because of CPU consumption
+ // it is n^2 algorithm - for lead-lead 20000x20000 combination are investigated
+ //
+ // Algorihm for curling tracks sign:
+ // 1 phase -makes a very rough fast cuts to minimize combinatorics
+ // a.) opposite sign
+ // b.) one of the tracks - not pointing to the primary vertex -
+ // c.) delta tan(theta)
+ // d.) delta phi
+ // 2 phase - calculates DCA between tracks - time consument
+
+ //
+ // fast cuts
+ //
+ // General cuts - for splitted tracks and for curling tracks
+ //
+ const Float_t kMaxdPhi = 0.2; // maximal distance in phi
+ //
+ // Curling tracks cuts
+ //
+ //
+ //
+ //
+ Int_t nentries = array->GetEntriesFast();
+ AliHelix *helixes = new AliHelix[nentries];
+ Float_t *xm = new Float_t[nentries];
+ Float_t *dz0 = new Float_t[nentries];
+ Float_t *dz1 = new Float_t[nentries];
+ //
+ //
+ TStopwatch timer;
+ timer.Start();
+ //
+ // Find track COG in x direction - point with best defined parameters
+ //
+ for (Int_t i=0;i<nentries;i++){
+ AliTPCseed* track = (AliTPCseed*)array->At(i);
+ if (!track) continue;
+ track->SetCircular(0);
+ new (&helixes[i]) AliHelix(*track);
+ Int_t ncl=0;
+ xm[i]=0;
+ Float_t dz[2];
+ track->GetDZ(GetX(),GetY(),GetZ(),GetBz(),dz);
+ dz0[i]=dz[0];
+ dz1[i]=dz[1];
+ for (Int_t icl=0; icl<160; icl++){
+ AliTPCclusterMI * cl = track->GetClusterPointer(icl);
+ if (cl) {
+ xm[i]+=cl->GetX();
+ ncl++;
+ }
+ }
+ if (ncl>0) xm[i]/=Float_t(ncl);
+ }
+ //
+ for (Int_t i0=0;i0<nentries;i0++){
+ AliTPCseed * track0 = (AliTPCseed*)array->At(i0);
+ if (!track0) continue;
+ Float_t xc0 = helixes[i0].GetHelix(6);
+ Float_t yc0 = helixes[i0].GetHelix(7);
+ Float_t r0 = helixes[i0].GetHelix(8);
+ Float_t rc0 = TMath::Sqrt(xc0*xc0+yc0*yc0);
+ Float_t fi0 = TMath::ATan2(yc0,xc0);
+
+ for (Int_t i1=i0+1;i1<nentries;i1++){
+ AliTPCseed * track1 = (AliTPCseed*)array->At(i1);
+ if (!track1) continue;
+ Int_t lab0=track0->GetLabel();
+ Int_t lab1=track1->GetLabel();
+ if (TMath::Abs(lab0)!=TMath::Abs(lab1)) continue;
+ //
+ Float_t xc1 = helixes[i1].GetHelix(6);
+ Float_t yc1 = helixes[i1].GetHelix(7);
+ Float_t r1 = helixes[i1].GetHelix(8);
+ Float_t rc1 = TMath::Sqrt(xc1*xc1+yc1*yc1);
+ Float_t fi1 = TMath::ATan2(yc1,xc1);
+ //
+ Float_t dfi = fi0-fi1;
+ //
+ //
+ if (dfi>1.5*TMath::Pi()) dfi-=TMath::Pi(); // take care about edge effect
+ if (dfi<-1.5*TMath::Pi()) dfi+=TMath::Pi(); //
+ if (TMath::Abs(dfi)>kMaxdPhi&&helixes[i0].GetHelix(4)*helixes[i1].GetHelix(4)<0){
+ //
+ // if short tracks with undefined sign
+ fi1 = -TMath::ATan2(yc1,-xc1);
+ dfi = fi0-fi1;
+ }
+ Float_t dtheta = TMath::Abs(track0->GetTgl()-track1->GetTgl())<TMath::Abs(track0->GetTgl()+track1->GetTgl())? track0->GetTgl()-track1->GetTgl():track0->GetTgl()+track1->GetTgl();
+
+ //
+ // debug stream to tune "fast cuts"
+ //
+ Double_t dist[3]; // distance at X
+ Double_t mdist[3]={0,0,0}; // mean distance X+-40cm
+ track0->GetDistance(track1,0.5*(xm[i0]+xm[i1])-40.,dist,AliTracker::GetBz());
+ for (Int_t i=0;i<3;i++) mdist[i]+=TMath::Abs(dist[i]);
+ track0->GetDistance(track1,0.5*(xm[i0]+xm[i1])+40.,dist,AliTracker::GetBz());
+ for (Int_t i=0;i<3;i++) mdist[i]+=TMath::Abs(dist[i]);
+ track0->GetDistance(track1,0.5*(xm[i0]+xm[i1]),dist,AliTracker::GetBz());
+ for (Int_t i=0;i<3;i++) mdist[i]+=TMath::Abs(dist[i]);
+ for (Int_t i=0;i<3;i++) mdist[i]*=0.33333;
+
+ Float_t sum =0;
+ Float_t sums=0;
+ for (Int_t icl=0; icl<160; icl++){
+ AliTPCclusterMI * cl0 = track0->GetClusterPointer(icl);
+ AliTPCclusterMI * cl1 = track1->GetClusterPointer(icl);
+ if (cl0&&cl1) {
+ sum++;
+ if (cl0==cl1) sums++;
+ }
+ }
+ //
+ if (AliTPCReconstructor::StreamLevel()>0) {
+ TTreeSRedirector &cstream = *fDebugStreamer;
+ cstream<<"Multi"<<
+ "iter="<<iter<<
+ "lab0="<<lab0<<
+ "lab1="<<lab1<<
+ "Tr0.="<<track0<< // seed0
+ "Tr1.="<<track1<< // seed1
+ "h0.="<<&helixes[i0]<<
+ "h1.="<<&helixes[i1]<<
+ //
+ "sum="<<sum<< //the sum of rows with cl in both
+ "sums="<<sums<< //the sum of shared clusters
+ "xm0="<<xm[i0]<< // the center of track
+ "xm1="<<xm[i1]<< // the x center of track
+ // General cut variables
+ "dfi="<<dfi<< // distance in fi angle
+ "dtheta="<<dtheta<< // distance int theta angle
+ //
+ "dz00="<<dz0[i0]<<
+ "dz01="<<dz0[i1]<<
+ "dz10="<<dz1[i1]<<
+ "dz11="<<dz1[i1]<<
+ "dist0="<<dist[0]<< //distance x
+ "dist1="<<dist[1]<< //distance y
+ "dist2="<<dist[2]<< //distance z
+ "mdist0="<<mdist[0]<< //distance x
+ "mdist1="<<mdist[1]<< //distance y
+ "mdist2="<<mdist[2]<< //distance z
+ //
+ "r0="<<r0<<
+ "rc0="<<rc0<<
+ "fi0="<<fi0<<
+ "fi1="<<fi1<<
+ "r1="<<r1<<
+ "rc1="<<rc1<<
+ "\n";
+ }
+ }
+ }
+ delete [] helixes;
+ delete [] xm;
+ if (AliTPCReconstructor::StreamLevel()>1) {
+ AliInfo("Time for curling tracks removal DEBUGGING MC");
+ timer.Print();
+ }
+}
+
+
+void AliTPCtrackerMI::FindSplitted(TObjArray * array, AliESDEvent */*esd*/, Int_t iter)
+{
+ //
+ //
+ // Two reasons to have multiple find tracks
+ // 1. Curling tracks can be find more than once
+ // 2. Splitted tracks
+ // a.) Multiple seeding to increase tracking efficiency - (~ 100% reached)
+ // b.) Edge effect on the sector boundaries
+ //
+ // This function tries to find tracks closed in the parametric space
+ //
+ // cut logic if distance is bigger than cut continue - Do Nothing
+ const Float_t kMaxdTheta = 0.05; // maximal distance in theta
+ const Float_t kMaxdPhi = 0.05; // maximal deistance in phi
+ const Float_t kdelta = 40.; //delta r to calculate track distance
+ //
+ // const Float_t kMaxDist0 = 1.; // maximal distance 0
+ //const Float_t kMaxDist1 = 0.3; // maximal distance 1 - cut if track in separate rows
+ //
+ /*
+ TCut csec("csec","abs(Tr0.fRelativeSector-Tr1.fRelativeSector)<2");
+ TCut cdtheta("cdtheta","abs(dtheta)<0.05");
+ */
+ //
+ //
+ //
+ Int_t nentries = array->GetEntriesFast();
+ AliHelix *helixes = new AliHelix[nentries];
+ Float_t *xm = new Float_t[nentries];
+ //
+ //
+ TStopwatch timer;
+ timer.Start();
+ //
+ //Sort tracks according quality
+ //
+ Int_t nseed = array->GetEntriesFast();
+ Float_t * quality = new Float_t[nseed];
+ Int_t * indexes = new Int_t[nseed];
+ for (Int_t i=0; i<nseed; i++) {
+ AliTPCseed *pt=(AliTPCseed*)array->UncheckedAt(i);
+ if (!pt){
+ quality[i]=-1;
+ continue;
+ }
+ pt->UpdatePoints(); //select first last max dens points
+ Float_t * points = pt->GetPoints();
+ if (points[3]<0.8) quality[i] =-1;
+ quality[i] = (points[2]-points[0])+pt->GetNumberOfClusters();
+ //prefer high momenta tracks if overlaps
+ quality[i] *= TMath::Sqrt(TMath::Abs(pt->Pt())+0.5);
+ }
+ TMath::Sort(nseed,quality,indexes);
+
+
+ //
+ // Find track COG in x direction - point with best defined parameters
+ //
+ for (Int_t i=0;i<nentries;i++){
+ AliTPCseed* track = (AliTPCseed*)array->At(i);
+ if (!track) continue;
+ track->SetCircular(0);
+ new (&helixes[i]) AliHelix(*track);
+ Int_t ncl=0;
+ xm[i]=0;
+ for (Int_t icl=0; icl<160; icl++){
+ AliTPCclusterMI * cl = track->GetClusterPointer(icl);
+ if (cl) {
+ xm[i]+=cl->GetX();
+ ncl++;
+ }
+ }
+ if (ncl>0) xm[i]/=Float_t(ncl);
+ }
+ //
+ for (Int_t is0=0;is0<nentries;is0++){
+ Int_t i0 = indexes[is0];
+ AliTPCseed * track0 = (AliTPCseed*)array->At(i0);
+ if (!track0) continue;
+ Float_t xc0 = helixes[i0].GetHelix(6);
+ Float_t yc0 = helixes[i0].GetHelix(7);
+ Float_t fi0 = TMath::ATan2(yc0,xc0);
+
+ for (Int_t is1=is0+1;is1<nentries;is1++){
+ Int_t i1 = indexes[is1];
+ AliTPCseed * track1 = (AliTPCseed*)array->At(i1);
+ if (!track1) continue;
+ //
+ Int_t dsec = TMath::Abs((track0->GetRelativeSector()%18)-(track1->GetRelativeSector()%18)); // sector distance
+ if (dsec>1 && dsec<17) continue;
+
+ if (track1->GetKinkIndexes()[0] == -track0->GetKinkIndexes()[0]) continue;
+
+ Float_t dtheta = TMath::Abs(track0->GetTgl()-track1->GetTgl())<TMath::Abs(track0->GetTgl()+track1->GetTgl())? track0->GetTgl()-track1->GetTgl():track0->GetTgl()+track1->GetTgl();
+ if (TMath::Abs(dtheta)>kMaxdTheta) continue;
+ //
+ Float_t xc1 = helixes[i1].GetHelix(6);
+ Float_t yc1 = helixes[i1].GetHelix(7);
+ Float_t fi1 = TMath::ATan2(yc1,xc1);
+ //
+ Float_t dfi = fi0-fi1;
+ if (dfi>TMath::Pi()) dfi-=TMath::TwoPi(); // take care about edge effect
+ if (dfi<-TMath::Pi()) dfi+=TMath::TwoPi(); //
+ if (TMath::Abs(dfi)>kMaxdPhi&&helixes[i0].GetHelix(4)*helixes[i1].GetHelix(4)<0){
+ //
+ // if short tracks with undefined sign
+ fi1 = -TMath::ATan2(yc1,-xc1);
+ dfi = fi0-fi1;
+ if (dfi>TMath::Pi()) dfi-=TMath::TwoPi(); // take care about edge effect
+ if (dfi<-TMath::Pi()) dfi+=TMath::TwoPi(); //
+ }
+ if (TMath::Abs(dfi)>kMaxdPhi) continue;
+ //
+ //
+ Float_t sum =0;
+ Float_t sums=0;
+ Float_t sum0=0;
+ Float_t sum1=0;
+ for (Int_t icl=0; icl<160; icl++){
+ Int_t index0=track0->GetClusterIndex2(icl);
+ Int_t index1=track1->GetClusterIndex2(icl);
+ Bool_t used0 = (index0>0 && !(index0&0x8000));
+ Bool_t used1 = (index1>0 && !(index1&0x8000));
+ //
+ if (used0) sum0++; // used cluster0
+ if (used1) sum1++; // used clusters1
+ if (used0&&used1) sum++;
+ if (index0==index1 && used0 && used1) sums++;
+ }
+
+ //
+ if (sums<10) continue;
+ if (sum<40) continue;
+ if (sums/Float_t(TMath::Min(sum0,sum1))<0.5) continue;
+ //
+ Double_t dist[5][4]; // distance at X
+ Double_t mdist[4]={0,0,0,0}; // mean distance on range +-delta
+
+ //
+ //
+ track0->GetDistance(track1,xm[i0],dist[0],AliTracker::GetBz());
+ for (Int_t i=0;i<3;i++) mdist[i]+=TMath::Abs(dist[0][i]);
+ track0->GetDistance(track1,xm[i1],dist[1],AliTracker::GetBz());
+ for (Int_t i=0;i<3;i++) mdist[i]+=TMath::Abs(dist[1][i]);
+ //
+ track0->GetDistance(track1,TMath::Min(xm[i1],xm[i0])-kdelta,dist[2],AliTracker::GetBz());
+ for (Int_t i=0;i<3;i++) mdist[i]+=TMath::Abs(dist[2][i]);
+ track0->GetDistance(track1,TMath::Max(xm[i1],xm[i0])+kdelta,dist[3],AliTracker::GetBz());
+ for (Int_t i=0;i<3;i++) mdist[i]+=TMath::Abs(dist[3][i]);
+ //
+ track0->GetDistance(track1,(xm[i1]+xm[i0])*0.5,dist[4],AliTracker::GetBz());
+ for (Int_t i=0;i<3;i++) mdist[i]+=TMath::Abs(dist[4][i]);
+ for (Int_t i=0;i<3;i++) mdist[i]*=0.2;
+ //
+ //
+ Int_t lab0=track0->GetLabel();
+ Int_t lab1=track1->GetLabel();
+ if( AliTPCReconstructor::StreamLevel()>5){
+ TTreeSRedirector &cstream = *fDebugStreamer;
+ cstream<<"Splitted"<<
+ "iter="<<iter<<
+ "lab0="<<lab0<<
+ "lab1="<<lab1<<
+ "Tr0.="<<track0<< // seed0
+ "Tr1.="<<track1<< // seed1
+ "h0.="<<&helixes[i0]<<
+ "h1.="<<&helixes[i1]<<
+ //
+ "sum="<<sum<< //the sum of rows with cl in both
+ "sum0="<<sum0<< //the sum of rows with cl in 0 track
+ "sum1="<<sum1<< //the sum of rows with cl in 1 track
+ "sums="<<sums<< //the sum of shared clusters
+ "xm0="<<xm[i0]<< // the center of track
+ "xm1="<<xm[i1]<< // the x center of track
+ // General cut variables
+ "dfi="<<dfi<< // distance in fi angle
+ "dtheta="<<dtheta<< // distance int theta angle
+ //
+ //
+ "dist0="<<dist[4][0]<< //distance x
+ "dist1="<<dist[4][1]<< //distance y
+ "dist2="<<dist[4][2]<< //distance z
+ "mdist0="<<mdist[0]<< //distance x
+ "mdist1="<<mdist[1]<< //distance y
+ "mdist2="<<mdist[2]<< //distance z
+ "\n";
+ }
+ delete array->RemoveAt(i1);
+ }
+ }
+ delete [] helixes;
+ delete [] xm;
+ delete [] quality;
+ delete [] indexes;
+ AliInfo("Time for splitted tracks removal");
+ timer.Print();
+}
+
+
+
+void AliTPCtrackerMI::FindCurling(const TObjArray * array, AliESDEvent */*esd*/, Int_t iter)
+{
+ //
+ // find Curling tracks
+ // Use AliTPCReconstructor::StreamLevel()>1 if you want to tune parameters - cuts
+ //
+ //
+ // Algorithm done in 2 phases - because of CPU consumption
+ // it is n^2 algorithm - for lead-lead 20000x20000 combination are investigated
+ // see detal in MC part what can be used to cut
+ //
+ //
+ //
+ const Float_t kMaxC = 400; // maximal curvature to of the track
+ const Float_t kMaxdTheta = 0.15; // maximal distance in theta
+ const Float_t kMaxdPhi = 0.15; // maximal distance in phi
+ const Float_t kPtRatio = 0.3; // ratio between pt
+ const Float_t kMinDCAR = 2.; // distance to the primary vertex in r - see cpipe cut
- Double_t alpha, cs,sn, xx2,yy2;
//
- alpha = (sec[1]-sec[2])*fSectors->GetAlpha();
- cs = TMath::Cos(alpha);
- sn = TMath::Sin(alpha);
- xx2= xyz[1][0]*cs-xyz[1][1]*sn;
- yy2= xyz[1][0]*sn+xyz[1][1]*cs;
- xyz[1][0] = xx2;
- xyz[1][1] = yy2;
+ // Curling tracks cuts
//
- alpha = (sec[0]-sec[2])*fSectors->GetAlpha();
- cs = TMath::Cos(alpha);
- sn = TMath::Sin(alpha);
- xx2= xyz[0][0]*cs-xyz[0][1]*sn;
- yy2= xyz[0][0]*sn+xyz[0][1]*cs;
- xyz[0][0] = xx2;
- xyz[0][1] = yy2;
//
+ const Float_t kMaxDeltaRMax = 40; // distance in outer radius
+ const Float_t kMaxDeltaRMin = 5.; // distance in lower radius - see cpipe cut
+ const Float_t kMinAngle = 2.9; // angle between tracks
+ const Float_t kMaxDist = 5; // biggest distance
//
+ // The cuts can be tuned using the "MC information stored in Multi tree ==> see FindMultiMC
+ /*
+ Fast cuts:
+ TCut csign("csign","Tr0.fP[4]*Tr1.fP[4]<0"); //opposite sign
+ TCut cmax("cmax","abs(Tr0.GetC())>1/400");
+ TCut cda("cda","sqrt(dtheta^2+dfi^2)<0.15");
+ TCut ccratio("ccratio","abs((Tr0.fP[4]+Tr1.fP[4])/(abs(Tr0.fP[4])+abs(Tr1.fP[4])))<0.3");
+ TCut cpipe("cpipe", "min(abs(r0-rc0),abs(r1-rc1))>5");
+ //
+ TCut cdrmax("cdrmax","abs(abs(rc0+r0)-abs(rc1+r1))<40")
+ TCut cdrmin("cdrmin","abs(abs(rc0+r0)-abs(rc1+r1))<10")
+ //
+ Multi->Draw("dfi","iter==0"+csign+cmax+cda+ccratio); ~94% of curling tracks fulfill
+ Multi->Draw("min(abs(r0-rc0),abs(r1-rc1))","iter==0&&abs(lab1)==abs(lab0)"+csign+cmax+cda+ccratio+cpipe+cdrmin+cdrmax); //80%
+ //
+ Curling2->Draw("dfi","iter==0&&abs(lab0)==abs(lab1)"+csign+cmax+cdtheta+cdfi+ccratio)
+
+ */
+ //
//
- Double_t x[5],c[15];
//
- x[0]=xyz[2][1];
- x[1]=xyz[2][2];
- x[4]=F1(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1]);
- x[2]=F2(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1]);
- x[3]=F3n(xyz[2][0],xyz[2][1],xyz[0][0],xyz[0][1],xyz[2][2],xyz[0][2],x[4]);
- //
- Double_t sy =0.1, sz =0.1;
+ Int_t nentries = array->GetEntriesFast();
+ AliHelix *helixes = new AliHelix[nentries];
+ for (Int_t i=0;i<nentries;i++){
+ AliTPCseed* track = (AliTPCseed*)array->At(i);
+ if (!track) continue;
+ track->SetCircular(0);
+ new (&helixes[i]) AliHelix(*track);
+ }
//
- Double_t sy1=0.2, sz1=0.2;
- Double_t sy2=0.2, sz2=0.2;
- Double_t sy3=0.2;
//
- Double_t f40=(F1(xyz[2][0],xyz[2][1]+sy,xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1])-x[4])/sy;
- Double_t f42=(F1(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1]+sy,xyz[0][0],xyz[0][1])-x[4])/sy;
- Double_t f43=(F1(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1]+sy)-x[4])/sy;
- Double_t f20=(F2(xyz[2][0],xyz[2][1]+sy,xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1])-x[2])/sy;
- Double_t f22=(F2(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1]+sy,xyz[0][0],xyz[0][1])-x[2])/sy;
- Double_t f23=(F2(xyz[2][0],xyz[2][1],xyz[1][0],xyz[1][1],xyz[0][0],xyz[0][1]+sy)-x[2])/sy;
+ TStopwatch timer;
+ timer.Start();
+ Double_t phase[2][2],radius[2];
//
- Double_t f30=(F3(xyz[2][0],xyz[2][1]+sy,xyz[0][0],xyz[0][1],xyz[2][2],xyz[0][2])-x[3])/sy;
- Double_t f31=(F3(xyz[2][0],xyz[2][1],xyz[0][0],xyz[0][1],xyz[2][2]+sz,xyz[0][2])-x[3])/sz;
- Double_t f32=(F3(xyz[2][0],xyz[2][1],xyz[0][0],xyz[0][1]+sy,xyz[2][2],xyz[0][2])-x[3])/sy;
- Double_t f34=(F3(xyz[2][0],xyz[2][1],xyz[0][0],xyz[0][1],xyz[2][2],xyz[0][2]+sz)-x[3])/sz;
-
-
- c[0]=sy1;
- c[1]=0.; c[2]=sz1;
- c[3]=f20*sy1; c[4]=0.; c[5]=f20*sy1*f20+f22*sy2*f22+f23*sy3*f23;
- c[6]=f30*sy1; c[7]=f31*sz1; c[8]=f30*sy1*f20+f32*sy2*f22;
- c[9]=f30*sy1*f30+f31*sz1*f31+f32*sy2*f32+f34*sz2*f34;
- c[10]=f40*sy1; c[11]=0.; c[12]=f40*sy1*f20+f42*sy2*f22+f43*sy3*f23;
- c[13]=f30*sy1*f40+f32*sy2*f42;
- c[14]=f40*sy1*f40+f42*sy2*f42+f43*sy3*f43;
-
- // Int_t row1 = fSectors->GetRowNumber(xyz[2][0]);
- AliTPCseed *seed=new AliTPCseed(xyz[2][0], sec[2]*fSectors->GetAlpha()+fSectors->GetAlphaShift(), x, c, 0);
- seed->SetLastPoint(row[2]);
- seed->SetFirstPoint(row[2]);
- for (Int_t i=row[0];i<row[2];i++){
- seed->SetClusterIndex(i, track->GetClusterIndex(i));
- }
+ // Find tracks
+ //
+ //
+ for (Int_t i0=0;i0<nentries;i0++){
+ AliTPCseed * track0 = (AliTPCseed*)array->At(i0);
+ if (!track0) continue;
+ if (TMath::Abs(track0->GetC())<1/kMaxC) continue;
+ Float_t xc0 = helixes[i0].GetHelix(6);
+ Float_t yc0 = helixes[i0].GetHelix(7);
+ Float_t r0 = helixes[i0].GetHelix(8);
+ Float_t rc0 = TMath::Sqrt(xc0*xc0+yc0*yc0);
+ Float_t fi0 = TMath::ATan2(yc0,xc0);
+
+ for (Int_t i1=i0+1;i1<nentries;i1++){
+ AliTPCseed * track1 = (AliTPCseed*)array->At(i1);
+ if (!track1) continue;
+ if (TMath::Abs(track1->GetC())<1/kMaxC) continue;
+ Float_t xc1 = helixes[i1].GetHelix(6);
+ Float_t yc1 = helixes[i1].GetHelix(7);
+ Float_t r1 = helixes[i1].GetHelix(8);
+ Float_t rc1 = TMath::Sqrt(xc1*xc1+yc1*yc1);
+ Float_t fi1 = TMath::ATan2(yc1,xc1);
+ //
+ Float_t dfi = fi0-fi1;
+ //
+ //
+ if (dfi>1.5*TMath::Pi()) dfi-=TMath::Pi(); // take care about edge effect
+ if (dfi<-1.5*TMath::Pi()) dfi+=TMath::Pi(); //
+ Float_t dtheta = TMath::Abs(track0->GetTgl()-track1->GetTgl())<TMath::Abs(track0->GetTgl()+track1->GetTgl())? track0->GetTgl()-track1->GetTgl():track0->GetTgl()+track1->GetTgl();
+ //
+ //
+ // FIRST fast cuts
+ if (track0->GetBConstrain()&&track1->GetBConstrain()) continue; // not constrained
+ if (track1->GetSigned1Pt()*track0->GetSigned1Pt()>0) continue; // not the same sign
+ if ( TMath::Abs(track1->GetTgl()+track0->GetTgl())>kMaxdTheta) continue; //distance in the Theta
+ if ( TMath::Abs(dfi)>kMaxdPhi) continue; //distance in phi
+ if ( TMath::Sqrt(dfi*dfi+dtheta*dtheta)>kMaxdPhi) continue; //common angular offset
+ //
+ Float_t pt0 = track0->GetSignedPt();
+ Float_t pt1 = track1->GetSignedPt();
+ if ((TMath::Abs(pt0+pt1)/(TMath::Abs(pt0)+TMath::Abs(pt1)))>kPtRatio) continue;
+ if ((iter==1) && TMath::Abs(TMath::Abs(rc0+r0)-TMath::Abs(rc1+r1))>kMaxDeltaRMax) continue;
+ if ((iter!=1) &&TMath::Abs(TMath::Abs(rc0-r0)-TMath::Abs(rc1-r1))>kMaxDeltaRMin) continue;
+ if (TMath::Min(TMath::Abs(rc0-r0),TMath::Abs(rc1-r1))<kMinDCAR) continue;
+ //
+ //
+ // Now find closest approach
+ //
+ //
+ //
+ Int_t npoints = helixes[i0].GetRPHIintersections(helixes[i1], phase, radius,10);
+ if (npoints==0) continue;
+ helixes[i0].GetClosestPhases(helixes[i1], phase);
+ //
+ Double_t xyz0[3];
+ Double_t xyz1[3];
+ Double_t hangles[3];
+ helixes[i0].Evaluate(phase[0][0],xyz0);
+ helixes[i1].Evaluate(phase[0][1],xyz1);
- return seed;
+ helixes[i0].GetAngle(phase[0][0],helixes[i1],phase[0][1],hangles);
+ Double_t deltah[2],deltabest;
+ if (TMath::Abs(hangles[2])<kMinAngle) continue;
+
+ if (npoints>0){
+ Int_t ibest=0;
+ helixes[i0].ParabolicDCA(helixes[i1],phase[0][0],phase[0][1],radius[0],deltah[0],2);
+ if (npoints==2){
+ helixes[i0].ParabolicDCA(helixes[i1],phase[1][0],phase[1][1],radius[1],deltah[1],2);
+ if (deltah[1]<deltah[0]) ibest=1;
+ }
+ deltabest = TMath::Sqrt(deltah[ibest]);
+ helixes[i0].Evaluate(phase[ibest][0],xyz0);
+ helixes[i1].Evaluate(phase[ibest][1],xyz1);
+ helixes[i0].GetAngle(phase[ibest][0],helixes[i1],phase[ibest][1],hangles);
+ Double_t radiusbest = TMath::Sqrt(radius[ibest]);
+ //
+ if (deltabest>kMaxDist) continue;
+ // if (mindcar+mindcaz<40 && (TMath::Abs(hangles[2])<kMinAngle ||deltabest>3)) continue;
+ Bool_t sign =kFALSE;
+ if (hangles[2]>kMinAngle) sign =kTRUE;
+ //
+ if (sign){
+ // circular[i0] = kTRUE;
+ // circular[i1] = kTRUE;
+ if (track0->OneOverPt()<track1->OneOverPt()){
+ track0->SetCircular(track0->GetCircular()+1);
+ track1->SetCircular(track1->GetCircular()+2);
+ }
+ else{
+ track1->SetCircular(track1->GetCircular()+1);
+ track0->SetCircular(track0->GetCircular()+2);
+ }
+ }
+ if (AliTPCReconstructor::StreamLevel()>1){
+ //
+ //debug stream to tune "fine" cuts
+ Int_t lab0=track0->GetLabel();
+ Int_t lab1=track1->GetLabel();
+ TTreeSRedirector &cstream = *fDebugStreamer;
+ cstream<<"Curling2"<<
+ "iter="<<iter<<
+ "lab0="<<lab0<<
+ "lab1="<<lab1<<
+ "Tr0.="<<track0<<
+ "Tr1.="<<track1<<
+ //
+ "r0="<<r0<<
+ "rc0="<<rc0<<
+ "fi0="<<fi0<<
+ "r1="<<r1<<
+ "rc1="<<rc1<<
+ "fi1="<<fi1<<
+ "dfi="<<dfi<<
+ "dtheta="<<dtheta<<
+ //
+ "npoints="<<npoints<<
+ "hangles0="<<hangles[0]<<
+ "hangles1="<<hangles[1]<<
+ "hangles2="<<hangles[2]<<
+ "xyz0="<<xyz0[2]<<
+ "xyzz1="<<xyz1[2]<<
+ "radius="<<radiusbest<<
+ "deltabest="<<deltabest<<
+ "phase0="<<phase[ibest][0]<<
+ "phase1="<<phase[ibest][1]<<
+ "\n";
+
+ }
+ }
+ }
+ }
+ delete [] helixes;
+ if (AliTPCReconstructor::StreamLevel()>1) {
+ AliInfo("Time for curling tracks removal");
+ timer.Print();
+ }
}
-void AliTPCtrackerMI::FindKinks(TObjArray * array, AliESD *esd)
+
+
+
+
+void AliTPCtrackerMI::FindKinks(TObjArray * array, AliESDEvent *esd)
{
//
// find kinks
//
// Find circling track
- TTreeSRedirector &cstream = *fDebugStreamer;
//
for (Int_t i0=0;i0<nentries;i0++){
AliTPCseed * track0 = (AliTPCseed*)array->At(i0);
if ( TMath::Abs(track1->GetTgl()+track0->GetTgl())>0.1) continue;
if (track0->GetBConstrain()&&track1->GetBConstrain()) continue;
if (TMath::Abs(1./track1->GetC())>200) continue;
- if (track1->Get1Pt()*track0->Get1Pt()>0) continue;
+ if (track1->GetSigned1Pt()*track0->GetSigned1Pt()>0) continue;
if (track1->GetTgl()*track0->GetTgl()>0) continue;
if (TMath::Max(TMath::Abs(1./track0->GetC()),TMath::Abs(1./track1->GetC()))>190) continue;
- if (track0->GetBConstrain()&&TMath::Abs(track1->Get1Pt())<TMath::Abs(track0->Get1Pt())) continue; //returning - lower momenta
- if (track1->GetBConstrain()&&TMath::Abs(track0->Get1Pt())<TMath::Abs(track1->Get1Pt())) continue; //returning - lower momenta
+ if (track0->GetBConstrain()&&track1->OneOverPt()<track0->OneOverPt()) continue; //returning - lower momenta
+ if (track1->GetBConstrain()&&track0->OneOverPt()<track1->OneOverPt()) continue; //returning - lower momenta
//
Float_t mindcar = TMath::Min(TMath::Abs(dca[i0]),TMath::Abs(dca[i1]));
if (mindcar<5) continue;
helixes[i0].GetAngle(phase[0][0],helixes[i1],phase[0][1],hangles);
Double_t deltah[2],deltabest;
if (hangles[2]<2.8) continue;
- /*
- cstream<<"C"<<track0->fLab<<track1->fLab<<
- track0->fP3<<track1->fP3<<
- track0->fP4<<track1->fP4<<
- delta<<rmean<<npoints<<
- hangles[0]<<hangles[2]<<
- xyz0[2]<<xyz1[2]<<radius[0]<<"\n";
- */
if (npoints>0){
Int_t ibest=0;
helixes[i0].ParabolicDCA(helixes[i1],phase[0][0],phase[0][1],radius[0],deltah[0],2);
//
if (deltabest>6) continue;
if (mindcar+mindcaz<40 && (hangles[2]<3.12||deltabest>3)) continue;
- Bool_t sign =kFALSE;
- if (hangles[2]>3.06) sign =kTRUE;
+ Bool_t lsign =kFALSE;
+ if (hangles[2]>3.06) lsign =kTRUE;
//
- if (sign){
+ if (lsign){
circular[i0] = kTRUE;
circular[i1] = kTRUE;
- if (TMath::Abs(track0->Get1Pt())<TMath::Abs(track1->Get1Pt())){
+ if (track0->OneOverPt()<track1->OneOverPt()){
track0->SetCircular(track0->GetCircular()+1);
track1->SetCircular(track1->GetCircular()+2);
}
track0->SetCircular(track0->GetCircular()+2);
}
}
- if (sign&&AliTPCReconstructor::StreamLevel()>1){
+ if (lsign&&AliTPCReconstructor::StreamLevel()>1){
//debug stream
Int_t lab0=track0->GetLabel();
Int_t lab1=track1->GetLabel();
+ TTreeSRedirector &cstream = *fDebugStreamer;
cstream<<"Curling"<<
"lab0="<<lab0<<
"lab1="<<lab1<<
for (Int_t i =0;i<nentries;i++){
if (sign[i]==0) continue;
AliTPCseed * track0 = (AliTPCseed*)array->At(i);
+ if (track0==0) {
+ AliInfo("seed==0");
+ continue;
+ }
ntracks++;
//
Double_t cradius0 = 40*40;
Float_t dens00=-1,dens01=-1;
Float_t dens10=-1,dens11=-1;
//
- Int_t found,foundable,shared;
- track0->GetClusterStatistic(0,row0-5, found, foundable,shared,kFALSE);
+ Int_t found,foundable,ishared;
+ track0->GetClusterStatistic(0,row0-5, found, foundable,ishared,kFALSE);
if (foundable>5) dens00 = Float_t(found)/Float_t(foundable);
- track0->GetClusterStatistic(row0+5,155, found, foundable,shared,kFALSE);
+ track0->GetClusterStatistic(row0+5,155, found, foundable,ishared,kFALSE);
if (foundable>5) dens01 = Float_t(found)/Float_t(foundable);
//
- track1->GetClusterStatistic(0,row0-5, found, foundable,shared,kFALSE);
+ track1->GetClusterStatistic(0,row0-5, found, foundable,ishared,kFALSE);
if (foundable>10) dens10 = Float_t(found)/Float_t(foundable);
- track1->GetClusterStatistic(row0+5,155, found, foundable,shared,kFALSE);
+ track1->GetClusterStatistic(row0+5,155, found, foundable,ishared,kFALSE);
if (foundable>10) dens11 = Float_t(found)/Float_t(foundable);
//
if (dens00<dens10 && dens01<dens11) continue;
if (TMath::Abs(ktrack0->GetC())>5) continue; // cut on the curvature for mother particle
AliExternalTrackParam paramm(*ktrack0);
AliExternalTrackParam paramd(*ktrack1);
- if (row0>60&&ktrack1->GetReference().GetX()>90.) new (¶md) AliExternalTrackParam(ktrack1->GetReference());
+ if (row0>60&&ktrack1->GetReference().GetX()>90.)new (¶md) AliExternalTrackParam(ktrack1->GetReference());
//
//
kink->SetMother(paramm);
Float_t x[3] = { kink->GetPosition()[0],kink->GetPosition()[1],kink->GetPosition()[2]};
Int_t index[4];
- fParam->Transform0to1(x,index);
- fParam->Transform1to2(x,index);
+ fkParam->Transform0to1(x,index);
+ fkParam->Transform1to2(x,index);
row0 = GetRowNumber(x[0]);
if (kink->GetR()<100) continue;
//
for (Int_t i=0;i<nkinks;i++){
quality[i] =100000;
- AliKink *kink = (AliKink*)kinks->At(i);
+ AliKink *kinkl = (AliKink*)kinks->At(i);
//
// refit kinks towards vertex
//
- Int_t index0 = kink->GetIndex(0);
- Int_t index1 = kink->GetIndex(1);
+ Int_t index0 = kinkl->GetIndex(0);
+ Int_t index1 = kinkl->GetIndex(1);
AliTPCseed * ktrack0 = (AliTPCseed*)array->At(index0);
AliTPCseed * ktrack1 = (AliTPCseed*)array->At(index1);
//
//
// Refit Kink under if too small angle
//
- if (kink->GetAngle(2)<0.05){
- kink->SetTPCRow0(GetRowNumber(kink->GetR()));
- Int_t row0 = kink->GetTPCRow0();
- Int_t drow = Int_t(2.+0.5/(0.05+kink->GetAngle(2)));
+ if (kinkl->GetAngle(2)<0.05){
+ kinkl->SetTPCRow0(GetRowNumber(kinkl->GetR()));
+ Int_t row0 = kinkl->GetTPCRow0();
+ Int_t drow = Int_t(2.+0.5/(0.05+kinkl->GetAngle(2)));
//
//
Int_t last = row0-drow;
AliTPCseed* seed1 = ReSeed(ktrack1,first,kTRUE);
//
if (seed0 && seed1){
- kink->SetStatus(1,8);
- if (RefitKink(*seed0,*seed1,*kink)) kink->SetStatus(1,9);
- row0 = GetRowNumber(kink->GetR());
+ kinkl->SetStatus(1,8);
+ if (RefitKink(*seed0,*seed1,*kinkl)) kinkl->SetStatus(1,9);
+ row0 = GetRowNumber(kinkl->GetR());
sumn = seed0->GetNumberOfClusters()+seed1->GetNumberOfClusters();
- new (&mothers[i]) AliTPCseed(*seed0);
- new (&daughters[i]) AliTPCseed(*seed1);
+ mothers[i] = *seed0;
+ daughters[i] = *seed1;
}
else{
delete kinks->RemoveAt(i);
if (seed1) delete seed1;
continue;
}
- if (kink->GetDistance()>0.5 || kink->GetR()<110 || kink->GetR()>240) {
+ if (kinkl->GetDistance()>0.5 || kinkl->GetR()<110 || kinkl->GetR()>240) {
delete kinks->RemoveAt(i);
if (seed0) delete seed0;
if (seed1) delete seed1;
delete seed1;
}
//
- if (kink) quality[i] = 160*((0.1+kink->GetDistance())*(2.-kink->GetTPCDensityFactor()))/(sumn+40.); //the longest -clossest will win
+ if (kinkl) quality[i] = 160*((0.1+kinkl->GetDistance())*(2.-kinkl->GetTPCDensityFactor()))/(sumn+40.); //the longest -clossest will win
}
TMath::Sort(nkinks,quality,indexes,kFALSE);
//
for (Int_t i=0;i<nkinks;i++){
- AliKink * kink = (AliKink*) kinks->At(indexes[i]);
- if (!kink) continue;
- kink->SetTPCRow0(GetRowNumber(kink->GetR()));
- Int_t index0 = kink->GetIndex(0);
- Int_t index1 = kink->GetIndex(1);
- if (circular[index0]||circular[index1]&&kink->GetDistance()>0.2) continue;
- kink->SetMultiple(usage[index0],0);
- kink->SetMultiple(usage[index1],1);
- if (kink->GetMultiple()[0]+kink->GetMultiple()[1]>2) continue;
- if (kink->GetMultiple()[0]+kink->GetMultiple()[1]>0 && quality[indexes[i]]>0.2) continue;
- if (kink->GetMultiple()[0]+kink->GetMultiple()[1]>0 && kink->GetDistance()>0.2) continue;
- if (circular[index0]||circular[index1]&&kink->GetDistance()>0.1) continue;
+ AliKink * kinkl = (AliKink*) kinks->At(indexes[i]);
+ if (!kinkl) continue;
+ kinkl->SetTPCRow0(GetRowNumber(kinkl->GetR()));
+ Int_t index0 = kinkl->GetIndex(0);
+ Int_t index1 = kinkl->GetIndex(1);
+ if (circular[index0]||(circular[index1]&&kinkl->GetDistance()>0.2)) continue;
+ kinkl->SetMultiple(usage[index0],0);
+ kinkl->SetMultiple(usage[index1],1);
+ if (kinkl->GetMultiple()[0]+kinkl->GetMultiple()[1]>2) continue;
+ if (kinkl->GetMultiple()[0]+kinkl->GetMultiple()[1]>0 && quality[indexes[i]]>0.2) continue;
+ if (kinkl->GetMultiple()[0]+kinkl->GetMultiple()[1]>0 && kinkl->GetDistance()>0.2) continue;
+ if (circular[index0]||(circular[index1]&&kinkl->GetDistance()>0.1)) continue;
AliTPCseed * ktrack0 = (AliTPCseed*)array->At(index0);
AliTPCseed * ktrack1 = (AliTPCseed*)array->At(index1);
if (!ktrack0 || !ktrack1) continue;
- Int_t index = esd->AddKink(kink);
+ Int_t index = esd->AddKink(kinkl);
//
//
if ( ktrack0->GetKinkIndex(0)==0 && ktrack1->GetKinkIndex(0)==0) { //best kink
if (mothers[indexes[i]].GetNumberOfClusters()>20 && daughters[indexes[i]].GetNumberOfClusters()>20 && (mothers[indexes[i]].GetNumberOfClusters()+daughters[indexes[i]].GetNumberOfClusters())>100){
- new (ktrack0) AliTPCseed(mothers[indexes[i]]);
- new (ktrack1) AliTPCseed(daughters[indexes[i]]);
+ *ktrack0 = mothers[indexes[i]];
+ *ktrack1 = daughters[indexes[i]];
}
}
//
for (Int_t i=0;i<nentries;i++){
AliTPCseed * track0 = (AliTPCseed*)array->At(i);
if (!track0) continue;
- if (track0->GetPt()<1.4) continue;
+ if (track0->Pt()<1.4) continue;
//remove double high momenta tracks - overlapped with kink candidates
- Int_t shared=0;
+ Int_t ishared=0;
Int_t all =0;
for (Int_t icl=track0->GetFirstPoint();icl<track0->GetLastPoint(); icl++){
if (track0->GetClusterPointer(icl)!=0){
all++;
- if (track0->GetClusterPointer(icl)->IsUsed(10)) shared++;
+ if (track0->GetClusterPointer(icl)->IsUsed(10)) ishared++;
}
}
- if (Float_t(shared+1)/Float_t(nall+1)>0.5) {
+ if (Float_t(ishared+1)/Float_t(all+1)>0.5) {
delete array->RemoveAt(i);
continue;
}
AliTPCseed & mother = *pmother;
AliTPCseed & daughter = *pdaughter;
- AliKink & kink = *pkink;
- if (CheckKinkPoint(track0,mother,daughter, kink)){
+ AliKink & kinkl = *pkink;
+ if (CheckKinkPoint(track0,mother,daughter, kinkl)){
if (mother.GetNumberOfClusters()<30||daughter.GetNumberOfClusters()<20) {
delete pmother;
delete pdaughter;
delete pkink;
continue; //too short tracks
}
- if (mother.GetPt()<1.4) {
+ if (mother.Pt()<1.4) {
delete pmother;
delete pdaughter;
delete pkink;
continue;
}
- Int_t row0= kink.GetTPCRow0();
- if (kink.GetDistance()>0.5 || kink.GetR()<110. || kink.GetR()>240.) {
+ Int_t row0= kinkl.GetTPCRow0();
+ if (kinkl.GetDistance()>0.5 || kinkl.GetR()<110. || kinkl.GetR()>240.) {
delete pmother;
delete pdaughter;
delete pkink;
continue;
}
//
- Int_t index = esd->AddKink(&kink);
+ Int_t index = esd->AddKink(&kinkl);
mother.SetKinkIndex(0,-(index+1));
daughter.SetKinkIndex(0,index+1);
if (mother.GetNumberOfClusters()>50) {
timer.Print();
}
-void AliTPCtrackerMI::FindV0s(TObjArray * array, AliESD *esd)
+void AliTPCtrackerMI::FindV0s(const TObjArray * array, AliESDEvent *esd)
{
//
// find V0s
//
//
// //
- TTreeSRedirector &cstream = *fDebugStreamer;
Float_t fprimvertex[3]={GetX(),GetY(),GetZ()};
AliV0 vertex;
Double_t cradius0 = 10*10;
if (sign[i]==0) continue;
AliTPCseed * track0 = (AliTPCseed*)array->At(i);
if (!track0) continue;
- if (AliTPCReconstructor::StreamLevel()>0){
+ if (AliTPCReconstructor::StreamLevel()>1){
+ TTreeSRedirector &cstream = *fDebugStreamer;
cstream<<"Tracks"<<
"Tr0.="<<track0<<
"dca="<<dca[i]<<
"\n";
}
//
- if (track0->Get1Pt()<0) continue;
+ if (track0->GetSigned1Pt()<0) continue;
if (track0->GetKinkIndex(0)>0||isPrim[i]) continue; //daughter kink
//
if (TMath::Abs(helixes[i].GetHelix(4))<0.000000001) continue;
Int_t eventNr = esd->GetEventNumberInFile(); // This is most likely NOT the event number you'd like to use. It has nothing to do with the 'real' event number
Double_t radiusm= (delta[0]<delta[1])? TMath::Sqrt(radius[0]):TMath::Sqrt(radius[1]);
Double_t deltam= (delta[0]<delta[1])? TMath::Sqrt(delta[0]):TMath::Sqrt(delta[1]);
- if (AliTPCReconstructor::StreamLevel()>0) {
+ if (AliTPCReconstructor::StreamLevel()>1) {
Int_t lab0=track0->GetLabel();
Int_t lab1=track1->GetLabel();
Char_t c0=track0->GetCircular();
Char_t c1=track1->GetCircular();
+ TTreeSRedirector &cstream = *fDebugStreamer;
cstream<<"V0"<<
"Event="<<eventNr<<
"vertex.="<<&vertex<<
AliTPCseed * track0 = (AliTPCseed*)array->At(index0);
AliTPCseed * track1 = (AliTPCseed*)array->At(index1);
if (track0->TPCrPID(0)>0.3&&track1->TPCrPID(0)>0.3&&v0->GetAnglep()[2]<0.15) quality[i]+=1000000; // gamma conversion candidate
- if (track0->TPCrPID(4)>0.9||track1->TPCrPID(4)>0.9&&minpulldca>4) quality[i]*=10; // lambda candidate candidate
+ if (track0->TPCrPID(4)>0.9||(track1->TPCrPID(4)>0.9&&minpulldca>4)) quality[i]*=10; // lambda candidate candidate
}
TMath::Sort(ncandidates,quality,indexes,kTRUE);
}
{
Int_t eventNr = esd->GetEventNumberInFile(); // This is most likely NOT the event number you'd like to use. It has nothing to do with the 'real' event number
- if (AliTPCReconstructor::StreamLevel()>0) {
+ if (AliTPCReconstructor::StreamLevel()>1) {
Int_t lab0=track0->GetLabel();
Int_t lab1=track1->GetLabel();
+ TTreeSRedirector &cstream = *fDebugStreamer;
cstream<<"V02"<<
"Event="<<eventNr<<
"vertex.="<<v0<<
for (Int_t irow=0;irow<kNdiv;irow++){
FollowBackProlongation(mother, rows[irow]);
FollowProlongation(daughter,rows[kNdiv-1-irow]);
- new(¶m0[irow]) AliTPCseed(mother);
- new(¶m1[kNdiv-1-irow]) AliTPCseed(daughter);
+ param0[irow] = mother;
+ param1[kNdiv-1-irow] = daughter;
}
//
// define kinks
param0[index].Reset(kTRUE);
FollowProlongation(param0[index],0);
//
- new (&mother) AliTPCseed(param0[index]);
- new (&daughter) AliTPCseed(param1[index]); // daughter in vertex
+ mother = param0[index];
+ daughter = param1[index]; // daughter in vertex
//
kink.SetMother(mother);
kink.SetDaughter(daughter);
if (index>=0) break;
index = TMath::Abs(index)-1;
AliESDkink * kink = fEvent->GetKink(index);
- //kink->fTPCdensity2[0][0]=-1;
- //kink->fTPCdensity2[0][1]=-1;
- kink->SetTPCDensity2(-1,0,0);
- kink->SetTPCDensity2(1,0,1);
+ kink->SetTPCDensity(-1,0,0);
+ kink->SetTPCDensity(1,0,1);
//
Int_t row0 = kink->GetTPCRow0() - 2 - Int_t( 0.5/ (0.05+kink->GetAngle(2)));
if (row0<15) row0=15;
//
Int_t found,foundable,shared;
seed->GetClusterStatistic(0,row0, found, foundable,shared,kFALSE);
- if (foundable>5) kink->SetTPCDensity2(Float_t(found)/Float_t(foundable),0,0);
+ if (foundable>5) kink->SetTPCDensity(Float_t(found)/Float_t(foundable),0,0);
seed->GetClusterStatistic(row1,155, found, foundable,shared,kFALSE);
- if (foundable>5) kink->SetTPCDensity2(Float_t(found)/Float_t(foundable),0,1);
+ if (foundable>5) kink->SetTPCDensity(Float_t(found)/Float_t(foundable),0,1);
}
}
if (index<=0) break;
index = TMath::Abs(index)-1;
AliESDkink * kink = fEvent->GetKink(index);
- kink->SetTPCDensity2(-1,1,0);
- kink->SetTPCDensity2(-1,1,1);
+ kink->SetTPCDensity(-1,1,0);
+ kink->SetTPCDensity(-1,1,1);
//
Int_t row0 = kink->GetTPCRow0() -2 - Int_t( 0.5/ (0.05+kink->GetAngle(2)));
if (row0<15) row0=15;
//
Int_t found,foundable,shared;
seed->GetClusterStatistic(0,row0, found, foundable,shared,kFALSE);
- if (foundable>5) kink->SetTPCDensity2(Float_t(found)/Float_t(foundable),1,0);
+ if (foundable>5) kink->SetTPCDensity(Float_t(found)/Float_t(foundable),1,0);
seed->GetClusterStatistic(row1,155, found, foundable,shared,kFALSE);
- if (foundable>5) kink->SetTPCDensity2(Float_t(found)/Float_t(foundable),1,1);
+ if (foundable>5) kink->SetTPCDensity(Float_t(found)/Float_t(foundable),1,1);
}
}
for (Int_t irow=0;irow<20;irow++){
FollowBackProlongation(*seed0, rows[irow]);
FollowProlongation(*seed1,rows[19-irow]);
- new(¶m0[irow]) AliTPCseed(*seed0);
- new(¶m1[19-irow]) AliTPCseed(*seed1);
+ param0[irow] = *seed0;
+ param1[19-irow] = *seed1;
}
//
// define kinks
seed1->ResetCovariance(10.);
FollowProlongation(*seed0,0);
FollowBackProlongation(*seed1,158);
- new (&mother) AliTPCseed(*seed0); // backup mother at position 0
+ mother = *seed0; // backup mother at position 0
seed0->Reset(kFALSE);
seed1->Reset(kFALSE);
seed0->ResetCovariance(10.);
for (Int_t irow=0;irow<20;irow++){
FollowBackProlongation(*seed0, rows[irow]);
FollowProlongation(*seed1,rows[19-irow]);
- new(¶m0[irow]) AliTPCseed(*seed0);
- new(¶m1[19-irow]) AliTPCseed(*seed1);
+ param0[irow] = *seed0;
+ param1[19-irow] = *seed1;
}
//
// define kinks
//
//
// new (&mother) AliTPCseed(param0[index]);
- new (&daughter) AliTPCseed(param1[index]);
+ daughter = param1[index];
daughter.SetLabel(kink.GetLabel(1));
param0[index].Reset(kTRUE);
FollowProlongation(param0[index],0);
- new (&mother) AliTPCseed(param0[index]);
+ mother = param0[index];
mother.SetLabel(kink.GetLabel(0));
delete seed0;
delete seed1;
return 0;
}
-Int_t AliTPCtrackerMI::Clusters2Tracks (AliESD *esd)
+Int_t AliTPCtrackerMI::Clusters2Tracks (AliESDEvent *esd)
{
//
if (fSeeds) DeleteSeeds();
//
RemoveUsed2(fSeeds,0.85,0.85,0);
if (AliTPCReconstructor::GetRecoParam()->GetDoKinks()) FindKinks(fSeeds,fEvent);
+ //FindCurling(fSeeds, fEvent,0);
+ if (AliTPCReconstructor::StreamLevel()>2) FindMultiMC(fSeeds, fEvent,-1); // find multi found tracks
RemoveUsed2(fSeeds,0.5,0.4,20);
+ FindSplitted(fSeeds, fEvent,0); // find multi found tracks
+ if (AliTPCReconstructor::StreamLevel()>2) FindMultiMC(fSeeds, fEvent,0); // find multi found tracks
+
// //
// // refit short tracks
// //
Int_t nseed=fSeeds->GetEntriesFast();
-// for (Int_t i=0; i<nseed; i++) {
-// AliTPCseed *pt=(AliTPCseed*)fSeeds->UncheckedAt(i), &t=*pt;
-// if (!pt) continue;
-// Int_t nc=t.GetNumberOfClusters();
-// if (nc<15) {
-// delete fSeeds->RemoveAt(i);
-// continue;
-// }
-// if (pt->GetKinkIndexes()[0]!=0) continue; // ignore kink candidates
-// if (nc>100) continue; // hopefully, enough for ITS
-// AliTPCseed *seed2 = new AliTPCseed(*pt);
-// //seed2->Reset(kFALSE);
-// //pt->ResetCovariance();
-// seed2->Modify(1);
-// FollowBackProlongation(*seed2,158);
-// //seed2->Reset(kFALSE);
-// seed2->Modify(10);
-// FollowProlongation(*seed2,0);
-// TTreeSRedirector &cstream = *fDebugStreamer;
-// cstream<<"Crefit"<<
-// "Tr0.="<<pt<<
-// "Tr1.="<<seed2<<
-// "\n";
-// if (seed2->fN>pt->fN){
-// delete fSeeds->RemoveAt(i);
-// fSeeds->AddAt(seed2,i);
-// }else{
-// delete seed2;
-// }
-// }
-// RemoveUsed2(fSeeds,0.6,0.6,50);
-
-// FindV0s(fSeeds,fEvent);
- //RemoveDouble(fSeeds,0.2,0.6,11);
-
//
Int_t found = 0;
for (Int_t i=0; i<nseed; i++) {
fdensity = 2.;
cuts[0]=0.0080;
+ Int_t fLastSeedRowSec=AliTPCReconstructor::GetRecoParam()->GetLastSeedRowSec();
+
// find secondaries
- for (Int_t delta = 30; delta<90; delta+=10){
+ for (Int_t delta = 30; delta<fLastSeedRowSec; delta+=10){
//
cuts[0] = 0.3;
cuts[1] = 3.5;
}
-void AliTPCtrackerMI::SumTracks(TObjArray *arr1,TObjArray *arr2) const
+void AliTPCtrackerMI::SumTracks(TObjArray *arr1,TObjArray *&arr2) const
{
//
//sum tracks to common container
}
}
}
- delete arr2;
+ delete arr2; arr2 = 0;
}
if (!pt) continue;
if (!t.IsActive()) continue;
// follow prolongation to the first layer
- if ( (fSectors ==fInnerSec) || (t.GetFirstPoint()-fParam->GetNRowLow()>rfirst+1) )
+ if ( (fSectors ==fInnerSec) || (t.GetFirstPoint()-fkParam->GetNRowLow()>rfirst+1) )
FollowProlongation(t, rfirst+1);
}
if (!pt) continue;
if (nr==80) pt->UpdateReference();
if (!pt->IsActive()) continue;
- // if ( (fSectors ==fOuterSec) && (pt->fFirstPoint-fParam->GetNRowLow())<nr) continue;
+ // if ( (fSectors ==fOuterSec) && (pt->fFirstPoint-fkParam->GetNRowLow())<nr) continue;
if (pt->GetRelativeSector()>17) {
continue;
}
AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
if (!pt) continue;
if (!pt->IsActive()) continue;
- // if ((fSectors ==fOuterSec) && (pt->fFirstPoint-fParam->GetNRowLow())<nr) continue;
+ // if ((fSectors ==fOuterSec) && (pt->fFirstPoint-fkParam->GetNRowLow())<nr) continue;
if (pt->GetRelativeSector()>17) {
continue;
}
fSectors = fInnerSec;
FollowBackProlongation(*pt,fInnerSec->GetNRows()+fOuterSec->GetNRows()-1);
}
-
+ CookLabel(pt,0.3);
}
return 0;
}
AliTPCseed *pt = (AliTPCseed*)arr->UncheckedAt(i);
if (pt) {
FollowProlongation(*pt,0);
+ CookLabel(pt,0.3);
}
+
}
return 0;
}
ParallelTracking(fSeeds,fOuterSec->GetNRows()+fInnerSec->GetNRows()-1,fInnerSec->GetNRows());
fSectors = fInnerSec;
ParallelTracking(fSeeds,fInnerSec->GetNRows()-1,0);
- //WriteTracks();
return 1;
}
void AliTPCtrackerMI::GetShape(AliTPCseed * seed, Int_t row)
{
//
- //
- Float_t sd2 = TMath::Abs((fParam->GetZLength()-TMath::Abs(seed->GetZ())))*fParam->GetDiffL()*fParam->GetDiffL();
- // Float_t padlength = fParam->GetPadPitchLength(seed->fSector);
- Float_t padlength = GetPadPitchLength(row);
- //
- Float_t sresy = (seed->GetSector() < fParam->GetNSector()/2) ? 0.2 :0.3;
+ //
+ AliTPCClusterParam * clparam = AliTPCcalibDB::Instance()->GetClusterParam();
+ Float_t zdrift = TMath::Abs((fkParam->GetZLength(0)-TMath::Abs(seed->GetZ())));
+ Int_t type = (seed->GetSector() < fkParam->GetNSector()/2) ? 0: (row>126) ? 1:2;
Double_t angulary = seed->GetSnp();
- angulary = angulary*angulary/(1-angulary*angulary);
- seed->SetCurrentSigmaY2(sd2+padlength*padlength*angulary/12.+sresy*sresy);
- //
- Float_t sresz = fParam->GetZSigma();
- Float_t angularz = seed->GetTgl();
- seed->SetCurrentSigmaZ2(sd2+padlength*padlength*angularz*angularz*(1+angulary)/12.+sresz*sresz);
+ angulary = angulary*angulary/(1.-angulary*angulary);
+ Double_t angularz = seed->GetTgl()*seed->GetTgl()*(1.+angulary);
+
+ Double_t sigmay = clparam->GetRMS0(0,type,zdrift,TMath::Sqrt(TMath::Abs(angulary)));
+ Double_t sigmaz = clparam->GetRMS0(1,type,zdrift,TMath::Sqrt(TMath::Abs(angularz)));
+ seed->SetCurrentSigmaY2(sigmay*sigmay);
+ seed->SetCurrentSigmaZ2(sigmaz*sigmaz);
+ // Float_t sd2 = TMath::Abs((fkParam->GetZLength(0)-TMath::Abs(seed->GetZ())))*fkParam->GetDiffL()*fkParam->GetDiffL();
+// // Float_t padlength = fkParam->GetPadPitchLength(seed->fSector);
+// Float_t padlength = GetPadPitchLength(row);
+// //
+// Float_t sresy = (seed->GetSector() < fkParam->GetNSector()/2) ? 0.2 :0.3;
+// seed->SetCurrentSigmaY2(sd2+padlength*padlength*angulary/12.+sresy*sresy);
+// //
+// Float_t sresz = fkParam->GetZSigma();
+// seed->SetCurrentSigmaZ2(sd2+padlength*padlength*angularz*angularz*(1+angulary)/12.+sresz*sresz);
/*
Float_t wy = GetSigmaY(seed);
Float_t wz = GetSigmaZ(seed);
}
-Float_t AliTPCtrackerMI::GetSigmaY(AliTPCseed * seed)
-{
- //
- //
- Float_t sd2 = TMath::Abs((fParam->GetZLength()-TMath::Abs(seed->GetZ())))*fParam->GetDiffL()*fParam->GetDiffL();
- Float_t padlength = fParam->GetPadPitchLength(seed->GetSector());
- Float_t sres = (seed->GetSector() < fParam->GetNSector()/2) ? 0.2 :0.3;
- Float_t angular = seed->GetSnp();
- angular = angular*angular/(1-angular*angular);
- // angular*=angular;
- //angular = TMath::Sqrt(angular/(1-angular));
- Float_t res = TMath::Sqrt(sd2+padlength*padlength*angular/12.+sres*sres);
- return res;
-}
-Float_t AliTPCtrackerMI::GetSigmaZ(AliTPCseed * seed)
-{
- //
- //
- Float_t sd2 = TMath::Abs((fParam->GetZLength()-TMath::Abs(seed->GetZ())))*fParam->GetDiffL()*fParam->GetDiffL();
- Float_t padlength = fParam->GetPadPitchLength(seed->GetSector());
- Float_t sres = fParam->GetZSigma();
- Float_t angular = seed->GetTgl();
- Float_t res = TMath::Sqrt(sd2+padlength*padlength*angular*angular/12.+sres*sres);
- return res;
-}
-
//__________________________________________________________________________
void AliTPCtrackerMI::CookLabel(AliKalmanTrack *tk, Float_t wrong) const {
//--------------------------------------------------------------------
//This function "cooks" a track label. If label<0, this track is fake.
//--------------------------------------------------------------------
- AliTPCseed * t = (AliTPCseed*)tk;
+ AliTPCseed * t = dynamic_cast<AliTPCseed*>(tk);
+ if(!t){
+ printf("%s:%d wrong type \n",(char*)__FILE__,__LINE__);
+ return;
+ }
+
Int_t noc=t->GetNumberOfClusters();
if (noc<10){
//printf("\nnot founded prolongation\n\n\n");
}
-Int_t AliTPCtrackerMI::AliTPCSector::GetRowNumber(Double_t x) const
-{
- //return pad row number for this x
- Double_t r;
- if (fN < 64){
- r=fRow[fN-1].GetX();
- if (x > r) return fN;
- r=fRow[0].GetX();
- if (x < r) return -1;
- return Int_t((x-r)/fPadPitchLength + 0.5);}
- else{
- r=fRow[fN-1].GetX();
- if (x > r) return fN;
- r=fRow[0].GetX();
- if (x < r) return -1;
- Double_t r1=fRow[64].GetX();
- if(x<r1){
- return Int_t((x-r)/f1PadPitchLength + 0.5);}
- else{
- return (Int_t((x-r1)/f2PadPitchLength + 0.5)+64);}
- }
-}
-
Int_t AliTPCtrackerMI::GetRowNumber(Double_t x[3]) const
{
//return pad row number for given x vector
return GetRowNumber(localx);
}
-//_________________________________________________________________________
-void AliTPCtrackerMI::AliTPCSector::Setup(const AliTPCParam *par, Int_t f) {
- //-----------------------------------------------------------------------
- // Setup inner sector
- //-----------------------------------------------------------------------
- if (f==0) {
- fAlpha=par->GetInnerAngle();
- fAlphaShift=par->GetInnerAngleShift();
- fPadPitchWidth=par->GetInnerPadPitchWidth();
- fPadPitchLength=par->GetInnerPadPitchLength();
- fN=par->GetNRowLow();
- fRow=new AliTPCRow[fN];
- for (Int_t i=0; i<fN; i++) {
- fRow[i].SetX(par->GetPadRowRadiiLow(i));
- fRow[i].SetDeadZone(1.5); //1.5 cm of dead zone
- }
- } else {
- fAlpha=par->GetOuterAngle();
- fAlphaShift=par->GetOuterAngleShift();
- fPadPitchWidth = par->GetOuterPadPitchWidth();
- fPadPitchLength = par->GetOuter1PadPitchLength();
- f1PadPitchLength = par->GetOuter1PadPitchLength();
- f2PadPitchLength = par->GetOuter2PadPitchLength();
-
- fN=par->GetNRowUp();
- fRow=new AliTPCRow[fN];
- for (Int_t i=0; i<fN; i++) {
- fRow[i].SetX(par->GetPadRowRadiiUp(i));
- fRow[i].SetDeadZone(1.5); // 1.5 cm of dead zone
- }
- }
-}
-
-AliTPCtrackerMI::AliTPCRow::AliTPCRow():
- fDeadZone(0.),
- fClusters1(0),
- fN1(0),
- fClusters2(0),
- fN2(0),
- fN(0),
- fX(0.)
-{
- //
- // default constructor
- //
-}
-
-AliTPCtrackerMI::AliTPCRow::~AliTPCRow(){
- //
-
-}
-
-
-
-//_________________________________________________________________________
-void
-AliTPCtrackerMI::AliTPCRow::InsertCluster(const AliTPCclusterMI* c, UInt_t index) {
- //-----------------------------------------------------------------------
- // Insert a cluster into this pad row in accordence with its y-coordinate
- //-----------------------------------------------------------------------
- if (fN==kMaxClusterPerRow) {
- cerr<<"AliTPCRow::InsertCluster(): Too many clusters !\n"; return;
- }
- if (fN==0) {fIndex[0]=index; fClusters[fN++]=c; return;}
- Int_t i=Find(c->GetZ());
- memmove(fClusters+i+1 ,fClusters+i,(fN-i)*sizeof(AliTPCclusterMI*));
- memmove(fIndex +i+1 ,fIndex +i,(fN-i)*sizeof(UInt_t));
- fIndex[i]=index; fClusters[i]=c; fN++;
-}
-
-void AliTPCtrackerMI::AliTPCRow::ResetClusters() {
- //
- // reset clusters
- fN = 0;
- fN1 = 0;
- fN2 = 0;
- //delete[] fClusterArray;
- if (fClusters1) delete []fClusters1;
- if (fClusters2) delete []fClusters2;
- //fClusterArray=0;
- fClusters1 = 0;
- fClusters2 = 0;
-}
-//___________________________________________________________________
-Int_t AliTPCtrackerMI::AliTPCRow::Find(Double_t z) const {
+void AliTPCtrackerMI::MakeBitmaps(AliTPCseed *t)
+{
//-----------------------------------------------------------------------
- // Return the index of the nearest cluster
+ // Fill the cluster and sharing bitmaps of the track
//-----------------------------------------------------------------------
- if (fN==0) return 0;
- if (z <= fClusters[0]->GetZ()) return 0;
- if (z > fClusters[fN-1]->GetZ()) return fN;
- Int_t b=0, e=fN-1, m=(b+e)/2;
- for (; b<e; m=(b+e)/2) {
- if (z > fClusters[m]->GetZ()) b=m+1;
- else e=m;
- }
- return m;
-}
-
-
-//___________________________________________________________________
-AliTPCclusterMI * AliTPCtrackerMI::AliTPCRow::FindNearest(Double_t y, Double_t z, Double_t roady, Double_t roadz) const {
- //-----------------------------------------------------------------------
- // Return the index of the nearest cluster in z y
- //-----------------------------------------------------------------------
- Float_t maxdistance = roady*roady + roadz*roadz;
-
- AliTPCclusterMI *cl =0;
- for (Int_t i=Find(z-roadz); i<fN; i++) {
- AliTPCclusterMI *c=(AliTPCclusterMI*)(fClusters[i]);
- if (c->GetZ() > z+roadz) break;
- if ( (c->GetY()-y) > roady ) continue;
- Float_t distance = (c->GetZ()-z)*(c->GetZ()-z)+(c->GetY()-y)*(c->GetY()-y);
- if (maxdistance>distance) {
- maxdistance = distance;
- cl=c;
- }
+ Int_t firstpoint = 0;
+ Int_t lastpoint = 159;
+ AliTPCTrackerPoint *point;
+ AliTPCclusterMI *cluster;
+
+ for (int iter=firstpoint; iter<lastpoint; iter++) {
+ // Change to cluster pointers to see if we have a cluster at given padrow
+ cluster = t->GetClusterPointer(iter);
+ if (cluster) {
+ t->SetClusterMapBit(iter, kTRUE);
+ point = t->GetTrackPoint(iter);
+ if (point->IsShared())
+ t->SetSharedMapBit(iter,kTRUE);
+ else
+ t->SetSharedMapBit(iter, kFALSE);
+ }
+ else {
+ t->SetClusterMapBit(iter, kFALSE);
+ t->SetSharedMapBit(iter, kFALSE);
+ }
}
- return cl;
}
-AliTPCclusterMI * AliTPCtrackerMI::AliTPCRow::FindNearest2(Double_t y, Double_t z, Double_t roady, Double_t roadz,UInt_t & index) const
-{
- //-----------------------------------------------------------------------
- // Return the index of the nearest cluster in z y
- //-----------------------------------------------------------------------
- Float_t maxdistance = roady*roady + roadz*roadz;
- AliTPCclusterMI *cl =0;
-
- //PH Check boundaries. 510 is the size of fFastCluster
- Int_t iz1 = Int_t(z-roadz+254.5);
- if (iz1<0 || iz1>=510) return cl;
- iz1 = TMath::Max(GetFastCluster(iz1)-1,0);
- Int_t iz2 = Int_t(z+roadz+255.5);
- if (iz2<0 || iz2>=510) return cl;
- iz2 = TMath::Min(GetFastCluster(iz2)+1,fN);
-
- //FindNearest3(y,z,roady,roadz,index);
- // for (Int_t i=Find(z-roadz); i<fN; i++) {
- for (Int_t i=iz1; i<iz2; i++) {
- AliTPCclusterMI *c=(AliTPCclusterMI*)(fClusters[i]);
- if (c->GetZ() > z+roadz) break;
- if ( c->GetY()-y > roady ) continue;
- if ( y-c->GetY() > roady ) continue;
- Float_t distance = (c->GetZ()-z)*(c->GetZ()-z)+(c->GetY()-y)*(c->GetY()-y);
- if (maxdistance>distance) {
- maxdistance = distance;
- cl=c;
- index =i;
- //roady = TMath::Sqrt(maxdistance);
- }
- }
- return cl;
+Bool_t AliTPCtrackerMI::IsFindable(AliTPCseed & track){
+ //
+ // return flag if there is findable cluster at given position
+ //
+ Float_t kDeltaZ=10;
+ Float_t z = track.GetZ();
+
+ if (TMath::Abs(z)<(AliTPCReconstructor::GetCtgRange()*track.GetX()+kDeltaZ) &&
+ TMath::Abs(z)<fkParam->GetZLength(0) &&
+ (TMath::Abs(track.GetSnp())<AliTPCReconstructor::GetMaxSnpTracker()))
+ return kTRUE;
+ return kFALSE;
}
-
-AliTPCclusterMI * AliTPCtrackerMI::AliTPCRow::FindNearest3(Double_t y, Double_t z, Double_t roady, Double_t roadz,UInt_t & index) const
-{
- //-----------------------------------------------------------------------
- // Return the index of the nearest cluster in z y
- //-----------------------------------------------------------------------
- Float_t maxdistance = roady*roady + roadz*roadz;
- // Int_t iz = Int_t(z+255.);
- AliTPCclusterMI *cl =0;
- for (Int_t i=Find(z-roadz); i<fN; i++) {
- //for (Int_t i=fFastCluster[iz-2]; i<fFastCluster[iz+2]; i++) {
- AliTPCclusterMI *c=(AliTPCclusterMI*)(fClusters[i]);
- if (c->GetZ() > z+roadz) break;
- if ( c->GetY()-y > roady ) continue;
- if ( y-c->GetY() > roady ) continue;
- Float_t distance = (c->GetZ()-z)*(c->GetZ()-z)+(c->GetY()-y)*(c->GetY()-y);
- if (maxdistance>distance) {
- maxdistance = distance;
- cl=c;
- index =i;
- //roady = TMath::Sqrt(maxdistance);
- }
- }
- return cl;
+void AliTPCtrackerMI::AddCovariance(AliTPCseed * seed){
+ //
+ // Adding systematic error
+ // !!!! the systematic error for element 4 is in 1/cm not in pt
+
+ const Double_t *param = AliTPCReconstructor::GetRecoParam()->GetSystematicError();
+ Double_t covar[15];
+ for (Int_t i=0;i<15;i++) covar[i]=0;
+ // 0
+ // 1 2
+ // 3 4 5
+ // 6 7 8 9
+ // 10 11 12 13 14
+ covar[0] = param[0]*param[0];
+ covar[2] = param[1]*param[1];
+ covar[5] = param[2]*param[2];
+ covar[9] = param[3]*param[3];
+ Double_t facC = AliTracker::GetBz()*kB2C;
+ covar[14]= param[4]*param[4]*facC*facC;
+ seed->AddCovariance(covar);
}
-
-
-
-
-
-// AliTPCTrackerPoint * AliTPCseed::GetTrackPoint(Int_t i)
-// {
-// //
-// //
-// return &fTrackPoints[i];
-// }
-
-