Logo ROOT   6.30.04
Reference Guide
 All Namespaces Files Pages
TGraphBentErrors.cxx
Go to the documentation of this file.
1 // @(#)root/hist:$Id$
2 // Author: Dave Morrison 30/06/2003
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2003, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 #include <string.h>
13 
14 #include "Riostream.h"
15 #include "TROOT.h"
16 #include "TGraphBentErrors.h"
17 #include "TStyle.h"
18 #include "TMath.h"
19 #include "TArrow.h"
20 #include "TBox.h"
21 #include "TVirtualPad.h"
22 #include "TH1.h"
23 #include "TF1.h"
24 #include "TClass.h"
25 
26 ClassImp(TGraphBentErrors);
27 
28 
29 ////////////////////////////////////////////////////////////////////////////////
30 
31 /** \class TGraphBentErrors
32  \ingroup Hist
33 A TGraphBentErrors is a TGraph with bent, asymmetric error bars.
34 
35 The TGraphBentErrors painting is performed thanks to the TGraphPainter
36 class. All details about the various painting options are given in this class.
37 
38 The picture below gives an example:
39 Begin_Macro(source)
40 {
41  auto c1 = new TCanvas("c1","A Simple Graph with bent error bars",200,10,700,500);
42  const Int_t n = 10;
43  Double_t x[n] = {-0.22, 0.05, 0.25, 0.35, 0.5, 0.61,0.7,0.85,0.89,0.95};
44  Double_t y[n] = {1,2.9,5.6,7.4,9,9.6,8.7,6.3,4.5,1};
45  Double_t exl[n] = {.05,.1,.07,.07,.04,.05,.06,.07,.08,.05};
46  Double_t eyl[n] = {.8,.7,.6,.5,.4,.4,.5,.6,.7,.8};
47  Double_t exh[n] = {.02,.08,.05,.05,.03,.03,.04,.05,.06,.03};
48  Double_t eyh[n] = {.6,.5,.4,.3,.2,.2,.3,.4,.5,.6};
49  Double_t exld[n] = {.0,.0,.0,.0,.0,.0,.0,.0,.0,.0};
50  Double_t eyld[n] = {.0,.0,.05,.0,.0,.0,.0,.0,.0,.0};
51  Double_t exhd[n] = {.0,.0,.0,.0,.0,.0,.0,.0,.0,.0};
52  Double_t eyhd[n] = {.0,.0,.0,.0,.0,.0,.0,.0,.05,.0};
53  auto gr = new TGraphBentErrors(n,x,y,exl,exh,eyl,eyh,exld,exhd,eyld,eyhd);
54  gr->SetTitle("TGraphBentErrors Example");
55  gr->SetMarkerColor(4);
56  gr->SetMarkerStyle(21);
57  gr->Draw("ALP");
58 }
59 End_Macro
60 */
61 
62 
63 ////////////////////////////////////////////////////////////////////////////////
64 /// TGraphBentErrors default constructor.
65 
66 TGraphBentErrors::TGraphBentErrors(): TGraph()
67 {
68  if (!CtorAllocate()) return;
69 }
70 
71 
72 ////////////////////////////////////////////////////////////////////////////////
73 /// TGraphBentErrors copy constructor
74 
75 TGraphBentErrors::TGraphBentErrors(const TGraphBentErrors &gr)
76  : TGraph(gr)
77 {
78  if (!CtorAllocate()) return;
79  Int_t n = fNpoints*sizeof(Double_t);
80  memcpy(fEXlow, gr.fEXlow, n);
81  memcpy(fEYlow, gr.fEYlow, n);
82  memcpy(fEXhigh, gr.fEXhigh, n);
83  memcpy(fEYhigh, gr.fEYhigh, n);
84  memcpy(fEXlowd, gr.fEXlowd, n);
85  memcpy(fEYlowd, gr.fEYlowd, n);
86  memcpy(fEXhighd, gr.fEXhighd, n);
87  memcpy(fEYhighd, gr.fEYhighd, n);
88 }
89 
90 
91 ////////////////////////////////////////////////////////////////////////////////
92 /// TGraphBentErrors normal constructor.
93 ///
94 /// the arrays are preset to zero
95 
96 TGraphBentErrors::TGraphBentErrors(Int_t n)
97  : TGraph(n)
98 {
99  if (!CtorAllocate()) return;
100  FillZero(0, fNpoints);
101 }
102 
103 
104 ////////////////////////////////////////////////////////////////////////////////
105 /// TGraphBentErrors normal constructor.
106 ///
107 /// if exl,h or eyl,h are null, the corresponding arrays are preset to zero
108 
109 TGraphBentErrors::TGraphBentErrors(Int_t n,
110  const Float_t *x, const Float_t *y,
111  const Float_t *exl, const Float_t *exh,
112  const Float_t *eyl, const Float_t *eyh,
113  const Float_t *exld, const Float_t *exhd,
114  const Float_t *eyld, const Float_t *eyhd)
115  : TGraph(n,x,y)
116 {
117  if (!CtorAllocate()) return;
118 
119  for (Int_t i=0;i<n;i++) {
120  if (exl) fEXlow[i] = exl[i];
121  else fEXlow[i] = 0;
122  if (exh) fEXhigh[i] = exh[i];
123  else fEXhigh[i] = 0;
124  if (eyl) fEYlow[i] = eyl[i];
125  else fEYlow[i] = 0;
126  if (eyh) fEYhigh[i] = eyh[i];
127  else fEYhigh[i] = 0;
128 
129  if (exld) fEXlowd[i] = exld[i];
130  else fEXlowd[i] = 0;
131  if (exhd) fEXhighd[i] = exhd[i];
132  else fEXhighd[i] = 0;
133  if (eyld) fEYlowd[i] = eyld[i];
134  else fEYlowd[i] = 0;
135  if (eyhd) fEYhighd[i] = eyhd[i];
136  else fEYhighd[i] = 0;
137  }
138 }
139 
140 
141 ////////////////////////////////////////////////////////////////////////////////
142 /// TGraphBentErrors normal constructor.
143 ///
144 /// if exl,h or eyl,h are null, the corresponding arrays are preset to zero
145 
146 TGraphBentErrors::TGraphBentErrors(Int_t n,
147  const Double_t *x, const Double_t *y,
148  const Double_t *exl, const Double_t *exh,
149  const Double_t *eyl, const Double_t *eyh,
150  const Double_t *exld, const Double_t *exhd,
151  const Double_t *eyld, const Double_t *eyhd)
152  : TGraph(n,x,y)
153 {
154  if (!CtorAllocate()) return;
155  n = sizeof(Double_t)*fNpoints;
156 
157  if (exl) memcpy(fEXlow, exl, n);
158  else memset(fEXlow, 0, n);
159  if (exh) memcpy(fEXhigh, exh, n);
160  else memset(fEXhigh, 0, n);
161  if (eyl) memcpy(fEYlow, eyl, n);
162  else memset(fEYlow, 0, n);
163  if (eyh) memcpy(fEYhigh, eyh, n);
164  else memset(fEYhigh, 0, n);
165 
166  if (exld) memcpy(fEXlowd, exld, n);
167  else memset(fEXlowd, 0, n);
168  if (exhd) memcpy(fEXhighd, exhd, n);
169  else memset(fEXhighd, 0, n);
170  if (eyld) memcpy(fEYlowd, eyld, n);
171  else memset(fEYlowd, 0, n);
172  if (eyhd) memcpy(fEYhighd, eyhd, n);
173  else memset(fEYhighd, 0, n);
174 }
175 
176 
177 ////////////////////////////////////////////////////////////////////////////////
178 /// TGraphBentErrors default destructor.
179 
180 TGraphBentErrors::~TGraphBentErrors()
181 {
182  delete [] fEXlow;
183  delete [] fEXhigh;
184  delete [] fEYlow;
185  delete [] fEYhigh;
186 
187  delete [] fEXlowd;
188  delete [] fEXhighd;
189  delete [] fEYlowd;
190  delete [] fEYhighd;
191 }
192 
193 
194 ////////////////////////////////////////////////////////////////////////////////
195 /// apply a function to all data points
196 /// y = f(x,y)
197 ///
198 /// Errors are calculated as eyh = f(x,y+eyh)-f(x,y) and
199 /// eyl = f(x,y)-f(x,y-eyl)
200 ///
201 /// Special treatment has to be applied for the functions where the
202 /// role of "up" and "down" is reversed.
203 /// function suggested/implemented by Miroslav Helbich <helbich@mail.desy.de>
204 
205 void TGraphBentErrors::Apply(TF1 *f)
206 {
207  Double_t x,y,exl,exh,eyl,eyh,eyl_new,eyh_new,fxy;
208 
209  if (fHistogram) {
210  delete fHistogram;
211  fHistogram = 0;
212  }
213  for (Int_t i=0;i<GetN();i++) {
214  GetPoint(i,x,y);
215  exl=GetErrorXlow(i);
216  exh=GetErrorXhigh(i);
217  eyl=GetErrorYlow(i);
218  eyh=GetErrorYhigh(i);
219 
220  fxy = f->Eval(x,y);
221  SetPoint(i,x,fxy);
222 
223  // in the case of the functions like y-> -1*y the roles of the
224  // upper and lower error bars is reversed
225  if (f->Eval(x,y-eyl)<f->Eval(x,y+eyh)) {
226  eyl_new = TMath::Abs(fxy - f->Eval(x,y-eyl));
227  eyh_new = TMath::Abs(f->Eval(x,y+eyh) - fxy);
228  }
229  else {
230  eyh_new = TMath::Abs(fxy - f->Eval(x,y-eyl));
231  eyl_new = TMath::Abs(f->Eval(x,y+eyh) - fxy);
232  }
233 
234  //error on x doesn't change
235  SetPointError(i,exl,exh,eyl_new,eyh_new);
236  }
237  if (gPad) gPad->Modified();
238 }
239 
240 
241 ////////////////////////////////////////////////////////////////////////////////
242 /// Compute range.
243 
244 void TGraphBentErrors::ComputeRange(Double_t &xmin, Double_t &ymin, Double_t &xmax, Double_t &ymax) const
245 {
246  TGraph::ComputeRange(xmin,ymin,xmax,ymax);
247 
248  for (Int_t i=0;i<fNpoints;i++) {
249  if (fX[i] -fEXlow[i] < xmin) {
250  if (gPad && gPad->GetLogx()) {
251  if (fEXlow[i] < fX[i]) xmin = fX[i]-fEXlow[i];
252  else xmin = TMath::Min(xmin,fX[i]/3);
253  } else {
254  xmin = fX[i]-fEXlow[i];
255  }
256  }
257  if (fX[i] +fEXhigh[i] > xmax) xmax = fX[i]+fEXhigh[i];
258  if (fY[i] -fEYlow[i] < ymin) {
259  if (gPad && gPad->GetLogy()) {
260  if (fEYlow[i] < fY[i]) ymin = fY[i]-fEYlow[i];
261  else ymin = TMath::Min(ymin,fY[i]/3);
262  } else {
263  ymin = fY[i]-fEYlow[i];
264  }
265  }
266  if (fY[i] +fEYhigh[i] > ymax) ymax = fY[i]+fEYhigh[i];
267  }
268 }
269 
270 
271 ////////////////////////////////////////////////////////////////////////////////
272 /// Copy and release.
273 
274 void TGraphBentErrors::CopyAndRelease(Double_t **newarrays,
275  Int_t ibegin, Int_t iend, Int_t obegin)
276 {
277  CopyPoints(newarrays, ibegin, iend, obegin);
278  if (newarrays) {
279  delete[] fEXlow;
280  fEXlow = newarrays[0];
281  delete[] fEXhigh;
282  fEXhigh = newarrays[1];
283  delete[] fEYlow;
284  fEYlow = newarrays[2];
285  delete[] fEYhigh;
286  fEYhigh = newarrays[3];
287  delete[] fEXlowd;
288  fEXlowd = newarrays[4];
289  delete[] fEXhighd;
290  fEXhighd = newarrays[5];
291  delete[] fEYlowd;
292  fEYlowd = newarrays[6];
293  delete[] fEYhighd;
294  fEYhighd = newarrays[7];
295  delete[] fX;
296  fX = newarrays[8];
297  delete[] fY;
298  fY = newarrays[9];
299  delete[] newarrays;
300  }
301 }
302 
303 
304 ////////////////////////////////////////////////////////////////////////////////
305 /// Copy errors from fE*** to arrays[***]
306 /// or to f*** Copy points.
307 
308 Bool_t TGraphBentErrors::CopyPoints(Double_t **arrays,
309  Int_t ibegin, Int_t iend, Int_t obegin)
310 {
311  if (TGraph::CopyPoints(arrays ? arrays+8 : 0, ibegin, iend, obegin)) {
312  Int_t n = (iend - ibegin)*sizeof(Double_t);
313  if (arrays) {
314  memmove(&arrays[0][obegin], &fEXlow[ibegin], n);
315  memmove(&arrays[1][obegin], &fEXhigh[ibegin], n);
316  memmove(&arrays[2][obegin], &fEYlow[ibegin], n);
317  memmove(&arrays[3][obegin], &fEYhigh[ibegin], n);
318  memmove(&arrays[4][obegin], &fEXlowd[ibegin], n);
319  memmove(&arrays[5][obegin], &fEXhighd[ibegin], n);
320  memmove(&arrays[6][obegin], &fEYlowd[ibegin], n);
321  memmove(&arrays[7][obegin], &fEYhighd[ibegin], n);
322  } else {
323  memmove(&fEXlow[obegin], &fEXlow[ibegin], n);
324  memmove(&fEXhigh[obegin], &fEXhigh[ibegin], n);
325  memmove(&fEYlow[obegin], &fEYlow[ibegin], n);
326  memmove(&fEYhigh[obegin], &fEYhigh[ibegin], n);
327  memmove(&fEXlowd[obegin], &fEXlowd[ibegin], n);
328  memmove(&fEXhighd[obegin], &fEXhighd[ibegin], n);
329  memmove(&fEYlowd[obegin], &fEYlowd[ibegin], n);
330  memmove(&fEYhighd[obegin], &fEYhighd[ibegin], n);
331  }
332  return kTRUE;
333  } else {
334  return kFALSE;
335  }
336 }
337 
338 
339 ////////////////////////////////////////////////////////////////////////////////
340 /// Should be called from ctors after fNpoints has been set
341 
342 Bool_t TGraphBentErrors::CtorAllocate(void)
343 {
344  if (!fNpoints) {
345  fEXlow = fEYlow = fEXhigh = fEYhigh = 0;
346  fEXlowd = fEYlowd = fEXhighd = fEYhighd = 0;
347  return kFALSE;
348  }
349  fEXlow = new Double_t[fMaxSize];
350  fEYlow = new Double_t[fMaxSize];
351  fEXhigh = new Double_t[fMaxSize];
352  fEYhigh = new Double_t[fMaxSize];
353  fEXlowd = new Double_t[fMaxSize];
354  fEYlowd = new Double_t[fMaxSize];
355  fEXhighd = new Double_t[fMaxSize];
356  fEYhighd = new Double_t[fMaxSize];
357  return kTRUE;
358 }
359 
360 ////////////////////////////////////////////////////////////////////////////////
361 /// protected function to perform the merge operation of a graph with asymmetric errors
362 
363 Bool_t TGraphBentErrors::DoMerge(const TGraph *g)
364 {
365  if (g->GetN() == 0) return kFALSE;
366 
367  Double_t * exl = g->GetEXlow();
368  Double_t * exh = g->GetEXhigh();
369  Double_t * eyl = g->GetEYlow();
370  Double_t * eyh = g->GetEYhigh();
371 
372  Double_t * exld = g->GetEXlowd();
373  Double_t * exhd = g->GetEXhighd();
374  Double_t * eyld = g->GetEYlowd();
375  Double_t * eyhd = g->GetEYhighd();
376 
377  if (exl == 0 || exh == 0 || eyl == 0 || eyh == 0 ||
378  exld == 0 || exhd == 0 || eyld == 0 || eyhd == 0) {
379  if (g->IsA() != TGraph::Class() )
380  Warning("DoMerge","Merging a %s is not compatible with a TGraphBentErrors - errors will be ignored",g->IsA()->GetName());
381  return TGraph::DoMerge(g);
382  }
383  for (Int_t i = 0 ; i < g->GetN(); i++) {
384  Int_t ipoint = GetN();
385  Double_t x = g->GetX()[i];
386  Double_t y = g->GetY()[i];
387  SetPoint(ipoint, x, y);
388  SetPointError(ipoint, exl[i], exh[i], eyl[i], eyh[i],
389  exld[i], exhd[i], eyld[i], eyhd[i] );
390  }
391 
392  return kTRUE;
393 
394 }
395 ////////////////////////////////////////////////////////////////////////////////
396 /// This function is called by GraphFitChisquare.
397 /// It returns the error along X at point i.
398 
399 Double_t TGraphBentErrors::GetErrorX(Int_t i) const
400 {
401  if (i < 0 || i >= fNpoints) return -1;
402  if (!fEXlow && !fEXhigh) return -1;
403  Double_t elow=0, ehigh=0;
404  if (fEXlow) elow = fEXlow[i];
405  if (fEXhigh) ehigh = fEXhigh[i];
406  return TMath::Sqrt(0.5*(elow*elow + ehigh*ehigh));
407 }
408 
409 
410 ////////////////////////////////////////////////////////////////////////////////
411 /// This function is called by GraphFitChisquare.
412 /// It returns the error along Y at point i.
413 
414 Double_t TGraphBentErrors::GetErrorY(Int_t i) const
415 {
416  if (i < 0 || i >= fNpoints) return -1;
417  if (!fEYlow && !fEYhigh) return -1;
418  Double_t elow=0, ehigh=0;
419  if (fEYlow) elow = fEYlow[i];
420  if (fEYhigh) ehigh = fEYhigh[i];
421  return TMath::Sqrt(0.5*(elow*elow + ehigh*ehigh));
422 }
423 
424 
425 ////////////////////////////////////////////////////////////////////////////////
426 /// Get high error on X[i].
427 
428 Double_t TGraphBentErrors::GetErrorXhigh(Int_t i) const
429 {
430  if (i<0 || i>fNpoints) return -1;
431  if (fEXhigh) return fEXhigh[i];
432  return -1;
433 }
434 
435 
436 ////////////////////////////////////////////////////////////////////////////////
437 /// Get low error on X[i].
438 
439 Double_t TGraphBentErrors::GetErrorXlow(Int_t i) const
440 {
441  if (i<0 || i>fNpoints) return -1;
442  if (fEXlow) return fEXlow[i];
443  return -1;
444 }
445 
446 
447 ////////////////////////////////////////////////////////////////////////////////
448 /// Get high error on Y[i].
449 
450 Double_t TGraphBentErrors::GetErrorYhigh(Int_t i) const
451 {
452  if (i<0 || i>fNpoints) return -1;
453  if (fEYhigh) return fEYhigh[i];
454  return -1;
455 }
456 
457 
458 ////////////////////////////////////////////////////////////////////////////////
459 /// Get low error on Y[i].
460 
461 Double_t TGraphBentErrors::GetErrorYlow(Int_t i) const
462 {
463  if (i<0 || i>fNpoints) return -1;
464  if (fEYlow) return fEYlow[i];
465  return -1;
466 }
467 
468 
469 ////////////////////////////////////////////////////////////////////////////////
470 /// Set zero values for point arrays in the range [begin, end)
471 
472 void TGraphBentErrors::FillZero(Int_t begin, Int_t end,
473  Bool_t from_ctor)
474 {
475  if (!from_ctor) {
476  TGraph::FillZero(begin, end, from_ctor);
477  }
478  Int_t n = (end - begin)*sizeof(Double_t);
479  memset(fEXlow + begin, 0, n);
480  memset(fEXhigh + begin, 0, n);
481  memset(fEYlow + begin, 0, n);
482  memset(fEYhigh + begin, 0, n);
483  memset(fEXlowd + begin, 0, n);
484  memset(fEXhighd + begin, 0, n);
485  memset(fEYlowd + begin, 0, n);
486  memset(fEYhighd + begin, 0, n);
487 }
488 
489 
490 ////////////////////////////////////////////////////////////////////////////////
491 /// Print graph and errors values.
492 
493 void TGraphBentErrors::Print(Option_t *) const
494 {
495  for (Int_t i=0;i<fNpoints;i++) {
496  printf("x[%d]=%g, y[%d]=%g, exl[%d]=%g, exh[%d]=%g, eyl[%d]=%g, eyh[%d]=%g\n"
497  ,i,fX[i],i,fY[i],i,fEXlow[i],i,fEXhigh[i],i,fEYlow[i],i,fEYhigh[i]);
498  }
499 }
500 
501 
502 ////////////////////////////////////////////////////////////////////////////////
503 /// Save primitive as a C++ statement(s) on output stream out
504 
505 void TGraphBentErrors::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/)
506 {
507  char quote = '"';
508  out << " " << std::endl;
509  static Int_t frameNumber = 2000;
510  frameNumber++;
511 
512  Int_t i;
513  TString fXName = TString(GetName()) + Form("_fx%d",frameNumber);
514  TString fYName = TString(GetName()) + Form("_fy%d",frameNumber);
515  TString fElXName = TString(GetName()) + Form("_felx%d",frameNumber);
516  TString fElYName = TString(GetName()) + Form("_fely%d",frameNumber);
517  TString fEhXName = TString(GetName()) + Form("_fehx%d",frameNumber);
518  TString fEhYName = TString(GetName()) + Form("_fehy%d",frameNumber);
519  TString fEldXName = TString(GetName()) + Form("_feldx%d",frameNumber);
520  TString fEldYName = TString(GetName()) + Form("_feldy%d",frameNumber);
521  TString fEhdXName = TString(GetName()) + Form("_fehdx%d",frameNumber);
522  TString fEhdYName = TString(GetName()) + Form("_fehdy%d",frameNumber);
523  out << " Double_t " << fXName << "[" << fNpoints << "] = {" << std::endl;
524  for (i = 0; i < fNpoints-1; i++) out << " " << fX[i] << "," << std::endl;
525  out << " " << fX[fNpoints-1] << "};" << std::endl;
526  out << " Double_t " << fYName << "[" << fNpoints << "] = {" << std::endl;
527  for (i = 0; i < fNpoints-1; i++) out << " " << fY[i] << "," << std::endl;
528  out << " " << fY[fNpoints-1] << "};" << std::endl;
529  out << " Double_t " << fElXName << "[" << fNpoints << "] = {" << std::endl;
530  for (i = 0; i < fNpoints-1; i++) out << " " << fEXlow[i] << "," << std::endl;
531  out << " " << fEXlow[fNpoints-1] << "};" << std::endl;
532  out << " Double_t " << fElYName << "[" << fNpoints << "] = {" << std::endl;
533  for (i = 0; i < fNpoints-1; i++) out << " " << fEYlow[i] << "," << std::endl;
534  out << " " << fEYlow[fNpoints-1] << "};" << std::endl;
535  out << " Double_t " << fEhXName << "[" << fNpoints << "] = {" << std::endl;
536  for (i = 0; i < fNpoints-1; i++) out << " " << fEXhigh[i] << "," << std::endl;
537  out << " " << fEXhigh[fNpoints-1] << "};" << std::endl;
538  out << " Double_t " << fEhYName << "[" << fNpoints << "] = {" << std::endl;
539  for (i = 0; i < fNpoints-1; i++) out << " " << fEYhigh[i] << "," << std::endl;
540  out << " " << fEYhigh[fNpoints-1] << "};" << std::endl;
541  out << " Double_t " << fEldXName << "[" << fNpoints << "] = {" << std::endl;
542  for (i = 0; i < fNpoints-1; i++) out << " " << fEXlowd[i] << "," << std::endl;
543  out << " " << fEXlowd[fNpoints-1] << "};" << std::endl;
544  out << " Double_t " << fEldYName << "[" << fNpoints << "] = {" << std::endl;
545  for (i = 0; i < fNpoints-1; i++) out << " " << fEYlowd[i] << "," << std::endl;
546  out << " " << fEYlowd[fNpoints-1] << "};" << std::endl;
547  out << " Double_t " << fEhdXName << "[" << fNpoints << "] = {" << std::endl;
548  for (i = 0; i < fNpoints-1; i++) out << " " << fEXhighd[i] << "," << std::endl;
549  out << " " << fEXhighd[fNpoints-1] << "};" << std::endl;
550  out << " Double_t " << fEhdYName << "[" << fNpoints << "] = {" << std::endl;
551  for (i = 0; i < fNpoints-1; i++) out << " " << fEYhighd[i] << "," << std::endl;
552  out << " " << fEYhighd[fNpoints-1] << "};" << std::endl;
553 
554  if (gROOT->ClassSaved(TGraphBentErrors::Class())) out << " ";
555  else out << " TGraphBentErrors *";
556  out << "grbe = new TGraphBentErrors("<< fNpoints << ","
557  << fXName << "," << fYName << ","
558  << fElXName << "," << fEhXName << ","
559  << fElYName << "," << fEhYName << ","
560  << fEldXName << "," << fEhdXName << ","
561  << fEldYName << "," << fEhdYName << ");"
562  << std::endl;
563 
564  out << " grbe->SetName(" << quote << GetName() << quote << ");" << std::endl;
565  out << " grbe->SetTitle(" << quote << GetTitle() << quote << ");" << std::endl;
566 
567  SaveFillAttributes(out,"grbe",0,1001);
568  SaveLineAttributes(out,"grbe",1,1,1);
569  SaveMarkerAttributes(out,"grbe",1,1,1);
570 
571  if (fHistogram) {
572  TString hname = fHistogram->GetName();
573  hname += frameNumber;
574  fHistogram->SetName(Form("Graph_%s",hname.Data()));
575  fHistogram->SavePrimitive(out,"nodraw");
576  out<<" grbe->SetHistogram("<<fHistogram->GetName()<<");"<<std::endl;
577  out<<" "<<std::endl;
578  }
579 
580  // save list of functions
581  TIter next(fFunctions);
582  TObject *obj;
583  while ((obj = next())) {
584  obj->SavePrimitive(out, Form("nodraw #%d\n",++frameNumber));
585  if (obj->InheritsFrom("TPaveStats")) {
586  out << " grbe->GetListOfFunctions()->Add(ptstats);" << std::endl;
587  out << " ptstats->SetParent(grbe->GetListOfFunctions());" << std::endl;
588  } else {
589  TString objname;
590  objname.Form("%s%d",obj->GetName(),frameNumber);
591  if (obj->InheritsFrom("TF1")) {
592  out << " " << objname << "->SetParent(grbe);\n";
593  }
594  out << " grbe->GetListOfFunctions()->Add("
595  << objname << ");" << std::endl;
596  }
597  }
598 
599  const char *l = strstr(option,"multigraph");
600  if (l) {
601  out<<" multigraph->Add(grbe,"<<quote<<l+10<<quote<<");"<<std::endl;
602  } else {
603  out<<" grbe->Draw("<<quote<<option<<quote<<");"<<std::endl;
604  }
605 }
606 
607 
608 ////////////////////////////////////////////////////////////////////////////////
609 /// Set ex and ey values for point pointed by the mouse.
610 
611 void TGraphBentErrors::SetPointError(Double_t exl, Double_t exh, Double_t eyl, Double_t eyh,
612  Double_t exld, Double_t exhd, Double_t eyld, Double_t eyhd)
613 {
614  Int_t px = gPad->GetEventX();
615  Int_t py = gPad->GetEventY();
616 
617  //localize point to be deleted
618  Int_t ipoint = -2;
619  Int_t i;
620  // start with a small window (in case the mouse is very close to one point)
621  for (i=0;i<fNpoints;i++) {
622  Int_t dpx = px - gPad->XtoAbsPixel(gPad->XtoPad(fX[i]));
623  Int_t dpy = py - gPad->YtoAbsPixel(gPad->YtoPad(fY[i]));
624  if (dpx*dpx+dpy*dpy < 25) {ipoint = i; break;}
625  }
626  if (ipoint == -2) return;
627 
628  fEXlow[ipoint] = exl;
629  fEYlow[ipoint] = eyl;
630  fEXhigh[ipoint] = exh;
631  fEYhigh[ipoint] = eyh;
632  fEXlowd[ipoint] = exld;
633  fEXhighd[ipoint] = exhd;
634  fEYlowd[ipoint] = eyld;
635  fEYhighd[ipoint] = eyhd;
636  gPad->Modified();
637 }
638 
639 
640 ////////////////////////////////////////////////////////////////////////////////
641 /// Set ex and ey values for point number i.
642 
643 void TGraphBentErrors::SetPointError(Int_t i, Double_t exl, Double_t exh, Double_t eyl, Double_t eyh,
644  Double_t exld, Double_t exhd, Double_t eyld, Double_t eyhd)
645 {
646  if (i < 0) return;
647  if (i >= fNpoints) {
648  // re-allocate the object
649  TGraphBentErrors::SetPoint(i,0,0);
650  }
651  fEXlow[i] = exl;
652  fEYlow[i] = eyl;
653  fEXhigh[i] = exh;
654  fEYhigh[i] = eyh;
655  fEXlowd[i] = exld;
656  fEXhighd[i] = exhd;
657  fEYlowd[i] = eyld;
658  fEYhighd[i] = eyhd;
659 }
660 
661 
662 ////////////////////////////////////////////////////////////////////////////////
663 /// Swap points.
664 
665 void TGraphBentErrors::SwapPoints(Int_t pos1, Int_t pos2)
666 {
667  SwapValues(fEXlow, pos1, pos2);
668  SwapValues(fEXhigh, pos1, pos2);
669  SwapValues(fEYlow, pos1, pos2);
670  SwapValues(fEYhigh, pos1, pos2);
671 
672  SwapValues(fEXlowd, pos1, pos2);
673  SwapValues(fEXhighd, pos1, pos2);
674  SwapValues(fEYlowd, pos1, pos2);
675  SwapValues(fEYhighd, pos1, pos2);
676 
677  TGraph::SwapPoints(pos1, pos2);
678 }