9 #ifndef LWH_Histogram1D_H
10 #define LWH_Histogram1D_H
15 #include "AIHistogram1D.h"
16 #include "ManagedObject.h"
42 : fax(new
Axis(n, lo, up)), vax(0),
43 sum(n + 2), sumw(n + 2), sumw2(n + 2), sumxw(n + 2), sumx2w(n + 2) {
52 sum(edges.size() + 1), sumw(edges.size() + 1), sumw2(edges.size() + 1),
53 sumxw(edges.size() + 1), sumx2w(edges.size() + 1) {
61 : IBaseHistogram(h), IHistogram(h), IHistogram1D(h),
ManagedObject(h),
62 fax(0), vax(0), sum(h.sum), sumw(h.sumw), sumw2(h.sumw2),
63 sumxw(h.sumxw), sumx2w(h.sumx2w) {
65 if ( vax ) ax = vax =
new VariAxis(*hvax);
66 else ax = fax =
new Axis(dynamic_cast<const Axis &>(*h.
ax));
104 throw std::runtime_error(
"LWH cannot handle annotations");
111 throw std::runtime_error(
"LWH cannot handle annotations");
127 sum = std::vector<int>(ax->bins() + 2);
128 sumw = std::vector<double>(ax->bins() + 2);
129 sumxw = std::vector<double>(ax->bins() + 2);
130 sumx2w = std::vector<double>(ax->bins() + 2);
131 sumw2 = std::vector<double>(ax->bins() + 2);
142 for (
int i = 2; i < ax->bins() + 2; ++i ) si += sum[i];
154 return entries() + extraEntries();
162 return sum[0] + sum[1];
173 for (
int i = 2; i < ax->bins() + 2; ++i ) {
188 for (
int i = 2; i < ax->bins() + 2; ++i ) sw += sumw[i];
198 return sumBinHeights() + sumExtraBinHeights();
206 return sumw[0] + sumw[1];
215 double minw = sumw[2];
216 for (
int i = 3; i < ax->bins() + 2; ++i ) minw = std::min(minw, sumw[i]);
226 double maxw = sumw[2];
227 for (
int i = 3; i < ax->bins() + 2; ++i ) maxw = std::max(maxw, sumw[i]);
238 bool fill(
double x,
double weight = 1.) {
239 int i = ax->coordToIndex(x) + 2;
242 sumxw[i] += x*weight;
243 sumx2w[i] += x*x*weight;
244 sumw2[i] += weight*weight;
245 return weight >= 0 && weight <= 1;
255 return sumw[i] != 0.0? sumxw[i]/sumw[i]:
256 ( vax? vax->binMidPoint(index): fax->binMidPoint(index) );
266 return sumw[i] == 0.0 || sum[i] < 2? ax->binWidth(index):
267 std::sqrt(std::max(sumw[i]*sumx2w[i] - sumxw[i]*sumxw[i], 0.0))/sumw[i];
277 return sum[index + 2];
287 return sumw[index + 2];
297 return std::sqrt(sumw2[index + 2]);
307 for (
int i = 2; i < ax->bins() + 2; ++i ) {
311 return s != 0.0? sx/s: 0.0;
322 for (
int i = 2; i < ax->bins() + 2; ++i ) {
327 return s != 0.0? std::sqrt(std::max(s*sx2 - sx*sx, 0.0))/s:
328 ax->upperEdge() - ax->lowerEdge();
347 return ax->coordToIndex(coord);
356 if ( ax->upperEdge() != h.
ax->upperEdge() ||
357 ax->lowerEdge() != h.
ax->lowerEdge() ||
358 ax->bins() != h.
ax->bins() )
return false;
359 for (
int i = 0; i < ax->bins() + 2; ++i ) {
361 sumw[i] += h.
sumw[i];
362 sumxw[i] += h.
sumxw[i];
364 sumw2[i] += h.
sumw2[i];
374 bool add(
const IHistogram1D & hist) {
375 return add(dynamic_cast<const Histogram1D &>(hist));
383 for (
int i = 0; i < ax->bins() + 2; ++i ) {
400 double oldintg = sumAllBinHeights();
401 if ( oldintg == 0.0 )
return;
402 for (
int i = 0; i < ax->bins() + 2; ++i ) {
403 double fac = intg/oldintg;
404 if ( i >= 2 ) fac /= (ax->binUpperEdge(i - 2) - ax->binLowerEdge(i - 2));
417 double intg = sumw[0] + sumw[1];
418 for (
int i = 2; i < ax->bins() + 2; ++i )
419 intg += sumw[i]*(ax->binUpperEdge(i - 2) - ax->binLowerEdge(i - 2));
427 void *
cast(
const std::string &)
const {
434 bool writeXML(std::ostream & os, std::string path, std::string name) {
435 os <<
" <histogram1d name=\"" << name
436 <<
"\"\n title=\"" << title()
437 <<
"\" path=\"" << path
438 <<
"\">\n <axis max=\"" << ax->upperEdge()
439 <<
"\" numberOfBins=\"" << ax->bins()
440 <<
"\" min=\"" << ax->lowerEdge()
441 <<
"\" direction=\"x\"";
444 for (
int i = 0, N = ax->bins() - 1; i < N; ++i )
445 os <<
" <binBorder value=\"" << ax->binUpperEdge(i) <<
"\"/>\n";
450 os <<
" <statistics entries=\"" << entries()
451 <<
"\">\n <statistic mean=\"" << mean()
452 <<
"\" direction=\"x\"\n rms=\"" << rms()
453 <<
"\"/>\n </statistics>\n <data1d>\n";
454 for (
int i = 0; i < ax->bins() + 2; ++i )
if ( sum[i] ) {
455 os <<
" <bin1d binNum=\"";
456 if ( i == 0 ) os <<
"UNDERFLOW";
457 else if ( i == 1 ) os <<
"OVERFLOW";
459 os <<
"\" entries=\"" << sum[i]
460 <<
"\" height=\"" << sumw[i]
461 <<
"\"\n error=\"" << std::sqrt(sumw2[i])
462 <<
"\" error2=\"" << sumw2[i]
463 <<
"\"\n weightedMean=\"" << binMean(i - 2)
464 <<
"\" weightedRms=\"" << binRms(i - 2)
467 os <<
" </data1d>\n </histogram1d>" << std::endl;
476 bool writeFLAT(std::ostream & os, std::string path, std::string name) {
477 os <<
"# " << path <<
"/" << name <<
" " << ax->lowerEdge()
478 <<
" " << ax->bins() <<
" " << ax->upperEdge()
479 <<
" \"" << title() <<
" \"" << std::endl;
480 for (
int i = 2; i < ax->bins() + 2; ++i )
481 os << 0.5*(ax->binLowerEdge(i - 2) + ax->binUpperEdge(i - 2)) <<
" "
482 << sumw[i] <<
" " << sqrt(sumw2[i]) <<
" " << sum[i] << std::endl;