(***********************************)
(***********************************)
(* The Partial Waves Decomposition *)
(* author: M. Gintner              *)
   version2 = "APR 30, 2002";
(***********************************)
(***********************************)
(* Formulas found in the book                            *)
(* Horejsi: Elektroslabe sjednoceni a stromova unitarita,*)
(* pp. 143-145, were used.                               *)
(* The error found in the eq.(F.2) was corrected.        *)
(***********************************)

(*************************************)
(* Complex-conjugated expression:    *)

Star[z_] := ReplaceAll[z,Complex[z1_,z2_]->Complex[z1,-z2]];

(*************************************)
(* The Wigner D-function definition: *)
(* j   ... a total angular momentum  *)
(* mf = h3 - h4, mi = h1 - h2        *)
(* h1,h2 ... initial state helicities*)
(* h3,h4 ... final state helicities  *)
(* theta, phi ... scattering angles  *)

Littled[j_,mf_,mi_,ksi_] :=
(-1)^(j-mf)*Sqrt[(j+mf)!/(2^(2j)*(j-mf)!*(j+mi)!*(j-mi)!)]*
(1+ksi)^(-(mf+mi)/2)*(1-ksi)^(-(mf-mi)/2)*
D[(1+ksi)^(j+mi)*(1-ksi)^(j-mi),{ksi,j-mf}];

WignerD[j_,mf_,mi_,theta_,phi_] :=
Exp[I*mi*phi]*Littled[j,mf,mi,Cos[theta]];

(*****************************************)
(* The Partial Wave Amplitude extracted: *)
(* j   ... a total angular momentum      *)
(* mf = h3 - h4, mi = h1 - h2            *)
(* h1,h2 ... initial state helicities    *)
(* h3,h4 ... final state helicities      *)
(* theta, phi ... scattering angles      *)
(* ampl ... scattering amplitude         *)

PWA[j_,h1_,h2_,h3_,h4_,ampl_] := Module[{WignerDstar,integrand,mi,mf},
mi = h1-h2; mf = h3-h4;
If[ Or[Abs[mi]>j,Abs[mf]>j],
    StringForm["Error: 'j' must be greater or equal to both,
'Abs[h1-h2]' and 'Abs[h3-h4]' !"],
(WignerDstar = Star[ WignerD[j,mf,mi,theta,phi]];
integrand = Simplify[ampl*WignerDstar*Sin[theta]];
1/(64*Pi^2)*Integrate[integrand,{theta,0,Pi},{phi,0,2Pi}])]  ];

(*****************************************)
(* The Partial Wave Amplitude extracted  *)
(* when 'phi = 0' :                      *)
(* j   ... a total angular momentum      *)
(* mf = h3 - h4, mi = h1 - h2            *)
(* h1,h2 ... initial state helicities    *)
(* h3,h4 ... final state helicities      *)
(* theta, phi ... scattering angles      *)
(* ampl ... scattering amplitude         *)

PWA0[j_,h1_,h2_,h3_,h4_,ampl_] := Module[{WignerDstar,integrand,mi,mf},
mi = h1-h2; mf = h3-h4;
If[ Or[Abs[mi]>j,Abs[mf]>j],
    StringForm["Error: 'j' must be greater or equal to both,
'Abs[h1-h2]' and 'Abs[h3-h4]' !"],
(WignerDstar = Star[ WignerD[j,mf,mi,theta,0]];
integrand = Simplify[ampl*WignerDstar*Sin[theta]];
1/(64*Pi^2)*(2Pi)*Integrate[integrand,{theta,0,Pi}])]  ];

(*****************************************)
(* Messages:                             *)

help1 = StringForm["The functions defined in this code:
Star, Littled, WignerD, PWA, PWA0. To get more information on
a specific function type 'help***' where the stars should be 
replaced by the name of the given function."];

helpStar = StringForm["Star[z]; gives complex conjugated of z"];

helpLittled = StringForm["Littled[j,mf,mi,ksi]; 'a little d' of the Wigner's 
D-function;
'j' is a total angular momentum;
'mf = h3 - h4';
'mi = h1 - h2';
'hi' is the helicity of the ith particle (+-1/2 for fermions);
'ksi' is typically the cos(theta)
"];

helpWignerD = StringForm["WignerD[j,mf,mi,theta,phi] :=
Exp[I*mi*phi]*Littled[j,mf,mi,Cos[theta]]"];

helpPWA = StringForm["PWA[j,h1,h2,h3,h4,ampl] is the 'j'th partial wave
amplitude of the 'ampl' of a process h1+h2 --> h3+h4 
with the helicities of particles 
given by 'hi'. The scattering angles in the 'ampl' must be denoted
as 'theta' and 'phi'."];

helpPWA0 = StringForm["PWA0[j,h1,h2,h3,h4,ampl] is
the 'PWA' function when 'phi = 0'."];

StringForm["The 'partial_waves.m' is a code for the partial wave
decomposition of an amplitude of a process 
h1+h2 --> h3+h4. To get a basic information about
the usage of the code type 'help1'. 
(Version: ``)",version1]

