| # Gaudi Framework, version v25r0 |

Home |  | Generated: Mon Feb 17 2014

- Main Page
- Related Pages
- Modules
- Namespaces
- Classes
- Files

- File List
- File Members

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
H1D.cpp
Go to the documentation of this file.
1 // $Id: H1D.cpp,v 1.14 2007/07/16 13:36:17 hmd Exp $
2 #ifdef __ICC
3 // disable icc remark #2259: non-pointer conversion from "X" to "Y" may lose significant bits
4 // TODO: To be removed, since it comes from ROOT TMathBase.h
5 #pragma warning(disable:2259)
6 #endif
7 #ifdef WIN32
8 // Disable warning
9 // warning C4996: 'sprintf': This function or variable may be unsafe.
10 // coming from TString.h
11 #pragma warning(disable:4996)
12 #endif
13 #include "H1D.h"
14 #include "GaudiPI.h"
15 #include "GaudiKernel/StreamBuffer.h"
16 #include "GaudiKernel/ObjectFactory.h"
17
18 std::pair<DataObject*,AIDA::IHistogram1D*> Gaudi::createH1D(const std::string& title,int nBins,double xlow, double xup) {
19 Histogram1D* p = new Histogram1D(new TH1D(title.c_str(),title.c_str(),nBins,xlow,xup));
20 return std::pair<DataObject*,AIDA::IHistogram1D*>(p,p);
21 }
22
23 std::pair<DataObject*,AIDA::IHistogram1D*> Gaudi::createH1D(const std::string& title, const Edges& e) {
24 Histogram1D* p = new Histogram1D(new TH1D(title.c_str(),title.c_str(),e.size()-1,&e.front()));
25 return std::pair<DataObject*,AIDA::IHistogram1D*>(p,p);
26 }
27
28 std::pair<DataObject*,AIDA::IHistogram1D*> Gaudi::createH1D(const AIDA::IHistogram1D& hist) {
29 TH1D *h = getRepresentation<AIDA::IHistogram1D,TH1D>(hist);
30 Histogram1D *n = h ? new Histogram1D(new TH1D(*h)) : 0;
31 return std::pair<DataObject*,AIDA::IHistogram1D*>(n,n);
32 }
33 namespace Gaudi {
34 template<> void *Generic1D<AIDA::IHistogram1D,TH1D>::cast(const std::string & className) const {
35 if (className == "AIDA::IHistogram1D")
36 return const_cast<AIDA::IHistogram1D*>((AIDA::IHistogram1D*)this);
37 if (className == "AIDA::IHistogram")
38 return const_cast<AIDA::IHistogram*>((AIDA::IHistogram*)this);
39 return 0;
40 }
41
42 template<> int Generic1D<AIDA::IHistogram1D,TH1D>::binEntries (int index) const {
43 if (binHeight(index)<=0) return 0;
44 double xx = binHeight(index)/binError(index);
45 return int(xx*xx+0.5);
46 }
47
48 template <>
49 void Generic1D<AIDA::IHistogram1D,TH1D>::adoptRepresentation(TObject* rep) {
50 TH1D* imp = dynamic_cast<TH1D*>(rep);
51 if ( imp ) {
52 if ( m_rep ) delete m_rep;
53 m_rep = imp;
54 return;
55 }
56 throw std::runtime_error("Cannot adopt native histogram representation.");
57 }
58 }
59
60 Gaudi::Histogram1D::Histogram1D() {
61 m_rep = new TH1D();
62 init("",false);
63 }
64
65 Gaudi::Histogram1D::Histogram1D(TH1D* rep) {
66 m_rep = rep;
67 init(m_rep->GetTitle());
68 initSums();
69 }
70
71 void Gaudi::Histogram1D::init(const std::string& title, bool initialize_axis) {
72 m_classType = "IHistogram1D";
73 if ( initialize_axis ) {
74 m_axis.initialize(m_rep->GetXaxis(),false);
75 }
76 const TArrayD* a = m_rep->GetSumw2();
77 if ( 0 == a || (a && a->GetSize()==0) ) m_rep->Sumw2();
78 setTitle(title);
79 m_rep->SetDirectory(0);
80 m_sumEntries = 0;
81 m_sumwx = 0;
82 }
83
84 void Gaudi::Histogram1D::initSums() {
85 m_sumwx = 0;
86 m_sumEntries = 0;
87 for(int i=1, n=m_rep->GetNbinsX(); i<=n; ++i) {
88 m_sumwx += m_rep->GetBinContent(i)*m_rep->GetBinCenter(i);
89 m_sumEntries += (int)m_rep->GetBinContent(i);
90 }
91 }
92
93 bool Gaudi::Histogram1D::reset() {
94 m_sumwx = 0;
95 m_sumEntries = 0;
96 return Base::reset();
97 }
98
100 void Gaudi::Histogram1D::adoptRepresentation(TObject* rep) {
101 Gaudi::Generic1D<AIDA::IHistogram1D,TH1D>::adoptRepresentation(rep);
102 if ( m_rep ) {
103 init(m_rep->GetTitle());
104 initSums();
105 }
106 }
107
108 bool Gaudi::Histogram1D::setBinContents(int i,int entries ,double height,double error,double centre) {
109 m_rep->SetBinContent(rIndex(i),height);
110 m_rep->SetBinError(rIndex(i),error);
111 // accumulate sumwx for in range bins
112 if (i != AIDA::IAxis::UNDERFLOW_BIN && i != AIDA::IAxis::OVERFLOW_BIN )
113 m_sumwx += centre*height;
114 m_sumEntries += entries;
115 return true;
116 }
117
118 #ifdef __ICC
119 // disable icc remark #1572: floating-point equality and inequality comparisons are unreliable
120 // The comparison are meant
121 #pragma warning(push)
122 #pragma warning(disable:1572)
123 #endif
124 bool Gaudi::Histogram1D::setRms(double rms) {
125 m_rep->SetEntries(m_sumEntries);
126 std::vector<double> stat(11);
127 // sum weights
128 stat[0] = sumBinHeights();
129 stat[1] = 0;
130 if (equivalentBinEntries() != 0)
131 stat[1] = ( sumBinHeights() * sumBinHeights() ) / equivalentBinEntries();
132 stat[2] = m_sumwx;
133 double mean = 0;
134 if ( sumBinHeights() != 0 ) mean = m_sumwx/ sumBinHeights();
135 stat[3] = ( mean*mean + rms*rms )* sumBinHeights();
136 m_rep->PutStats(&stat.front());
137 return true;
138 }
139
140 // set histogram statistics
141 bool Gaudi::Histogram1D::setStatistics(int allEntries,double eqBinEntries,double mean,double rms) {
142 m_rep->SetEntries(allEntries);
143 // fill statistcal vector for Root
144 std::vector<double> stat(11);
145 // sum weights
146 stat[0] = sumBinHeights();
147 // sum weights **2
148 stat[1] = 0;
149 if (eqBinEntries != 0)
150 stat[1] = ( sumBinHeights() * sumBinHeights() ) / eqBinEntries;
151 // sum weights * x
152 stat[2] = mean*sumBinHeights();
153 // sum weight * x **2
154 stat[3] = ( mean*mean + rms*rms )* sumBinHeights();
155 m_rep->PutStats(&stat.front());
156 return true;
157 }
158
159 bool Gaudi::Histogram1D::fill ( double x,double weight ) {
160 (weight == 1.) ? m_rep->Fill(x) : m_rep->Fill(x,weight);
161 return true;
162 }
163
164 void Gaudi::Histogram1D::copyFromAida(const AIDA::IHistogram1D & h) {
165 // implement here the copy
166 std::string tit = h.title()+"Copy";
167 delete m_rep;
168 if (h.axis().isFixedBinning() ) {
169 m_rep = new TH1D(tit.c_str(),tit.c_str(),h.axis().bins(),h.axis().lowerEdge(),h.axis().upperEdge());
170 }
171 else {
172 Edges e;
173 for (int i =0; i < h.axis().bins(); ++i) {
174 e.push_back(h.axis().binLowerEdge(i));
175 }
176 // add also upperedges at the end
177 e.push_back(h.axis().upperEdge() );
178 m_rep = new TH1D(tit.c_str(),tit.c_str(),e.size()-1,&e.front());
179 }
180 m_axis.initialize(m_rep->GetXaxis(),false);
181 m_rep->Sumw2();
182 m_sumEntries = 0;
183 m_sumwx = 0;
184 // sumw
185 double sumw = h.sumBinHeights();
186 // sumw2
187 double sumw2 = 0;
188 if (h.equivalentBinEntries() != 0)
189 sumw2 = ( sumw * sumw ) /h.equivalentBinEntries();
190
191 double sumwx = h.mean()*h.sumBinHeights();
192 double sumwx2 = (h.mean()*h.mean() + h.rms()*h.rms() )*h.sumBinHeights();
193
194 // copy the contents in
195 for (int i=-2; i < axis().bins(); ++i) {
196 // root binning starts from one !
197 m_rep->SetBinContent(rIndex(i),h.binHeight(i) );
198 m_rep->SetBinError(rIndex(i),h.binError(i) );
199 }
200 // need to do set entries after setting contents otherwise root will recalulate them
201 // taking into account how many time SetBinContents() has been called
202 m_rep->SetEntries(h.allEntries());
203 // stat vector
204 std::vector<double> stat(11);
205 stat[0] = sumw;
206 stat[1] = sumw2;
207 stat[2] = sumwx;
208 stat[3] = sumwx2;
209 m_rep->PutStats(&stat.front());
210 }
211
212 #ifdef __ICC
213 // re-enable icc remark #1572
214 #pragma warning(pop)
215 #endif
216
217 StreamBuffer& Gaudi::Histogram1D::serialize(StreamBuffer& s) {
218 //DataObject::serialize(s);
219 std::string title;
220 int size;
221 s >> size;
222 for (int j = 0; j < size; j++) {
223 std::string key, value;
224 s >> key >> value;
225 annotation().addItem (key, value);
226 if ("Title" == key) {
227 title = value;
228 }
229 }
230 double lowerEdge, upperEdge, binHeight, binError;
231 int isFixedBinning, bins;
232 s >> isFixedBinning >> bins;
233
234 if ( m_rep ) delete m_rep;
235 if ( isFixedBinning ) {
236 s >> lowerEdge >> upperEdge;
237 m_rep = new TH1D(title.c_str(),title.c_str(),bins,lowerEdge,upperEdge);
238 } else {
239 Edges edges;
240 edges.resize(bins);
241 for ( int i = 0; i <= bins; ++i )
242 s >> *(double*)&edges[i];
243 m_rep = new TH1D(title.c_str(),title.c_str(),edges.size()-1,&edges.front());
244 }
245 m_axis.initialize(m_rep->GetXaxis(),true);
246 m_rep->Sumw2();
247 m_sumEntries = 0;
248 m_sumwx = 0;
249
250 for ( int i = 0; i <= bins + 1; ++i ) {
251 s >> binHeight >> binError;
252 m_rep->SetBinContent( i, binHeight );
253 m_rep->SetBinError( i, binError );
254 }
255 Stat_t allEntries;
256 s >> allEntries;
257 m_rep->SetEntries( allEntries );
258 Stat_t stats[4]; // stats array
259 s >> stats[0] >> stats[1] >> stats[2] >> stats[3];
260 m_rep->PutStats( stats );
261 return s;
262 }
263
264 StreamBuffer& Gaudi::Histogram1D::serialize(StreamBuffer& s) const {
265 //DataObject::serialize(s);
266 s << static_cast<int>( annotation().size() );
267 for (int i = 0; i < annotation().size(); i++) {
268 s << annotation().key(i);
269 s << annotation().value(i);
270 }
271 const AIDA::IAxis & axis( this->axis() );
272 const int isFixedBinning = axis.isFixedBinning();
273 const int bins = axis.bins();
274 s << isFixedBinning << bins;
275 if ( isFixedBinning ) {
276 s << axis.lowerEdge();
277 }
278 else {
279 for ( int i = 0; i < bins; ++i )
280 s << axis.binLowerEdge(i);
281 }
282 s << axis.upperEdge();
283 for ( int i = 0; i <= bins + 1; ++i )
284 s << m_rep->GetBinContent(i) << m_rep->GetBinError( i );
285
286 s << m_rep->GetEntries();
287 Stat_t stats[4]; // stats array
288 m_rep->GetStats( stats );
289 s << stats[0] << stats[1] << stats[2] << stats[3];
290 return s;
291 }
292
293 typedef Gaudi::Histogram1D H1D;
294 DECLARE_DATAOBJECT_FACTORY(H1D)