Warning, /include/Geant4/tools/histo/b2 is written in an unsupported language. File is not indexed.
0001 // Copyright (C) 2010, Guy Barrand. All rights reserved.
0002 // See the file tools.license for terms.
0003
0004 #ifndef tools_histo_b2
0005 #define tools_histo_b2
0006
0007 #include "base_histo"
0008
0009 #include <ostream>
0010
0011 namespace tools {
0012 namespace histo {
0013
0014 template <class TC,class TO,class TN,class TW,class TH>
0015 class b2 : public base_histo<TC,TO,TN,TW,TH> {
0016 typedef base_histo<TC,TO,TN,TW,TH> parent;
0017 protected:
0018 enum {AxisX=0,AxisY=1};
0019 public:
0020 typedef base_histo<TC,TO,TN,TW,TH> base_histo_t;
0021 typedef typename parent::axis_t axis_t;
0022 typedef typename parent::bn_t bn_t;
0023 public:
0024 virtual TH bin_error(int,int) const = 0; //for print
0025 public:
0026 // Partition :
0027 TC mean_x() const {
0028 if(parent::m_in_range_Sw==0) return 0;
0029 return parent::m_in_range_Sxw[0]/parent::m_in_range_Sw;
0030 }
0031
0032 TC mean_y() const {
0033 if(parent::m_in_range_Sw==0) return 0;
0034 return parent::m_in_range_Sxw[1]/parent::m_in_range_Sw;
0035 }
0036
0037 TC rms_x() const {
0038 if(parent::m_in_range_Sw==0) return 0;
0039 TC mean = parent::m_in_range_Sxw[0]/parent::m_in_range_Sw;
0040 return ::sqrt(::fabs((parent::m_in_range_Sx2w[0] / parent::m_in_range_Sw) - mean * mean));
0041 }
0042
0043 TC rms_y() const {
0044 if(parent::m_in_range_Sw==0) return 0;
0045 TC mean = parent::m_in_range_Sxw[1]/parent::m_in_range_Sw;
0046 return ::sqrt(::fabs((parent::m_in_range_Sx2w[1] / parent::m_in_range_Sw) - mean * mean));
0047 }
0048
0049 int coord_to_index_x(TC aCoord) const {
0050 return axis_x().coord_to_index(aCoord);
0051 }
0052 int coord_to_index_y(TC aCoord) const {
0053 return axis_y().coord_to_index(aCoord);
0054 }
0055
0056 // bins :
0057 TN bin_entries(int aI,int aJ) const {
0058 TO offset;
0059 if(!_find_offset(aI,aJ,offset)) return 0;
0060 return parent::m_bin_entries[offset];
0061 }
0062
0063 TW bin_Sw(int aI,int aJ) const {
0064 TO offset;
0065 if(!_find_offset(aI,aJ,offset)) return 0;
0066 return parent::m_bin_Sw[offset];
0067 }
0068
0069 TW bin_Sw2(int aI,int aJ) const {
0070 TO offset;
0071 if(!_find_offset(aI,aJ,offset)) return 0;
0072 return parent::m_bin_Sw2[offset];
0073 }
0074 TC bin_Sxw(int aI,int aJ) const {
0075 TO offset;
0076 if(!_find_offset(aI,aJ,offset)) return 0;
0077 return parent::m_bin_Sxw[offset][AxisX];
0078 }
0079 TC bin_Sx2w(int aI,int aJ) const {
0080 TO offset;
0081 if(!_find_offset(aI,aJ,offset)) return 0;
0082 return parent::m_bin_Sx2w[offset][AxisX];
0083 }
0084 TC bin_Syw(int aI,int aJ) const {
0085 TO offset;
0086 if(!_find_offset(aI,aJ,offset)) return 0;
0087 return parent::m_bin_Sxw[offset][AxisY];
0088 }
0089 TC bin_Sy2w(int aI,int aJ) const {
0090 TO offset;
0091 if(!_find_offset(aI,aJ,offset)) return 0;
0092 return parent::m_bin_Sx2w[offset][AxisY];
0093 }
0094
0095 TH bin_height(int aI,int aJ) const {
0096 TO offset;
0097 if(!_find_offset(aI,aJ,offset)) return 0;
0098 return this->get_bin_height(offset);
0099 }
0100
0101 TC bin_center_x(int aI) const {
0102 return parent::m_axes[0].bin_center(aI);
0103 }
0104 TC bin_center_y(int aJ) const {
0105 return parent::m_axes[1].bin_center(aJ);
0106 }
0107
0108 TC bin_mean_x(int aI,int aJ) const {
0109 TO offset;
0110 if(!_find_offset(aI,aJ,offset)) return 0;
0111 TW sw = parent::m_bin_Sw[offset];
0112 if(sw==0) return 0;
0113 return parent::m_bin_Sxw[offset][AxisX]/sw;
0114 }
0115
0116 TC bin_mean_y(int aI,int aJ) const {
0117 TO offset;
0118 if(!_find_offset(aI,aJ,offset)) return 0;
0119 TW sw = parent::m_bin_Sw[offset];
0120 if(sw==0) return 0;
0121 return parent::m_bin_Sxw[offset][AxisY]/sw;
0122 }
0123
0124 TC bin_rms_x(int aI,int aJ) const {
0125 TO offset;
0126 if(!_find_offset(aI,aJ,offset)) return 0;
0127 TW sw = parent::m_bin_Sw[offset];
0128 if(sw==0) return 0;
0129 TC sxw = parent::m_bin_Sxw[offset][AxisX];
0130 TC sx2w = parent::m_bin_Sx2w[offset][AxisX];
0131 TC mean = sxw/sw;
0132 return ::sqrt(::fabs((sx2w / sw) - mean * mean));
0133 }
0134
0135 TC bin_rms_y(int aI,int aJ) const {
0136 TO offset;
0137 if(!_find_offset(aI,aJ,offset)) return 0;
0138 TW sw = parent::m_bin_Sw[offset];
0139 if(sw==0) return 0;
0140 TC sxw = parent::m_bin_Sxw[offset][AxisY];
0141 TC sx2w = parent::m_bin_Sx2w[offset][AxisY];
0142 TC mean = sxw/sw;
0143 return ::sqrt(::fabs((sx2w / sw) - mean * mean));
0144 }
0145
0146 // Axes :
0147 const axis_t& axis_x() const {return parent::m_axes[0];}
0148 const axis_t& axis_y() const {return parent::m_axes[1];}
0149 axis_t& axis_x() {return parent::m_axes[0];} //touchy
0150 axis_t& axis_y() {return parent::m_axes[1];} //touchy
0151 // Projection :
0152
0153 TN bin_entries_x(int aI) const {
0154 if(!parent::m_dimension) return 0;
0155 bn_t ibin;
0156 if(!parent::m_axes[0].in_range_to_absolute_index(aI,ibin)) return 0;
0157 bn_t ybins = parent::m_axes[1].bins()+2;
0158 TO offset;
0159 TN _entries = 0;
0160 for(bn_t jbin=0;jbin<ybins;jbin++) {
0161 offset = ibin + jbin * parent::m_axes[1].m_offset;
0162 _entries += parent::m_bin_entries[offset];
0163 }
0164 return _entries;
0165 }
0166
0167 TW bin_height_x(int aI) const {
0168 if(!parent::m_dimension) return 0;
0169 bn_t ibin;
0170 if(!parent::m_axes[0].in_range_to_absolute_index(aI,ibin)) return 0;
0171 bn_t ybins = parent::m_axes[1].bins()+2;
0172 TO offset;
0173 TW sw = 0;
0174 for(bn_t jbin=0;jbin<ybins;jbin++) {
0175 offset = ibin + jbin * parent::m_axes[1].m_offset;
0176 sw += this->get_bin_height(offset);
0177 }
0178 return sw;
0179 }
0180
0181 TN bin_entries_y(int aJ) const {
0182 if(!parent::m_dimension) return 0;
0183 bn_t jbin;
0184 if(!parent::m_axes[1].in_range_to_absolute_index(aJ,jbin)) return 0;
0185 bn_t xbins = parent::m_axes[0].bins()+2;
0186 TO offset;
0187 TN _entries = 0;
0188 for(bn_t ibin=0;ibin<xbins;ibin++) {
0189 offset = ibin + jbin * parent::m_axes[1].m_offset;
0190 _entries += parent::m_bin_entries[offset];
0191 }
0192 return _entries;
0193 }
0194
0195 TW bin_height_y(int aJ) const {
0196 if(!parent::m_dimension) return 0;
0197 bn_t jbin;
0198 if(!parent::m_axes[1].in_range_to_absolute_index(aJ,jbin)) return 0;
0199 bn_t xbins = parent::m_axes[0].bins()+2;
0200 TO offset;
0201 TW sw = 0;
0202 for(bn_t ibin=0;ibin<xbins;ibin++) {
0203 offset = ibin + jbin * parent::m_axes[1].m_offset;
0204 sw += this->get_bin_height(offset);
0205 }
0206 return sw;
0207 }
0208
0209 TC Sxyw() const {return parent::m_in_range_plane_Sxyw[0];}
0210 public:
0211 //NOTE : print is a Python keyword.
0212 void hprint(std::ostream& a_out) {
0213 // A la HPRINT.
0214 a_out << parent::dimension() << parent::title() << std::endl;
0215 a_out
0216 << " * ENTRIES = " << parent::all_entries() << std::endl;
0217
0218 // 6 | 7 | 8
0219 // -----------
0220 // 3 | 4 | 5
0221 // -----------
0222 // 0 | 1 | 2
0223
0224 TW height_0 = bin_height(axis_t::UNDERFLOW_BIN,
0225 axis_t::UNDERFLOW_BIN);
0226 TW height_2 = bin_height(axis_t::OVERFLOW_BIN,
0227 axis_t::UNDERFLOW_BIN);
0228 TW height_6 = bin_height(axis_t::UNDERFLOW_BIN,
0229 axis_t::OVERFLOW_BIN);
0230 TW height_8 = bin_height(axis_t::OVERFLOW_BIN,
0231 axis_t::OVERFLOW_BIN);
0232
0233 bn_t i,j;
0234
0235 TW height_1 = 0;
0236 TW height_7 = 0;
0237 for(i=0;i<axis_x().bins();i++){
0238 height_1 += bin_height(i,axis_t::UNDERFLOW_BIN);
0239 height_7 += bin_height(i,axis_t::OVERFLOW_BIN);
0240 }
0241
0242 TW height_3 = 0;
0243 TW height_5 = 0;
0244 for(j=0;j<axis_y().bins();j++){
0245 height_3 += bin_height(axis_t::UNDERFLOW_BIN,j);
0246 height_5 += bin_height(axis_t::OVERFLOW_BIN,j);
0247 }
0248
0249 TW height_4 = 0;
0250 for(i=0;i<axis_x().bins();i++){
0251 for(j=0;j<axis_y().bins();j++){
0252 height_4 += bin_height(i,j);
0253 }
0254 }
0255
0256 a_out
0257 << " " << height_6 << " " << height_7 << " " << height_8 << std::endl;
0258 a_out
0259 << " " << height_3 << " " << height_4 << " " << height_5 << std::endl;
0260 a_out
0261 << " " << height_0 << " " << height_1 << " " << height_2 << std::endl;
0262
0263 // Some bins :
0264 bn_t xbins = axis_x().bins();
0265 bn_t ybins = axis_y().bins();
0266 a_out
0267 << " * ENTRIES[0,0] = "
0268 << bin_entries(0,0)
0269 << " * HEIGHT[0,0] = "
0270 << bin_height(0,0)
0271 << " * ERROR[0,0] = "
0272 << bin_error(0,0)
0273 << std::endl;
0274 a_out
0275 << " * ENTRIES[N/2,N/2] = "
0276 << bin_entries(xbins/2,ybins/2)
0277 << " * HEIGHT[N/2,N/2] = "
0278 << bin_height(xbins/2,ybins/2)
0279 << " * ERROR[N/2,N/2] = "
0280 << bin_error(xbins/2,ybins/2)
0281 << std::endl;
0282 a_out
0283 << " * ENTRIES[N-1,N-1] = "
0284 << bin_entries(xbins-1,ybins-1)
0285 << " * HEIGHT[N-1,N-1] = "
0286 << bin_height(xbins-1,ybins-1)
0287 << " * ERROR[N-1,N-1] = "
0288 << bin_error(xbins-1,ybins-1)
0289 << std::endl;
0290 }
0291 protected:
0292 b2(const std::string& a_title,bn_t aXnumber,TC aXmin,TC aXmax,bn_t aYnumber,TC aYmin,TC aYmax) {
0293 parent::m_title = a_title;
0294 std::vector<bn_t> nbins;
0295 nbins.push_back(aXnumber);
0296 nbins.push_back(aYnumber);
0297 std::vector<TC> mins;
0298 mins.push_back(aXmin);
0299 mins.push_back(aYmin);
0300 std::vector<TC> maxs;
0301 maxs.push_back(aXmax);
0302 maxs.push_back(aYmax);
0303 parent::configure(2,nbins,mins,maxs);
0304 }
0305
0306 b2(const std::string& a_title,const std::vector<TC>& a_edges_x,const std::vector<TC>& a_edges_y) {
0307 parent::m_title = a_title;
0308 std::vector< std::vector<TC> > edges(2);
0309 edges[0] = a_edges_x;
0310 edges[1] = a_edges_y;
0311 parent::configure(2,edges);
0312 }
0313
0314 virtual ~b2(){}
0315 protected:
0316 b2(const b2& a_from):parent(a_from) {}
0317 b2& operator=(const b2& a_from) {parent::operator=(a_from);return *this;}
0318 public:
0319 bool configure(bn_t aXnumber,TC aXmin,TC aXmax,bn_t aYnumber,TC aYmin,TC aYmax){
0320 std::vector<bn_t> nbins;
0321 nbins.push_back(aXnumber);
0322 nbins.push_back(aYnumber);
0323 std::vector<TC> mins;
0324 mins.push_back(aXmin);
0325 mins.push_back(aYmin);
0326 std::vector<TC> maxs;
0327 maxs.push_back(aXmax);
0328 maxs.push_back(aYmax);
0329 return parent::configure(2,nbins,mins,maxs);
0330 }
0331
0332 bool configure(const std::vector<TC>& a_edges_x,const std::vector<TC>& a_edges_y){
0333 std::vector< std::vector<TC> > edges(2);
0334 edges[0] = a_edges_x;
0335 edges[1] = a_edges_y;
0336 return parent::configure(2,edges);
0337 }
0338
0339 protected:
0340 bool _find_offset(int aI,int aJ,TO& a_offset) const {
0341 if(parent::m_dimension!=2) {a_offset=0;return false;}
0342 bn_t ibin,jbin;
0343 if(!parent::m_axes[0].in_range_to_absolute_index(aI,ibin)) {a_offset=0;return false;}
0344 if(!parent::m_axes[1].in_range_to_absolute_index(aJ,jbin)) {a_offset=0;return false;}
0345 a_offset = ibin + jbin * parent::m_axes[1].m_offset;
0346 return true;
0347 }
0348 };
0349
0350 }}
0351
0352 #endif
0353
0354
0355
0356
