Int_t AliITSUClusterPix::Compare(const TObject* obj) const
{
// compare clusters accodring to specific mode
- const AliITSUClusterPix* px = dynamic_cast<const AliITSUClusterPix*>(obj);
+ const AliITSUClusterPix* px = (const AliITSUClusterPix*)obj;
float xyz[3],xyz1[3];
if (fgMode & kSortIdLocXZ) { // sorting in local frame
if (GetVolumeId()==px->GetVolumeId()) {
Bool_t AliITSUClusterPix::IsEqual(const TObject* obj) const
{
// compare clusters accodring to specific mode
- const AliITSUClusterPix* px = dynamic_cast<const AliITSUClusterPix*>(obj);
+ const AliITSUClusterPix* px = (const AliITSUClusterPix*)obj;
const Float_t kTol = 1e-5;
float xyz[3],xyz1[3];
if (fgMode & kSortIdLocXZ) { // sorting in local frame
Bool_t AliITSURecoLayer::IsEqual(const TObject* obj) const
{
// check if layers are equal in R
- const AliITSURecoLayer* lr = dynamic_cast<const AliITSURecoLayer*>(obj);
+ const AliITSURecoLayer* lr = (const AliITSURecoLayer*)obj;
return Abs(lr->GetR()-GetR())<1e-6 ? kTRUE : kFALSE;
}
Int_t AliITSURecoLayer::Compare(const TObject* obj) const
{
// compare two layers
- const AliITSURecoLayer* lr = dynamic_cast<const AliITSURecoLayer*>(obj);
+ const AliITSURecoLayer* lr = (const AliITSURecoLayer*)obj;
double dr = GetR() - lr->GetR();
if (Abs(dr)<1e-6) return 0;
return dr>0 ? 1:-1;
delete clFinder;
}
//
- delete[] fClusters;
+ delete fClusters[i];
}
+ delete[] fClusters;
//
delete fGeom;
}
TBranch *lrBranch[fGeom->GetNLayers()];
//
for (int ilr=0;ilr<fGeom->GetNLayers();ilr++) {
+ lrBranch[ilr] = 0;
if (clustersTree) { // do we write clusters tree?
int tp = fGeom->GetLayerDetTypeID(ilr)/AliITSUGeomTGeo::kMaxSegmPerDetType;
if (tp==AliITSUGeomTGeo::kDetTypePix) {
,fSignalAfterElect(0.0)
{
// Default constructor
+ for (int i=kBuffSize;i--;) {
+ fTrack[i] = -2;
+ fHits[i] = -1;
+ fSignal[i] = 0;
+ }
}
//______________________________________________________________________
,fParent(0)
{
// def c-tor
+ ResetFMatrix();
}
//_________________________________________________________________________
,fParent(src.fParent)
{
// def c-tor
+ for (int i=kNFElem;i--;) fFMatrix[i] = src.fFMatrix[i];
}
//_________________________________________________________________________
fChi2Glo = src.fChi2Glo;
fChi2Cl = src.fChi2Cl;
fParent = src.fParent;
+ for (int i=kNFElem;i--;) fFMatrix[i] = src.fFMatrix[i];
AliExternalTrackParam::operator=(src);
return *this;
}
if (IsKilled() != sd->IsKilled()) return kFALSE;
return Abs(GetChi2Glo() - sd->GetChi2Glo())<kTol;
}
+
+//______________________________________________________________________________
+Bool_t AliITSUSeed::PropagateToX(Double_t xk, Double_t b)
+{
+ // Propagate this track to the plane X=xk (cm) in the field "b" (kG)
+ Double_t dx=xk-fX;
+ if (TMath::Abs(dx)<=kAlmost0) return kTRUE;
+
+ Double_t crv=GetC(b);
+ if (TMath::Abs(b) < kAlmost0Field) crv=0.;
+ Double_t x2r = crv*dx;
+ Double_t f1=fP[2], f2=f1 + x2r;
+ if (TMath::Abs(f1) >= kAlmost1) return kFALSE;
+ if (TMath::Abs(f2) >= kAlmost1) return kFALSE;
+ if (TMath::Abs(fP[4])< kAlmost0) return kFALSE;
+
+ Double_t &fP0=fP[0], &fP1=fP[1], &fP2=fP[2], &fP3=fP[3], &fP4=fP[4];
+ Double_t
+ &fC00=fC[0],
+ &fC10=fC[1], &fC11=fC[2],
+ &fC20=fC[3], &fC21=fC[4], &fC22=fC[5],
+ &fC30=fC[6], &fC31=fC[7], &fC32=fC[8], &fC33=fC[9],
+ &fC40=fC[10], &fC41=fC[11], &fC42=fC[12], &fC43=fC[13], &fC44=fC[14];
+
+ Double_t r1=TMath::Sqrt((1.-f1)*(1.+f1)), r2=TMath::Sqrt((1.-f2)*(1.+f2));
+ if (TMath::Abs(r1)<kAlmost0) return kFALSE;
+ if (TMath::Abs(r2)<kAlmost0) return kFALSE;
+
+ fX=xk;
+ double dy2dx = (f1+f2)/(r1+r2);
+ fP0 += dx*dy2dx;
+ if (TMath::Abs(x2r)<0.05) {
+ fP1 += dx*(r2 + f2*dy2dx)*fP3; // Many thanks to P.Hristov !
+ fP2 += x2r;
+ }
+ else {
+ // for small dx/R the linear apporximation of the arc by the segment is OK,
+ // but at large dx/R the error is very large and leads to incorrect Z propagation
+ // angle traversed delta = 2*asin(dist_start_end / R / 2), hence the arc is: R*deltaPhi
+ // The dist_start_end is obtained from sqrt(dx^2+dy^2) = x/(r1+r2)*sqrt(2+f1*f2+r1*r2)
+ // Similarly, the rotation angle in linear in dx only for dx<<R
+ double chord = dx*TMath::Sqrt(1+dy2dx*dy2dx); // distance from old position to new one
+ double rot = 2*TMath::ASin(0.5*chord*crv); // angular difference seen from the circle center
+ fP1 += rot/crv*fP3;
+ fP2 = TMath::Sin(rot + TMath::ASin(fP2));
+ }
+
+ //f = F - 1
+ double r1i = 1./r1;
+ double r2i = 1./r2;
+ double tg1 = f1*r1i;
+ double tg2 = f2*r2i;
+ double v0 = 1. + dy2dx*tg2;
+ double v1 = (r1i+r2i)*(dy2dx*(tg1+tg2)+2);
+ double v2 = (r1i+r2i)*v0;
+ //
+ double f24 = dx*crv/fP4;
+ double f02 = dx*v1;
+ double f04 = dx*v2*f24;
+ double f12 = dx*fP3* (f2*v1+dy2dx-tg2);
+ double f13 = dx*r2*v0;
+ double f14 = dx*f24*fP3*(f2*v2+dy2dx-tg2);
+ //
+ //b = C*ft
+ Double_t b00=f02*fC20 + f04*fC40, b01=f12*fC20 + f14*fC40 + f13*fC30;
+ Double_t b02=f24*fC40;
+ Double_t b10=f02*fC21 + f04*fC41, b11=f12*fC21 + f14*fC41 + f13*fC31;
+ Double_t b12=f24*fC41;
+ Double_t b20=f02*fC22 + f04*fC42, b21=f12*fC22 + f14*fC42 + f13*fC32;
+ Double_t b22=f24*fC42;
+ Double_t b40=f02*fC42 + f04*fC44, b41=f12*fC42 + f14*fC44 + f13*fC43;
+ Double_t b42=f24*fC44;
+ Double_t b30=f02*fC32 + f04*fC43, b31=f12*fC32 + f14*fC43 + f13*fC33;
+ Double_t b32=f24*fC43;
+
+ //a = f*b = f*C*ft
+ Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a02=f02*b22+f04*b42;
+ Double_t a11=f12*b21+f14*b41+f13*b31,a12=f12*b22+f14*b42+f13*b32;
+ Double_t a22=f24*b42;
+
+ //F*C*Ft = C + (b + bt + a)
+ fC00 += b00 + b00 + a00;
+ fC10 += b10 + b01 + a01;
+ fC20 += b20 + b02 + a02;
+ fC30 += b30;
+ fC40 += b40;
+ fC11 += b11 + b11 + a11;
+ fC21 += b21 + b12 + a12;
+ fC31 += b31;
+ fC41 += b41;
+ fC22 += b22 + b22 + a22;
+ fC32 += b32;
+ fC42 += b42;
+ //
+ // update stored transformation matrix F = Fnew*Fold
+ fFMatrix[kF04] += f04 + f24*fFMatrix[kF02];
+ fFMatrix[kF14] += f14 + f24*fFMatrix[kF12];
+ fFMatrix[kF02] += f02;
+ fFMatrix[kF12] += f12;
+ fFMatrix[kF13] += f13;
+ fFMatrix[kF24] += f24;
+ //
+ CheckCovariance();
+
+ return kTRUE;
+}
{
public:
enum {kKilled=BIT(14)};
+ enum {kF02,kF04,kF12,kF13,kF14,kF24, kF44,kNFElem}; // non-trivial elems of propagation matrix
+ enum {kB00,kB01,kB02,kB03,kB04,kB10,kB11,kB12,kB13,kB14, kNBElem}; // non-trivial elems of B matrix (I - K*H)
//
AliITSUSeed();
AliITSUSeed(const AliITSUSeed& src);
virtual Bool_t IsEqual(const TObject* obj) const;
virtual Int_t Compare(const TObject* obj) const;
//
+ // test
+ void ResetFMatrix();
+ void ApplyELoss2FMatrix(Double_t frac, Bool_t beforeProp);
+ Bool_t ApplyMaterialCorrection(Double_t xOverX0, Double_t xTimesRho, Double_t mass, Bool_t beforeProp);
+ Bool_t PropagateToX(Double_t xk, Double_t b);
+ //
protected:
+ //
+ Double_t fFMatrix[kNFElem]; // matxif of propagation from prev layer (non-trivial elements)
+ Double_t fResid[2]; // residuals vector
+ Double_t fCombErrI[3]; // inverse combined error matrix
+ Double_t fBMatix[kNBElem]; // I - K*H matix non-trivial elements
UShort_t fHitsPattern; // bit pattern of hits
UInt_t fClID; // packed cluster info (see AliITSUAux::PackCluster)
Float_t fChi2Glo; // current chi2 global
if (cl>=0) fHitsPattern |= 0x1<<lr;
}
+//_________________________________________________________________________
+inline void AliITSUSeed::ResetFMatrix()
+{
+ // reset transport matrix
+ fFMatrix[kF02] = fFMatrix[kF04] = fFMatrix[kF12] = fFMatrix[kF13] = fFMatrix[kF14] = fFMatrix[kF24] = 0;
+ fFMatrix[kF44] = 1.0; // this element accumulates eloss
+}
+
+//_________________________________________________________________________
+inline Bool_t AliITSUSeed::ApplyMaterialCorrection(Double_t xOverX0, Double_t xTimesRho, Double_t mass, Bool_t beforeProp)
+{
+ // apply material correction and modify transport matrix
+ double pold = Get1P();
+ if (!CorrectForMeanMaterial(xOverX0,xTimesRho,mass)) return kFALSE;
+ ApplyELoss2FMatrix( Get1P()/pold, beforeProp);
+ return kTRUE;
+}
+
+
+//_________________________________________________________________________
+inline void AliITSUSeed::ApplyELoss2FMatrix(Double_t frac, Bool_t beforeProp)
+{
+ // Accounts for the energy loss in the transport matrix
+ // equivalent to multiplying Fmatix by E=diag{1,1,1,1,P4new/P4old}, where P4 is the 1/pt param.
+ // If beforeProp is true, then it is assumed that the eloss was applied before the transport,
+ // i.e. F' = F * E, otherwise, after transport, F' = E * F
+ fFMatrix[kF44] *= frac;
+ if (beforeProp) {
+ fFMatrix[kF04] *= frac;
+ fFMatrix[kF14] *= frac;
+ fFMatrix[kF24] *= frac;
+ }
+}
+
#endif
{
if (source.fBTree) {
fBTree = new TBtree();
- for (int i=fItems->GetEntriesFast();i--;) {
- TObject* obj = fItems->At(i);
- if (obj && ! IsDisabled(obj)) continue;
- RegisterItem(obj);
+ if (fItems) {
+ for (int i=fItems->GetEntriesFast();i--;) {
+ TObject* obj = fItems->At(i);
+ if (obj && ! IsDisabled(obj)) continue;
+ RegisterItem(obj);
+ }
}
}
}
if (doLayer) {
if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
// go via layer to its boundary, applying material correction.
- if (!PropagateTrackTo(seed,xToGo,fCurrMass, lrFr->GetMaxStep(), kFALSE, -1, 0, kTRUE)) return kFALSE;
+ if (!PropagateSeed(seed,xToGo,fCurrMass, lrFr->GetMaxStep())) return kFALSE;
}
}
AliITSURecoLayer* lrTo = fITS->GetLayer( (lFrom+=dir) );
// go the entrance of the layer, assuming no materials in between
double xToGo = dir>0 ? lrTo->GetRMin() : lrTo->GetRMax();
if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
- if (!seed->PropagateTo(xToGo, GetBz())) return kFALSE; // RS: do we need BxByBz?
+ if (!PropagateSeed(seed,xToGo,fCurrMass,100, kFALSE )) return kFALSE;
lrFr = lrTo;
}
return kTRUE;
}
return kTRUE;
}
+
+//______________________________________________________________________________
+Bool_t AliITSUTrackerGlo::PropagateSeed(AliITSUSeed *seed, Double_t xToGo, Double_t mass, Double_t maxStep, Bool_t matCorr)
+{
+ // propagate seed to given x applying material correction if requested
+ const Double_t kEpsilon = 1e-5;
+ Double_t xpos = seed->GetX();
+ Int_t dir = (xpos<xToGo) ? 1:-1;
+ Double_t xyz0[3],xyz1[3],param[7];
+ //
+ if (matCorr) seed->GetXYZ(xyz1); //starting global position
+ while ( (xToGo-xpos)*dir > kEpsilon){
+ Double_t step = dir*TMath::Min(TMath::Abs(xToGo-xpos), maxStep);
+ Double_t x = xpos+step;
+ Double_t bz=GetBz(); // getting the local Bz
+ if (!seed->PropagateToX(x,bz)) return kFALSE;
+ if (matCorr) {
+ xyz0[0]=xyz1[0]; // global pos at the beginning of step
+ xyz0[1]=xyz1[1];
+ xyz0[2]=xyz1[2];
+ seed->GetXYZ(xyz1); // // global pos at the end of step
+ MeanMaterialBudget(xyz0,xyz1,param);
+ Double_t xrho=param[0]*param[4], xx0=param[1];
+ if (dir>0) xrho = -xrho; // outward should be negative
+ if (!seed->ApplyMaterialCorrection(xx0,xrho,mass,kFALSE)) return kFALSE;
+ }
+ xpos = seed->GetX();
+ }
+ return kTRUE;
+}
Int_t GetNSeeds(Int_t lr) const {return fSeedsLr[lr].GetEntriesFast();} //RS TOCHECK
AliITSUSeed* GetSeed(Int_t lr, Int_t sID) const {return (AliITSUSeed*)fSeedsLr[lr].UncheckedAt(sID);} //RS TOCHECK
Bool_t TransportToLayer(AliITSUSeed* seed, Int_t lFrom, Int_t lTo);
+ Bool_t PropagateSeed(AliITSUSeed *seed, Double_t xToGo, Double_t mass, Double_t maxStep=1.0, Bool_t matCorr=kTRUE);
+ //
Bool_t NeedToKill(AliITSUSeed* seed, Int_t flag);
Bool_t GetRoadWidth(AliITSUSeed* seed, int ilrA);
Int_t CheckCluster(AliITSUSeed* seed, Int_t lr, Int_t clID);
void DeleteLastSeedFromPool() {fSeedsPool.RemoveLast();}
void ResetSeedTree(); // RS TOCHECK
//
+
private:
AliITSUTrackerGlo(const AliITSUTrackerGlo&);
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff