00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #include "EvtGenBase/EvtPatches.hh"
00022 #include <iostream>
00023 #include <math.h>
00024 #include <assert.h>
00025 #include "EvtGenBase/EvtComplex.hh"
00026 #include "EvtGenBase/EvtVector4C.hh"
00027 using std::ostream;
00028
00029
00030 EvtVector4C::EvtVector4C(){
00031
00032 v[0]=EvtComplex(0.0); v[1]=EvtComplex(0.0);
00033 v[2]=EvtComplex(0.0); v[3]=EvtComplex(0.0);
00034
00035 }
00036
00037 EvtVector4C::~EvtVector4C(){}
00038
00039 EvtVector4C::EvtVector4C(const EvtComplex& e0,const EvtComplex& e1,
00040 const EvtComplex& e2,const EvtComplex& e3){
00041
00042 v[0]=e0; v[1]=e1; v[2]=e2; v[3]=e3;
00043 }
00044
00045
00046 EvtVector4C rotateEuler(const EvtVector4C& rs,
00047 double alpha,double beta,double gamma){
00048
00049 EvtVector4C tmp(rs);
00050 tmp.applyRotateEuler(alpha,beta,gamma);
00051 return tmp;
00052
00053 }
00054
00055 EvtVector4C boostTo(const EvtVector4C& rs,
00056 const EvtVector4R p4){
00057
00058 EvtVector4C tmp(rs);
00059 tmp.applyBoostTo(p4);
00060 return tmp;
00061
00062 }
00063
00064 EvtVector4C boostTo(const EvtVector4C& rs,
00065 const EvtVector3R boost){
00066
00067 EvtVector4C tmp(rs);
00068 tmp.applyBoostTo(boost);
00069 return tmp;
00070
00071 }
00072
00073 void EvtVector4C::applyBoostTo(const EvtVector4R& p4){
00074
00075 double e=p4.get(0);
00076
00077 EvtVector3R boost(p4.get(1)/e,p4.get(2)/e,p4.get(3)/e);
00078
00079 applyBoostTo(boost);
00080
00081 return;
00082
00083 }
00084
00085 void EvtVector4C::applyBoostTo(const EvtVector3R& boost){
00086
00087 double bx,by,bz,gamma,b2;
00088
00089 bx=boost.get(0);
00090 by=boost.get(1);
00091 bz=boost.get(2);
00092
00093 double bxx=bx*bx;
00094 double byy=by*by;
00095 double bzz=bz*bz;
00096
00097 b2=bxx+byy+bzz;
00098
00099
00100 if (b2==0.0){
00101 return;
00102 }
00103
00104 assert(b2<1.0);
00105
00106 gamma=1.0/sqrt(1-b2);
00107
00108
00109 double gb2=(gamma-1.0)/b2;
00110
00111 double gb2xy=gb2*bx*by;
00112 double gb2xz=gb2*bx*bz;
00113 double gb2yz=gb2*by*bz;
00114
00115 double gbx=gamma*bx;
00116 double gby=gamma*by;
00117 double gbz=gamma*bz;
00118
00119 EvtComplex e2=v[0];
00120 EvtComplex px2=v[1];
00121 EvtComplex py2=v[2];
00122 EvtComplex pz2=v[3];
00123
00124 v[0]=gamma*e2+gbx*px2+gby*py2+gbz*pz2;
00125
00126 v[1]=gbx*e2+gb2*bxx*px2+px2+gb2xy*py2+gb2xz*pz2;
00127
00128 v[2]=gby*e2+gb2*byy*py2+py2+gb2xy*px2+gb2yz*pz2;
00129
00130 v[3]=gbz*e2+gb2*bzz*pz2+pz2+gb2yz*py2+gb2xz*px2;
00131
00132 return;
00133
00134 }
00135
00136 void EvtVector4C::applyRotateEuler(double phi,double theta,double ksi){
00137
00138 double sp=sin(phi);
00139 double st=sin(theta);
00140 double sk=sin(ksi);
00141 double cp=cos(phi);
00142 double ct=cos(theta);
00143 double ck=cos(ksi);
00144
00145 EvtComplex x=( ck*ct*cp-sk*sp)*v[1]+( -sk*ct*cp-ck*sp)*v[2]+st*cp*v[3];
00146 EvtComplex y=( ck*ct*sp+sk*cp)*v[1]+(-sk*ct*sp+ck*cp)*v[2]+st*sp*v[3];
00147 EvtComplex z=-ck*st*v[1]+sk*st*v[2]+ct*v[3];
00148
00149 v[1]=x;
00150 v[2]=y;
00151 v[3]=z;
00152
00153 }
00154
00155
00156 ostream& operator<<(ostream& s, const EvtVector4C& v){
00157
00158 s<<"("<<v.v[0]<<","<<v.v[1]<<","<<v.v[2]<<","<<v.v[3]<<")";
00159
00160 return s;
00161
00162 }
00163