// 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 "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 "TRandom.h"
#include "AliTPCcalibDB.h"
#include "AliTPCTransform.h"
+#include "AliTPCClusterParam.h"
//
ClassImp(AliTPCtrackerMI)
+
class AliTPCFastMath {
public:
AliTPCFastMath();
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;
//---------------------------------------------------------------------
// 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);
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):
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, AliESDEvent * 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)
{
AliTPCseed *pt=(AliTPCseed*)arr->UncheckedAt(i);
if (!pt) continue;
pt->UpdatePoints();
+ AddCovariance(pt);
+ if (AliTPCReconstructor::StreamLevel()>1) {
+ (*fDebugStreamer)<<"Track0"<<
+ "Tr0.="<<pt<<
+ "\n";
+ }
// pt->PropagateTo(fParam->GetInnerRadiusLow());
if (pt->GetKinkIndex(0)<=0){ //don't propagate daughter tracks
pt->PropagateTo(fParam->GetInnerRadiusLow());
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){
//
//
- //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[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);
- }
-
- //
- 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;
- }
- //
+ // Use calibrated cluster error from OCDB
//
+ AliTPCClusterParam * clparam = AliTPCcalibDB::Instance()->GetClusterParam();
//
- 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(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];
- }
- }
- }
+ 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 (ctype<0){
- res+=0.005;
- res*=2.4; // overestimate error 2 times
- }
- res+= snoise2;
+// //
+// 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(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 (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){
//
//
- //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 (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;
- }
- //
+ // Use calibrated cluster error from OCDB
//
+ AliTPCClusterParam * clparam = AliTPCcalibDB::Instance()->GetClusterParam();
//
- 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(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];
- }
- }
- 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;
+ 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
}
- 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(0)-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;
-
- if ((ctype==0) && (fSectors ==fOuterSec))
- res *= 0.81 +TMath::Exp(6.8*(rsigmaz-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;
+// //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];
+// //
- seed->SetErrorZ2(res);
- return res;
+// //
+// 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(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];
+// }
+// }
+// }
+
+// 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;
}
-*/
+
return LoadClusters();
}
+
+Int_t AliTPCtrackerMI::LoadClusters(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;
+ fParam->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 i=0;i<tpcrow->GetN1();i++)
+ tpcrow->SetCluster1(i, *(AliTPCclusterMI*)(clrow->GetArray()->At(i)));
+ }
+ 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)));
+ }
+ }
+ //
+ delete clrow;
+ LoadOuterSectors();
+ LoadInnerSectors();
+ return 0;
+}
+
+Int_t AliTPCtrackerMI::LoadClusters(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], fParam);
+ }
+ else{
+ left = (sec-fkNIS*2)/fkNOS;
+ fOuterSec[(sec-fkNIS*2)%fkNOS].InsertCluster(clust, count[sec][row], fParam);
+ }
+ }
+
+ // Load functions must be called behind LoadCluster(TClonesArray*)
+ // needed by HLT
+ //LoadOuterSectors();
+ //LoadInnerSectors();
+
+ return 0;
+}
+
+
Int_t AliTPCtrackerMI::LoadClusters()
{
//
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]);
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){
//
//
Double_t x[3]={cluster->GetRow(),cluster->GetPad(),cluster->GetTimeBin()};
Int_t i[1]={cluster->GetDetector()};
transform->Transform(x,i,0,1);
- if (!AliTPCReconstructor::GetRecoParam()->GetBYMirror()){
- if (cluster->GetDetector()%36>17){
- x[1]*=-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";
}
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) {
//
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 (!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.;
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);
t.GetProlongation(x,y,z);
}
//
- // update current shape info every 3 pad-row
- if ( (nr%6==0) || t.GetNumberOfClusters()<2 || (t.GetCurrentSigmaY2()<0.0001) ){
+ // update current shape info every 2 pad-row
+ if ( (nr%2==0) || t.GetNumberOfClusters()<2 || (t.GetCurrentSigmaY2()<0.0001) ){
// t.fCurrentSigmaY = GetSigmaY(&t);
//t.fCurrentSigmaZ = GetSigmaZ(&t);
GetShape(&t,nr);
//Int_t nr2 = nr;
- 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.;
if (cl) {
t.SetCurrentCluster(cl);
- // Int_t accept = AcceptCluster(&t,t.fCurrentCluster,1.);
+ // Int_t accept = AcceptCluster(&t,t.fCurrentCluster);
//if (accept<3){
t.SetClusterIndex2(row,index);
t.SetClusterPointer(row, cl);
}
//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);
}
// 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)){
//
seed->PropagateTo(fParam->GetInnerRadiusLow());
seed->UpdatePoints();
+ AddCovariance(seed);
MakeBitmaps(seed);
AliESDtrack *esd=event->GetTrack(i);
seed->CookdEdx(0.02,0.6);
}
//FindKinks(fSeeds,event);
Info("RefitInward","Number of refitted tracks %d",ntracks);
- fEvent =0;
- //WriteTracks();
return 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);
seed->CookdEdx(0.02,0.6);
CookLabel(seed,0.1); //For comparison only
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
- if (AliTPCReconstructor::StreamLevel()>1) {
+ if (AliTPCReconstructor::StreamLevel()>1 && esd) {
(*fDebugStreamer)<<"Cback"<<
"Tr0.="<<seed<<
+ "esd.="<<esd<<
"EventNrInFile="<<eventnumber<<
- "\n"; // patch 28 fev 06
+ "\n";
}
}
}
//FindKinks(fSeeds,event);
Info("PropagateBack","Number of back propagated tracks %d",ntracks);
fEvent =0;
- //WriteTracks();
return 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;
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);
// 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;
//
//
// 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)) {
}
if (ncl>0) xm[i]/=Float_t(ncl);
}
- TTreeSRedirector &cstream = *fDebugStreamer;
//
for (Int_t i0=0;i0<nentries;i0++){
AliTPCseed * track0 = (AliTPCseed*)array->At(i0);
}
}
//
+ if (AliTPCReconstructor::StreamLevel()>0) {
+ TTreeSRedirector &cstream = *fDebugStreamer;
cstream<<"Multi"<<
"iter="<<iter<<
"lab0="<<lab0<<
"r1="<<r1<<
"rc1="<<rc1<<
"\n";
+ }
}
}
delete [] helixes;
}
if (ncl>0) xm[i]/=Float_t(ncl);
}
- TTreeSRedirector &cstream = *fDebugStreamer;
//
for (Int_t is0=0;is0<nentries;is0++){
Int_t i0 = indexes[is0];
//
Int_t lab0=track0->GetLabel();
Int_t lab1=track1->GetLabel();
- 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";
+ 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(TObjArray * array, AliESDEvent *esd, Int_t iter)
+void AliTPCtrackerMI::FindCurling(TObjArray * array, AliESDEvent */*esd*/, Int_t iter)
{
//
// find Curling tracks
//
// Find tracks
//
- TTreeSRedirector &cstream = *fDebugStreamer;
//
for (Int_t i0=0;i0<nentries;i0++){
AliTPCseed * track0 = (AliTPCseed*)array->At(i0);
//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<<
//
// Find circling track
- TTreeSRedirector &cstream = *fDebugStreamer;
//
for (Int_t i0=0;i0<nentries;i0++){
AliTPCseed * track0 = (AliTPCseed*)array->At(i0);
//debug stream
Int_t lab0=track0->GetLabel();
Int_t lab1=track1->GetLabel();
+ TTreeSRedirector &cstream = *fDebugStreamer;
cstream<<"Curling"<<
"lab0="<<lab0<<
"lab1="<<lab1<<
//
//
// //
- TTreeSRedirector &cstream = *fDebugStreamer;
Float_t fprimvertex[3]={GetX(),GetY(),GetZ()};
AliV0 vertex;
Double_t cradius0 = 10*10;
AliTPCseed * track0 = (AliTPCseed*)array->At(i);
if (!track0) continue;
if (AliTPCReconstructor::StreamLevel()>1){
+ TTreeSRedirector &cstream = *fDebugStreamer;
cstream<<"Tracks"<<
"Tr0.="<<track0<<
"dca="<<dca[i]<<
Int_t lab1=track1->GetLabel();
Char_t c0=track0->GetCircular();
Char_t c1=track1->GetCircular();
+ TTreeSRedirector &cstream = *fDebugStreamer;
cstream<<"V0"<<
"Event="<<eventNr<<
"vertex.="<<&vertex<<
if (AliTPCReconstructor::StreamLevel()>1) {
Int_t lab0=track0->GetLabel();
Int_t lab1=track1->GetLabel();
+ TTreeSRedirector &cstream = *fDebugStreamer;
cstream<<"V02"<<
"Event="<<eventNr<<
"vertex.="<<v0<<
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(0)-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((fParam->GetZLength(0)-TMath::Abs(seed->GetZ())));
+ Int_t type = (seed->GetSector() < fParam->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((fParam->GetZLength(0)-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;
+// seed->SetCurrentSigmaY2(sd2+padlength*padlength*angulary/12.+sresy*sresy);
+// //
+// Float_t sresz = fParam->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(0)-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(0)-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 {
}
-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();
- if(fRow)delete [] fRow;fRow = 0;
- 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();
- if(fRow)delete [] fRow;fRow = 0;
- 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(){
- //
- for (Int_t i = 0; i < fN1; i++)
- fClusters1[i].~AliTPCclusterMI();
- delete [] fClusters1;
- for (Int_t i = 0; i < fN2; i++)
- fClusters2[i].~AliTPCclusterMI();
- delete [] fClusters2;
-}
-
-
-
-//_________________________________________________________________________
-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
- // MvL: Need to call destructors for AliTPCclusterMI, to delete fInfo
- for (Int_t i = 0; i < fN1; i++)
- fClusters1[i].~AliTPCclusterMI();
- delete [] fClusters1; fClusters1=0;
- for (Int_t i = 0; i < fN2; i++)
- fClusters2[i].~AliTPCclusterMI();
- delete [] fClusters2; fClusters2=0;
-
- fN = 0;
- fN1 = 0;
- fN2 = 0;
- //delete[] fClusterArray;
-
- //fClusterArray=0;
-}
-
-
-//___________________________________________________________________
-Int_t AliTPCtrackerMI::AliTPCRow::Find(Double_t z) const {
- //-----------------------------------------------------------------------
- // Return the index of the nearest cluster
- //-----------------------------------------------------------------------
- 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;
- }
- }
- 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);
- Bool_t skipUsed = !(AliTPCReconstructor::GetRecoParam()->GetClusterSharing());
- //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;
- if (skipUsed && c->IsUsed(11)) 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::MakeBitmaps(AliTPCseed *t)
}
}
}
+
+
+
+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);
+}