TEvtGen/EvtGenBase/EvtdFunctionSingle.cxx

   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) 2000      Caltech, UCSB
  10 //
  11 // Module: EvtdFunctionSingle.cc
  12 //
  13 // Description: Evaluates one Wigner d-Functions.
  14 //
  15 // Modification history:
  16 //
  17 //    fkw           February 2, 2001     changes to satisfy KCC
  18 //    RYD            August 10, 2000         Module created
  19 //
  20 //------------------------------------------------------------------------
  21 //
  22 #include "EvtGenBase/EvtPatches.hh"
  23 #include <stdlib.h>
  24 #include <math.h>
  25 #include <iostream>
  26 #include <assert.h>
  27 #include "EvtGenBase/EvtdFunctionSingle.hh"
  28
  29 EvtdFunctionSingle::EvtdFunctionSingle(){
  30   _j=0;
  31   _m1=0;
  32   _m2=0;
  33   _coef=0;
  34   _kmin=0;
  35   _kmax=0;
  36 }
  37
  38
  39 EvtdFunctionSingle::~EvtdFunctionSingle(){
  40   if (_coef!=0) delete [] _coef;
  41 }
  42
  43
  44 void EvtdFunctionSingle::init(int j,int m1,int m2){
  45
  46   assert(abs(m2)>=abs(m1));
  47   assert(m2>=0);
  48
  49   _j=j;
  50   _m1=m1;
  51   _m2=m2;
  52
  53   _kmin=_m2-_m1;
  54   _kmax=_j-_m1;
  55
  56   assert(_kmin<=_kmax);
  57
  58   _coef=new double[(_kmax-_kmin)/2+1];
  59
  60   int k;
  61
  62   for(k=_kmin;k<=_kmax;k+=2){
  63     int sign=1;
  64     if ((k-_m2+_m1)%4!=0) sign=-sign;
  65     double tmp = fact((_j+_m2)/2)*fact((_j-_m2)/2)
  66                    *fact((_j+_m1)/2)*fact((_j-_m1)/2);
  67     _coef[(k-_kmin)/2]=sign*sqrt(tmp)/
  68       (fact((_j+_m2-k)/2)*fact(k/2)*fact((_j-_m1-k)/2)*fact((k-_m2+_m1)/2));
  69
  70   }
  71
  72 }
  73
  74
  75 double EvtdFunctionSingle::d(int j,int m1,int m2, double theta){
  76
  77   assert(j==_j);
  78   assert(m1==_m1);
  79   assert(m2==_m2);
  80
  81   double c2=cos(0.5*theta);
  82   double s2=sin(0.5*theta);
  83
  84   double d=0.0;
  85
  86   int k;
  87   for(k=_kmin;k<=_kmax;k+=2){
  88     d+=_coef[(k-_kmin)/2]*pow(c2,(2*_j-2*k+m2-m1)/2)*pow(s2,(2*k-m2+m1)/2);
  89   }
  90
  91   return d;
  92
  93 }
  94
  95
  96 int EvtdFunctionSingle::fact(int n){
  97
  98   assert(n>=0);
  99
 100   int f=1;
 101
 102   int k;
 103   for(k=2;k<=n;k++) f*=k;
 104
 105   return f;
 106
 107 }
 108
 109
 110