00001
00018 #ifdef _MSC_VER
00019
00020 #include "msdevstudio/MSconfig.h"
00021 #endif
00022
00023 #include "Bins2DProfile.h"
00024
00025 #include "datasrcs/DataPointTuple.h"
00026 #include "datasrcs/NTuple.h"
00027
00028 #include <cmath>
00029
00030 #include <cassert>
00031
00032 #ifdef ITERATOR_MEMBER_DEFECT
00033 using namespace std;
00034 #else
00035 using std::fill;
00036 using std::list;
00037 using std::sqrt;
00038 using std::vector;
00039 #endif
00040
00041 using namespace hippodraw;
00042
00043 Bins2DProfile::Bins2DProfile ( )
00044 : Bins2DBase ( "Bins2DProfile" )
00045 {
00046 }
00047
00048 Bins2DProfile::Bins2DProfile( const Bins2DProfile & binner )
00049 : Bins2DBase( binner ),
00050 m_sumZ( binner.m_sumZ ),
00051 m_variance( binner.m_variance ),
00052 m_z_range ( binner.m_z_range )
00053 {
00054 m_num_bins = binner.m_num_bins;
00055 }
00056
00057 Bins2DProfile::~Bins2DProfile ()
00058 {
00059 }
00060
00061 BinsBase *
00062 Bins2DProfile::clone () const
00063 {
00064 return new Bins2DProfile ( *this );
00065 }
00066
00067 void
00068 Bins2DProfile::
00069 setNumberOfBins ( hippodraw::Axes::Type axis, int nb )
00070 {
00071 assert ( axis == Axes::X || axis == Axes::Y );
00072
00073 Bins2DBase::setNumberOfBins ( axis, nb );
00074 if ( axis == Axes::X ) {
00075 int number = numberOfBins ( Axes::X );
00076 m_data.resize( number + 2 );
00077 m_sumZ.resize( number + 2 );
00078 m_variance.resize( number );
00079 }
00080 else {
00081 int number_y = numberOfBins ( Axes::Y );
00082 unsigned int i = 0;
00083 for ( ; i < m_data.size(); i++ ) {
00084 m_data[i].resize( number_y + 2 );
00085 m_sumZ[i].resize( number_y + 2 );
00086 }
00087 for ( i = 0; i < m_variance.size(); i++ ) {
00088 m_variance[i].resize( number_y );
00089 }
00090 }
00091
00092 reset ();
00093
00094 m_num_bins = numberOfBins ( Axes::X ) * numberOfBins ( Axes::Y );
00095 }
00096
00097 void Bins2DProfile::reset()
00098 {
00099 unsigned int i = 0;
00100 for ( ; i < m_data.size(); i++ ) {
00101 fill ( m_data[i].begin(), m_data[i].end(), 0.0 );
00102 fill ( m_sumZ[i].begin(), m_sumZ[i].end(), 0.0 );
00103 }
00104
00105 for ( i = 0; i < m_variance.size(); i++ ) {
00106 fill ( m_variance[i].begin(), m_variance[i].end(), 0.0 );
00107 }
00108
00109 m_empty = true;
00110 }
00111
00112 void Bins2DProfile::accumulate( double x, double y, double z, double wt )
00113 {
00114 int i = binNumberX( x );
00115 int j = binNumberY( y );
00116
00117 if( i > 0 && i <= numberOfBins ( Axes::X ) &&
00118 j > 0 && j <= numberOfBins ( Axes::Y ) ) {
00119 m_variance[i-1][j-1] += wt * wt;
00120 }
00121 m_data[i][j] += wt;
00122 m_sumZ[i][j] += z;
00123
00124 m_empty = false;
00125 }
00126
00127 double Bins2DProfile::getZValue ( double x, double y ) const
00128 {
00129
00130 int i = binNumberX( x );
00131 int j = binNumberY( y );
00132
00133
00134 if (m_data[i][j] == 0) return 0.0;
00135
00136 double v = m_sumZ[i][j] / m_data[i][j];
00137 return v;
00138 }
00139
00140 NTuple *
00141 Bins2DProfile::
00142 createNTuple () const
00143 {
00144 unsigned int size = numberOfBins ( Axes::X ) * numberOfBins ( Axes::Y );
00145 NTuple * ntuple = prepareNTuple ( size );
00146
00147 fillDataSource ( ntuple );
00148
00149 return ntuple;
00150 }
00151
00152 namespace dp = hippodraw::DataPoint3DTuple;
00153
00157 void
00158 Bins2DProfile::
00159 fillDataSource ( DataSource * ntuple ) const
00160 {
00161 ntuple -> clear ();
00162 vector < double > row ( dp::SIZE );
00163
00164 int v_inc = 0;
00165 vector<double>::const_iterator vity;
00166 double next_x = getLow ( Axes::X );
00167
00168 int num_x = numberOfBins ( Axes::X );
00169 int num_y = numberOfBins ( Axes::Y );
00170
00171 for ( int x_inc = 1; x_inc <= num_x; x_inc++ ) {
00172 vity = m_variance[v_inc++].begin();
00173 double widthX = binWidthX ( x_inc - 1 );
00174 double half_widthX = 0.5 * widthX;
00175 double x = next_x + half_widthX;
00176 next_x += widthX;
00177 double next_y = getLow ( Axes::Y );
00178
00179 for ( int y_inc = 1; y_inc <= num_y; y_inc++ ) {
00180 double widthY = binWidthY ( y_inc - 1 );
00181 double half_widthY = 0.5 * widthY;
00182 double y = next_y;
00183 y += half_widthY;
00184 next_y += widthY;
00185
00186 double var2 = *vity++;
00187 double v = 0.0;
00188 double verr = 0.0;
00189 if ( m_data[x_inc][y_inc] != 0 ) {
00190 double num = m_data[x_inc][y_inc];
00191 v = m_sumZ[x_inc][y_inc] / num;
00192 verr = sqrt ( ( var2 / ( num - 1. ) - v * v ) );
00193 }
00194 row[dp::X] = x;
00195 row[dp::Y] = y;
00196 row[dp::Z] = v;
00197 row [dp::XERR] = half_widthX;
00198 row [dp::YERR] = half_widthY;
00199 row [dp::ZERR] = verr;
00200
00201 ntuple -> addRow ( row );
00202 }
00203 }
00204 vector < unsigned int > shape ( 3 );
00205 shape[0] = num_x;
00206 shape[1] = num_y;
00207 shape[2] = dp::SIZE;
00208
00209 ntuple -> setShape ( shape );
00210 }
00211
00215 void
00216 Bins2DProfile::
00217 setBinContents ( const DataSource * )
00218 {
00219 }