/*****************************************************************************
* Project: RooFit                                                           *
* Package: RooFitModels                                                     *
*    File: $Id$
* Authors:                                                                  *
*   GR, Gerhard Raven,   Nikhef & VU, Gerhard.Raven@nikhef.nl
*                                                                           *
*                                                                           *
* Redistribution and use in source and binary forms,                        *
* with or without modification, are permitted according to the terms        *
*****************************************************************************/

//////////////////////////////////////////////////////////////////////////////
//
// BEGIN_HTML
//   Implementation of the so-called real spherical harmonics, using the orthonormal normalization,
// which are related to spherical harmonics as:
//
//  Y_{l0} = Y_l^0   (m=0)
//  Y_{lm} = \frac{1}{\sqrt{2}}  \left( Y_l^m     + (-1)^m     Y_l^{-m}   \right) (m>0)
//  Y_{lm} = \frac{1}{i\sqrt{2}} \left( Y_l^{|m|} - (-1)^{|m|} Y_l^{-|m|} \right) (m<0)
//
// which implies:
//
// Y_{l0}(\cos\theta,\phi) =          N_{l0}   P_l^0    (\cos\theta)              (m=0)
// Y_{lm}(\cos\theta,\phi) = \sqrt{2} N_{lm}   P_l^m    (\cos\theta) cos(|m|\phi) (m>0)
// Y_{lm}(\cos\theta,\phi) = \sqrt{2} N_{l|m|} P_l^{|m|}(\cos\theta) sin(|m|\phi) (m<0)
//
// where
//  N_{lm} = \sqrt{ \frac{2l+1}{4\pi} \frac{ (l-m)! }{ (l+m)! } }
//
// Note that the normalization corresponds to the orthonormal case,
// and thus we have Int d\cos\theta d\phi Y_{lm} Y_{l'm'} = \delta_{ll'} \delta{mm'}
//
// Note that in addition, this code can also represent the product of two
// (real) spherical harmonics -- it actually uses the fact that Y_{00} = \sqrt{\frac{1}{4\pi}}
// in order to represent a single spherical harmonics by multiplying it
// by \sqrt{4\pi} Y_00, as this makes it trivial to compute the analytical
// integrals, using the orthogonality properties of Y_l^m...
//
// END_HTML
//

#include "RooFit.h"
#include "Riostream.h"
#include <math.h>

#include "RooSpHarmonic.h"
#include "Math/SpecFunc.h"
#include "TMath.h"

using namespace std;

ClassImp(RooSpHarmonic)
;

//_____________________________________________________________________________
namespace {
inline double N(int l, int m=0) {
double n = sqrt( double(2*l+1)/(4*TMath::Pi())*TMath::Factorial(l-m)/TMath::Factorial(l+m) );
return m==0 ? n : TMath::Sqrt2() * n;
}
}

//_____________________________________________________________________________
RooSpHarmonic::RooSpHarmonic() :
_n(0),
_sgn1(0),
_sgn2(0)
{
}

//_____________________________________________________________________________
RooSpHarmonic::RooSpHarmonic(const char* name, const char* title, RooAbsReal& ctheta, RooAbsReal& phi, int l, int m)
: RooLegendre(name, title,ctheta,l,m<0?-m:m)
, _phi("phi", "phi", this, phi)
, _n( 2*sqrt(TMath::Pi()))
, _sgn1( m==0 ? 0 : m<0 ? -1 : +1 )
, _sgn2( 0 )
{
}

//_____________________________________________________________________________
RooSpHarmonic::RooSpHarmonic(const char* name, const char* title, RooAbsReal& ctheta, RooAbsReal& phi, int l1, int m1, int l2, int m2)
: RooLegendre(name, title,ctheta,l1, m1<0?-m1:m1,l2,m2<0?-m2:m2)
, _phi("phi", "phi", this, phi)
, _n(1)
, _sgn1( m1==0 ? 0 : m1<0 ? -1 : +1 )
, _sgn2( m2==0 ? 0 : m2<0 ? -1 : +1 )
{
}

//_____________________________________________________________________________
RooSpHarmonic::RooSpHarmonic(const RooSpHarmonic& other, const char* name)
: RooLegendre(other, name)
, _phi("phi", this,other._phi)
, _n(other._n)
, _sgn1( other._sgn1 )
, _sgn2( other._sgn2 )
{
}

//_____________________________________________________________________________
Double_t RooSpHarmonic::evaluate() const
{
double n = _n*N(_l1,_m1)*N(_l2,_m2)*RooLegendre::evaluate();
if (_sgn1!=0) n *= (_sgn1<0 ? sin(_m1*_phi) : cos(_m1*_phi) );
if (_sgn2!=0) n *= (_sgn2<0 ? sin(_m2*_phi) : cos(_m2*_phi) );
return n;
}

namespace {
bool fullRange(const RooRealProxy& x ,const char* range)  {
return ( x.min(range) == x.min() && x.max(range) == x.max() ) ;
}
}

//_____________________________________________________________________________
Int_t RooSpHarmonic::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
{
// TODO: check that phi.max - phi.min = 2 pi... ctheta.max = +1, and ctheta.min = -1
// we don't support indefinite integrals... maybe one day, when there is a use for it.....
bool noRange  = ( rangeName == 0 || strlen(rangeName)==0 );
bool phiOK    = ( noRange || fullRange(_phi,rangeName) );
bool cthetaOK = ( noRange || fullRange(_ctheta,rangeName) );
if (cthetaOK && phiOK && matchArgs(allVars, analVars, _ctheta,_phi)) return 3;
if (            phiOK && matchArgs(allVars, analVars,         _phi)) return 2;
return RooLegendre::getAnalyticalIntegral(allVars,analVars,rangeName);
}

//_____________________________________________________________________________
Double_t RooSpHarmonic::analyticalIntegral(Int_t code, const char* range) const
{
if (code==3) {
return (_l1==_l2 && _sgn1*_m1==_sgn2*_m2 ) ? _n : 0 ;
} else if (code == 2) {
if ( _sgn1*_m1 != _sgn2*_m2) return 0;
return ( _m1==0 ? 2 : 1 ) * TMath::Pi()*_n*N(_l1,_m1)*N(_l2,_m2)*RooLegendre::evaluate();
} else {
double n = _n*N(_l1,_m1)*N(_l2,_m2)*RooLegendre::analyticalIntegral(code,range);
if (_sgn1!=0) n *= (_sgn1<0 ? sin(_m1*_phi) : cos(_m1*_phi) );
if (_sgn2!=0) n *= (_sgn2<0 ? sin(_m2*_phi) : cos(_m2*_phi) );
return n;
}
}

Int_t RooSpHarmonic::getMaxVal( const RooArgSet& vars) const {
return RooLegendre::getMaxVal(vars);
}

Double_t RooSpHarmonic::maxVal( Int_t code) const {
double n = _n*N(_l1,_m1)*N(_l2,_m2);
return n*RooLegendre::maxVal(code);
}

RooSpHarmonic.cxx:1
RooSpHarmonic.cxx:2
RooSpHarmonic.cxx:3
RooSpHarmonic.cxx:4
RooSpHarmonic.cxx:5
RooSpHarmonic.cxx:6
RooSpHarmonic.cxx:7
RooSpHarmonic.cxx:8
RooSpHarmonic.cxx:9
RooSpHarmonic.cxx:10
RooSpHarmonic.cxx:11
RooSpHarmonic.cxx:12
RooSpHarmonic.cxx:13
RooSpHarmonic.cxx:14
RooSpHarmonic.cxx:15
RooSpHarmonic.cxx:16
RooSpHarmonic.cxx:17
RooSpHarmonic.cxx:18
RooSpHarmonic.cxx:19
RooSpHarmonic.cxx:20
RooSpHarmonic.cxx:21
RooSpHarmonic.cxx:22
RooSpHarmonic.cxx:23
RooSpHarmonic.cxx:24
RooSpHarmonic.cxx:25
RooSpHarmonic.cxx:26
RooSpHarmonic.cxx:27
RooSpHarmonic.cxx:28
RooSpHarmonic.cxx:29
RooSpHarmonic.cxx:30
RooSpHarmonic.cxx:31
RooSpHarmonic.cxx:32
RooSpHarmonic.cxx:33
RooSpHarmonic.cxx:34
RooSpHarmonic.cxx:35
RooSpHarmonic.cxx:36
RooSpHarmonic.cxx:37
RooSpHarmonic.cxx:38
RooSpHarmonic.cxx:39
RooSpHarmonic.cxx:40
RooSpHarmonic.cxx:41
RooSpHarmonic.cxx:42
RooSpHarmonic.cxx:43
RooSpHarmonic.cxx:44
RooSpHarmonic.cxx:45
RooSpHarmonic.cxx:46
RooSpHarmonic.cxx:47
RooSpHarmonic.cxx:48
RooSpHarmonic.cxx:49
RooSpHarmonic.cxx:50
RooSpHarmonic.cxx:51
RooSpHarmonic.cxx:52
RooSpHarmonic.cxx:53
RooSpHarmonic.cxx:54
RooSpHarmonic.cxx:55
RooSpHarmonic.cxx:56
RooSpHarmonic.cxx:57
RooSpHarmonic.cxx:58
RooSpHarmonic.cxx:59
RooSpHarmonic.cxx:60
RooSpHarmonic.cxx:61
RooSpHarmonic.cxx:62
RooSpHarmonic.cxx:63
RooSpHarmonic.cxx:64
RooSpHarmonic.cxx:65
RooSpHarmonic.cxx:66
RooSpHarmonic.cxx:67
RooSpHarmonic.cxx:68
RooSpHarmonic.cxx:69
RooSpHarmonic.cxx:70
RooSpHarmonic.cxx:71
RooSpHarmonic.cxx:72
RooSpHarmonic.cxx:73
RooSpHarmonic.cxx:74
RooSpHarmonic.cxx:75
RooSpHarmonic.cxx:76
RooSpHarmonic.cxx:77
RooSpHarmonic.cxx:78
RooSpHarmonic.cxx:79
RooSpHarmonic.cxx:80
RooSpHarmonic.cxx:81
RooSpHarmonic.cxx:82
RooSpHarmonic.cxx:83
RooSpHarmonic.cxx:84
RooSpHarmonic.cxx:85
RooSpHarmonic.cxx:86
RooSpHarmonic.cxx:87
RooSpHarmonic.cxx:88
RooSpHarmonic.cxx:89
RooSpHarmonic.cxx:90
RooSpHarmonic.cxx:91
RooSpHarmonic.cxx:92
RooSpHarmonic.cxx:93
RooSpHarmonic.cxx:94
RooSpHarmonic.cxx:95
RooSpHarmonic.cxx:96
RooSpHarmonic.cxx:97
RooSpHarmonic.cxx:98
RooSpHarmonic.cxx:99
RooSpHarmonic.cxx:100
RooSpHarmonic.cxx:101
RooSpHarmonic.cxx:102
RooSpHarmonic.cxx:103
RooSpHarmonic.cxx:104
RooSpHarmonic.cxx:105
RooSpHarmonic.cxx:106
RooSpHarmonic.cxx:107
RooSpHarmonic.cxx:108
RooSpHarmonic.cxx:109
RooSpHarmonic.cxx:110
RooSpHarmonic.cxx:111
RooSpHarmonic.cxx:112
RooSpHarmonic.cxx:113
RooSpHarmonic.cxx:114
RooSpHarmonic.cxx:115
RooSpHarmonic.cxx:116
RooSpHarmonic.cxx:117
RooSpHarmonic.cxx:118
RooSpHarmonic.cxx:119
RooSpHarmonic.cxx:120
RooSpHarmonic.cxx:121
RooSpHarmonic.cxx:122
RooSpHarmonic.cxx:123
RooSpHarmonic.cxx:124
RooSpHarmonic.cxx:125
RooSpHarmonic.cxx:126
RooSpHarmonic.cxx:127
RooSpHarmonic.cxx:128
RooSpHarmonic.cxx:129
RooSpHarmonic.cxx:130
RooSpHarmonic.cxx:131
RooSpHarmonic.cxx:132
RooSpHarmonic.cxx:133
RooSpHarmonic.cxx:134
RooSpHarmonic.cxx:135
RooSpHarmonic.cxx:136
RooSpHarmonic.cxx:137
RooSpHarmonic.cxx:138
RooSpHarmonic.cxx:139
RooSpHarmonic.cxx:140
RooSpHarmonic.cxx:141
RooSpHarmonic.cxx:142
RooSpHarmonic.cxx:143
RooSpHarmonic.cxx:144
RooSpHarmonic.cxx:145
RooSpHarmonic.cxx:146
RooSpHarmonic.cxx:147
RooSpHarmonic.cxx:148
RooSpHarmonic.cxx:149
RooSpHarmonic.cxx:150
RooSpHarmonic.cxx:151
RooSpHarmonic.cxx:152
RooSpHarmonic.cxx:153
RooSpHarmonic.cxx:154
RooSpHarmonic.cxx:155
RooSpHarmonic.cxx:156
RooSpHarmonic.cxx:157