+AliTPCseed *AliTPCtrackerMI::MakeSeed(AliTPCseed *track, Float_t r0, Float_t r1, Float_t r2)
+{
+ //
+ //
+ //reseed using track points
+ Int_t p0 = int(r0*track->GetNumberOfClusters()); // point 0
+ Int_t p1 = int(r1*track->GetNumberOfClusters());
+ Int_t p2 = int(r2*track->GetNumberOfClusters()); // last point
+ Int_t pp2=0;
+ Double_t x0[3],x1[3],x2[3];
+ x0[0]=-1;
+ x0[0]=-1;
+ x0[0]=-1;
+
+ // find track position at given ratio of the length
+ Int_t sec0, sec1, sec2;
+ sec0=0;
+ sec1=0;
+ sec2=0;
+ Int_t index=-1;
+ Int_t clindex;
+ for (Int_t i=0;i<160;i++){
+ if (track->fClusterPointer[i]){
+ index++;
+ AliTPCTrackerPoint *trpoint =track->GetTrackPoint(i);
+ if ( (index<p0) || x0[0]<0 ){
+ if (trpoint->GetX()>1){
+ clindex = track->GetClusterIndex2(i);
+ if (clindex>0){
+ x0[0] = trpoint->GetX();
+ x0[1] = trpoint->GetY();
+ x0[2] = trpoint->GetZ();
+ sec0 = ((clindex&0xff000000)>>24)%18;
+ }
+ }
+ }
+
+ if ( (index<p1) &&(trpoint->GetX()>1)){
+ clindex = track->GetClusterIndex2(i);
+ if (clindex>0){
+ x1[0] = trpoint->GetX();
+ x1[1] = trpoint->GetY();
+ x1[2] = trpoint->GetZ();
+ sec1 = ((clindex&0xff000000)>>24)%18;
+ }
+ }
+ if ( (index<p2) &&(trpoint->GetX()>1)){
+ clindex = track->GetClusterIndex2(i);
+ if (clindex>0){
+ x2[0] = trpoint->GetX();
+ x2[1] = trpoint->GetY();
+ x2[2] = trpoint->GetZ();
+ sec2 = ((clindex&0xff000000)>>24)%18;
+ pp2 = i;
+ }
+ }
+ }
+ }
+
+ Double_t alpha, cs,sn, xx2,yy2;
+ //
+ alpha = (sec1-sec2)*fSectors->GetAlpha();
+ cs = TMath::Cos(alpha);
+ sn = TMath::Sin(alpha);
+ xx2= x1[0]*cs-x1[1]*sn;
+ yy2= x1[0]*sn+x1[1]*cs;
+ x1[0] = xx2;
+ x1[1] = yy2;
+ //
+ alpha = (sec0-sec2)*fSectors->GetAlpha();
+ cs = TMath::Cos(alpha);
+ sn = TMath::Sin(alpha);
+ xx2= x0[0]*cs-x0[1]*sn;
+ yy2= x0[0]*sn+x0[1]*cs;
+ x0[0] = xx2;
+ x0[1] = yy2;
+ //
+ //
+ //
+ Double_t x[5],c[15];
+ //
+ x[0]=x2[1];
+ x[1]=x2[2];
+ x[4]=F1(x2[0],x2[1],x1[0],x1[1],x0[0],x0[1]);
+ // if (x[4]>1) return 0;
+ x[2]=F2(x2[0],x2[1],x1[0],x1[1],x0[0],x0[1]);
+ x[3]=F3n(x2[0],x2[1],x0[0],x0[1],x2[2],x0[2],x[4]);
+ //if (TMath::Abs(x[3]) > 2.2) return 0;
+ //if (TMath::Abs(x[2]) > 1.99) return 0;
+ //
+ Double_t sy =0.1, sz =0.1;
+ //
+ Double_t sy1=0.02+track->GetSigmaY2(), sz1=0.02+track->GetSigmaZ2();
+ Double_t sy2=0.01+track->GetSigmaY2(), sz2=0.01+track->GetSigmaZ2();
+ Double_t sy3=0.01+track->GetSigmaY2();
+ //
+ Double_t f40=(F1(x2[0],x2[1]+sy,x1[0],x1[1],x0[0],x0[1])-x[4])/sy;
+ Double_t f42=(F1(x2[0],x2[1],x1[0],x1[1]+sy,x0[0],x0[1])-x[4])/sy;
+ Double_t f43=(F1(x2[0],x2[1],x1[0],x1[1],x0[0],x0[1]+sy)-x[4])/sy;
+ Double_t f20=(F2(x2[0],x2[1]+sy,x1[0],x1[1],x0[0],x0[1])-x[2])/sy;
+ Double_t f22=(F2(x2[0],x2[1],x1[0],x1[1]+sy,x0[0],x0[1])-x[2])/sy;
+ Double_t f23=(F2(x2[0],x2[1],x1[0],x1[1],x0[0],x0[1]+sy)-x[2])/sy;
+ //
+ Double_t f30=(F3(x2[0],x2[1]+sy,x0[0],x0[1],x2[2],x0[2])-x[3])/sy;
+ Double_t f31=(F3(x2[0],x2[1],x0[0],x0[1],x2[2]+sz,x0[2])-x[3])/sz;
+ Double_t f32=(F3(x2[0],x2[1],x0[0],x0[1]+sy,x2[2],x0[2])-x[3])/sy;
+ Double_t f34=(F3(x2[0],x2[1],x0[0],x0[1],x2[2],x0[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(x2[0]);
+ AliTPCseed *seed=new AliTPCseed(0, x, c, x2[0], sec2*fSectors->GetAlpha()+fSectors->GetAlphaShift());
+ // Double_t y0,z0,y1,z1, y2,z2;
+ //seed->GetProlongation(x0[0],y0,z0);
+ // seed->GetProlongation(x1[0],y1,z1);
+ //seed->GetProlongation(x2[0],y2,z2);
+ // seed =0;
+ seed->fLastPoint = pp2;
+ seed->fFirstPoint = pp2;
+
+
+ return seed;
+}
+
+
+AliTPCseed *AliTPCtrackerMI::ReSeed(AliTPCseed *track, Float_t r0, Float_t r1, Float_t r2)
+{
+ //
+ //
+ //reseed using founded clusters
+ //
+ // Find the number of clusters
+ Int_t nclusters = 0;
+ for (Int_t irow=0;irow<160;irow++){
+ if (track->GetClusterIndex(irow)>0) nclusters++;
+ }
+ //
+ Int_t ipos[3];
+ ipos[0] = TMath::Max(int(r0*nclusters),0); // point 0 cluster
+ ipos[1] = TMath::Min(int(r1*nclusters),nclusters-1); //
+ ipos[2] = TMath::Min(int(r2*nclusters),nclusters-1); // last point
+ //
+ //
+ Double_t xyz[3][3];
+ Int_t row[3],sec[3]={0,0,0};
+ //
+ // find track row position at given ratio of the length
+ Int_t index=-1;
+ for (Int_t irow=0;irow<160;irow++){
+ if (track->GetClusterIndex2(irow)<0) continue;
+ index++;
+ for (Int_t ipoint=0;ipoint<3;ipoint++){
+ if (index<=ipos[ipoint]) row[ipoint] = irow;
+ }
+ }
+ //
+ //Get cluster and sector position
+ for (Int_t ipoint=0;ipoint<3;ipoint++){
+ Int_t clindex = track->GetClusterIndex2(row[ipoint]);
+ AliTPCclusterMI * cl = GetClusterMI(clindex);
+ if (cl==0) {
+ //Error("Bug\n");
+ // AliTPCclusterMI * cl = GetClusterMI(clindex);
+ return 0;
+ }
+ sec[ipoint] = ((clindex&0xff000000)>>24)%18;
+ xyz[ipoint][0] = GetXrow(row[ipoint]);
+ xyz[ipoint][1] = cl->GetY();
+ xyz[ipoint][2] = cl->GetZ();
+ }
+ //
+ //
+ // 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(0, x, c, xyz[2][0], sec[2]*fSectors->GetAlpha()+fSectors->GetAlphaShift());
+ seed->fLastPoint = row[2];
+ seed->fFirstPoint = row[2];
+ return seed;
+}
+
+Int_t AliTPCtrackerMI::CheckKinkPoint(AliTPCseed*seed, Float_t th)
+{
+ //
+ //
+ //
+ for (Int_t i=0;i<12;i++) seed->fKinkPoint[i]=0;
+ //
+ if (TMath::Abs(seed->GetC())>0.01) return 0;
+ //
+
+ Float_t x[160], y[160], erry[160], z[160], errz[160];
+ Int_t sec[160];
+ Float_t xt[160], yt[160], zt[160];
+ Int_t i1 = 200;
+ Int_t i2 = 0;
+ Int_t secm = -1;
+ Int_t padm = -1;
+ Int_t middle = seed->GetNumberOfClusters()/2;
+ //
+ //
+ // find central sector, get local cooordinates
+ Int_t count = 0;
+ for (Int_t i=seed->fFirstPoint;i<=seed->fLastPoint;i++) {
+ sec[i]= seed->GetClusterSector(i)%18;
+ x[i] = GetXrow(i);
+ if (sec[i]>=0) {
+ AliTPCclusterMI * cl = seed->fClusterPointer[i];
+ // if (cl==0) cl = GetClusterMI(seed->GetClusterIndex2(i));
+ if (cl==0) {
+ sec[i] = -1;
+ continue;
+ }
+ //
+ //
+ if (i>i2) i2 = i; //last point with cluster
+ if (i2<i1) i1 = i; //first point with cluster
+ y[i] = cl->GetY();
+ z[i] = cl->GetZ();
+ AliTPCTrackerPoint * point = seed->GetTrackPoint(i);
+ xt[i] = x[i];
+ yt[i] = point->GetY();
+ zt[i] = point->GetZ();
+
+ if (point->GetX()>0){
+ erry[i] = point->GetErrY();
+ errz[i] = point->GetErrZ();
+ }
+
+ count++;
+ if (count<middle) {
+ secm = sec[i]; //central sector
+ padm = i; //middle point with cluster
+ }
+ }
+ }
+ //
+ // rotate position to global coordinate system connected to sector at last the point
+ //
+ for (Int_t i=i1;i<=i2;i++){
+ //
+ if (sec[i]<0) continue;
+ Double_t alpha = (sec[i2]-sec[i])*fSectors->GetAlpha();
+ Double_t cs = TMath::Cos(alpha);
+ Double_t sn = TMath::Sin(alpha);
+ Float_t xx2= x[i]*cs+y[i]*sn;
+ Float_t yy2= -x[i]*sn+y[i]*cs;
+ x[i] = xx2;
+ y[i] = yy2;
+ //
+ xx2= xt[i]*cs+yt[i]*sn;
+ yy2= -xt[i]*sn+yt[i]*cs;
+ xt[i] = xx2;
+ yt[i] = yy2;
+
+ }
+ //get "state" vector
+ Double_t xh[5],xm = x[padm];
+ xh[0]=yt[i2];
+ xh[1]=zt[i2];
+ xh[4]=F1(xt[i2],yt[i2],xt[padm],yt[padm],xt[i1],yt[i1]);
+ xh[2]=F2(xt[i2],yt[i2],xt[padm],yt[padm],xt[i1],yt[i1]);
+ xh[3]=F3n(xt[i2],yt[i2],xt[i1],yt[i1],zt[i2],zt[i1],xh[4]);
+ //
+ //
+ for (Int_t i=i1;i<=i2;i++){
+ Double_t yy,zz;
+ if (sec[i]<0) continue;
+ GetProlongation(x[i2], x[i],xh,yy,zz);
+ if (TMath::Abs(y[i]-yy)>4||TMath::Abs(z[i]-zz)>4){
+ //Double_t xxh[5];
+ //xxh[4]=F1old(x[i2],y[i2],x[padm],y[padm],x[i1],y[i1]);
+ //xxh[2]=F2old(x[i2],y[i2],x[padm],y[padm],x[i1],y[i1]);
+ Error("AliTPCtrackerMI::CheckKinkPoint","problem\n");
+ }
+ y[i] = y[i] - yy;
+ z[i] = z[i] - zz;
+ }
+ Float_t dyup[160],dydown[160], dzup[160], dzdown[160];
+ Float_t yup[160], ydown[160], zup[160], zdown[160];
+
+ AliTPCpolyTrack ptrack1,ptrack2;
+ //
+ // derivation up
+ for (Int_t i=i1;i<=i2;i++){
+ AliTPCclusterMI * cl = seed->fClusterPointer[i];
+ if (!cl) continue;
+ if (cl->GetType()<0) continue;
+ if (cl->GetType()>10) continue;
+
+ if (sec[i]>=0){
+ ptrack1.AddPoint(x[i]-xm,y[i],z[i],0.1,0.1);
+ }
+ if (ptrack1.GetN()>4.){
+ ptrack1.UpdateParameters();
+ Double_t ddy,ddz;
+ ptrack1.GetFitDerivation(x[i]-xm,ddy,ddz);
+ Double_t yy,zz;
+ ptrack1.GetFitPoint(x[i]-xm,yy,zz);
+
+ dyup[i] = ddy;
+ dzup[i] = ddz;
+ yup[i] = yy;
+ zup[i] = zz;
+
+ }
+ else{
+ dyup[i]=0.; //not enough points
+ }
+ }
+ //
+ // derivation down
+ for (Int_t i=i2;i>=i1;i--){
+ AliTPCclusterMI * cl = seed->fClusterPointer[i];
+ if (!cl) continue;
+ if (cl->GetType()<0) continue;
+ if (cl->GetType()>10) continue;
+ if (sec[i]>=0){
+ ptrack2.AddPoint(x[i]-xm,y[i],z[i],0.1,0.1);
+ }
+ if (ptrack2.GetN()>4){
+ ptrack2.UpdateParameters();
+ Double_t ddy,ddz;
+ ptrack2.GetFitDerivation(x[i]-xm,ddy,ddz);
+ Double_t yy,zz;
+ ptrack2.GetFitPoint(x[i]-xm,yy,zz);
+
+ dydown[i] = ddy;
+ dzdown[i] = ddz;
+ ydown[i] = yy;
+ zdown[i] = zz;
+ }
+ else{
+ dydown[i]=0.; //not enough points
+ }
+ }
+ //
+ //
+ // find maximal difference of the derivation
+ for (Int_t i=0;i<12;i++) seed->fKinkPoint[i]=0;
+
+
+ for (Int_t i=i1+10;i<i2-10;i++){
+ if ( (TMath::Abs(dydown[i])<0.00000001) || (TMath::Abs(dyup[i])<0.00000001) ||i<30)continue;
+ // printf("%f\t%f\t%f\t%f\t%f\n",x[i],dydown[i],dyup[i],dzdown[i],dzup[i]);
+ //
+ Float_t ddy = TMath::Abs(dydown[i]-dyup[i]);
+ Float_t ddz = TMath::Abs(dzdown[i]-dzup[i]);
+ if ( (ddy+ddz)> th){
+ seed->fKinkPoint[0] = i;
+ seed->fKinkPoint[1] = ddy;
+ seed->fKinkPoint[2] = ddz;
+ th = ddy+ddz;
+ }
+ }
+
+ if (fTreeDebug){
+ //
+ //write information to the debug tree
+ TBranch * br = fTreeDebug->GetBranch("debug");
+ TClonesArray * arr = new TClonesArray("AliTPCTrackPoint2");
+ arr->ExpandCreateFast(i2-i1);
+ br->SetAddress(&arr);
+ //
+ AliTPCclusterMI cldummy;
+ cldummy.SetQ(0);
+ AliTPCTrackPoint2 pdummy;
+ pdummy.GetTPoint().fIsShared = 10;
+ //
+ Double_t alpha = sec[i2]*fSectors->GetAlpha();
+ Double_t cs = TMath::Cos(alpha);
+ Double_t sn = TMath::Sin(alpha);
+
+ for (Int_t i=i1;i<i2;i++){
+ AliTPCTrackPoint2 *trpoint = (AliTPCTrackPoint2*)arr->UncheckedAt(i-i1);
+ //cluster info
+ AliTPCclusterMI * cl0 = seed->fClusterPointer[i];
+ //
+ AliTPCTrackerPoint * point = seed->GetTrackPoint(i);