inlib
1.2.0
|
00001 // Copyright (C) 2010, Guy Barrand. All rights reserved. 00002 // See the file inlib.license for terms. 00003 00004 #ifndef inlib_histo_h3 00005 #define inlib_histo_h3 00006 00007 #include "b3" 00008 00009 namespace inlib { 00010 namespace histo { 00011 00012 template <class TC,class TN,class TW,class TH> 00013 class h3 : public b3<TC,TN,TW,TH> { 00014 typedef b3<TC,TN,TW,TH> parent; 00015 public: 00016 typedef typename b3<TC,TN,TW,TH>::bn_t bn_t; 00017 protected: 00018 virtual TH get_bin_height(int a_offset) const { //TH should be the same as TW 00019 return parent::m_bin_Sw[a_offset]; 00020 } 00021 public: 00022 00023 virtual TH bin_error(int aI,int aJ,int aK) const { 00024 if(parent::m_bin_number==0) return 0; 00025 bn_t ibin; 00026 if(!parent::m_axes[0].in_range_to_absolute_index(aI,ibin)) return 0; 00027 bn_t jbin; 00028 if(!parent::m_axes[1].in_range_to_absolute_index(aJ,jbin)) return 0; 00029 bn_t kbin; 00030 if(!parent::m_axes[2].in_range_to_absolute_index(aK,kbin)) return 0; 00031 bn_t offset = ibin + jbin * parent::m_axes[1].m_offset + 00032 kbin * parent::m_axes[2].m_offset; 00033 00034 return ::sqrt(parent::m_bin_Sw2[offset]); 00035 } 00036 00037 public: 00038 bool multiply(TW aFactor){ 00039 if(!parent::base_multiply(aFactor)) return false; 00040 this->update_fast_getters(); 00041 return true; 00042 } 00043 bool scale(TW aFactor) {return multiply(aFactor);} 00044 00045 void copy_from_data(const histo_data<TC,TN,TW>& a_from) { 00046 parent::base_from_data(a_from); 00047 } 00048 histo_data<TC,TN,TW> get_histo_data() const { 00049 return parent::base_get_data(); 00050 } 00051 00052 bool reset() { 00053 parent::base_reset(); 00054 this->update_fast_getters(); 00055 return true; 00056 } 00057 00058 bool fill(TC aX,TC aY,TC aZ,TW aWeight = 1) { 00059 //m_coords[0] = aX; 00060 //m_coords[1] = aY; 00061 //m_coords[2] = aZ; 00062 //return fill_bin(m_coords,aWeight); 00063 00064 if(parent::m_dimension<=0) return false; 00065 00066 bn_t ibin,jbin,kbin; 00067 if(!parent::m_axes[0].coord_to_absolute_index(aX,ibin)) return false; 00068 if(!parent::m_axes[1].coord_to_absolute_index(aY,jbin)) return false; 00069 if(!parent::m_axes[3].coord_to_absolute_index(aZ,kbin)) return false; 00070 00071 bn_t offset = ibin + jbin * parent::m_axes[1].m_offset + 00072 kbin * parent::m_axes[2].m_offset; 00073 00074 parent::m_bin_entries[offset]++; 00075 parent::m_bin_Sw[offset] += aWeight; 00076 parent::m_bin_Sw2[offset] += aWeight * aWeight; 00077 00078 TC xw = aX * aWeight; 00079 TC x2w = aX * xw; 00080 parent::m_bin_Sxw[offset][0] += xw; 00081 parent::m_bin_Sx2w[offset][0] += x2w; 00082 00083 TC yw = aY * aWeight; 00084 TC y2w = aY * yw; 00085 parent::m_bin_Sxw[offset][1] += yw; 00086 parent::m_bin_Sx2w[offset][1] += y2w; 00087 00088 TC zw = aZ * aWeight; 00089 TC z2w = aZ * zw; 00090 parent::m_bin_Sxw[offset][2] += zw; 00091 parent::m_bin_Sx2w[offset][2] += z2w; 00092 00093 bool inRange = true; 00094 if(ibin==0) inRange = false; 00095 else if(ibin==(parent::m_axes[0].m_number_of_bins+1)) inRange = false; 00096 00097 if(jbin==0) inRange = false; 00098 else if(jbin==(parent::m_axes[1].m_number_of_bins+1)) inRange = false; 00099 00100 if(kbin==0) inRange = false; 00101 else if(kbin==(parent::m_axes[2].m_number_of_bins+1)) inRange = false; 00102 00103 if(inRange) { 00104 parent::m_in_range_entries++; 00105 parent::m_in_range_Sw += aWeight; 00106 00107 parent::m_in_range_Sxw += xw; 00108 parent::m_in_range_Sx2w += x2w; 00109 00110 parent::m_in_range_Syw += yw; 00111 parent::m_in_range_Sy2w += y2w; 00112 00113 parent::m_in_range_Szw += zw; 00114 parent::m_in_range_Sz2w += z2w; 00115 } 00116 00117 return true; 00118 } 00119 00120 bool add(const h3& a_histo){ 00121 parent::base_add(a_histo); 00122 this->update_fast_getters(); 00123 return true; 00124 } 00125 bool subtract(const h3& a_histo){ 00126 parent::base_subtract(a_histo); 00127 this->update_fast_getters(); 00128 return true; 00129 } 00130 00131 bool multiply(const h3& a_histo) { 00132 if(!parent::base_multiply(a_histo)) return false; 00133 this->update_fast_getters(); 00134 return true; 00135 } 00136 00137 bool divide(const h3& a_histo) { 00138 if(!parent::base_divide(a_histo)) return false; 00139 this->update_fast_getters(); 00140 return true; 00141 } 00142 00143 public: 00144 /* 00145 // Slices : 00146 h2d* slice_xy(int aKbeg,int aKend) const; 00147 h2d* slice_yz(int aIbeg,int aIend) const; 00148 h2d* slice_xz(int aJbeg,int aJend) const; 00149 */ 00150 public: 00151 h3(const std::string& a_title, 00152 bn_t aXnumber,TC aXmin,TC aXmax, 00153 bn_t aYnumber,TC aYmin,TC aYmax, 00154 bn_t aZnumber,TC aZmin,TC aZmax) 00155 : parent(a_title,aXnumber,aXmin,aXmax, 00156 aYnumber,aYmin,aYmax, 00157 aZnumber,aZmin,aZmax) 00158 {} 00159 00160 h3(const std::string& a_title, 00161 const std::vector<TC>& aEdgesX, 00162 const std::vector<TC>& aEdgesY, 00163 const std::vector<TC>& aEdgesZ) 00164 : parent(a_title,aEdgesX,aEdgesY,aEdgesZ) 00165 {} 00166 00167 virtual ~h3(){} 00168 public: 00169 h3(const h3& a_from): parent(a_from){} 00170 h3& operator=(const h3& a_from){ 00171 parent::operator=(a_from); 00172 return *this; 00173 } 00174 }; 00175 00176 }} 00177 00178 #endif 00179 00180 00181 00182