e3b3e3ed66f0f1bc664e6857bba6107c7e48132b
[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  
81 ClassImp(Ali3Vector) // Class implementation to enable ROOT I/O
82  
83 Ali3Vector::Ali3Vector()
84 {
85 // Creation of an Ali3Vector object and initialisation of parameters
86 // All attributes initialised to 0
87  fV=0;
88  fTheta=0;
89  fPhi=0;
90  fDx=0;
91  fDy=0;
92  fDz=0;
93  fDresult=0;
94 }
95 ///////////////////////////////////////////////////////////////////////////
96 Ali3Vector::~Ali3Vector()
97 {
98 // Destructor to delete dynamically allocated memory
99 }
100 ///////////////////////////////////////////////////////////////////////////
101 void Ali3Vector::SetVector(Double_t* v,TString f)
102 {
103 // Store vector according to reference frame f
104 // All errors will be reset to 0
105  fDx=0;
106  fDy=0;
107  fDz=0;
108  fDresult=0;
109
110  Double_t pi=acos(-1.);
111
112  Int_t frame=0;
113  if (f == "car") frame=1;
114  if (f == "sph") frame=2;
115  if (f == "cyl") frame=3;
116
117  Double_t x,y,z,rho,phi;
118
119  switch (frame)
120  {
121   case 1: // Cartesian coordinates
122    x=v[0];
123    y=v[1];
124    z=v[2];
125    fV=sqrt(x*x+y*y+z*z);
126    fTheta=0;
127    if (fV && fabs(z/fV)<=1.)
128    {
129     fTheta=acos(z/fV);
130    }
131    else
132    {
133     if (z<0.) fTheta=pi;
134    }
135    if (fTheta<0.) fTheta+=2.*pi;
136    fPhi=0;
137    if (x || y) fPhi=atan2(y,x);
138    if (fPhi<0.) fPhi+=2.*pi;
139    break;
140
141   case 2: // Spherical coordinates
142    fV=v[0];
143    fTheta=v[1];
144    fPhi=v[2];
145    break;
146
147   case 3: // Cylindrical coordinates
148    rho=v[0];
149    phi=v[1];
150    z=v[2];
151    fV=sqrt(rho*rho+z*z);
152    fPhi=phi;
153    if (fPhi<0.) fPhi+=2.*pi;
154    fTheta=0;
155    if (fV && fabs(z/fV)<=1.)
156    {
157     fTheta=acos(z/fV);
158    }
159    else
160    {
161     if (z<0.) fTheta=pi;
162    }
163    if (fTheta<0.) fTheta+=2.*pi;
164    break;
165
166   default: // Unsupported reference frame
167    cout << "*Ali3Vector::SetVector* Unsupported frame : " << f << endl
168         << " Possible frames are 'car', 'sph' and 'cyl'." << endl; 
169    fV=0;
170    fTheta=0;
171    fPhi=0;
172    break;
173  }
174 }
175 ///////////////////////////////////////////////////////////////////////////
176 void Ali3Vector::GetVector(Double_t* v,TString f)
177 {
178 // Provide vector according to reference frame f
179  Int_t frame=0;
180  if (f == "car") frame=1;
181  if (f == "sph") frame=2;
182  if (f == "cyl") frame=3;
183
184  switch (frame)
185  {
186   case 1: // Cartesian coordinates
187    v[0]=fV*sin(fTheta)*cos(fPhi);
188    v[1]=fV*sin(fTheta)*sin(fPhi);
189    v[2]=fV*cos(fTheta);
190    break;
191
192   case 2: // Spherical coordinates
193    v[0]=fV;
194    v[1]=fTheta;
195    v[2]=fPhi;
196    break;
197
198   case 3: // Cylindrical coordinates
199    v[0]=fV*sin(fTheta);
200    v[1]=fPhi;
201    v[2]=fV*cos(fTheta);
202    break;
203
204   default: // Unsupported reference frame
205    cout << "*Ali3Vector::GetVector* Unsupported frame : " << f << endl
206         << " Possible frames are 'car', 'sph' and 'cyl'." << endl; 
207    for (Int_t i=0; i<3; i++)
208    {
209     v[i]=0;
210    }
211    break;
212  }
213 }
214 ///////////////////////////////////////////////////////////////////////////
215 void Ali3Vector::SetVector(Float_t* v,TString f)
216 {
217 // Store vector according to reference frame f
218 // All errors will be reset to 0
219  Double_t vec[3];
220  for (Int_t i=0; i<3; i++)
221  {
222   vec[i]=v[i];
223  }
224  SetVector(vec,f);
225 }
226 ///////////////////////////////////////////////////////////////////////////
227 void Ali3Vector::GetVector(Float_t* v,TString f)
228 {
229 // Provide vector according to reference frame f
230  Double_t vec[3];
231  GetVector(vec,f);
232  for (Int_t i=0; i<3; i++)
233  {
234   v[i]=vec[i];
235  }
236 }
237 ///////////////////////////////////////////////////////////////////////////
238 void Ali3Vector::SetErrors(Double_t* e,TString f)
239 {
240 // Store errors according to reference frame f
241 // The error on scalar results is reset to 0
242  fDresult=0;
243
244  Int_t frame=0;
245  if (f == "car") frame=1;
246  if (f == "sph") frame=2;
247  if (f == "cyl") frame=3;
248
249  Double_t dx2,dy2,dz2,rho;
250
251  switch (frame)
252  {
253   case 1: // Cartesian coordinates
254    fDx=fabs(e[0]);
255    fDy=fabs(e[1]);
256    fDz=fabs(e[2]);
257    break;
258
259   case 2: // Spherical coordinates
260    dx2=pow((cos(fPhi)*sin(fTheta)*e[0]),2)+pow((fV*cos(fTheta)*cos(fPhi)*e[1]),2)
261        +pow((fV*sin(fTheta)*sin(fPhi)*e[2]),2);
262    dy2=pow((sin(fPhi)*sin(fTheta)*e[0]),2)+pow((fV*cos(fTheta)*sin(fPhi)*e[1]),2)
263        +pow((fV*sin(fTheta)*cos(fPhi)*e[2]),2);
264    dz2=pow((cos(fTheta)*e[0]),2)+pow((fV*sin(fTheta)*e[1]),2);
265    fDx=sqrt(dx2);
266    fDy=sqrt(dy2);
267    fDz=sqrt(dz2);
268    break;
269
270   case 3: // Cylindrical coordinates
271    rho=fV*sin(fTheta);
272    dx2=pow((cos(fPhi)*e[0]),2)+pow((rho*sin(fPhi)*e[1]),2);
273    dy2=pow((sin(fPhi)*e[0]),2)+pow((rho*cos(fPhi)*e[1]),2);
274    fDx=sqrt(dx2);
275    fDy=sqrt(dy2);
276    fDz=fabs(e[2]);
277    break;
278
279   default: // Unsupported reference frame
280    cout << "*Ali3Vector::SetErrors* Unsupported frame : " << f << endl
281         << " Possible frames are 'car', 'sph' and 'cyl'." << endl; 
282    fDx=0;
283    fDy=0;
284    fDz=0;
285    break;
286  }
287 }
288 ///////////////////////////////////////////////////////////////////////////
289 void Ali3Vector::GetErrors(Double_t* e,TString f)
290 {
291 // Provide errors according to reference frame f
292  Int_t frame=0;
293  if (f == "car") frame=1;
294  if (f == "sph") frame=2;
295  if (f == "cyl") frame=3;
296
297  Double_t pi=acos(-1.);
298
299  Double_t dr2,dtheta2,dphi2,rho,drho2;
300  Double_t v[3]; 
301  Double_t rxy2; // Shorthand for (x*x+y*y)
302
303  switch (frame)
304  {
305   case 1: // Cartesian coordinates
306    e[0]=fDx;
307    e[1]=fDy;
308    e[2]=fDz;
309    break;
310
311   case 2: // Spherical coordinates
312    GetVector(v,"car");
313    rxy2=pow(v[0],2)+pow(v[1],2);
314    if (sqrt(rxy2)<(fV*1e-10)) rxy2=0;
315    if (fV) 
316    {
317     dr2=(pow((v[0]*fDx),2)+pow((v[1]*fDy),2)+pow((v[2]*fDz),2))/(fV*fV);
318    }
319    else
320    {
321     dr2=0;
322    }
323    if (fV)
324    {
325     dtheta2=rxy2*pow(fDz,2)/pow(fV,4);
326     if (v[2] && rxy2)
327     {
328      dtheta2+=rxy2*pow(v[2],2)*(pow((v[0]*fDx),2)+pow((v[1]*fDy),2)) /
329               pow(((pow(v[2],2)*rxy2)+pow(rxy2,2)),2);
330     }
331    }
332    else
333    {
334     dtheta2=0;
335    }
336    if (rxy2)
337    {
338     dphi2=(pow((v[1]*fDx),2)+pow((v[0]*fDy),2))/(pow(rxy2,2));
339    }
340    else
341    {
342     dphi2=0;
343    }
344    e[0]=sqrt(dr2);
345    e[1]=sqrt(dtheta2);
346    if (e[1]>pi) e[1]=pi;
347    e[2]=sqrt(dphi2);
348    if (e[2]>(2.*pi)) e[2]=2.*pi;
349    break;
350
351   case 3: // Cylindrical coordinates
352    GetVector(v,"car");
353    rho=fabs(fV*sin(fTheta));
354    if (rho<(fV*1e-10)) rho=0;
355    if (rho) 
356    {
357     drho2=(pow((v[0]*fDx),2)+pow((v[1]*fDy),2))/(rho*rho);
358    }
359    else
360    {
361     drho2=0;
362    }
363    if (rho)
364    {
365     dphi2=(pow((v[1]*fDx),2)+pow((v[0]*fDy),2))/(pow(rho,4));
366    }
367    else
368    {
369     dphi2=0;
370    }
371    e[0]=sqrt(drho2);
372    e[1]=sqrt(dphi2);
373    if (e[1]>(2.*pi)) e[1]=2.*pi;
374    e[2]=fDz;
375    break;
376
377   default: // Unsupported reference frame
378    cout << "*Ali3Vector::GetErrors* Unsupported frame : " << f << endl
379         << " Possible frames are 'car', 'sph' and 'cyl'." << endl; 
380    for (Int_t i=0; i<3; i++)
381    {
382     e[i]=0;
383    }
384    break;
385  }
386 }
387 ///////////////////////////////////////////////////////////////////////////
388 void Ali3Vector::SetErrors(Float_t* e,TString f)
389 {
390 // Store errors according to reference frame f
391 // The error on scalar results is reset to 0
392  Double_t vec[3];
393  for (Int_t i=0; i<3; i++)
394  {
395   vec[i]=e[i];
396  }
397  SetErrors(vec,f);
398 }
399 ///////////////////////////////////////////////////////////////////////////
400 void Ali3Vector::GetErrors(Float_t* e,TString f)
401 {
402 // Provide errors according to reference frame f
403  Double_t vec[3];
404  GetErrors(vec,f);
405  for (Int_t i=0; i<3; i++)
406  {
407   e[i]=vec[i];
408  }
409 }
410 ///////////////////////////////////////////////////////////////////////////
411 void Ali3Vector::Data(TString f)
412 {
413 // Print vector components according to reference frame f
414  if (f=="car" || f=="sph" || f=="cyl")
415  {
416   Double_t vec[3],err[3];
417   GetVector(vec,f);
418   GetErrors(err,f);
419   cout << " Vector in " << f << " coordinates : "
420        << vec[0] << " " << vec[1] << " " << vec[2] << endl; 
421   cout << "   Err. in " << f << " coordinates : "
422        << err[0] << " " << err[1] << " " << err[2] << endl; 
423  }
424  else
425  {
426   cout << " *Ali3Vector::Data* Unsupported frame : " << f << endl
427        << "  Possible frames are 'car', 'sph' and 'cyl'." << endl; 
428  }
429 }
430 ///////////////////////////////////////////////////////////////////////////
431 Double_t Ali3Vector::GetNorm()
432 {
433 // Provide the norm of the current vector
434 // The error on the scalar result (norm) is updated accordingly
435  Double_t e[3];
436  GetErrors(e,"sph");
437  fDresult=e[0]; 
438  return fV;
439 }
440 ///////////////////////////////////////////////////////////////////////////
441 Double_t Ali3Vector::GetPseudoRapidity()
442 {
443 // Provide the pseudo-rapidity w.r.t. the z-axis.
444 // In other words : eta=-log(tan(theta/2))
445 // The error on the scalar result (pseudo-rap.) is updated accordingly
446  Double_t pi=acos(-1.);
447  Double_t v[3];
448  GetVector(v,"sph");
449  Double_t thetahalf=v[1]/2.;
450  Double_t arg=0;
451  if (v[1]<pi) arg=tan(thetahalf);
452  Double_t eta=9999;
453  if (arg>0) eta=-log(arg);
454  Double_t e[3];
455  GetErrors(e,"sph");
456  Double_t prod=cos(thetahalf)*sin(thetahalf);
457  fDresult=0;
458  if (prod) fDresult=fabs(e[1]/2.*prod);
459  return eta;
460 }
461 ///////////////////////////////////////////////////////////////////////////
462 Double_t Ali3Vector::Dot(Ali3Vector& q)
463 {
464 // Provide the dot product of the current vector with vector q
465 // The error on the scalar result (dotproduct) is updated accordingly
466
467  Double_t dotpro=0;
468
469  if ((this) == &q) // Check for special case v.Dot(v)
470  {
471   Double_t norm=GetNorm();
472   Double_t dnorm=GetResultError();
473   dotpro=pow(norm,2);
474   fDresult=2.*norm*dnorm;
475  }
476  else
477  {
478   Double_t a[3],b[3];
479   Double_t ea[3],eb[3];
480   Double_t d2=0;
481
482   GetVector(a,"car");
483   GetErrors(ea,"car");
484   q.GetVector(b,"car");
485   q.GetErrors(eb,"car");
486   for (Int_t i=0; i<3; i++)
487   {
488    dotpro+=a[i]*b[i];
489    d2+=pow(b[i]*ea[i],2)+pow(a[i]*eb[i],2);
490   }
491   fDresult=sqrt(d2);
492  }
493
494  return dotpro;
495 }
496 ///////////////////////////////////////////////////////////////////////////
497 Double_t Ali3Vector::GetResultError()
498 {
499 // Provide the error on the result of an operation yielding a scalar
500 // E.g. GetNorm() or Dot()
501  return fDresult;
502 }
503 ///////////////////////////////////////////////////////////////////////////
504 Ali3Vector Ali3Vector::Cross(Ali3Vector& q)
505 {
506 // Provide the cross product of the current vector with vector q
507 // Error propagation is performed automatically
508  Double_t a[3],b[3],c[3];
509  Double_t ea[3],eb[3],ec[3],d2;
510
511  GetVector(a,"car");
512  GetErrors(ea,"car");
513  q.GetVector(b,"car");
514  q.GetErrors(eb,"car");
515
516  c[0]=a[1]*b[2]-a[2]*b[1];
517  c[1]=a[2]*b[0]-a[0]*b[2];
518  c[2]=a[0]*b[1]-a[1]*b[0];
519
520  d2=pow(b[2]*ea[1],2)+pow(a[1]*eb[2],2)
521    +pow(b[1]*ea[2],2)+pow(a[2]*eb[1],2);
522  ec[0]=sqrt(d2);
523
524  d2=pow(b[0]*ea[2],2)+pow(a[2]*eb[0],2)
525    +pow(b[2]*ea[0],2)+pow(a[0]*eb[2],2);
526  ec[1]=sqrt(d2);
527
528  d2=pow(b[1]*ea[0],2)+pow(a[0]*eb[1],2)
529    +pow(b[0]*ea[1],2)+pow(a[1]*eb[0],2);
530  ec[2]=sqrt(d2);
531
532  Ali3Vector v;
533  v.SetVector(c,"car");
534  v.SetErrors(ec,"car");
535   
536  return v;
537 }
538 ///////////////////////////////////////////////////////////////////////////
539 Ali3Vector Ali3Vector::operator+(Ali3Vector& q)
540 {
541 // Add vector q to the current vector
542 // Error propagation is performed automatically
543  Double_t a[3],b[3],ea[3],eb[3];
544
545  GetVector(a,"car");
546  GetErrors(ea,"car");
547  q.GetVector(b,"car");
548  q.GetErrors(eb,"car");
549
550  for (Int_t i=0; i<3; i++)
551  {
552   a[i]+=b[i];
553   ea[i]=sqrt(pow(ea[i],2)+pow(eb[i],2));
554  }
555
556  Ali3Vector v;
557  v.SetVector(a,"car");
558  v.SetErrors(ea,"car");
559   
560  return v;
561 }
562 ///////////////////////////////////////////////////////////////////////////
563 Ali3Vector Ali3Vector::operator-(Ali3Vector& q)
564 {
565 // Subtract vector q from the current vector
566 // Error propagation is performed automatically
567  Double_t a[3],b[3],ea[3],eb[3];
568
569  GetVector(a,"car");
570  GetErrors(ea,"car");
571  q.GetVector(b,"car");
572  q.GetErrors(eb,"car");
573
574  for (Int_t i=0; i<3; i++)
575  {
576   a[i]-=b[i];
577   ea[i]=sqrt(pow(ea[i],2)+pow(eb[i],2));
578  }
579
580  Ali3Vector v;
581  v.SetVector(a,"car");
582  v.SetErrors(ea,"car");
583   
584  return v;
585 }
586 ///////////////////////////////////////////////////////////////////////////
587 Ali3Vector Ali3Vector::operator*(Double_t s)
588 {
589 // Multiply the current vector with a scalar s.
590 // Error propagation is performed automatically.
591  Double_t a[3],ea[3];
592
593  GetVector(a,"car");
594  GetErrors(ea,"car");
595
596  for (Int_t i=0; i<3; i++)
597  {
598   a[i]*=s;
599   ea[i]*=s;
600  }
601
602  Ali3Vector v;
603  v.SetVector(a,"car");
604  v.SetErrors(ea,"car");
605   
606  return v;
607 }
608 ///////////////////////////////////////////////////////////////////////////
609 Ali3Vector Ali3Vector::operator/(Double_t s)
610 {
611 // Divide the current vector by a scalar s
612 // Error propagation is performed automatically
613
614  if (fabs(s)<1.e-20) // Protect against division by 0
615  {
616   cout << " *Ali3Vector::/* Division by 0 detected. No action taken." << endl;
617   return *this;
618  }
619  else
620  {
621   Double_t a[3],ea[3];
622
623   GetVector(a,"car");
624   GetErrors(ea,"car");
625
626   for (Int_t i=0; i<3; i++)
627   {
628    a[i]/=s;
629    ea[i]/=s;
630   }
631
632   Ali3Vector v;
633   v.SetVector(a,"car");
634   v.SetErrors(ea,"car");
635   
636   return v;
637  }
638 }
639 ///////////////////////////////////////////////////////////////////////////
640 Ali3Vector& Ali3Vector::operator+=(Ali3Vector& q)
641 {
642 // Add vector q to the current vector
643 // Error propagation is performed automatically
644  Double_t a[3],b[3],ea[3],eb[3];
645
646  GetVector(a,"car");
647  GetErrors(ea,"car");
648  q.GetVector(b,"car");
649  q.GetErrors(eb,"car");
650
651  for (Int_t i=0; i<3; i++)
652  {
653   a[i]+=b[i];
654   ea[i]=sqrt(pow(ea[i],2)+pow(eb[i],2));
655  }
656
657  SetVector(a,"car");
658  SetErrors(ea,"car");
659   
660  return *this;
661 }
662 ///////////////////////////////////////////////////////////////////////////
663 Ali3Vector& Ali3Vector::operator-=(Ali3Vector& q)
664 {
665 // Subtract vector q from the current vector
666 // Error propagation is performed automatically
667  Double_t a[3],b[3],ea[3],eb[3];
668
669  GetVector(a,"car");
670  GetErrors(ea,"car");
671  q.GetVector(b,"car");
672  q.GetErrors(eb,"car");
673
674  for (Int_t i=0; i<3; i++)
675  {
676   a[i]-=b[i];
677   ea[i]=sqrt(pow(ea[i],2)+pow(eb[i],2));
678  }
679
680  SetVector(a,"car");
681  SetErrors(ea,"car");
682   
683  return *this;
684 }
685 ///////////////////////////////////////////////////////////////////////////
686 Ali3Vector& Ali3Vector::operator*=(Double_t s)
687 {
688 // Multiply the current vector with a scalar s
689 // Error propagation is performed automatically
690  Double_t a[3],ea[3];
691
692  GetVector(a,"car");
693  GetErrors(ea,"car");
694
695  for (Int_t i=0; i<3; i++)
696  {
697   a[i]*=s;
698   ea[i]*=s;
699  }
700
701  SetVector(a,"car");
702  SetErrors(ea,"car");
703   
704  return *this;
705 }
706 ///////////////////////////////////////////////////////////////////////////
707 Ali3Vector& Ali3Vector::operator/=(Double_t s)
708 {
709 // Divide the current vector by a scalar s
710 // Error propagation is performed automatically
711
712  if (fabs(s)<1.e-20) // Protect against division by 0
713  {
714   cout << " *Ali3Vector::/=* Division by 0 detected. No action taken." << endl;
715   return *this;
716  }
717  else
718  {
719   Double_t a[3],ea[3];
720
721   GetVector(a,"car");
722   GetErrors(ea,"car");
723
724   for (Int_t i=0; i<3; i++)
725   {
726    a[i]/=s;
727    ea[i]/=s;
728   }
729
730   SetVector(a,"car");
731   SetErrors(ea,"car");
732   
733   return *this;
734  }
735 }
736 ///////////////////////////////////////////////////////////////////////////
737 Ali3Vector Ali3Vector::GetVecTrans()
738 {
739 // Provide the transverse vector w.r.t. z-axis.
740 // Error propagation is performed automatically
741  Double_t pi=acos(-1.);
742  Double_t a[3],ea[3];
743
744  GetVector(a,"sph");
745  GetErrors(ea,"sph");
746
747  Double_t vt,dvt2;
748  vt=a[0]*sin(a[1]);
749  dvt2=pow((sin(a[1])*ea[0]),2)+pow((a[0]*cos(a[1])*ea[1]),2);
750
751  a[0]=fabs(vt);
752  a[1]=pi/2.;
753
754  ea[0]=sqrt(dvt2);
755  ea[1]=0;
756
757  Ali3Vector v;
758  v.SetVector(a,"sph");
759  v.SetErrors(ea,"sph");
760   
761  return v;
762 }
763 ///////////////////////////////////////////////////////////////////////////
764 Ali3Vector Ali3Vector::GetVecLong()
765 {
766 // Provide the longitudinal vector w.r.t. z-axis.
767 // Error propagation is performed automatically
768  Double_t pi=acos(-1.);
769  Double_t a[3],ea[3];
770
771  GetVector(a,"sph");
772  GetErrors(ea,"sph");
773
774  Double_t vl,dvl2;
775  vl=a[0]*cos(a[1]);
776  dvl2=pow((cos(a[1])*ea[0]),2)+pow((a[0]*sin(a[1])*ea[1]),2);
777
778  a[0]=fabs(vl);
779  a[1]=0;
780  if (vl<0) a[1]=pi;
781  a[2]=0;
782
783  ea[0]=sqrt(dvl2);
784  ea[1]=0;
785  ea[2]=0;
786
787  Ali3Vector v;
788  v.SetVector(a,"sph");
789  v.SetErrors(ea,"sph");
790   
791  return v;
792 }
793 ///////////////////////////////////////////////////////////////////////////