[u/mrichter/AliRoot.git] / TEvtGen / EvtGenBase / EvtPhotonParticle.cpp

//--------------------------------------------------------------------------
//
// 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: EvtPhotonParticle.cc
//
// Description: Class to describe massless vectors
//
// Modification history:
//
//    DJL/RYD   September 25, 1996           Module created
//
//------------------------------------------------------------------------
// 
#include "EvtGenBase/EvtPatches.hh"
#include <stdlib.h>
#include <iostream>
#include <math.h>
#include "EvtGenBase/EvtComplex.hh"
#include "EvtGenBase/EvtPhotonParticle.hh"
#include "EvtGenBase/EvtVector4C.hh"
#include "EvtGenBase/EvtReport.hh"
using std::endl;

EvtPhotonParticle::EvtPhotonParticle() {
}

EvtPhotonParticle::~EvtPhotonParticle(){}

void EvtPhotonParticle::init(EvtId part_n,const EvtVector4R& p4){

  init(part_n,p4.get(0),p4.get(1),
       p4.get(2),p4.get(3));
}

void EvtPhotonParticle::init(EvtId part_n,double e,double px,double py,double pz){

  _validP4=true;
  setp(e,px,py,pz);
  setpart_num(part_n);

  setLifetime();
  
  //defere calculation of basis vectors untill they are needed!
  _evalBasis=0;

}


EvtVector4C EvtPhotonParticle::epsParentPhoton(int i){

  if (!_evalBasis){

    _evalBasis=1;
    eps1.set(EvtComplex(0.0,0.0), EvtComplex(-1.0/sqrt(2.0),0.0),
	     EvtComplex(0.0,-1.0/sqrt(2.0)), EvtComplex(0.0,0.0));
    eps2.set(EvtComplex(0.0,0.0), EvtComplex(1.0/sqrt(2.0),0.0),
	     EvtComplex(0.0,-1.0/sqrt(2.0)), EvtComplex(0.0,0.0));
    
    // These are for photon along z axis.  Rotate to get 
    // correct direction...
    
    double phi,theta;

    EvtVector4R p=this->getP4();

    double px=p.get(1);
    double py=p.get(2);
    double pz=p.get(3);

    phi = atan2(py,px);
    theta = acos(pz/sqrt(px*px+py*py+pz*pz));
    eps1.applyRotateEuler(phi,theta,-phi);
    eps2.applyRotateEuler(phi,theta,-phi);
    
  }


  EvtVector4C temp;
  
  switch(i) {

  case 0:
    temp=eps1;
    break;
  case 1:
    temp=eps2;
    break;
  default:
    report(ERROR,"EvtGen") << "EvtPhotonParticle.cc: Asked "
			   << "for state:"<<i<<endl;
    ::abort();
    break;
  }

  return temp;
}

EvtVector4C EvtPhotonParticle::epsPhoton(int ){

  report(ERROR,"EvtGen") << "EvtPhotonParticle.cc: Can not get "
			 << "state in photons restframe."<<endl;;
  ::abort();
  return EvtVector4C();

}


EvtSpinDensity EvtPhotonParticle::rotateToHelicityBasis() const {

  EvtVector4C eplus(0.0,-1.0/sqrt(2.0),
		    EvtComplex(0.0,-1.0/sqrt(2.0)),0.0);
  EvtVector4C eminus(0.0,1.0/sqrt(2.0),
		     EvtComplex(0.0,-1.0/sqrt(2.0)),0.0);

  //Really uggly have to cast away constness because the
  //function epsParentPhoton caches the state vectors...
  EvtVector4C e1=((EvtParticle*)this)->epsParentPhoton(0);
  EvtVector4C e2=((EvtParticle*)this)->epsParentPhoton(1);


  EvtSpinDensity R;
  R.setDim(2);

  R.set(0,0,(eplus.conj())*e1);
  R.set(0,1,(eplus.conj())*e2);
  
  R.set(1,0,(eminus.conj())*e1);
  R.set(1,1,(eminus.conj())*e2);
  
  return R;
	
}


EvtSpinDensity EvtPhotonParticle::rotateToHelicityBasis(double alpha,
							double beta,
							double gamma) const{

  EvtVector4C eplus(0.0,-1.0/sqrt(2.0),
		    EvtComplex(0.0,-1.0/sqrt(2.0)),0.0);
  EvtVector4C eminus(0.0,1.0/sqrt(2.0),
		     EvtComplex(0.0,-1.0/sqrt(2.0)),0.0);
  
  eplus.applyRotateEuler(alpha,beta,gamma);
  eminus.applyRotateEuler(alpha,beta,gamma);
  
  
  //Really uggly have to cast away constness because the
  //function epsParentPhoton caches the state vectors...
  EvtVector4C e1=((EvtParticle*)this)->epsParentPhoton(0);
  EvtVector4C e2=((EvtParticle*)this)->epsParentPhoton(1);
  
  EvtSpinDensity R;
  R.setDim(2);
  
  R.set(0,0,(eplus.conj())*e1);
  R.set(0,1,(eplus.conj())*e2);
  
  R.set(1,0,(eminus.conj())*e1);
  R.set(1,1,(eminus.conj())*e2);
  
  return R;

}
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: EvtPhotonParticle.cc
	12	//
	13	// Description: Class to describe massless vectors
	14	//
	15	// Modification history:
	16	//
	17	// DJL/RYD September 25, 1996 Module created
	18	//
	19	//------------------------------------------------------------------------
	20	//
	21	#include "EvtGenBase/EvtPatches.hh"
	22	#include <stdlib.h>
	23	#include <iostream>
	24	#include <math.h>
	25	#include "EvtGenBase/EvtComplex.hh"
	26	#include "EvtGenBase/EvtPhotonParticle.hh"
	27	#include "EvtGenBase/EvtVector4C.hh"
	28	#include "EvtGenBase/EvtReport.hh"
	29	using std::endl;
	30
0ca57c2f	31	EvtPhotonParticle::EvtPhotonParticle() {
	32	}
	33
da0e9ce3	34	EvtPhotonParticle::~EvtPhotonParticle(){}
	35
	36	void EvtPhotonParticle::init(EvtId part_n,const EvtVector4R& p4){
	37
	38	init(part_n,p4.get(0),p4.get(1),
	39	p4.get(2),p4.get(3));
	40	}
	41
	42	void EvtPhotonParticle::init(EvtId part_n,double e,double px,double py,double pz){
	43
	44	_validP4=true;
	45	setp(e,px,py,pz);
	46	setpart_num(part_n);
	47
	48	setLifetime();
	49
	50	//defere calculation of basis vectors untill they are needed!
	51	_evalBasis=0;
	52
	53	}
	54
	55
	56	EvtVector4C EvtPhotonParticle::epsParentPhoton(int i){
	57
	58	if (!_evalBasis){
	59
	60	_evalBasis=1;
	61	eps1.set(EvtComplex(0.0,0.0), EvtComplex(-1.0/sqrt(2.0),0.0),
	62	EvtComplex(0.0,-1.0/sqrt(2.0)), EvtComplex(0.0,0.0));
	63	eps2.set(EvtComplex(0.0,0.0), EvtComplex(1.0/sqrt(2.0),0.0),
	64	EvtComplex(0.0,-1.0/sqrt(2.0)), EvtComplex(0.0,0.0));
	65
	66	// These are for photon along z axis. Rotate to get
	67	// correct direction...
	68
	69	double phi,theta;
	70
	71	EvtVector4R p=this->getP4();
	72
	73	double px=p.get(1);
	74	double py=p.get(2);
	75	double pz=p.get(3);
	76
	77	phi = atan2(py,px);
	78	theta = acos(pz/sqrt(pxpx+pypy+pz*pz));
	79	eps1.applyRotateEuler(phi,theta,-phi);
	80	eps2.applyRotateEuler(phi,theta,-phi);
	81
	82	}
	83
	84
	85	EvtVector4C temp;
	86
	87	switch(i) {
	88
	89	case 0:
	90	temp=eps1;
	91	break;
	92	case 1:
	93	temp=eps2;
	94	break;
	95	default:
	96	report(ERROR,"EvtGen") << "EvtPhotonParticle.cc: Asked "
	97	<< "for state:"<<i<<endl;
98	::abort();
99	break;
100	}
101
102	return temp;
103	}
104
105	EvtVector4C EvtPhotonParticle::epsPhoton(int ){
106
107	report(ERROR,"EvtGen") << "EvtPhotonParticle.cc: Can not get "
108	<< "state in photons restframe."<<endl;;
109	::abort();
110	return EvtVector4C();
111
112	}
113
114
115	EvtSpinDensity EvtPhotonParticle::rotateToHelicityBasis() const {
116
117	EvtVector4C eplus(0.0,-1.0/sqrt(2.0),
118	EvtComplex(0.0,-1.0/sqrt(2.0)),0.0);
119	EvtVector4C eminus(0.0,1.0/sqrt(2.0),
120	EvtComplex(0.0,-1.0/sqrt(2.0)),0.0);
121
122	//Really uggly have to cast away constness because the
123	//function epsParentPhoton caches the state vectors...
124	EvtVector4C e1=((EvtParticle*)this)->epsParentPhoton(0);
125	EvtVector4C e2=((EvtParticle*)this)->epsParentPhoton(1);
126
127
128	EvtSpinDensity R;
129	R.setDim(2);
130
131	R.set(0,0,(eplus.conj())*e1);
132	R.set(0,1,(eplus.conj())*e2);
133
134	R.set(1,0,(eminus.conj())*e1);
135	R.set(1,1,(eminus.conj())*e2);
136
137	return R;
138
139	}
140
141
142	EvtSpinDensity EvtPhotonParticle::rotateToHelicityBasis(double alpha,
143	double beta,
144	double gamma) const{
145
146	EvtVector4C eplus(0.0,-1.0/sqrt(2.0),
147	EvtComplex(0.0,-1.0/sqrt(2.0)),0.0);
148	EvtVector4C eminus(0.0,1.0/sqrt(2.0),
149	EvtComplex(0.0,-1.0/sqrt(2.0)),0.0);
150
151	eplus.applyRotateEuler(alpha,beta,gamma);
152	eminus.applyRotateEuler(alpha,beta,gamma);
153
154
155	//Really uggly have to cast away constness because the
156	//function epsParentPhoton caches the state vectors...
157	EvtVector4C e1=((EvtParticle*)this)->epsParentPhoton(0);
158	EvtVector4C e2=((EvtParticle*)this)->epsParentPhoton(1);
159
160	EvtSpinDensity R;
161	R.setDim(2);
162
163	R.set(0,0,(eplus.conj())*e1);
164	R.set(0,1,(eplus.conj())*e2);
165
166	R.set(1,0,(eminus.conj())*e1);
167	R.set(1,1,(eminus.conj())*e2);
168
169	return R;
170
171	}
172
173