00001 00012 #ifdef _MSC_VER 00013 #include "msdevstudio/MSconfig.h" 00014 #endif 00015 00016 #include "NTupleLikeliHoodFCN.h" 00017 00018 #include "datasrcs/DataPointTuple.h" 00019 #include "datasrcs/DataSource.h" 00020 #include "functions/FunctionBase.h" 00021 00022 #include <cmath> 00023 #include <cassert> 00024 00025 using std::vector; 00026 using std::log; 00027 using std::pow; 00028 using std::exp; 00029 00030 using namespace hippodraw; 00031 00032 NTupleLikeliHoodFCN:: 00033 NTupleLikeliHoodFCN () 00034 { 00035 } 00036 00037 NTupleLikeliHoodFCN:: 00038 NTupleLikeliHoodFCN ( const NTupleLikeliHoodFCN & fcn ) 00039 : NTupleFCN ( fcn ) 00040 { 00041 } 00042 00043 StatedFCN * 00044 NTupleLikeliHoodFCN:: 00045 clone () const 00046 { 00047 return new NTupleLikeliHoodFCN ( *this ); 00048 } 00049 00050 namespace dp = hippodraw::DataPoint2DTuple; 00051 00052 double 00053 NTupleLikeliHoodFCN:: 00054 objectiveValue () const 00055 { 00056 double result = 0.0; 00057 double sum1 = 0.0; 00058 double sum2 = 0.0; 00059 00060 int ix = m_indices [ dp::X ]; 00061 int iy = m_indices [ dp::Y ]; 00062 int ie = m_indices [ dp::XERR ]; 00063 00064 int num_bins = m_ntuple -> rows (); 00065 for( int i = 0; i < num_bins; i++ ) 00066 { 00067 if ( acceptRow ( i ) ) { 00068 const vector < double > & row = m_ntuple -> getRow ( i ); 00069 00070 // Get the edges of the bin. 00071 double x = row [ ix ]; 00072 double hwidth = 0.5 * row [ ie ]; 00073 double edgeL = x - hwidth; 00074 double edgeR = x + hwidth; 00075 00076 // Integrate function over the entire i-th bin 00077 double theta = m_function-> integrate( edgeL, edgeR ); 00078 00079 // Get number of enteries in the ith bin 00080 double n = 2.0 * hwidth * row [ iy ]; 00081 00082 // Computing the factorial of n and taking log of it. 00083 // Exact computation is done when N is small ( < 100) else 00084 // we use sterliing's approximations for it. 00085 // CAUTION: Exactly calculating factorial itself usually leads to 00086 // a very large and unrepresentable number so we avoid it. 00087 //double logfactn = 0.0; 00088 //if( n < 100 ) 00089 //{ 00090 // 00091 // for( int j = 2; j <= n; j++) 00092 // logfactn += log( (double) j ); 00093 // 00094 // result += ( n * log ( theta ) - theta - logfactn ); 00095 //} 00096 //else 00097 //{ 00098 // logfactn = n * log( (double) n ) - n; 00099 // result +=( n * ( log ( theta ) - log( (double) n ) + 1 ) - theta ); 00100 //} 00101 00102 //result += ( n * log ( theta ) - theta ); 00103 00104 sum1 += n * log ( theta ); 00105 sum2 += theta; 00106 00107 } 00108 } 00109 result = sum1 - sum2; 00110 // The quantity that is (asymptotically) equivalent to 00111 // chi-square is -2 log( L ) i.e we return -2 * result 00112 00113 return ( -2 * result ); 00114 } 00115 00116 bool 00117 NTupleLikeliHoodFCN:: 00118 needsIntegrated () const 00119 { 00120 return true; 00121 }