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_c3d 00005 #define inlib_histo_c3d 00006 00007 #include "base_cloud" 00008 00009 #include <inlib/mnmx> 00010 00011 #include "h3d" 00012 00013 namespace inlib { 00014 namespace histo { 00015 00016 class c3d : public base_cloud { 00017 public: 00018 bool set_title(const std::string&); 00019 unsigned int dimension() const {return 3;} 00020 bool reset(); 00021 unsigned int entries() const; 00022 public: 00023 double sum_of_weights() const; 00024 bool convert_to_histogram(); 00025 bool is_converted() const; 00026 bool scale(double); 00027 public: 00028 bool fill(double,double,double,double = 1); 00029 double lower_edge_x() const; 00030 double upper_edge_x() const; 00031 double lower_edge_y() const; 00032 double upper_edge_y() const; 00033 double lower_edge_z() const; 00034 double upper_edge_z() const; 00035 double value_x(unsigned int) const; 00036 double value_y(unsigned int) const; 00037 double value_z(unsigned int) const; 00038 double weight(unsigned int) const; 00039 double mean_x() const; 00040 double mean_y() const; 00041 double mean_z() const; 00042 double rms_x() const; 00043 double rms_y() const; 00044 double rms_z() const; 00045 bool convert(unsigned int,double,double, 00046 unsigned int,double,double, 00047 unsigned int,double,double); 00048 bool convert(const std::vector<double>&, 00049 const std::vector<double>&, 00050 const std::vector<double>&); 00051 const histo::h3d& histogram() const; 00052 bool fill_histogram(histo::h3d& a_histo) const { 00053 unsigned int number = m_xs.size(); 00054 for(unsigned int index=0;index<number;index++) { 00055 if(!a_histo.fill(m_xs[index],m_ys[index],m_zs[index],m_ws[index])) 00056 return false; 00057 } 00058 return true; 00059 } 00060 bool set_conversion_parameters(unsigned int,double,double, 00061 unsigned int,double,double, 00062 unsigned int,double,double); 00063 public: 00064 c3d(); 00065 c3d(const std::string&,int=-1); 00066 virtual ~c3d(){delete m_histo;} 00067 public: 00068 c3d(const c3d& a_from) 00069 :base_cloud(a_from) 00070 ,m_xs(a_from.m_xs) 00071 ,m_ys(a_from.m_ys) 00072 ,m_zs(a_from.m_zs) 00073 ,m_lower_x(a_from.m_lower_x) 00074 ,m_upper_x(a_from.m_upper_x) 00075 ,m_lower_y(a_from.m_lower_y) 00076 ,m_upper_y(a_from.m_upper_y) 00077 ,m_lower_z(a_from.m_lower_z) 00078 ,m_upper_z(a_from.m_upper_z) 00079 ,m_Sxw(a_from.m_Sxw) 00080 ,m_Sx2w(a_from.m_Sx2w) 00081 ,m_Syw(a_from.m_Syw) 00082 ,m_Sy2w(a_from.m_Sy2w) 00083 ,m_Szw(a_from.m_Szw) 00084 ,m_Sz2w(a_from.m_Sz2w) 00085 ,m_cnv_x_num(a_from.m_cnv_x_num) 00086 ,m_cnv_x_min(a_from.m_cnv_x_min) 00087 ,m_cnv_x_max(a_from.m_cnv_x_max) 00088 ,m_cnv_y_num(a_from.m_cnv_y_num) 00089 ,m_cnv_y_min(a_from.m_cnv_y_min) 00090 ,m_cnv_y_max(a_from.m_cnv_y_max) 00091 ,m_cnv_z_num(a_from.m_cnv_z_num) 00092 ,m_cnv_z_min(a_from.m_cnv_z_min) 00093 ,m_cnv_z_max(a_from.m_cnv_z_max) 00094 ,m_histo(0) 00095 { 00096 if(a_from.m_histo) { 00097 m_histo = new histo::h3d(*a_from.m_histo); 00098 } 00099 } 00100 00101 c3d& operator=(const c3d& a_from) { 00102 base_cloud::operator=(a_from); 00103 m_xs = a_from.m_xs; 00104 m_ys = a_from.m_ys; 00105 m_zs = a_from.m_zs; 00106 m_lower_x = a_from.m_lower_x; 00107 m_upper_x = a_from.m_upper_x; 00108 m_lower_y = a_from.m_lower_y; 00109 m_upper_y = a_from.m_upper_y; 00110 m_lower_z = a_from.m_lower_z; 00111 m_upper_z = a_from.m_upper_z; 00112 m_Sxw = a_from.m_Sxw; 00113 m_Sx2w = a_from.m_Sx2w; 00114 m_Syw = a_from.m_Syw; 00115 m_Sy2w = a_from.m_Sy2w; 00116 m_Szw = a_from.m_Szw; 00117 m_Sz2w = a_from.m_Sz2w; 00118 m_cnv_x_num = a_from.m_cnv_x_num; 00119 m_cnv_x_min = a_from.m_cnv_x_min; 00120 m_cnv_x_max = a_from.m_cnv_x_max; 00121 m_cnv_y_num = a_from.m_cnv_y_num; 00122 m_cnv_y_min = a_from.m_cnv_y_min; 00123 m_cnv_y_max = a_from.m_cnv_y_max; 00124 m_cnv_z_num = a_from.m_cnv_z_num; 00125 m_cnv_z_min = a_from.m_cnv_z_min; 00126 m_cnv_z_max = a_from.m_cnv_z_max; 00127 delete m_histo; 00128 m_histo = 0; 00129 if(a_from.m_histo) { 00130 m_histo = new histo::h3d(*a_from.m_histo); 00131 } 00132 return *this; 00133 } 00134 00135 protected: 00136 void clear(); 00137 protected: 00138 std::vector<double> m_xs; 00139 std::vector<double> m_ys; 00140 std::vector<double> m_zs; 00141 double m_lower_x; 00142 double m_upper_x; 00143 double m_lower_y; 00144 double m_upper_y; 00145 double m_lower_z; 00146 double m_upper_z; 00147 double m_Sxw; 00148 double m_Sx2w; 00149 double m_Syw; 00150 double m_Sy2w; 00151 double m_Szw; 00152 double m_Sz2w; 00153 // 00154 unsigned int m_cnv_x_num; 00155 double m_cnv_x_min; 00156 double m_cnv_x_max; 00157 unsigned int m_cnv_y_num; 00158 double m_cnv_y_min; 00159 double m_cnv_y_max; 00160 unsigned int m_cnv_z_num; 00161 double m_cnv_z_min; 00162 double m_cnv_z_max; 00163 histo::h3d* m_histo; 00164 }; 00165 00166 }} 00167 00168 namespace inlib { 00169 namespace histo { 00170 00171 inline 00172 c3d::c3d() 00173 :base_cloud(UNLIMITED()) 00174 ,m_lower_x(0) 00175 ,m_upper_x(0) 00176 ,m_lower_y(0) 00177 ,m_upper_y(0) 00178 ,m_lower_z(0) 00179 ,m_upper_z(0) 00180 ,m_Sxw(0) 00181 ,m_Sx2w(0) 00182 ,m_Syw(0) 00183 ,m_Sy2w(0) 00184 ,m_Szw(0) 00185 ,m_Sz2w(0) 00186 ,m_cnv_x_num(0) 00187 ,m_cnv_x_min(0) 00188 ,m_cnv_x_max(0) 00189 ,m_cnv_y_num(0) 00190 ,m_cnv_y_min(0) 00191 ,m_cnv_y_max(0) 00192 ,m_cnv_z_num(0) 00193 ,m_cnv_z_min(0) 00194 ,m_cnv_z_max(0) 00195 ,m_histo(0) 00196 {} 00197 00198 inline 00199 c3d::c3d(const std::string& a_title,int aLimit) 00200 :base_cloud(aLimit) 00201 ,m_lower_x(0) 00202 ,m_upper_x(0) 00203 ,m_lower_y(0) 00204 ,m_upper_y(0) 00205 ,m_lower_z(0) 00206 ,m_upper_z(0) 00207 ,m_Sxw(0) 00208 ,m_Sx2w(0) 00209 ,m_Syw(0) 00210 ,m_Sy2w(0) 00211 ,m_Szw(0) 00212 ,m_Sz2w(0) 00213 ,m_cnv_x_num(0) 00214 ,m_cnv_x_min(0) 00215 ,m_cnv_x_max(0) 00216 ,m_cnv_y_num(0) 00217 ,m_cnv_y_min(0) 00218 ,m_cnv_y_max(0) 00219 ,m_cnv_z_num(0) 00220 ,m_cnv_z_min(0) 00221 ,m_cnv_z_max(0) 00222 ,m_histo(0) 00223 { 00224 set_title(a_title); 00225 } 00226 00227 inline 00228 bool c3d::is_converted() const {return m_histo ? true : false;} 00229 00230 inline 00231 void c3d::clear(){ 00232 m_lower_x = 0; 00233 m_upper_x = 0; 00234 m_lower_y = 0; 00235 m_upper_y = 0; 00236 m_lower_z = 0; 00237 m_upper_z = 0; 00238 m_Sw = 0; 00239 m_Sxw = 0; 00240 m_Sx2w = 0; 00241 m_Syw = 0; 00242 m_Sy2w = 0; 00243 m_Szw = 0; 00244 m_Sz2w = 0; 00245 m_xs.clear(); 00246 m_ys.clear(); 00247 m_zs.clear(); 00248 m_ws.clear(); 00249 } 00250 00251 inline 00252 bool c3d::convert( 00253 unsigned int aBinsX,double aLowerEdgeX,double aUpperEdgeX 00254 ,unsigned int aBinsY,double aLowerEdgeY,double aUpperEdgeY 00255 ,unsigned int aBinsZ,double aLowerEdgeZ,double aUpperEdgeZ 00256 ) { 00257 if(m_histo) return true; // Done. 00258 m_histo = new histo::h3d(this->title(), 00259 aBinsX,aLowerEdgeX,aUpperEdgeX, 00260 aBinsY,aLowerEdgeY,aUpperEdgeY, 00261 aBinsZ,aLowerEdgeZ,aUpperEdgeZ); 00262 if(!m_histo) return false; 00263 bool status = fill_histogram(*m_histo); 00264 clear(); 00265 return status; 00266 } 00267 00268 inline 00269 bool c3d::convert_to_histogram(){ 00270 if( (m_cnv_x_num<=0) || (m_cnv_x_max<=m_cnv_x_min) || 00271 (m_cnv_y_num<=0) || (m_cnv_y_max<=m_cnv_y_min) || 00272 (m_cnv_z_num<=0) || (m_cnv_z_max<=m_cnv_z_min) ) { 00273 double dx = 0.01 * (upper_edge_x() - lower_edge_x())/BINS(); 00274 double dy = 0.01 * (upper_edge_y() - lower_edge_y())/BINS(); 00275 double dz = 0.01 * (upper_edge_z() - lower_edge_z())/BINS(); 00276 return convert(BINS(),lower_edge_x(),upper_edge_x()+dx, 00277 BINS(),lower_edge_y(),upper_edge_y()+dy, 00278 BINS(),lower_edge_z(),upper_edge_z()+dz); 00279 } else { 00280 return convert(m_cnv_x_num,m_cnv_x_min,m_cnv_x_max, 00281 m_cnv_y_num,m_cnv_y_min,m_cnv_y_max, 00282 m_cnv_z_num,m_cnv_z_min,m_cnv_z_max); 00283 } 00284 } 00285 00286 inline 00287 bool c3d::set_title(const std::string& a_title){ 00288 m_title = a_title; 00289 if(m_histo) m_histo->set_title(a_title); 00290 return true; 00291 } 00292 00293 inline 00294 bool c3d::scale(double a_scale) { 00295 if(m_histo) { 00296 return m_histo->scale(a_scale); 00297 } else { 00298 unsigned int number = m_ws.size(); 00299 for(unsigned int index=0;index<number;index++) m_ws[index] *= a_scale; 00300 m_Sw *= a_scale; 00301 m_Sxw *= a_scale; 00302 m_Sx2w *= a_scale; 00303 m_Syw *= a_scale; 00304 m_Sy2w *= a_scale; 00305 m_Szw *= a_scale; 00306 m_Sz2w *= a_scale; 00307 return true; 00308 } 00309 } 00310 00311 inline 00312 bool c3d::set_conversion_parameters( 00313 unsigned int aCnvXnumber,double aCnvXmin,double aCnvXmax 00314 ,unsigned int aCnvYnumber,double aCnvYmin,double aCnvYmax 00315 ,unsigned int aCnvZnumber,double aCnvZmin,double aCnvZmax 00316 ){ 00317 m_cnv_x_num = aCnvXnumber; 00318 m_cnv_x_min = aCnvXmin; 00319 m_cnv_x_max = aCnvXmax; 00320 m_cnv_y_num = aCnvYnumber; 00321 m_cnv_y_min = aCnvYmin; 00322 m_cnv_y_max = aCnvYmax; 00323 m_cnv_z_num = aCnvZnumber; 00324 m_cnv_z_min = aCnvZmin; 00325 m_cnv_z_max = aCnvZmax; 00326 return true; 00327 } 00328 00329 00330 inline 00331 const h3d& c3d::histogram() const { 00332 if(!m_histo) const_cast<c3d&>(*this).convert_to_histogram(); 00333 return *m_histo; 00334 } 00335 00336 inline 00337 bool c3d::reset() { 00338 clear(); 00339 delete m_histo; 00340 m_histo = 0; 00341 return true; 00342 } 00343 00344 inline 00345 bool c3d::fill(double aX,double aY,double aZ,double aW){ 00346 if(!m_histo && (m_limit!=UNLIMITED()) && ((int)m_xs.size()>=m_limit)){ 00347 convert_to_histogram(); 00348 } 00349 00350 if(m_histo) { 00351 return m_histo->fill(aX,aY,aZ,aW); 00352 } else { 00353 if(m_xs.size()) { 00354 m_lower_x = inlib::mn<double>(aX,m_lower_x); 00355 m_upper_x = inlib::mx<double>(aX,m_upper_x); 00356 } else { 00357 m_lower_x = aX; 00358 m_upper_x = aX; 00359 } 00360 if(m_ys.size()) { 00361 m_lower_y = inlib::mn<double>(aY,m_lower_y); 00362 m_upper_y = inlib::mx<double>(aY,m_upper_y); 00363 } else { 00364 m_lower_y = aY; 00365 m_upper_y = aY; 00366 } 00367 if(m_zs.size()) { 00368 m_lower_z = inlib::mn<double>(aZ,m_lower_z); 00369 m_upper_z = inlib::mx<double>(aZ,m_upper_z); 00370 } else { 00371 m_lower_z = aZ; 00372 m_upper_z = aZ; 00373 } 00374 m_xs.push_back(aX); 00375 m_ys.push_back(aY); 00376 m_zs.push_back(aZ); 00377 m_ws.push_back(aW); 00378 m_Sw += aW; 00379 double xw = aX * aW; 00380 m_Sxw += xw; 00381 m_Sx2w += aX * xw; 00382 double yw = aY * aW; 00383 m_Syw += yw; 00384 m_Sy2w += aY * yw; 00385 double zw = aZ * aW; 00386 m_Szw += zw; 00387 m_Sz2w += aZ * zw; 00388 return true; 00389 } 00390 } 00391 00392 inline 00393 bool c3d::convert( 00394 const std::vector<double>& aEdgesX 00395 ,const std::vector<double>& aEdgesY 00396 ,const std::vector<double>& aEdgesZ 00397 ) { 00398 if(m_histo) return true; 00399 m_histo = new histo::h3d(this->title(), 00400 aEdgesX,aEdgesY,aEdgesZ); 00401 if(!m_histo) return false; 00402 bool status = fill_histogram(*m_histo); 00403 clear(); 00404 return status; 00405 } 00406 00407 inline 00408 double c3d::sum_of_weights() const { 00409 return (m_histo ? m_histo->sum_bin_heights() : m_Sw); 00410 } 00411 00412 inline 00413 unsigned int c3d::entries() const { 00414 return m_histo ? m_histo->all_entries() : m_ws.size(); 00415 } 00416 00417 inline 00418 double c3d::lower_edge_x() const { 00419 return m_histo ? m_histo->axis_x().lower_edge() : m_lower_x; 00420 } 00421 inline 00422 double c3d::lower_edge_y() const { 00423 return m_histo ? m_histo->axis_y().lower_edge() : m_lower_y; 00424 } 00425 inline 00426 double c3d::lower_edge_z() const { 00427 return m_histo ? m_histo->axis_z().lower_edge() : m_lower_z; 00428 } 00429 inline 00430 double c3d::upper_edge_x() const { 00431 return m_histo ? m_histo->axis_x().upper_edge() : m_upper_x; 00432 } 00433 inline 00434 double c3d::upper_edge_y() const { 00435 return m_histo ? m_histo->axis_y().upper_edge() : m_upper_y; 00436 } 00437 inline 00438 double c3d::upper_edge_z() const { 00439 return m_histo ? m_histo->axis_z().upper_edge() : m_upper_z; 00440 } 00441 inline 00442 double c3d::value_x(unsigned int aIndex) const { 00443 return m_histo ? 0 : m_xs[aIndex]; 00444 } 00445 inline 00446 double c3d::value_y(unsigned int aIndex) const { 00447 return m_histo ? 0 : m_ys[aIndex]; 00448 } 00449 inline 00450 double c3d::value_z(unsigned int aIndex) const { 00451 return m_histo ? 0 : m_zs[aIndex]; 00452 } 00453 inline 00454 double c3d::weight(unsigned int aIndex) const { 00455 return m_histo ? 0 : m_ws[aIndex]; 00456 } 00457 inline 00458 double c3d::mean_x() const { 00459 return m_histo ? m_histo->mean_x() : (m_Sw?m_Sxw/m_Sw:0); 00460 } 00461 inline 00462 double c3d::mean_y() const { 00463 return m_histo ? m_histo->mean_y() : (m_Sw?m_Syw/m_Sw:0); 00464 } 00465 inline 00466 double c3d::mean_z() const { 00467 return m_histo ? m_histo->mean_z() : (m_Sw?m_Szw/m_Sw:0); 00468 } 00469 inline 00470 double c3d::rms_x() const { 00471 double rms = 0; //FIXME nan. 00472 if(m_histo) { 00473 rms = m_histo->rms_x(); 00474 } else { 00475 if(m_Sw==0) { 00476 } else { 00477 double mean = m_Sxw / m_Sw; 00478 rms = ::sqrt(::fabs( (m_Sx2w / m_Sw) - mean * mean)); 00479 } 00480 } 00481 return rms; 00482 } 00483 inline 00484 double c3d::rms_y() const { 00485 double rms = 0; //FIXME nan. 00486 if(m_histo) { 00487 rms = m_histo->rms_y(); 00488 } else { 00489 if(m_Sw==0) { 00490 } else { 00491 double mean = m_Syw / m_Sw; 00492 rms = ::sqrt(::fabs( (m_Sy2w / m_Sw) - mean * mean)); 00493 } 00494 } 00495 return rms; 00496 } 00497 inline 00498 double c3d::rms_z() const { 00499 double rms = 0; //FIXME nan. 00500 if(m_histo) { 00501 rms = m_histo->rms_z(); 00502 } else { 00503 if(m_Sw==0) { 00504 } else { 00505 double mean = m_Szw / m_Sw; 00506 rms = ::sqrt(::fabs( (m_Sz2w / m_Sw) - mean * mean)); 00507 } 00508 } 00509 return rms; 00510 } 00511 00512 }} 00513 00514 #endif