]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TRD/AliTRDtrack.cxx
Add tracking without tilted pads as option
[u/mrichter/AliRoot.git] / TRD / AliTRDtrack.cxx
index 983add1394417bf3106a855d6fb38a163e65291e..1cb894dbfc5e2b6263f933c72c127da70059c53b 100644 (file)
@@ -15,6 +15,9 @@
 
 /*
 $Log$
+Revision 1.15  2003/01/27 16:34:49  cblume
+Update of tracking by Sergei and Chuncheng
+
 Revision 1.14  2002/11/07 15:52:09  cblume
 Update of tracking code for tilted pads
 
@@ -88,8 +91,8 @@ AliTRDtrack::AliTRDtrack(const AliTRDcluster *c, UInt_t index,
   fCzy=cc[1];  fCzz=cc[2];
   fCey=cc[3];  fCez=cc[4];  fCee=cc[5];
   fCty=cc[6];  fCtz=cc[7];  fCte=cc[8];  fCtt=cc[9];
-  fCcy=cc[10]; fCcz=cc[11]; fCce=cc[12]; fCct=cc[13]; fCcc=cc[14];
-
+  fCcy=cc[10]; fCcz=cc[11]; fCce=cc[12]; fCct=cc[13]; fCcc=cc[14];  
+  
   fIndex[0]=index;
   SetNumberOfClusters(1);
 
@@ -102,6 +105,7 @@ AliTRDtrack::AliTRDtrack(const AliTRDcluster *c, UInt_t index,
   if(s*s < 1) q *= TMath::Sqrt((1-s*s)/(1+t*t));
 
   fdQdl[0] = q;
+
 }                              
            
 //_____________________________________________________________________________
@@ -109,7 +113,7 @@ AliTRDtrack::AliTRDtrack(const AliTRDtrack& t) : AliKalmanTrack(t) {
   //
   // Copy constructor.
   //
-
+  
   SetLabel(t.GetLabel());
   fSeedLab=t.GetSeedLabel();
 
@@ -121,13 +125,13 @@ AliTRDtrack::AliTRDtrack(const AliTRDtrack& t) : AliKalmanTrack(t) {
   fAlpha=t.fAlpha;
   fX=t.fX;
 
-  fY=t.fY; fZ=t.fZ; fC=t.fC; fE=t.fE; fT=t.fT;
+  fY=t.fY; fZ=t.fZ; fE=t.fE; fT=t.fT; fC=t.fC;
 
   fCyy=t.fCyy;
   fCzy=t.fCzy;  fCzz=t.fCzz;
   fCey=t.fCey;  fCez=t.fCez;  fCee=t.fCee;
   fCty=t.fCty;  fCtz=t.fCtz;  fCte=t.fCte;  fCtt=t.fCtt;
-  fCcy=t.fCcy;  fCcz=t.fCcz;  fCce=t.fCce;  fCct=t.fCct;  fCcc=t.fCcc;
+  fCcy=t.fCcy;  fCcz=t.fCcz;  fCce=t.fCce;  fCct=t.fCct;  fCcc=t.fCcc;  
 
   Int_t n=t.GetNumberOfClusters(); 
   SetNumberOfClusters(n);
@@ -135,6 +139,7 @@ AliTRDtrack::AliTRDtrack(const AliTRDtrack& t) : AliKalmanTrack(t) {
     fIndex[i]=t.fIndex[i];
     fdQdl[i]=t.fdQdl[i];
   }
+
 }                                
 
 //_____________________________________________________________________________
@@ -146,6 +151,7 @@ AliTRDtrack::AliTRDtrack(const AliKalmanTrack& t, Double_t alpha)
 
   SetLabel(t.GetLabel());
   SetChi2(0.);
+  SetMass(t.GetMass());
   SetNumberOfClusters(0);
 
   fdEdx=0;
@@ -162,11 +168,11 @@ AliTRDtrack::AliTRDtrack(const AliKalmanTrack& t, Double_t alpha)
 
   x = GetConvConst();  
 
-  fY=p[0]; 
-  fZ=p[1]; 
+  fY=p[0];
+  fZ=p[1];
   fT=p[3];
   fC=p[4]/x;
-  fE=fX*fC-p[2]; 
+  fE=fC*fX - p[2];   
 
   //Conversion of the covariance matrix
   Double_t c[15]; t.GetExternalCovariance(c);
@@ -177,11 +183,11 @@ AliTRDtrack::AliTRDtrack(const AliKalmanTrack& t, Double_t alpha)
   Double_t c32=fX*c[13] - c[8];
   Double_t c20=fX*c[10] - c[3], c21=fX*c[11] - c[4], c42=fX*c[14] - c[12];
 
-  fCyy=c[0];
-  fCzy=c[1];   fCzz=c[2];
-  fCey=c20;    fCez=c21;    fCee=c22;
-  fCty=c[6];   fCtz=c[7];   fCte=c32;   fCtt=c[9];
-  fCcy=c[10];  fCcz=c[11];  fCce=c42;   fCct=c[13];   fCcc=c[14];
+  fCyy=c[0 ];
+  fCzy=c[1 ];   fCzz=c[2 ];
+  fCey=c20;     fCez=c21;     fCee=c22;
+  fCty=c[6 ];   fCtz=c[7 ];   fCte=c32;   fCtt=c[9 ];
+  fCcy=c[10];   fCcz=c[11];   fCce=c42;   fCct=c[13]; fCcc=c[14];  
 
 }              
 
@@ -195,7 +201,7 @@ void AliTRDtrack::GetExternalParameters(Double_t& xr, Double_t x[5]) const {
      x[1]=GetZ();
      x[2]=GetSnp();
      x[3]=GetTgl();
-     x[4]=fC*GetConvConst();
+     x[4]=Get1Pt();
 }           
 
 //_____________________________________________________________________________
@@ -213,17 +219,20 @@ void AliTRDtrack::GetExternalCovariance(Double_t cc[15]) const {
   cc[1 ]=fCzy;   cc[2 ]=fCzz;
   cc[3 ]=c20;    cc[4 ]=c21;    cc[5 ]=c22;
   cc[6 ]=fCty;   cc[7 ]=fCtz;   cc[8 ]=c32;   cc[9 ]=fCtt;
-  cc[10]=fCcy*a; cc[11]=fCcz*a; cc[12]=c42*a; cc[13]=fCct*a; cc[14]=fCcc*a*a;
+  cc[10]=fCcy*a; cc[11]=fCcz*a; cc[12]=c42*a; cc[13]=fCct*a; cc[14]=fCcc*a*a; 
+  
 }               
                        
 
 //_____________________________________________________________________________
 void AliTRDtrack::GetCovariance(Double_t cc[15]) const {
+
   cc[0]=fCyy;
   cc[1]=fCzy;  cc[2]=fCzz;
   cc[3]=fCey;  cc[4]=fCez;  cc[5]=fCee;
   cc[6]=fCcy;  cc[7]=fCcz;  cc[8]=fCce;  cc[9]=fCcc;
   cc[10]=fCty; cc[11]=fCtz; cc[12]=fCte; cc[13]=fCct; cc[14]=fCtt;
+  
 }    
 
 //_____________________________________________________________________________
@@ -278,17 +287,14 @@ void AliTRDtrack::CookdEdx(Double_t low, Double_t up) {
 }                     
 
 
-
-
 //_____________________________________________________________________________
 Int_t AliTRDtrack::PropagateTo(Double_t xk,Double_t x0,Double_t rho)
 {
   // Propagates a track of particle with mass=pm to a reference plane 
   // defined by x=xk through media of density=rho and radiationLength=x0
 
-  
   if (TMath::Abs(fC*xk - fE) >= 0.99999) {
-    Int_t n=GetNumberOfClusters(); 
+    Int_t n=GetNumberOfClusters();
     if (n>4) cerr<<n<<" AliTRDtrack warning: Propagation failed !\n";
     return 0;
   }
@@ -297,7 +303,7 @@ Int_t AliTRDtrack::PropagateTo(Double_t xk,Double_t x0,Double_t rho)
   Double_t c1=fC*x1 - fE;
   if((c1*c1) > 1) return 0;
   Double_t r1=sqrt(1.- c1*c1);
-  Double_t c2=fC*x2 - fE;
+  Double_t c2=fC*x2 - fE; 
   if((c2*c2) > 1) return 0;
   Double_t r2=sqrt(1.- c2*c2);
 
@@ -311,7 +317,7 @@ Int_t AliTRDtrack::PropagateTo(Double_t xk,Double_t x0,Double_t rho)
   Double_t cr=c1*r2+c2*r1;
   Double_t f12=-dx*fT*(2*cr + cc*(c2*c1/r1-r1 + c1*c2/r2-r2))/(cr*cr);
   Double_t f13= dx*cc/cr;
-  Double_t f14= dx*fT*(cr*xx-cc*(r1*x2-c2*c1*x1/r1+r2*x1-c1*c2*x2/r2))/(cr*cr);
+  Double_t f14=dx*fT*(cr*xx-cc*(r1*x2-c2*c1*x1/r1+r2*x1-c1*c2*x2/r2))/(cr*cr);
 
   //b = C*ft
   Double_t b00=f02*fCey + f04*fCcy, b01=f12*fCey + f14*fCcy + f13*fCty;
@@ -332,23 +338,19 @@ Int_t AliTRDtrack::PropagateTo(Double_t xk,Double_t x0,Double_t rho)
   fCzz += a11 + 2*b11;
   fCez += b21;
   fCtz += b31;
-  fCcz += b41;                  
-
-  fX=x2;
+  fCcz += b41;
 
+  fX=x2;                                                     
 
   //Multiple scattering  ******************
-
   Double_t d=sqrt((x1-fX)*(x1-fX)+(y1-fY)*(y1-fY)+(z1-fZ)*(z1-fZ));
-  Double_t p2=GetPt()*GetPt()*(1.+fT*fT);
-  p2 = TMath::Min(p2,1e+08);  // to avoid division by (1-1) for stiff tracks
+  Double_t p2=(1.+ GetTgl()*GetTgl())/(Get1Pt()*Get1Pt());
   Double_t beta2=p2/(p2 + GetMass()*GetMass());
+  Double_t theta2=14.1*14.1/(beta2*p2*1e6)*d/x0*rho;
 
   Double_t ey=fC*fX - fE, ez=fT;
   Double_t xz=fC*ez, zz1=ez*ez+1, xy=fE+ey;
 
-  Double_t theta2=14.1*14.1/(beta2*p2*1e6)*d/x0*rho;
-
   fCee += (2*ey*ez*ez*fE+1-ey*ey+ez*ez+fE*fE*ez*ez)*theta2;
   fCte += ez*zz1*xy*theta2;
   fCtt += zz1*zz1*theta2;
@@ -356,18 +358,16 @@ Int_t AliTRDtrack::PropagateTo(Double_t xk,Double_t x0,Double_t rho)
   fCct += xz*zz1*theta2;
   fCcc += xz*xz*theta2;
 
-
   //Energy losses************************
-  if (x1 < x2) d=-d;
-  Double_t dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2)) - beta2)*d*rho;
-  //PH  SetLength(GetLength()+d);
+  if((5940*beta2/(1-beta2+1e-10) - beta2) < 0) return 0;
 
-  cc = fC;
+  Double_t dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2+1e-10)) - beta2)*d*rho;
+  if (x1 < x2) dE=-dE;
+  cc=fC;
   fC*=(1.- sqrt(p2+GetMass()*GetMass())/p2*dE);
-  fE+=fX*(fC-cc);
-
-  return 1;        
+  fE+=fX*(fC-cc);    
 
+  return 1;            
 }     
 
 
@@ -376,81 +376,99 @@ Int_t AliTRDtrack::Update(const AliTRDcluster *c, Double_t chisq, UInt_t index,
 {
   // Assignes found cluster to the track and updates track information
 
+  Bool_t fNoTilt = kTRUE;
+  if(TMath::Abs(h01) > 0.003) fNoTilt = kFALSE;
 
-  Double_t r00=c->GetSigmaY2(), r01=0., r11=c->GetSigmaZ2()*100.0;
-
-  r00+=(fCyy+2.0*h01*fCzy+h01*h01*fCzz);
-  r01+=(fCzy+h01*fCzz);
-  r11+=fCzz;
-
+  Double_t r00=c->GetSigmaY2(), r01=0., r11=c->GetSigmaZ2();
+  r00+=fCyy; r01+=fCzy; r11+=fCzz;
   Double_t det=r00*r11 - r01*r01;
   Double_t tmp=r00; r00=r11/det; r11=tmp/det; r01=-r01/det;
 
-//  Double_t k00=fCyy*r00+fCzy*r01, k01=fCyy*r01+fCzy*r11;
-//  Double_t k10=fCzy*r00+fCzz*r01, k11=fCzy*r01+fCzz*r11;
-//  Double_t k20=fCey*r00+fCez*r01, k21=fCey*r01+fCez*r11;
-//  Double_t k30=fCty*r00+fCtz*r01, k31=fCty*r01+fCtz*r11;
-//  Double_t k40=fCcy*r00+fCcz*r01, k41=fCcy*r01+fCcz*r11;
-
-  Double_t k00=fCyy*r00+fCzy*(r01+h01*r00),k01=fCyy*r01+fCzy*(r11+h01*r01);
-  Double_t k10=fCzy*r00+fCzz*(r01+h01*r00),k11=fCzy*r01+fCzz*(r11+h01*r01);
-  Double_t k20=fCey*r00+fCez*(r01+h01*r00),k21=fCey*r01+fCez*(r11+h01*r01);
-  Double_t k30=fCty*r00+fCtz*(r01+h01*r00),k31=fCty*r01+fCtz*(r11+h01*r01);
-  Double_t k40=fCcy*r00+fCcz*(r01+h01*r00),k41=fCcy*r01+fCcz*(r11+h01*r01);
+  Double_t k00=fCyy*r00+fCzy*r01, k01=fCyy*r01+fCzy*r11;
+  Double_t k10=fCzy*r00+fCzz*r01, k11=fCzy*r01+fCzz*r11;
+  Double_t k20=fCey*r00+fCez*r01, k21=fCey*r01+fCez*r11;
+  Double_t k30=fCty*r00+fCtz*r01, k31=fCty*r01+fCtz*r11;
+  Double_t k40=fCcy*r00+fCcz*r01, k41=fCcy*r01+fCcz*r11;
 
   Double_t dy=c->GetY() - fY, dz=c->GetZ() - fZ;
-
-  dy=dy+h01*dz;
-
   Double_t cur=fC + k40*dy + k41*dz, eta=fE + k20*dy + k21*dz;
-  if (TMath::Abs(cur*fX-eta) >= 0.99999) {
-    Int_t n=GetNumberOfClusters(); 
-    if (n>4) cerr<<n<<" AliTRDtrack warning: Filtering failed !\n";
-    return 0;
-  }
-
-  fY += k00*dy + k01*dz;
-  fZ += k10*dy + k11*dz;
-  fE  = eta;
-  fT += k30*dy + k31*dz;
-  fC  = cur;
-
-
-  k01+=h01*k00;
-  k11+=h01*k10;
-  k21+=h01*k20;
-  k31+=h01*k30;
-  k41+=h01*k40;
 
   Double_t c01=fCzy, c02=fCey, c03=fCty, c04=fCcy;
   Double_t c12=fCez, c13=fCtz, c14=fCcz;
 
+  if(fNoTilt) {
+    if (TMath::Abs(cur*fX-eta) >= 0.99999) {
+      Int_t n=GetNumberOfClusters();
+      if (n>4) cerr<<n<<" AliTRDtrack warning: Filtering failed !\n";
+      return 0;
+    }
+    fY += k00*dy + k01*dz;
+    fZ += k10*dy + k11*dz;
+    fE  = eta;
+    fT += k30*dy + k31*dz;
+    fC  = cur;
+  }
+  else {
+    Double_t xu_factor = 100.;  // empirical factor set by C.Xu
+                                // in the first tilt version      
+    r00=c->GetSigmaY2(), r01=0., r11=c->GetSigmaZ2()*xu_factor; 
+    r00+=(fCyy+2.0*h01*fCzy+h01*h01*fCzz);
+    r01+=(fCzy+h01*fCzz);  
+    det=r00*r11 - r01*r01;
+    tmp=r00; r00=r11/det; r11=tmp/det; r01=-r01/det;
+
+    k00=fCyy*r00+fCzy*(r01+h01*r00),k01=fCyy*r01+fCzy*(r11+h01*r01);
+    k10=fCzy*r00+fCzz*(r01+h01*r00),k11=fCzy*r01+fCzz*(r11+h01*r01);
+    k20=fCey*r00+fCez*(r01+h01*r00),k21=fCey*r01+fCez*(r11+h01*r01);
+    k30=fCty*r00+fCtz*(r01+h01*r00),k31=fCty*r01+fCtz*(r11+h01*r01);
+    k40=fCcy*r00+fCcz*(r01+h01*r00),k41=fCcy*r01+fCcz*(r11+h01*r01);  
+
+    dy=c->GetY() - fY; dz=c->GetZ() - fZ; 
+    dy=dy+h01*dz;
+
+    cur=fC + k40*dy + k41*dz; eta=fE + k20*dy + k21*dz;
+    if (TMath::Abs(cur*fX-eta) >= 0.99999) {
+      Int_t n=GetNumberOfClusters();
+      if (n>4) cerr<<n<<" AliTRDtrack warning: Filtering failed !\n";
+      return 0;
+    }                           
+    fY += k00*dy + k01*dz;
+    fZ += k10*dy + k11*dz;
+    fE  = eta;
+    fT += k30*dy + k31*dz;
+    fC  = cur;
+    
+    k01+=h01*k00;
+    k11+=h01*k10;
+    k21+=h01*k20;
+    k31+=h01*k30;
+    k41+=h01*k40;  
+  }
+
   fCyy-=k00*fCyy+k01*fCzy; fCzy-=k00*c01+k01*fCzz;
   fCey-=k00*c02+k01*c12;   fCty-=k00*c03+k01*c13;
   fCcy-=k00*c04+k01*c14;
-
+  
   fCzz-=k10*c01+k11*fCzz;
   fCez-=k10*c02+k11*c12;   fCtz-=k10*c03+k11*c13;
   fCcz-=k10*c04+k11*c14;
-
+  
   fCee-=k20*c02+k21*c12;   fCte-=k20*c03+k21*c13;
   fCce-=k20*c04+k21*c14;
-
+  
   fCtt-=k30*c03+k31*c13;
-  fCct-=k40*c03+k41*c13;        
-
-  fCcc-=k40*c04+k41*c14;
+  fCct-=k40*c03+k41*c13;
+  
+  fCcc-=k40*c04+k41*c14;                 
 
-  Int_t n=GetNumberOfClusters();  
+  Int_t n=GetNumberOfClusters();
   fIndex[n]=index;
-  SetNumberOfClusters(n+1);  
-
-  SetChi2(GetChi2()+chisq); 
+  SetNumberOfClusters(n+1);
 
+  SetChi2(GetChi2()+chisq);
   //  cerr<<"in update: fIndex["<<fN<<"] = "<<index<<endl;
 
-  return 1;
-
+  return 1;     
 }                     
 
 
@@ -461,39 +479,35 @@ Int_t AliTRDtrack::Rotate(Double_t alpha)
 
   fAlpha += alpha;
   if (fAlpha<-TMath::Pi()) fAlpha += 2*TMath::Pi();
-  if (fAlpha>=TMath::Pi()) fAlpha -= 2*TMath::Pi();    
+  if (fAlpha>=TMath::Pi()) fAlpha -= 2*TMath::Pi();
 
   Double_t x1=fX, y1=fY;
   Double_t ca=cos(alpha), sa=sin(alpha);
   Double_t r1=fC*fX - fE;
 
-  if (TMath::Abs(r1) >= 0.99999) {
-    Int_t n=GetNumberOfClusters(); 
-    if (n>4) cerr<<n<<" AliTRDtrack warning: Rotation failed !\n";
-    return 0;
-  }
-
   fX = x1*ca + y1*sa;
-  fY=-x1*sa + y1*ca;
+  fY =-x1*sa + y1*ca;
+  if((r1*r1) > 1) return 0;
   fE=fE*ca + (fC*y1 + sqrt(1.- r1*r1))*sa;
 
   Double_t r2=fC*fX - fE;
   if (TMath::Abs(r2) >= 0.99999) {
-    Int_t n=GetNumberOfClusters(); 
+    Int_t n=GetNumberOfClusters();
     if (n>4) cerr<<n<<" AliTRDtrack warning: Rotation failed !\n";
     return 0;
   }
 
+  if((r2*r2) > 1) return 0;
   Double_t y0=fY + sqrt(1.- r2*r2)/fC;
   if ((fY-y0)*fC >= 0.) {
-    Int_t n=GetNumberOfClusters(); 
+    Int_t n=GetNumberOfClusters();
     if (n>4) cerr<<n<<" AliTRDtrack warning: Rotation failed !!!\n";
     return 0;
   }
 
   //f = F - 1
   Double_t f00=ca-1,    f24=(y1 - r1*x1/sqrt(1.- r1*r1))*sa,
-           f20=fC*sa, f22=(ca + sa*r1/sqrt(1.- r1*r1))-1;
+           f20=fC*sa,  f22=(ca + sa*r1/sqrt(1.- r1*r1))-1;
 
   //b = C*ft
   Double_t b00=fCyy*f00, b02=fCyy*f20+fCcy*f24+fCey*f22;
@@ -505,8 +519,6 @@ Int_t AliTRDtrack::Rotate(Double_t alpha)
   //a = f*b = f*C*ft
   Double_t a00=f00*b00, a02=f00*b02, a22=f20*b02+f24*b42+f22*b22;
 
-  // *** Double_t dy2=fCyy;  
-          
   //F*C*Ft = C + (a + b + bt)
   fCyy += a00 + 2*b00;
   fCzy += b10;
@@ -518,10 +530,7 @@ Int_t AliTRDtrack::Rotate(Double_t alpha)
   fCee += a22 + 2*b22;
   fCce += b42;
 
-  // *** fCyy+=dy2*sa*sa*r1*r1/(1.- r1*r1);
-  // *** fCzz+=d2y*sa*sa*fT*fT/(1.- r1*r1);   
-
-  return 1;
+  return 1;                            
 }                         
 
 
@@ -529,24 +538,32 @@ Int_t AliTRDtrack::Rotate(Double_t alpha)
 Double_t AliTRDtrack::GetPredictedChi2(const AliTRDcluster *c, Double_t h01) const
 {
   
-  Double_t r00=c->GetSigmaY2(), r01=0., r11=c->GetSigmaZ2();
-
-  r00+=fCyy; r01+=fCzy; r11+=fCzz;
-
-  Double_t det=r00*r11 - r01*r01;
-  if (TMath::Abs(det) < 1.e-10) {
-    Int_t n=GetNumberOfClusters(); 
-    if (n>4) cerr<<n<<" AliTRDtrack warning: Singular matrix !\n";
-    return 1e10;
+  Bool_t fNoTilt = kTRUE;
+  if(TMath::Abs(h01) > 0.003) fNoTilt = kFALSE;
+  Double_t chi2, dy, r00, r01, r11;
+
+  if(fNoTilt) {
+    dy=c->GetY() - fY;
+    r00=c->GetSigmaY2();    
+    chi2 = (dy*dy)/r00;    
   }
-  Double_t tmp=r00; r00=r11; r11=tmp; r01=-r01;
-
-  Double_t dy=c->GetY() - fY, dz=c->GetZ() - fZ;
-
-  dy=dy+h01*dz;
-
-  return (dy*r00*dy + 2*r01*dy*dz + dz*r11*dz)/det; 
+  else {
+    r00=c->GetSigmaY2(); r01=0.; r11=c->GetSigmaZ2();
+    r00+=fCyy; r01+=fCzy; r11+=fCzz;
+
+    Double_t det=r00*r11 - r01*r01;
+    if (TMath::Abs(det) < 1.e-10) {
+      Int_t n=GetNumberOfClusters(); 
+      if (n>4) cerr<<n<<" AliTRDtrack warning: Singular matrix !\n";
+      return 1e10;
+    }
+    Double_t tmp=r00; r00=r11; r11=tmp; r01=-r01;
+    Double_t dy=c->GetY() - fY, dz=c->GetZ() - fZ;
+    dy=dy+h01*dz;
 
+    chi2 = (dy*r00*dy + 2*r01*dy*dz + dz*r11*dz)/det; 
+  }
+  return chi2;
 }      
 
 
@@ -592,9 +609,9 @@ void AliTRDtrack::ResetCovariance() {
   //
 
   fCyy*=10.;
-  fCzy=0.;   fCzz*=10.;
-  fCey=0.;   fCez=0.;   fCee*=10.;
-  fCty=0.;   fCtz=0.;   fCte=0.;   fCtt*=10.;
-  fCcy=0.;   fCcz=0.;   fCce=0.;   fCct=0.;   fCcc*=10.;
-
+  fCzy=0.;  fCzz*=10.;
+  fCey=0.;  fCez=0.;  fCee*=10.;
+  fCty=0.;  fCtz=0.;  fCte=0.;  fCtt*=10.;
+  fCcy=0.;  fCcz=0.;  fCce=0.;  fCct=0.;  fCcc*=10.;  
 }                                                         
+