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 <stdlib.h>
00025 #include "EvtGenBase/EvtComplex.hh"
00026 #include "EvtGenBase/EvtNeutrinoParticle.hh"
00027 #include "EvtGenBase/EvtVector4R.hh"
00028 #include "EvtGenBase/EvtPDL.hh"
00029 #include "EvtGenBase/EvtReport.hh"
00030
00031 using std::endl;
00032
00033
00034 EvtNeutrinoParticle::~EvtNeutrinoParticle(){}
00035
00036 EvtNeutrinoParticle::EvtNeutrinoParticle(){
00037
00038 return;
00039 }
00040
00041 void EvtNeutrinoParticle::init(EvtId part_n,const EvtVector4R& p4){
00042
00043 _validP4=true;
00044 setp(p4);
00045 setpart_num(part_n);
00046
00047 double e,px,py,pz;
00048 e=p4.get(0);
00049 px=p4.get(1);
00050 py=p4.get(2);
00051 pz=p4.get(3);
00052
00053 if (EvtPDL::getStdHep(part_n)==0){
00054 report(ERROR,"EvtGen") << "Error in EvtNeutrinoParticle::init, part_n="
00055 << part_n.getId()<<endl;
00056 }
00057
00058 if (EvtPDL::getStdHep(part_n)>0){
00059
00060 double beta,alpha,p2,norm;
00061
00062
00063
00064
00065 p2=px*px+py*py+pz*pz;
00066
00067 beta=acos(pz/sqrt(p2));
00068 alpha=atan2(py,px);
00069
00070 norm=sqrt(2*e);
00071
00072 double cosb,sinb,cosa,sina;
00073
00074 cosb=cos(0.5*beta);
00075 sinb=sin(0.5*beta);
00076
00077 cosa=cos(0.5*alpha);
00078 sina=sin(0.5*alpha);
00079
00080 spinor_parent.set(-norm*sinb*EvtComplex(cosa,-sina),
00081 norm*cosb*EvtComplex(cosa,sina),
00082 norm*sinb*EvtComplex(cosa,-sina),
00083 -norm*cosb*EvtComplex(cosa,sina));
00084
00085 }
00086 else{
00087
00088 px=-p4.get(1);
00089 py=-p4.get(2);
00090 pz=-p4.get(3);
00091
00092 double pn,sqrpn;
00093
00094 pn=e;
00095 sqrpn=sqrt(pn-pz);
00096
00097 spinor_parent.set((1.0/sqrpn)*EvtComplex(px,-py),
00098 EvtComplex(sqrpn,0.0),
00099 (-1.0/sqrpn)*EvtComplex(px,-py),
00100 -EvtComplex(sqrpn,0.0));
00101
00102
00103 }
00104
00105 setLifetime();
00106
00107 }
00108
00109
00110 EvtDiracSpinor EvtNeutrinoParticle::spParentNeutrino() const {
00111
00112 return spinor_parent;
00113 }
00114
00115 EvtDiracSpinor EvtNeutrinoParticle::spNeutrino() const {
00116
00117 report(ERROR,"EvtGen") << "Tried to get neutrino spinor in restframe";
00118 report(ERROR,"EvtGen") << "Will terminate execution.";
00119
00120 ::abort();
00121
00122 return spinor_rest;
00123 }
00124
00125
00126 EvtSpinDensity EvtNeutrinoParticle::rotateToHelicityBasis() const{
00127
00128 report(ERROR,"EvtGen") << "rotateToHelicityBasis not implemented for neutrino.";
00129 report(ERROR,"EvtGen") << "Will terminate execution.";
00130
00131 ::abort();
00132
00133 EvtSpinDensity rho;
00134 return rho;
00135
00136 }
00137
00138 EvtSpinDensity EvtNeutrinoParticle::rotateToHelicityBasis(double alpha,
00139 double beta,
00140 double gamma) const{
00141
00142 report(ERROR,"EvtGen") << "rotateToHelicityBasis(alpha,beta,gamma) not implemented for neutrino.";
00143 report(ERROR,"EvtGen") << "Will terminate execution.";
00144
00145 ::abort();
00146
00147 EvtSpinDensity rho;
00148 return rho;
00149
00150
00151 }