04-apr-2004 NvE SetMass() invoked from AliTrack::Set3Momentum to get also the energy...
[u/mrichter/AliRoot.git] / RALICE / Ali3Vector.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 // $Id$
17
18 ///////////////////////////////////////////////////////////////////////////
19 // Class Ali3Vector
20 // Handling of 3-vectors in various reference frames.
21 //
22 // This class is meant to serve as a base class for ALICE objects
23 // that have 3-dimensional vector characteristics.
24 // Error propagation is performed automatically. 
25 //
26 // Note :
27 // ------
28 // Vectors (v), Errors (e) and reference frames (f) are specified via
29 // SetVector(Float_t* v,TString f)
30 // SetErrors(Float_t* e,TString f)
31 // under the following conventions :
32 //
33 // f="car" ==> v in Cartesian coordinates   (x,y,z)
34 // f="sph" ==> v in Spherical coordinates   (r,theta,phi)
35 // f="cyl" ==> v in Cylindrical coordinates (rho,phi,z)
36 //
37 // All angles are in radians.
38 //
39 // Example :
40 // ---------
41 //
42 // Ali3Vector a;
43 // Float_t v[3]={-1,25,7};
44 // Float_t e[3]={0.03,0.5,0.21};
45 // a.SetVector(v,"car");
46 // a.SetErrors(e,"car");
47 // a.Data();
48 //
49 // Float_t vec[3];
50 // Float_t err[3];
51 // a.GetVector(vec,"sph");
52 // a.GetErrors(vec,"sph");
53 //
54 // Ali3Vector b;
55 // Float_t v2[3]={6,-18,33};
56 // Float_t e2[3]={0.19,0.45,0.93};
57 // b.SetVector(v2,"car");
58 // b.SetErrors(e2,"car");
59 //
60 // Float_t dotpro=a.Dot(b);
61 // Float_t doterror=a.GetResultError();
62 //
63 // Ali3Vector c=a.Cross(b);
64 // c.Data("sph");
65 // c.GetVector(vec,"cyl");
66 // c.GetErrors(err,"cyl");
67 //
68 // Float_t norm=c.GetNorm();
69 // Float_t normerror=c.GetResultError();
70 //
71 // c=a+b;
72 // c=a-b;
73 // c=a*5;
74 //
75 //--- Author: Nick van Eijndhoven 30-mar-1999 UU-SAP Utrecht
76 //- Modified: NvE $Date$ UU-SAP Utrecht
77 ///////////////////////////////////////////////////////////////////////////
78
79 #include "Ali3Vector.h"
80 #include "Riostream.h"
81  
82 ClassImp(Ali3Vector) // Class implementation to enable ROOT I/O
83  
84 Ali3Vector::Ali3Vector()
85 {
86 // Creation of an Ali3Vector object and initialisation of parameters
87 // All attributes initialised to 0
88  SetZero();
89 }
90 ///////////////////////////////////////////////////////////////////////////
91 Ali3Vector::~Ali3Vector()
92 {
93 // Destructor to delete dynamically allocated memory
94 }
95 ///////////////////////////////////////////////////////////////////////////
96 Ali3Vector::Ali3Vector(const Ali3Vector& v)
97 {
98 // Copy constructor
99  fV=v.fV;
100  fTheta=v.fTheta;
101  fPhi=v.fPhi;
102  fDx=v.fDx;
103  fDy=v.fDy;
104  fDz=v.fDz;
105  fDresult=v.fDresult;
106 }
107 ///////////////////////////////////////////////////////////////////////////
108 void Ali3Vector::Load(Ali3Vector& q)
109 {
110 // Load all attributes of the input Ali3Vector into this Ali3Vector object.
111  Double_t temp=q.GetResultError();
112  Double_t a[3];
113  q.GetVector(a,"sph");
114  SetVector(a,"sph");
115  q.GetErrors(a,"car");
116  SetErrors(a,"car");
117  fDresult=temp;
118 }
119 ///////////////////////////////////////////////////////////////////////////
120 void Ali3Vector::SetZero()
121 {
122 // (Re)set all attributes to zero.
123  fV=0;
124  fTheta=0;
125  fPhi=0;
126  fDx=0;
127  fDy=0;
128  fDz=0;
129  fDresult=0;
130 }
131 ///////////////////////////////////////////////////////////////////////////
132 void Ali3Vector::SetVector(Double_t* v,TString f)
133 {
134 // Store vector according to reference frame f
135 // All errors will be reset to 0
136  fDx=0;
137  fDy=0;
138  fDz=0;
139  fDresult=0;
140
141  Double_t pi=acos(-1.);
142
143  Int_t frame=0;
144  if (f == "car") frame=1;
145  if (f == "sph") frame=2;
146  if (f == "cyl") frame=3;
147
148  Double_t x,y,z,rho,phi;
149
150  switch (frame)
151  {
152   case 1: // Cartesian coordinates
153    x=v[0];
154    y=v[1];
155    z=v[2];
156    fV=sqrt(x*x+y*y+z*z);
157    fTheta=0;
158    if (fV && fabs(z/fV)<=1.)
159    {
160     fTheta=acos(z/fV);
161    }
162    else
163    {
164     if (z<0.) fTheta=pi;
165    }
166    if (fTheta<0.) fTheta+=2.*pi;
167    fPhi=0;
168    if (x || y) fPhi=atan2(y,x);
169    if (fPhi<0.) fPhi+=2.*pi;
170    break;
171
172   case 2: // Spherical coordinates
173    fV=v[0];
174    fTheta=v[1];
175    fPhi=v[2];
176    break;
177
178   case 3: // Cylindrical coordinates
179    rho=v[0];
180    phi=v[1];
181    z=v[2];
182    fV=sqrt(rho*rho+z*z);
183    fPhi=phi;
184    if (fPhi<0.) fPhi+=2.*pi;
185    fTheta=0;
186    if (fV && fabs(z/fV)<=1.)
187    {
188     fTheta=acos(z/fV);
189    }
190    else
191    {
192     if (z<0.) fTheta=pi;
193    }
194    if (fTheta<0.) fTheta+=2.*pi;
195    break;
196
197   default: // Unsupported reference frame
198    cout << "*Ali3Vector::SetVector* Unsupported frame : " << f.Data() << endl
199         << " Possible frames are 'car', 'sph' and 'cyl'." << endl; 
200    fV=0;
201    fTheta=0;
202    fPhi=0;
203    break;
204  }
205 }
206 ///////////////////////////////////////////////////////////////////////////
207 void Ali3Vector::GetVector(Double_t* v,TString f) const
208 {
209 // Provide vector according to reference frame f
210  Int_t frame=0;
211  if (f == "car") frame=1;
212  if (f == "sph") frame=2;
213  if (f == "cyl") frame=3;
214
215  switch (frame)
216  {
217   case 1: // Cartesian coordinates
218    v[0]=fV*sin(fTheta)*cos(fPhi);
219    v[1]=fV*sin(fTheta)*sin(fPhi);
220    v[2]=fV*cos(fTheta);
221    break;
222
223   case 2: // Spherical coordinates
224    v[0]=fV;
225    v[1]=fTheta;
226    v[2]=fPhi;
227    break;
228
229   case 3: // Cylindrical coordinates
230    v[0]=fV*sin(fTheta);
231    v[1]=fPhi;
232    v[2]=fV*cos(fTheta);
233    break;
234
235   default: // Unsupported reference frame
236    cout << "*Ali3Vector::GetVector* Unsupported frame : " << f.Data() << endl
237         << " Possible frames are 'car', 'sph' and 'cyl'." << endl; 
238    for (Int_t i=0; i<3; i++)
239    {
240     v[i]=0;
241    }
242    break;
243  }
244 }
245 ///////////////////////////////////////////////////////////////////////////
246 void Ali3Vector::SetVector(Float_t* v,TString f)
247 {
248 // Store vector according to reference frame f
249 // All errors will be reset to 0
250  Double_t vec[3];
251  for (Int_t i=0; i<3; i++)
252  {
253   vec[i]=v[i];
254  }
255  SetVector(vec,f);
256 }
257 ///////////////////////////////////////////////////////////////////////////
258 void Ali3Vector::GetVector(Float_t* v,TString f) const
259 {
260 // Provide vector according to reference frame f
261  Double_t vec[3];
262  GetVector(vec,f);
263  for (Int_t i=0; i<3; i++)
264  {
265   v[i]=vec[i];
266  }
267 }
268 ///////////////////////////////////////////////////////////////////////////
269 void Ali3Vector::SetErrors(Double_t* e,TString f)
270 {
271 // Store errors according to reference frame f
272 // The error on scalar results is reset to 0
273  fDresult=0;
274
275  Int_t frame=0;
276  if (f == "car") frame=1;
277  if (f == "sph") frame=2;
278  if (f == "cyl") frame=3;
279
280  Double_t dx2,dy2,dz2,rho;
281
282  switch (frame)
283  {
284   case 1: // Cartesian coordinates
285    fDx=fabs(e[0]);
286    fDy=fabs(e[1]);
287    fDz=fabs(e[2]);
288    break;
289
290   case 2: // Spherical coordinates
291    dx2=pow((cos(fPhi)*sin(fTheta)*e[0]),2)+pow((fV*cos(fTheta)*cos(fPhi)*e[1]),2)
292        +pow((fV*sin(fTheta)*sin(fPhi)*e[2]),2);
293    dy2=pow((sin(fPhi)*sin(fTheta)*e[0]),2)+pow((fV*cos(fTheta)*sin(fPhi)*e[1]),2)
294        +pow((fV*sin(fTheta)*cos(fPhi)*e[2]),2);
295    dz2=pow((cos(fTheta)*e[0]),2)+pow((fV*sin(fTheta)*e[1]),2);
296    fDx=sqrt(dx2);
297    fDy=sqrt(dy2);
298    fDz=sqrt(dz2);
299    break;
300
301   case 3: // Cylindrical coordinates
302    rho=fV*sin(fTheta);
303    dx2=pow((cos(fPhi)*e[0]),2)+pow((rho*sin(fPhi)*e[1]),2);
304    dy2=pow((sin(fPhi)*e[0]),2)+pow((rho*cos(fPhi)*e[1]),2);
305    fDx=sqrt(dx2);
306    fDy=sqrt(dy2);
307    fDz=fabs(e[2]);
308    break;
309
310   default: // Unsupported reference frame
311    cout << "*Ali3Vector::SetErrors* Unsupported frame : " << f.Data() << endl
312         << " Possible frames are 'car', 'sph' and 'cyl'." << endl; 
313    fDx=0;
314    fDy=0;
315    fDz=0;
316    break;
317  }
318 }
319 ///////////////////////////////////////////////////////////////////////////
320 void Ali3Vector::GetErrors(Double_t* e,TString f) const
321 {
322 // Provide errors according to reference frame f
323  Int_t frame=0;
324  if (f == "car") frame=1;
325  if (f == "sph") frame=2;
326  if (f == "cyl") frame=3;
327
328  Double_t pi=acos(-1.);
329
330  Double_t dr2,dtheta2,dphi2,rho,drho2;
331  Double_t v[3]; 
332  Double_t rxy2; // Shorthand for (x*x+y*y)
333
334  switch (frame)
335  {
336   case 1: // Cartesian coordinates
337    e[0]=fDx;
338    e[1]=fDy;
339    e[2]=fDz;
340    break;
341
342   case 2: // Spherical coordinates
343    GetVector(v,"car");
344    rxy2=pow(v[0],2)+pow(v[1],2);
345    if (sqrt(rxy2)<(fV*1e-10)) rxy2=0;
346    if (fV) 
347    {
348     dr2=(pow((v[0]*fDx),2)+pow((v[1]*fDy),2)+pow((v[2]*fDz),2))/(fV*fV);
349    }
350    else
351    {
352     dr2=0;
353    }
354    if (fV)
355    {
356     dtheta2=rxy2*pow(fDz,2)/pow(fV,4);
357     if (v[2] && rxy2)
358     {
359      dtheta2+=rxy2*pow(v[2],2)*(pow((v[0]*fDx),2)+pow((v[1]*fDy),2)) /
360               pow(((pow(v[2],2)*rxy2)+pow(rxy2,2)),2);
361     }
362    }
363    else
364    {
365     dtheta2=0;
366    }
367    if (rxy2)
368    {
369     dphi2=(pow((v[1]*fDx),2)+pow((v[0]*fDy),2))/(pow(rxy2,2));
370    }
371    else
372    {
373     dphi2=0;
374    }
375    e[0]=sqrt(dr2);
376    e[1]=sqrt(dtheta2);
377    if (e[1]>pi) e[1]=pi;
378    e[2]=sqrt(dphi2);
379    if (e[2]>(2.*pi)) e[2]=2.*pi;
380    break;
381
382   case 3: // Cylindrical coordinates
383    GetVector(v,"car");
384    rho=fabs(fV*sin(fTheta));
385    if (rho<(fV*1e-10)) rho=0;
386    if (rho) 
387    {
388     drho2=(pow((v[0]*fDx),2)+pow((v[1]*fDy),2))/(rho*rho);
389    }
390    else
391    {
392     drho2=0;
393    }
394    if (rho)
395    {
396     dphi2=(pow((v[1]*fDx),2)+pow((v[0]*fDy),2))/(pow(rho,4));
397    }
398    else
399    {
400     dphi2=0;
401    }
402    e[0]=sqrt(drho2);
403    e[1]=sqrt(dphi2);
404    if (e[1]>(2.*pi)) e[1]=2.*pi;
405    e[2]=fDz;
406    break;
407
408   default: // Unsupported reference frame
409    cout << "*Ali3Vector::GetErrors* Unsupported frame : " << f.Data() << endl
410         << " Possible frames are 'car', 'sph' and 'cyl'." << endl; 
411    for (Int_t i=0; i<3; i++)
412    {
413     e[i]=0;
414    }
415    break;
416  }
417 }
418 ///////////////////////////////////////////////////////////////////////////
419 void Ali3Vector::SetErrors(Float_t* e,TString f)
420 {
421 // Store errors according to reference frame f
422 // The error on scalar results is reset to 0
423  Double_t vec[3];
424  for (Int_t i=0; i<3; i++)
425  {
426   vec[i]=e[i];
427  }
428  SetErrors(vec,f);
429 }
430 ///////////////////////////////////////////////////////////////////////////
431 void Ali3Vector::GetErrors(Float_t* e,TString f) const
432 {
433 // Provide errors according to reference frame f
434  Double_t vec[3];
435  GetErrors(vec,f);
436  for (Int_t i=0; i<3; i++)
437  {
438   e[i]=vec[i];
439  }
440 }
441 ///////////////////////////////////////////////////////////////////////////
442 void Ali3Vector::Data(TString f) const
443 {
444 // Print vector components according to reference frame f
445  if (f=="car" || f=="sph" || f=="cyl")
446  {
447   Double_t vec[3],err[3];
448   GetVector(vec,f);
449   GetErrors(err,f);
450   cout << " Vector in " << f.Data() << " coordinates : "
451        << vec[0] << " " << vec[1] << " " << vec[2] << endl; 
452   cout << "   Err. in " << f.Data() << " coordinates : "
453        << err[0] << " " << err[1] << " " << err[2] << endl; 
454  }
455  else
456  {
457   cout << " *Ali3Vector::Data* Unsupported frame : " << f.Data() << endl
458        << "  Possible frames are 'car', 'sph' and 'cyl'." << endl; 
459  }
460 }
461 ///////////////////////////////////////////////////////////////////////////
462 Double_t Ali3Vector::GetNorm()
463 {
464 // Provide the norm of the current vector
465 // The error on the scalar result (norm) is updated accordingly
466  Double_t e[3];
467  GetErrors(e,"sph");
468  fDresult=e[0]; 
469  return fV;
470 }
471 ///////////////////////////////////////////////////////////////////////////
472 Double_t Ali3Vector::GetPseudoRapidity()
473 {
474 // Provide the pseudo-rapidity w.r.t. the z-axis.
475 // In other words : eta=-log(tan(theta/2))
476 // The error on the scalar result (pseudo-rap.) is updated accordingly
477  Double_t pi=acos(-1.);
478  Double_t v[3];
479  GetVector(v,"sph");
480  Double_t thetahalf=v[1]/2.;
481  Double_t arg=0;
482  if (v[1]<pi) arg=tan(thetahalf);
483  Double_t eta=9999;
484  if (arg>0) eta=-log(arg);
485  Double_t e[3];
486  GetErrors(e,"sph");
487  Double_t prod=cos(thetahalf)*sin(thetahalf);
488  fDresult=0;
489  if (prod) fDresult=fabs(e[1]/2.*prod);
490  return eta;
491 }
492 ///////////////////////////////////////////////////////////////////////////
493 Double_t Ali3Vector::Dot(Ali3Vector& q)
494 {
495 // Provide the dot product of the current vector with vector q
496 // The error on the scalar result (dotproduct) is updated accordingly
497
498  Double_t dotpro=0;
499
500  if ((this) == &q) // Check for special case v.Dot(v)
501  {
502   Double_t norm=GetNorm();
503   Double_t dnorm=GetResultError();
504   dotpro=pow(norm,2);
505   fDresult=2.*norm*dnorm;
506  }
507  else
508  {
509   Double_t a[3],b[3];
510   Double_t ea[3],eb[3];
511   Double_t d2=0;
512
513   GetVector(a,"car");
514   GetErrors(ea,"car");
515   q.GetVector(b,"car");
516   q.GetErrors(eb,"car");
517   for (Int_t i=0; i<3; i++)
518   {
519    dotpro+=a[i]*b[i];
520    d2+=pow(b[i]*ea[i],2)+pow(a[i]*eb[i],2);
521   }
522   fDresult=sqrt(d2);
523  }
524
525  return dotpro;
526 }
527 ///////////////////////////////////////////////////////////////////////////
528 Double_t Ali3Vector::GetResultError() const
529 {
530 // Provide the error on the result of an operation yielding a scalar
531 // E.g. GetNorm() or Dot()
532  return fDresult;
533 }
534 ///////////////////////////////////////////////////////////////////////////
535 Ali3Vector Ali3Vector::Cross(Ali3Vector& q) const
536 {
537 // Provide the cross product of the current vector with vector q
538 // Error propagation is performed automatically
539  Double_t a[3],b[3],c[3];
540  Double_t ea[3],eb[3],ec[3],d2;
541
542  GetVector(a,"car");
543  GetErrors(ea,"car");
544  q.GetVector(b,"car");
545  q.GetErrors(eb,"car");
546
547  c[0]=a[1]*b[2]-a[2]*b[1];
548  c[1]=a[2]*b[0]-a[0]*b[2];
549  c[2]=a[0]*b[1]-a[1]*b[0];
550
551  d2=pow(b[2]*ea[1],2)+pow(a[1]*eb[2],2)
552    +pow(b[1]*ea[2],2)+pow(a[2]*eb[1],2);
553  ec[0]=sqrt(d2);
554
555  d2=pow(b[0]*ea[2],2)+pow(a[2]*eb[0],2)
556    +pow(b[2]*ea[0],2)+pow(a[0]*eb[2],2);
557  ec[1]=sqrt(d2);
558
559  d2=pow(b[1]*ea[0],2)+pow(a[0]*eb[1],2)
560    +pow(b[0]*ea[1],2)+pow(a[1]*eb[0],2);
561  ec[2]=sqrt(d2);
562
563  Ali3Vector v;
564  v.SetVector(c,"car");
565  v.SetErrors(ec,"car");
566   
567  return v;
568 }
569 ///////////////////////////////////////////////////////////////////////////
570 Ali3Vector Ali3Vector::operator+(Ali3Vector& q) const
571 {
572 // Add vector q to the current vector
573 // Error propagation is performed automatically
574  Double_t a[3],b[3],ea[3],eb[3];
575
576  GetVector(a,"car");
577  GetErrors(ea,"car");
578  q.GetVector(b,"car");
579  q.GetErrors(eb,"car");
580
581  for (Int_t i=0; i<3; i++)
582  {
583   a[i]+=b[i];
584   ea[i]=sqrt(pow(ea[i],2)+pow(eb[i],2));
585  }
586
587  Ali3Vector v;
588  v.SetVector(a,"car");
589  v.SetErrors(ea,"car");
590   
591  return v;
592 }
593 ///////////////////////////////////////////////////////////////////////////
594 Ali3Vector Ali3Vector::operator-(Ali3Vector& q) const
595 {
596 // Subtract vector q from the current vector
597 // Error propagation is performed automatically
598  Double_t a[3],b[3],ea[3],eb[3];
599
600  GetVector(a,"car");
601  GetErrors(ea,"car");
602  q.GetVector(b,"car");
603  q.GetErrors(eb,"car");
604
605  for (Int_t i=0; i<3; i++)
606  {
607   a[i]-=b[i];
608   ea[i]=sqrt(pow(ea[i],2)+pow(eb[i],2));
609  }
610
611  Ali3Vector v;
612  v.SetVector(a,"car");
613  v.SetErrors(ea,"car");
614   
615  return v;
616 }
617 ///////////////////////////////////////////////////////////////////////////
618 Ali3Vector Ali3Vector::operator*(Double_t s) const
619 {
620 // Multiply the current vector with a scalar s.
621 // Error propagation is performed automatically.
622  Double_t a[3],ea[3];
623
624  GetVector(a,"car");
625  GetErrors(ea,"car");
626
627  for (Int_t i=0; i<3; i++)
628  {
629   a[i]*=s;
630   ea[i]*=s;
631  }
632
633  Ali3Vector v;
634  v.SetVector(a,"car");
635  v.SetErrors(ea,"car");
636   
637  return v;
638 }
639 ///////////////////////////////////////////////////////////////////////////
640 Ali3Vector Ali3Vector::operator/(Double_t s) const
641 {
642 // Divide the current vector by a scalar s
643 // Error propagation is performed automatically
644
645  if (fabs(s)<1.e-20) // Protect against division by 0
646  {
647   cout << " *Ali3Vector::/* Division by 0 detected. No action taken." << endl;
648   return *this;
649  }
650  else
651  {
652   Double_t a[3],ea[3];
653
654   GetVector(a,"car");
655   GetErrors(ea,"car");
656
657   for (Int_t i=0; i<3; i++)
658   {
659    a[i]/=s;
660    ea[i]/=s;
661   }
662
663   Ali3Vector v;
664   v.SetVector(a,"car");
665   v.SetErrors(ea,"car");
666   
667   return v;
668  }
669 }
670 ///////////////////////////////////////////////////////////////////////////
671 Ali3Vector& Ali3Vector::operator+=(Ali3Vector& q)
672 {
673 // Add vector q to the current vector
674 // Error propagation is performed automatically
675  Double_t a[3],b[3],ea[3],eb[3];
676
677  GetVector(a,"car");
678  GetErrors(ea,"car");
679  q.GetVector(b,"car");
680  q.GetErrors(eb,"car");
681
682  for (Int_t i=0; i<3; i++)
683  {
684   a[i]+=b[i];
685   ea[i]=sqrt(pow(ea[i],2)+pow(eb[i],2));
686  }
687
688  SetVector(a,"car");
689  SetErrors(ea,"car");
690   
691  return *this;
692 }
693 ///////////////////////////////////////////////////////////////////////////
694 Ali3Vector& Ali3Vector::operator-=(Ali3Vector& q)
695 {
696 // Subtract vector q from the current vector
697 // Error propagation is performed automatically
698  Double_t a[3],b[3],ea[3],eb[3];
699
700  GetVector(a,"car");
701  GetErrors(ea,"car");
702  q.GetVector(b,"car");
703  q.GetErrors(eb,"car");
704
705  for (Int_t i=0; i<3; i++)
706  {
707   a[i]-=b[i];
708   ea[i]=sqrt(pow(ea[i],2)+pow(eb[i],2));
709  }
710
711  SetVector(a,"car");
712  SetErrors(ea,"car");
713   
714  return *this;
715 }
716 ///////////////////////////////////////////////////////////////////////////
717 Ali3Vector& Ali3Vector::operator*=(Double_t s)
718 {
719 // Multiply the current vector with a scalar s
720 // Error propagation is performed automatically
721  Double_t a[3],ea[3];
722
723  GetVector(a,"car");
724  GetErrors(ea,"car");
725
726  for (Int_t i=0; i<3; i++)
727  {
728   a[i]*=s;
729   ea[i]*=s;
730  }
731
732  SetVector(a,"car");
733  SetErrors(ea,"car");
734   
735  return *this;
736 }
737 ///////////////////////////////////////////////////////////////////////////
738 Ali3Vector& Ali3Vector::operator/=(Double_t s)
739 {
740 // Divide the current vector by a scalar s
741 // Error propagation is performed automatically
742
743  if (fabs(s)<1.e-20) // Protect against division by 0
744  {
745   cout << " *Ali3Vector::/=* Division by 0 detected. No action taken." << endl;
746   return *this;
747  }
748  else
749  {
750   Double_t a[3],ea[3];
751
752   GetVector(a,"car");
753   GetErrors(ea,"car");
754
755   for (Int_t i=0; i<3; i++)
756   {
757    a[i]/=s;
758    ea[i]/=s;
759   }
760
761   SetVector(a,"car");
762   SetErrors(ea,"car");
763   
764   return *this;
765  }
766 }
767 ///////////////////////////////////////////////////////////////////////////
768 Ali3Vector Ali3Vector::GetVecTrans() const
769 {
770 // Provide the transverse vector w.r.t. z-axis.
771 // Error propagation is performed automatically
772  Double_t pi=acos(-1.);
773  Double_t a[3],ea[3];
774
775  GetVector(a,"sph");
776  GetErrors(ea,"sph");
777
778  Double_t vt,dvt2;
779  vt=a[0]*sin(a[1]);
780  dvt2=pow((sin(a[1])*ea[0]),2)+pow((a[0]*cos(a[1])*ea[1]),2);
781
782  a[0]=fabs(vt);
783  a[1]=pi/2.;
784
785  ea[0]=sqrt(dvt2);
786  ea[1]=0;
787
788  Ali3Vector v;
789  v.SetVector(a,"sph");
790  v.SetErrors(ea,"sph");
791   
792  return v;
793 }
794 ///////////////////////////////////////////////////////////////////////////
795 Ali3Vector Ali3Vector::GetVecLong() const
796 {
797 // Provide the longitudinal vector w.r.t. z-axis.
798 // Error propagation is performed automatically
799  Double_t pi=acos(-1.);
800  Double_t a[3],ea[3];
801
802  GetVector(a,"sph");
803  GetErrors(ea,"sph");
804
805  Double_t vl,dvl2;
806  vl=a[0]*cos(a[1]);
807  dvl2=pow((cos(a[1])*ea[0]),2)+pow((a[0]*sin(a[1])*ea[1]),2);
808
809  a[0]=fabs(vl);
810  a[1]=0;
811  if (vl<0) a[1]=pi;
812  a[2]=0;
813
814  ea[0]=sqrt(dvl2);
815  ea[1]=0;
816  ea[2]=0;
817
818  Ali3Vector v;
819  v.SetVector(a,"sph");
820  v.SetErrors(ea,"sph");
821   
822  return v;
823 }
824 ///////////////////////////////////////////////////////////////////////////
825 Ali3Vector Ali3Vector::GetPrimed(TRotMatrix* m) const
826 {
827 // Provide vector components (and errors) in a rotated frame.
828 // The orientation of the rotated frame is described by the TRotMatrix
829 // input argument.
830  Ali3Vector v=*this;
831  if (!m) return v;
832
833  Double_t* mat=m->GetMatrix();
834
835  Double_t a[3],aprim[3];
836
837  GetVector(a,"car");
838  aprim[0]=a[0]*mat[0]+a[1]*mat[1]+a[2]*mat[2];
839  aprim[1]=a[0]*mat[3]+a[1]*mat[4]+a[2]*mat[5];
840  aprim[2]=a[0]*mat[6]+a[1]*mat[7]+a[2]*mat[8];
841  v.SetVector(aprim,"car");
842
843  GetErrors(a,"car");
844  aprim[0]=sqrt(pow(a[0]*mat[0],2)+pow(a[1]*mat[1],2)+pow(a[2]*mat[2],2));
845  aprim[1]=sqrt(pow(a[0]*mat[3],2)+pow(a[1]*mat[4],2)+pow(a[2]*mat[5],2));
846  aprim[2]=sqrt(pow(a[0]*mat[6],2)+pow(a[1]*mat[7],2)+pow(a[2]*mat[8],2));
847  v.SetErrors(aprim,"car");
848
849  return v;
850 }
851 ///////////////////////////////////////////////////////////////////////////