00001
00012 #ifdef _MSC_VER
00013 #include "msdevstudio/MSconfig.h"
00014 #endif
00015
00016 #include "binners/Bins1DProfile.h"
00017
00018 #include "datasrcs/DataPointTuple.h"
00019 #include "datasrcs/NTuple.h"
00020
00021 #include <algorithm>
00022 #include <numeric>
00023
00024 #include <cmath>
00025
00026 #include <cassert>
00027
00028
00029 #ifdef ITERATOR_MEMBER_DEFECT
00030 using namespace std;
00031 #else
00032 using std::fill;
00033 using std::list;
00034 using std::sqrt;
00035 using std::string;
00036 using std::vector;
00037 using std::min_element;
00038 using std::max_element;
00039 #endif
00040
00041 using namespace hippodraw;
00042
00043 Bins1DProfile::Bins1DProfile ( )
00044 : Bins1DBase ( "Bins1DProfile" )
00045 {
00046 }
00047
00048 Bins1DProfile::Bins1DProfile ( const Bins1DProfile & binner )
00049 : Bins1DBase ( binner ),
00050 m_sum( binner.m_sum ),
00051 m_sumsq( binner.m_sumsq ),
00052 m_num( binner.m_num ),
00053 m_sumwtsq (binner.m_sumwtsq )
00054 {
00055 }
00056
00057 BinsBase *
00058 Bins1DProfile::clone () const
00059 {
00060 return new Bins1DProfile(*this);
00061 }
00062
00063 Bins1DProfile::~Bins1DProfile ()
00064 {
00065 }
00066
00067 double Bins1DProfile::minBin()
00068 {
00069 return *min_element ( m_sum.begin() + 1, m_sum.end ( ) - 1 );
00070 }
00071
00072 double Bins1DProfile::maxBin()
00073 {
00074 return *max_element ( m_sum.begin () + 1, m_sum.end () - 1 );
00075 }
00076
00077 void Bins1DProfile::resize ( int number )
00078 {
00079 m_sum.resize( number + 2 );
00080 m_num.resize( number + 2 );
00081 m_sumsq.resize( number + 2 );
00082 m_sumwtsq.resize( number + 2);
00083
00084 reset();
00085
00086 m_values_dirty = true;
00087 }
00088
00089 void Bins1DProfile::reset()
00090 {
00091 fill( m_sum.begin(), m_sum.end(), 0.0 );
00092 fill( m_num.begin(), m_num.end(), 0 );
00093 fill( m_sumsq.begin(), m_sumsq.end(), 0.0 );
00094 fill( m_sumwtsq.begin(), m_sumwtsq.end(), 0.0 );
00095
00096 m_empty = true;
00097 }
00098
00099 int Bins1DProfile::getNumberOfEntries () const
00100 {
00101 return static_cast < int > ( std::accumulate ( m_num.begin() + 1,
00102 m_num.end() - 2, 0.0 ) );
00103 }
00104
00105 int Bins1DProfile::getNumberOfEntries ( int i ) const
00106 {
00107
00108 return static_cast < int > (m_num [ i+1] );
00109 }
00110
00111 int Bins1DProfile::getUnderflow () const
00112 {
00113 return -1;
00114 }
00115
00116 int Bins1DProfile::getOverflow () const
00117 {
00118 return -1;
00119 }
00120
00121
00122 void Bins1DProfile::accumulate( double x, double y, double wt, double )
00123 {
00124 int i = binNumber ( x );
00125
00126 m_sum[i] += y * wt;
00127 m_sumsq[i] += y * y * wt;
00128 m_num[i] += wt;
00129
00130 m_sumwtsq[i] += wt * wt;
00131
00132 m_empty = false;
00133 }
00134
00135
00136 NTuple *
00137 Bins1DProfile::
00138 createNTuple () const
00139 {
00140 unsigned int size = m_sum.size ();
00141 NTuple * ntuple = prepareNTuple ( size );
00142
00143 fillDataSource ( ntuple );
00144
00145 return ntuple;
00146 }
00147
00148 namespace dp = hippodraw::DataPoint2DTuple;
00149
00150 void
00151 Bins1DProfile::
00152 fillDataSource ( DataSource * ntuple ) const
00153 {
00154 ntuple -> clear();
00155 vector < double > row ( dp::SIZE );
00156
00157 double width;
00158 double x = getLow ( Axes::X );
00159 int x_number = numberOfBins ( Axes::X );
00160 for ( int i = 0; i < x_number; i++ ) {
00161 width = binWidth ( i );
00162 double half_width = 0.5 * width;
00163 x += half_width;
00164 double y = 0;
00165 double num = m_num [ i+1 ];
00166 if ( num <= 0 ) {
00167 x += half_width;
00168 continue;
00169 }
00170 if ( m_sum[i+1] != 0 ) {
00171 y = m_sum[i+1] / m_num[i+1];
00172 }
00173 else y = 0.;
00174
00175 double yerr = 0.0;
00176
00177 if ( m_num[i+1] > 1. ) {
00178 double num = m_num[i+1];
00179
00180
00181
00182
00183
00184
00185
00186
00187 yerr = sqrt ( (num/((num*num)-m_sumwtsq[i+1]))*
00188 (m_sumsq[i+1]-2.0*y*m_sum[i+1]+y*y*num) );
00189
00190
00191
00192
00193 }
00194 else {
00195 yerr = 0.0;
00196 }
00197
00198 row[dp::X] = x;
00199 row[dp::Y] = y;
00200 row[dp::XERR] = half_width;
00201 row[dp::YERR] = yerr;
00202
00203 ntuple -> addRow ( row );
00204
00205 x += half_width;
00206 }
00207 }
00208
00212 void
00213 Bins1DProfile::
00214 setBinContents ( const DataSource * )
00215 {
00216 }