[u/mrichter/AliRoot.git] / TEvtGen / EvtGenBase / EvtDiracSpinor.cxx

//--------------------------------------------------------------------------
//
// Environment:
//      This software is part of the EvtGen package developed jointly
//      for the BaBar and CLEO collaborations.  If you use all or part
//      of it, please give an appropriate acknowledgement.
//
// Copyright Information: See EvtGen/COPYRIGHT
//      Copyright (C) 1998      Caltech, UCSB
//
// Module: EvtDiracSpinor.cc
//
// Description:  Class to describe (EvtDiracParticle) spinors.
//
// Modification history:
//
//    DJL/RYD     September 25, 1996         Module created
//
//------------------------------------------------------------------------
// 
#include "EvtGenBase/EvtPatches.hh"
#include <math.h>
#include <assert.h>
#include "EvtGenBase/EvtDiracSpinor.hh"
#include "EvtGenBase/EvtGammaMatrix.hh"
#include "EvtGenBase/EvtComplex.hh"
#include "EvtGenBase/EvtReport.hh"
#include "EvtGenBase/EvtVector4C.hh"
#include "EvtGenBase/EvtTensor4C.hh"
using std::ostream;


EvtDiracSpinor::~EvtDiracSpinor(){}

EvtDiracSpinor::EvtDiracSpinor(const EvtComplex& sp0,const EvtComplex& sp1,
				    const EvtComplex& sp2,const EvtComplex& sp3){
  set(sp0,sp1,sp2,sp3);
}

void EvtDiracSpinor::set(const EvtComplex& sp0,const EvtComplex& sp1,
			 const EvtComplex& sp2,const EvtComplex& sp3){

  spinor[0]=sp0;spinor[1]=sp1;spinor[2]=sp2;spinor[3]=sp3;
}

void EvtDiracSpinor::set_spinor(int i,const EvtComplex& sp){

  spinor[i]=sp;
}

ostream& operator<<(ostream& s, const EvtDiracSpinor& sp){

  s <<"["<<sp.spinor[0]<<","<<sp.spinor[1]<<","
    <<sp.spinor[2]<<","<<sp.spinor[3]<<"]";
  return s;

}


const EvtComplex& EvtDiracSpinor::get_spinor(int i) const { 
   
  return spinor[i];

}

EvtDiracSpinor rotateEuler(const EvtDiracSpinor& sp,
			   double alpha,double beta,double gamma){

  EvtDiracSpinor tmp(sp);
  tmp.applyRotateEuler(alpha,beta,gamma);
  return tmp;

}

EvtDiracSpinor boostTo(const EvtDiracSpinor& sp,
		       const EvtVector4R p4){

  EvtDiracSpinor tmp(sp);
  tmp.applyBoostTo(p4);
  return tmp;

}

EvtDiracSpinor boostTo(const EvtDiracSpinor& sp,
		       const EvtVector3R boost){

  EvtDiracSpinor tmp(sp);
  tmp.applyBoostTo(boost);
  return tmp;

}

void EvtDiracSpinor::applyBoostTo(const EvtVector4R& p4){

  double e=p4.get(0);

  EvtVector3R boost(p4.get(1)/e,p4.get(2)/e,p4.get(3)/e);

  applyBoostTo(boost);

  return;

}


void EvtDiracSpinor::applyBoostTo(const EvtVector3R& boost) {

  double bx,by,bz,gamma,b2,f1,f2;
  EvtComplex spinorp[4];

  bx=boost.get(0);
  by=boost.get(1);
  bz=boost.get(2);
  b2=bx*bx+by*by+bz*bz;

  if (b2==0.0){
    return;
  }

  assert(b2<1.0);

  gamma=1.0/sqrt(1-b2);
  
  f1=sqrt((gamma+1.0)/2.0);
  f2=f1*gamma/(gamma+1.0);

  spinorp[0]=f1*spinor[0]+f2*bz*spinor[2]+
    f2*EvtComplex(bx,-by)*spinor[3];
  spinorp[1]=f1*spinor[1]+f2*EvtComplex(bx,by)*spinor[2]-
    f2*bz*spinor[3];
  spinorp[2]=f2*bz*spinor[0]+f2*EvtComplex(bx,-by)*spinor[1]+
    f1*spinor[2];
  spinorp[3]=f2*EvtComplex(bx,by)*spinor[0]-
    f2*bz*spinor[1]+f1*spinor[3];
  
  spinor[0]=spinorp[0];
  spinor[1]=spinorp[1];
  spinor[2]=spinorp[2];
  spinor[3]=spinorp[3];

  return;
}

void EvtDiracSpinor::applyRotateEuler(double alpha,double beta,
				      double gamma) {

  EvtComplex retVal[4];
  
  double cb2=cos(0.5*beta);
  double sb2=sin(0.5*beta);
  double capg2=cos(0.5*(alpha+gamma));
  double camg2=cos(0.5*(alpha-gamma));
  double sapg2=sin(0.5*(alpha+gamma));
  double samg2=sin(0.5*(alpha-gamma));

  EvtComplex m11(cb2*capg2,-cb2*sapg2);
  EvtComplex m12(-sb2*camg2,sb2*samg2);
  EvtComplex m21(sb2*camg2,sb2*samg2);
  EvtComplex m22(cb2*capg2,cb2*sapg2);

  retVal[0]=m11*spinor[0]+m12*spinor[1];
  retVal[1]=m21*spinor[0]+m22*spinor[1];
  retVal[2]=m11*spinor[2]+m12*spinor[3];
  retVal[3]=m21*spinor[2]+m22*spinor[3];

  spinor[0]=retVal[0];
  spinor[1]=retVal[1];
  spinor[2]=retVal[2];
  spinor[3]=retVal[3];

  return;
}


EvtDiracSpinor EvtDiracSpinor::conj() const {

  EvtDiracSpinor sp;

  for ( int i=0; i<4; i++)
    sp.set_spinor(i,::conj(spinor[i]));
  
  return sp;
}

EvtVector4C EvtLeptonVACurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){

  //Old code; below is a new specialized code that does it more efficiently.
  //EvtGammaMatrix mat;
  //EvtVector4C temp;
  //mat.va0();
  //temp.set(0,d*(mat*dp));
  //mat.va1();
  //temp.set(1,d*(mat*dp));
  //mat.va2();
  //temp.set(2,d*(mat*dp));
  //mat.va3();
  //temp.set(3,d*(mat*dp));
  //return temp;
 

  EvtComplex u02=::conj(d.spinor[0]-d.spinor[2]);  
  EvtComplex u13=::conj(d.spinor[1]-d.spinor[3]);  

  EvtComplex v02=dp.spinor[0]-dp.spinor[2];
  EvtComplex v13=dp.spinor[1]-dp.spinor[3];

  EvtComplex a=u02*v02;
  EvtComplex b=u13*v13;

  EvtComplex c=u02*v13;
  EvtComplex e=u13*v02;

  return EvtVector4C(a+b,-(c+e),EvtComplex(0,1)*(c-e),b-a);

  
}

EvtVector4C EvtLeptonVCurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){

  EvtVector4C temp;

  // no conjugate here; done in the multiplication
  // yes this is stupid and fooled me to for a long time (ryd)

  temp.set(0,d*(EvtGammaMatrix::v0()*dp));
  temp.set(1,d*(EvtGammaMatrix::v1()*dp));
  temp.set(2,d*(EvtGammaMatrix::v2()*dp));
  temp.set(3,d*(EvtGammaMatrix::v3()*dp));
  
  return temp;
}


EvtVector4C EvtLeptonACurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){

  EvtVector4C temp;

  EvtGammaMatrix mat;

  // no conjugate here; done in the multiplication
  // yes this is stupid and fooled me to for a long time (ryd)

  mat = EvtGammaMatrix::v0()-EvtGammaMatrix::va0();
  temp.set(0,d*(mat*dp));

  mat = EvtGammaMatrix::v1()-EvtGammaMatrix::va1();
  temp.set(1,d*(mat*dp));

  mat = EvtGammaMatrix::v2()-EvtGammaMatrix::va2();
  temp.set(2,d*(mat*dp));

  mat = EvtGammaMatrix::v3()-EvtGammaMatrix::va3();
  temp.set(3,d*(mat*dp));
  
  return temp;
}

EvtComplex EvtLeptonSCurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){

  EvtComplex temp;

  // no conjugate here; done in the multiplication
  // yes this is stupid and fooled me to for a long time (ryd)

  temp=d*(EvtGammaMatrix::g0()*dp);
  
  return temp;
}

EvtComplex EvtLeptonPCurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){

  EvtComplex temp;

  // no conjugate here; done in the multiplication
  // yes this is stupid and fooled me to for a long time (ryd)
  static EvtGammaMatrix m=EvtGammaMatrix::g0()*EvtGammaMatrix::g5();
  temp=d*(m*dp);
  
  return temp;
}

EvtTensor4C EvtLeptonTCurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){

  EvtTensor4C temp;
  temp.zero();
  EvtComplex i2(0,0.5);

  static EvtGammaMatrix mat01=EvtGammaMatrix::g0()*
    (EvtGammaMatrix::g0()*EvtGammaMatrix::g1()-
     EvtGammaMatrix::g1()*EvtGammaMatrix::g0());
  static EvtGammaMatrix mat02=EvtGammaMatrix::g0()*
    (EvtGammaMatrix::g0()*EvtGammaMatrix::g2()-
     EvtGammaMatrix::g2()*EvtGammaMatrix::g0());
  static EvtGammaMatrix mat03=EvtGammaMatrix::g0()*
    (EvtGammaMatrix::g0()*EvtGammaMatrix::g3()-
    EvtGammaMatrix::g3()*EvtGammaMatrix::g0());
  static EvtGammaMatrix mat12=EvtGammaMatrix::g0()*
    (EvtGammaMatrix::g1()*EvtGammaMatrix::g2()-
    EvtGammaMatrix::g2()*EvtGammaMatrix::g1());
  static EvtGammaMatrix mat13=EvtGammaMatrix::g0()*
    (EvtGammaMatrix::g1()*EvtGammaMatrix::g3()-
     EvtGammaMatrix::g3()*EvtGammaMatrix::g1());
  static EvtGammaMatrix mat23=EvtGammaMatrix::g0()*
    (EvtGammaMatrix::g2()*EvtGammaMatrix::g3()-
     EvtGammaMatrix::g3()*EvtGammaMatrix::g2());

 
  temp.set(0,1,i2*(d*(mat01*dp)));
  temp.set(1,0,-temp.get(0,1));

  temp.set(0,2,i2*(d*(mat02*dp)));
  temp.set(2,0,-temp.get(0,2));

  temp.set(0,3,i2*(d*(mat03*dp)));
  temp.set(3,0,-temp.get(0,3));

  temp.set(1,2,i2*(d*(mat12*dp)));
  temp.set(2,1,-temp.get(1,2));

  temp.set(1,3,i2*(d*(mat13*dp)));
  temp.set(3,1,-temp.get(1,3));

  temp.set(2,3,i2*(d*(mat23*dp)));
  temp.set(3,2,-temp.get(2,3));
  
  return temp;
}


EvtDiracSpinor operator*(const EvtComplex& c, const EvtDiracSpinor& d) {
     EvtDiracSpinor result;
     result.spinor[0] = c*d.spinor[0];
     result.spinor[1] = c*d.spinor[1];
     result.spinor[2] = c*d.spinor[2];
     result.spinor[3] = c*d.spinor[3];
     return result;
 }


EvtComplex operator*(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){

  int i;
  EvtComplex temp;

  temp=EvtComplex(0.0,0.0);

  for(i=0;i<4;i++){
    temp+=::conj(d.get_spinor(i))*dp.get_spinor(i);
  }
  return temp;
}


EvtDiracSpinor EvtDiracSpinor::adjoint() const
{
    EvtDiracSpinor d = this->conj(); // first conjugate, then multiply with gamma0
    EvtGammaMatrix g0 = EvtGammaMatrix::g0();
    EvtDiracSpinor result; // automatically initialized to 0

    for (int i=0; i<4; ++i)
        for (int j=0; j<4; ++j)
            result.spinor[i] += d.spinor[j] * g0._gamma[i][j];

    return result;
}
Commit	Line	Data
da0e9ce3	1	//--------------------------------------------------------------------------
	2	//
	3	// Environment:
	4	// This software is part of the EvtGen package developed jointly
	5	// for the BaBar and CLEO collaborations. If you use all or part
	6	// of it, please give an appropriate acknowledgement.
	7	//
	8	// Copyright Information: See EvtGen/COPYRIGHT
	9	// Copyright (C) 1998 Caltech, UCSB
	10	//
	11	// Module: EvtDiracSpinor.cc
	12	//
	13	// Description: Class to describe (EvtDiracParticle) spinors.
	14	//
	15	// Modification history:
	16	//
	17	// DJL/RYD September 25, 1996 Module created
	18	//
	19	//------------------------------------------------------------------------
	20	//
	21	#include "EvtGenBase/EvtPatches.hh"
	22	#include <math.h>
	23	#include <assert.h>
	24	#include "EvtGenBase/EvtDiracSpinor.hh"
	25	#include "EvtGenBase/EvtGammaMatrix.hh"
	26	#include "EvtGenBase/EvtComplex.hh"
	27	#include "EvtGenBase/EvtReport.hh"
	28	#include "EvtGenBase/EvtVector4C.hh"
	29	#include "EvtGenBase/EvtTensor4C.hh"
	30	using std::ostream;
	31
	32
	33	EvtDiracSpinor::~EvtDiracSpinor(){}
	34
	35	EvtDiracSpinor::EvtDiracSpinor(const EvtComplex& sp0,const EvtComplex& sp1,
	36	const EvtComplex& sp2,const EvtComplex& sp3){
	37	set(sp0,sp1,sp2,sp3);
	38	}
	39
	40	void EvtDiracSpinor::set(const EvtComplex& sp0,const EvtComplex& sp1,
	41	const EvtComplex& sp2,const EvtComplex& sp3){
	42
	43	spinor[0]=sp0;spinor[1]=sp1;spinor[2]=sp2;spinor[3]=sp3;
	44	}
	45
	46	void EvtDiracSpinor::set_spinor(int i,const EvtComplex& sp){
	47
	48	spinor[i]=sp;
	49	}
	50
	51	ostream& operator<<(ostream& s, const EvtDiracSpinor& sp){
	52
	53	s <<"["<<sp.spinor[0]<<","<<sp.spinor[1]<<","
	54	<<sp.spinor[2]<<","<<sp.spinor[3]<<"]";
	55	return s;
	56
	57	}
	58
	59
	60	const EvtComplex& EvtDiracSpinor::get_spinor(int i) const {
	61
	62	return spinor[i];
	63
	64	}
65
66	EvtDiracSpinor rotateEuler(const EvtDiracSpinor& sp,
67	double alpha,double beta,double gamma){
68
69	EvtDiracSpinor tmp(sp);
70	tmp.applyRotateEuler(alpha,beta,gamma);
71	return tmp;
72
73	}
74
75	EvtDiracSpinor boostTo(const EvtDiracSpinor& sp,
76	const EvtVector4R p4){
77
78	EvtDiracSpinor tmp(sp);
79	tmp.applyBoostTo(p4);
80	return tmp;
81
82	}
83
84	EvtDiracSpinor boostTo(const EvtDiracSpinor& sp,
85	const EvtVector3R boost){
86
87	EvtDiracSpinor tmp(sp);
88	tmp.applyBoostTo(boost);
89	return tmp;
90
91	}
92
93	void EvtDiracSpinor::applyBoostTo(const EvtVector4R& p4){
94
95	double e=p4.get(0);
96
97	EvtVector3R boost(p4.get(1)/e,p4.get(2)/e,p4.get(3)/e);
98
99	applyBoostTo(boost);
100
101	return;
102
103	}
104
105
106
107	void EvtDiracSpinor::applyBoostTo(const EvtVector3R& boost) {
108
109	double bx,by,bz,gamma,b2,f1,f2;
110	EvtComplex spinorp[4];
111
112	bx=boost.get(0);
113	by=boost.get(1);
114	bz=boost.get(2);
115	b2=bxbx+byby+bz*bz;
116
117	if (b2==0.0){
118	return;
119	}
120
121	assert(b2<1.0);
122
123	gamma=1.0/sqrt(1-b2);
124
125	f1=sqrt((gamma+1.0)/2.0);
126	f2=f1*gamma/(gamma+1.0);
127
128	spinorp[0]=f1spinor[0]+f2bz*spinor[2]+
129	f2EvtComplex(bx,-by)spinor[3];
130	spinorp[1]=f1spinor[1]+f2EvtComplex(bx,by)*spinor[2]-
131	f2bzspinor[3];
132	spinorp[2]=f2bzspinor[0]+f2EvtComplex(bx,-by)spinor[1]+
133	f1*spinor[2];
134	spinorp[3]=f2EvtComplex(bx,by)spinor[0]-
135	f2bzspinor[1]+f1*spinor[3];
136
137	spinor[0]=spinorp[0];
138	spinor[1]=spinorp[1];
139	spinor[2]=spinorp[2];
140	spinor[3]=spinorp[3];
141
142	return;
143	}
144
145	void EvtDiracSpinor::applyRotateEuler(double alpha,double beta,
146	double gamma) {
147
148	EvtComplex retVal[4];
149
150	double cb2=cos(0.5*beta);
151	double sb2=sin(0.5*beta);
152	double capg2=cos(0.5*(alpha+gamma));
153	double camg2=cos(0.5*(alpha-gamma));
154	double sapg2=sin(0.5*(alpha+gamma));
155	double samg2=sin(0.5*(alpha-gamma));
156
157	EvtComplex m11(cb2capg2,-cb2sapg2);
158	EvtComplex m12(-sb2camg2,sb2samg2);
159	EvtComplex m21(sb2camg2,sb2samg2);
160	EvtComplex m22(cb2capg2,cb2sapg2);
161
162	retVal[0]=m11spinor[0]+m12spinor[1];
163	retVal[1]=m21spinor[0]+m22spinor[1];
164	retVal[2]=m11spinor[2]+m12spinor[3];
165	retVal[3]=m21spinor[2]+m22spinor[3];
166
167	spinor[0]=retVal[0];
168	spinor[1]=retVal[1];
169	spinor[2]=retVal[2];
170	spinor[3]=retVal[3];
171
172	return;
173	}
174
175
176
177	EvtDiracSpinor EvtDiracSpinor::conj() const {
178
179	EvtDiracSpinor sp;
180
181	for ( int i=0; i<4; i++)
182	sp.set_spinor(i,::conj(spinor[i]));
183
184	return sp;
185	}
186
187	EvtVector4C EvtLeptonVACurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){
188
189	//Old code; below is a new specialized code that does it more efficiently.
190	//EvtGammaMatrix mat;
191	//EvtVector4C temp;
192	//mat.va0();
193	//temp.set(0,d(matdp));
194	//mat.va1();
195	//temp.set(1,d(matdp));
196	//mat.va2();
197	//temp.set(2,d(matdp));
198	//mat.va3();
199	//temp.set(3,d(matdp));
200	//return temp;
201
202
203	EvtComplex u02=::conj(d.spinor[0]-d.spinor[2]);
204	EvtComplex u13=::conj(d.spinor[1]-d.spinor[3]);
205
206	EvtComplex v02=dp.spinor[0]-dp.spinor[2];
207	EvtComplex v13=dp.spinor[1]-dp.spinor[3];
208
209	EvtComplex a=u02*v02;
210	EvtComplex b=u13*v13;
211
212	EvtComplex c=u02*v13;
213	EvtComplex e=u13*v02;
214
215	return EvtVector4C(a+b,-(c+e),EvtComplex(0,1)*(c-e),b-a);
216
217
218	}
219
220	EvtVector4C EvtLeptonVCurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){
221
222	EvtVector4C temp;
223
224	// no conjugate here; done in the multiplication
225	// yes this is stupid and fooled me to for a long time (ryd)
226
227	temp.set(0,d(EvtGammaMatrix::v0()dp));
228	temp.set(1,d(EvtGammaMatrix::v1()dp));
229	temp.set(2,d(EvtGammaMatrix::v2()dp));
230	temp.set(3,d(EvtGammaMatrix::v3()dp));
231
232	return temp;
233	}
234
235
236	EvtVector4C EvtLeptonACurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){
237
238	EvtVector4C temp;
239
240	EvtGammaMatrix mat;
241
242	// no conjugate here; done in the multiplication
243	// yes this is stupid and fooled me to for a long time (ryd)
244
245	mat = EvtGammaMatrix::v0()-EvtGammaMatrix::va0();
246	temp.set(0,d(matdp));
247
248	mat = EvtGammaMatrix::v1()-EvtGammaMatrix::va1();
249	temp.set(1,d(matdp));
250
251	mat = EvtGammaMatrix::v2()-EvtGammaMatrix::va2();
252	temp.set(2,d(matdp));
253
254	mat = EvtGammaMatrix::v3()-EvtGammaMatrix::va3();
255	temp.set(3,d(matdp));
256
257	return temp;
258	}
259
260	EvtComplex EvtLeptonSCurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){
261
262	EvtComplex temp;
263
264	// no conjugate here; done in the multiplication
265	// yes this is stupid and fooled me to for a long time (ryd)
266
267	temp=d(EvtGammaMatrix::g0()dp);
268
269	return temp;
270	}
271
272	EvtComplex EvtLeptonPCurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){
273
274	EvtComplex temp;
275
276	// no conjugate here; done in the multiplication
277	// yes this is stupid and fooled me to for a long time (ryd)
278	static EvtGammaMatrix m=EvtGammaMatrix::g0()*EvtGammaMatrix::g5();
279	temp=d(mdp);
280
281	return temp;
282	}
283
284	EvtTensor4C EvtLeptonTCurrent(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){
285
286	EvtTensor4C temp;
287	temp.zero();
288	EvtComplex i2(0,0.5);
289
290	static EvtGammaMatrix mat01=EvtGammaMatrix::g0()*
291	(EvtGammaMatrix::g0()*EvtGammaMatrix::g1()-
292	EvtGammaMatrix::g1()*EvtGammaMatrix::g0());
293	static EvtGammaMatrix mat02=EvtGammaMatrix::g0()*
294	(EvtGammaMatrix::g0()*EvtGammaMatrix::g2()-
295	EvtGammaMatrix::g2()*EvtGammaMatrix::g0());
296	static EvtGammaMatrix mat03=EvtGammaMatrix::g0()*
297	(EvtGammaMatrix::g0()*EvtGammaMatrix::g3()-
298	EvtGammaMatrix::g3()*EvtGammaMatrix::g0());
299	static EvtGammaMatrix mat12=EvtGammaMatrix::g0()*
300	(EvtGammaMatrix::g1()*EvtGammaMatrix::g2()-
301	EvtGammaMatrix::g2()*EvtGammaMatrix::g1());
302	static EvtGammaMatrix mat13=EvtGammaMatrix::g0()*
303	(EvtGammaMatrix::g1()*EvtGammaMatrix::g3()-
304	EvtGammaMatrix::g3()*EvtGammaMatrix::g1());
305	static EvtGammaMatrix mat23=EvtGammaMatrix::g0()*
306	(EvtGammaMatrix::g2()*EvtGammaMatrix::g3()-
307	EvtGammaMatrix::g3()*EvtGammaMatrix::g2());
308
309
310	temp.set(0,1,i2(d(mat01*dp)));
311	temp.set(1,0,-temp.get(0,1));
312
313	temp.set(0,2,i2(d(mat02*dp)));
314	temp.set(2,0,-temp.get(0,2));
315
316	temp.set(0,3,i2(d(mat03*dp)));
317	temp.set(3,0,-temp.get(0,3));
318
319	temp.set(1,2,i2(d(mat12*dp)));
320	temp.set(2,1,-temp.get(1,2));
321
322	temp.set(1,3,i2(d(mat13*dp)));
323	temp.set(3,1,-temp.get(1,3));
324
325	temp.set(2,3,i2(d(mat23*dp)));
326	temp.set(3,2,-temp.get(2,3));
327
328	return temp;
329	}
330
331
332	EvtDiracSpinor operator*(const EvtComplex& c, const EvtDiracSpinor& d) {
333	EvtDiracSpinor result;
334	result.spinor[0] = c*d.spinor[0];
335	result.spinor[1] = c*d.spinor[1];
336	result.spinor[2] = c*d.spinor[2];
337	result.spinor[3] = c*d.spinor[3];
338	return result;
339	}
340
341
342	EvtComplex operator*(const EvtDiracSpinor& d,const EvtDiracSpinor& dp){
343
344	int i;
345	EvtComplex temp;
346
347	temp=EvtComplex(0.0,0.0);
348
349	for(i=0;i<4;i++){
350	temp+=::conj(d.get_spinor(i))*dp.get_spinor(i);
351	}
352	return temp;
353	}
354
355
356	EvtDiracSpinor EvtDiracSpinor::adjoint() const
357	{
358	EvtDiracSpinor d = this->conj(); // first conjugate, then multiply with gamma0
359	EvtGammaMatrix g0 = EvtGammaMatrix::g0();
360	EvtDiracSpinor result; // automatically initialized to 0
361
362	for (int i=0; i<4; ++i)
363	for (int j=0; j<4; ++j)
364	result.spinor[i] += d.spinor[j] * g0._gamma[i][j];
365
366	return result;
367	}