Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TF3.cxx
Go to the documentation of this file.
1// @(#)root/hist:$Id$
2// Author: Rene Brun 27/10/95
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, 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 "TROOT.h"
13#include "TF3.h"
14#include "TBuffer.h"
15#include "TMath.h"
16#include "TH3.h"
17#include "TVirtualPad.h"
18#include "TRandom.h"
19#include "TVectorD.h"
20#include "Riostream.h"
21#include "TColor.h"
22#include "TVirtualFitter.h"
23#include "TVirtualHistPainter.h"
25#include <cassert>
26
28
29/** \class TF3
30 \ingroup Functions
31A 3-Dim function with parameters
32*/
33
34////////////////////////////////////////////////////////////////////////////////
35/// F3 default constructor
36
38{
39 fNpz = 0;
40 fZmin = 0;
41 fZmax = 1;
42}
43
44
45////////////////////////////////////////////////////////////////////////////////
46/// F3 constructor using a formula definition
47///
48/// See TFormula constructor for explanation of the formula syntax.
49
50TF3::TF3(const char *name,const char *formula, Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax, Option_t * opt)
51 :TF2(name,formula,xmin,xmax,ymax,ymin,opt)
52{
53 fZmin = zmin;
54 fZmax = zmax;
55 fNpz = 30;
56 Int_t ndim = GetNdim();
57 // accept 1-d or 2-d formula
58 if (ndim < 3) fNdim = 3;
59 if (ndim > 3 && xmin < xmax && ymin < ymax && zmin < zmax) {
60 Error("TF3","function: %s/%s has dimension %d instead of 3",name,formula,ndim);
61 MakeZombie();
62 }
63}
64
65////////////////////////////////////////////////////////////////////////////////
66/// F3 constructor using a pointer to real function
67///
68/// \param[in] name object name
69/// \param[in] fcn pointer to real function
70/// \param[in] xmin,xmax x axis limits
71/// \param[in] ymin,ymax y axis limits
72/// \param[in] zmin,zmax z axis limits
73/// \param[in] npar is the number of free parameters used by the function
74/// \param[in] ndim number of dimensions
75///
76/// For example, for a 3-dim function with 3 parameters, the user function
77/// looks like:
78///
79/// Double_t fun1(Double_t *x, Double_t *par)
80/// return par[0]*x[2] + par[1]*exp(par[2]*x[0]*x[1]);
81///
82/// \warning A function created with this constructor cannot be Cloned.
83
85 :TF2(name,fcn,xmin,xmax,ymin,ymax,npar,ndim)
86{
87 fZmin = zmin;
88 fZmax = zmax;
89 fNpz = 30;
90}
91
92////////////////////////////////////////////////////////////////////////////////
93/// F3 constructor using a pointer to real function---
94///
95/// \param[in] name object name
96/// \param[in] fcn pointer to real function
97/// \param[in] xmin,xmax x axis limits
98/// \param[in] ymin,ymax y axis limits
99/// \param[in] zmin,zmax z axis limits
100/// \param[in] npar is the number of free parameters used by the function
101/// \param[in] ndim number of dimensions
102///
103/// For example, for a 3-dim function with 3 parameters, the user function
104/// looks like:
105///
106/// Double_t fun1(Double_t *x, Double_t *par)
107/// return par[0]*x[2] + par[1]*exp(par[2]*x[0]*x[1]);
108///
109/// WARNING! A function created with this constructor cannot be Cloned.
110
111TF3::TF3(const char *name,Double_t (*fcn)(const Double_t *, const Double_t *), Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax, Int_t npar, Int_t ndim)
112 : TF2(name,fcn,xmin,xmax,ymin,ymax,npar,ndim),
113 fZmin(zmin),
114 fZmax(zmax),
115 fNpz(30)
116{
117}
118
119////////////////////////////////////////////////////////////////////////////////
120/// F3 constructor using a ParamFunctor
121///
122/// a functor class implementing operator() (double *, double *)
123///
124/// \param[in] name object name
125/// \param[in] f parameter functor
126/// \param[in] xmin,xmax x axis limits
127/// \param[in] ymin,ymax y axis limits
128/// \param[in] zmin,zmax z axis limits
129/// \param[in] npar is the number of free parameters used by the function
130/// \param[in] ndim number of dimensions
131///
132/// \warning A function created with this constructor cannot be Cloned.
133
135 : TF2(name, f, xmin, xmax, ymin, ymax, npar, ndim),
136 fZmin(zmin),
137 fZmax(zmax),
138 fNpz(30)
139{
140}
141
142////////////////////////////////////////////////////////////////////////////////
143/// Operator =
144
146{
147 if (this != &rhs)
148 rhs.TF3::Copy(*this);
149 return *this;
150}
151
152////////////////////////////////////////////////////////////////////////////////
153/// F3 default destructor
154
156{
157}
158
159////////////////////////////////////////////////////////////////////////////////
160/// Copy constructor.
161
162TF3::TF3(const TF3 &f3) : TF2()
163{
164 f3.TF3::Copy(*this);
165}
166
167////////////////////////////////////////////////////////////////////////////////
168/// Copy this F3 to a new F3
169
170void TF3::Copy(TObject &obj) const
171{
172 TF2::Copy(obj);
173 ((TF3&)obj).fZmin = fZmin;
174 ((TF3&)obj).fZmax = fZmax;
175 ((TF3&)obj).fNpz = fNpz;
176}
177
178////////////////////////////////////////////////////////////////////////////////
179/// Compute distance from point px,py to a function
180///
181/// Compute the closest distance of approach from point px,py to this function.
182/// The distance is computed in pixels units.
183
184
186{
187 return TF1::DistancetoPrimitive(px, py);
188
189}
190
191////////////////////////////////////////////////////////////////////////////////
192/// Draw this function with its current attributes
193
195{
196 TString opt = option;
197 opt.ToLower();
198 if (gPad && !opt.Contains("same")) gPad->Clear();
199
201
202}
203
204////////////////////////////////////////////////////////////////////////////////
205/// Execute action corresponding to one event
206///
207/// This member function is called when a F3 is clicked with the locator
208
210{
211 TF1::ExecuteEvent(event, px, py);
212}
213
214////////////////////////////////////////////////////////////////////////////////
215/// Return minimum/maximum value of the function
216///
217/// To find the minimum on a range, first set this range via the SetRange function
218/// If a vector x of coordinate is passed it will be used as starting point for the minimum.
219/// In addition on exit x will contain the coordinate values at the minimuma
220/// If x is NULL or x is inifinity or NaN, first, a grid search is performed to find the initial estimate of the
221/// minimum location. The range of the function is divided into fNpx and fNpy
222/// sub-ranges. If the function is "good" (or "bad"), these values can be changed
223/// by SetNpx and SetNpy functions
224///
225/// Then, a minimization is used with starting values found by the grid search
226/// The minimizer algorithm used (by default Minuit) can be changed by callinga
227/// ROOT::Math::Minimizer::SetDefaultMinimizerType("..")
228/// Other option for the minimizer can be set using the static method of the MinimizerOptions class
229
231{
232 //First do a grid search with step size fNpx and fNpy
233
234 Double_t xx[3];
235 Double_t rsign = (findmax) ? -1. : 1.;
236 TF3 & function = const_cast<TF3&>(*this); // needed since EvalPar is not const
237 Double_t xxmin = 0, yymin = 0, zzmin = 0, ttmin = 0;
238 if (x == nullptr || ( (x!= nullptr) && ( !TMath::Finite(x[0]) || !TMath::Finite(x[1]) || !TMath::Finite(x[2]) ) ) ){
239 Double_t dx = (fXmax - fXmin)/fNpx;
240 Double_t dy = (fYmax - fYmin)/fNpy;
241 Double_t dz = (fZmax - fZmin)/fNpz;
242 xxmin = fXmin;
243 yymin = fYmin;
244 zzmin = fZmin;
245 ttmin = rsign * TMath::Infinity();
246 for (Int_t i=0; i<fNpx; i++){
247 xx[0]=fXmin + (i+0.5)*dx;
248 for (Int_t j=0; j<fNpy; j++){
249 xx[1]=fYmin+(j+0.5)*dy;
250 for (Int_t k=0; k<fNpz; k++){
251 xx[2] = fZmin+(k+0.5)*dz;
252 Double_t tt = function(xx);
253 if (rsign*tt < rsign*ttmin) {xxmin = xx[0], yymin = xx[1]; zzmin = xx[2]; ttmin=tt;}
254 }
255 }
256 }
257
258 xxmin = TMath::Min(fXmax, xxmin);
259 yymin = TMath::Min(fYmax, yymin);
260 zzmin = TMath::Min(fZmax, zzmin);
261 }
262 else {
263 xxmin = x[0];
264 yymin = x[1];
265 zzmin = x[2];
266 zzmin = function(x);
267 }
268 xx[0] = xxmin;
269 xx[1] = yymin;
270 xx[2] = zzmin;
271
272 double fmin = GetMinMaxNDim(xx,findmax);
273 if (rsign*fmin < rsign*zzmin) {
274 if (x) {x[0] = xx[0]; x[1] = xx[1]; x[2] = xx[2];}
275 return fmin;
276 }
277 // here if minimization failed
278 if (x) { x[0] = xxmin; x[1] = yymin; x[2] = zzmin; }
279 return ttmin;
280}
281
282////////////////////////////////////////////////////////////////////////////////
283/// Compute the X, Y and Z values corresponding to the minimum value of the function
284/// on its range.
285///
286/// Returns the function value at the minimum.
287/// To find the minimum on a subrange, use the SetRange() function first.
288///
289/// Method:
290/// First, a grid search is performed to find the initial estimate of the
291/// minimum location. The range of the function is divided
292/// into fNpx,fNpy and fNpz sub-ranges. If the function is "good" (or "bad"),
293/// these values can be changed by SetNpx(), SetNpy() and SetNpz() functions.
294/// Then, Minuit minimization is used with starting values found by the grid search
295///
296/// Note that this method will always do first a grid search in contrast to GetMinimum
297
299{
300 double xx[3] = { 0,0,0 };
301 xx[0] = TMath::QuietNaN(); // to force to do grid search in TF3::FindMinMax
302 double fmin = FindMinMax(xx, false);
303 x = xx[0]; y = xx[1]; z = xx[2];
304 return fmin;
305
306}
307
308////////////////////////////////////////////////////////////////////////////////
309/// Compute the X, Y and Z values corresponding to the maximum value of the function
310/// on its range.
311///
312/// Return the function value at the maximum. See TF3::GetMinimumXYZ
313
315{
316 double xx[3] = { 0,0,0 };
317 xx[0] = TMath::QuietNaN(); // to force to do grid search in TF3::FindMinMax
318 double fmax = FindMinMax(xx, true);
319 x = xx[0]; y = xx[1]; z = xx[2];
320 return fmax;
321
322}
323
324////////////////////////////////////////////////////////////////////////////////
325/// Return 3 random numbers following this function shape
326///
327/// The distribution contained in this TF3 function is integrated
328/// over the cell contents.
329/// It is normalized to 1.
330/// Getting the three random numbers implies:
331/// - Generating a random number between 0 and 1 (say r1)
332/// - Look in which cell in the normalized integral r1 corresponds to
333/// - make a linear interpolation in the returned cell
334///
335/// IMPORTANT NOTE
336///
337/// The integral of the function is computed at fNpx * fNpy * fNpz points.
338/// If the function has sharp peaks, you should increase the number of
339/// points (SetNpx, SetNpy, SetNpz) such that the peak is correctly tabulated
340/// at several points.
341
342void TF3::GetRandom3(Double_t &xrandom, Double_t &yrandom, Double_t &zrandom, TRandom * rng)
343{
344 // Check if integral array must be built
345 Int_t i,j,k,cell;
346 Double_t dx = (fXmax-fXmin)/fNpx;
347 Double_t dy = (fYmax-fYmin)/fNpy;
348 Double_t dz = (fZmax-fZmin)/fNpz;
349 Int_t ncells = fNpx*fNpy*fNpz;
350 Double_t xx[3];
351 Double_t *parameters = GetParameters();
352 InitArgs(xx,parameters);
353 if (fIntegral.empty() ) {
354 fIntegral.resize(ncells+1);
355 //fIntegral = new Double_t[ncells+1];
356 fIntegral[0] = 0;
357 Double_t integ;
358 Int_t intNegative = 0;
359 cell = 0;
360 for (k=0;k<fNpz;k++) {
361 xx[2] = fZmin+(k+0.5)*dz;
362 for (j=0;j<fNpy;j++) {
363 xx[1] = fYmin+(j+0.5)*dy;
364 for (i=0;i<fNpx;i++) {
365 xx[0] = fXmin+(i+0.5)*dx;
366 integ = EvalPar(xx,parameters);
367 if (integ < 0) {intNegative++; integ = -integ;}
368 fIntegral[cell+1] = fIntegral[cell] + integ;
369 cell++;
370 }
371 }
372 }
373 if (intNegative > 0) {
374 Warning("GetRandom3","function:%s has %d negative values: abs assumed",GetName(),intNegative);
375 }
376 if (fIntegral[ncells] == 0) {
377 Error("GetRandom3","Integral of function is zero");
378 return;
379 }
380 for (i=1;i<=ncells;i++) { // normalize integral to 1
381 fIntegral[i] /= fIntegral[ncells];
382 }
383 }
384
385// return random numbers
386 Double_t r;
387 if (!rng) rng = gRandom;
388 r = rng->Rndm();
389 cell = TMath::BinarySearch(ncells,fIntegral.data(),r);
390 k = cell/(fNpx*fNpy);
391 j = (cell -k*fNpx*fNpy)/fNpx;
392 i = cell -fNpx*(j +fNpy*k);
393 xrandom = fXmin +dx*i +dx*rng->Rndm();
394 yrandom = fYmin +dy*j +dy*rng->Rndm();
395 zrandom = fZmin +dz*k +dz*rng->Rndm();
396}
397
398////////////////////////////////////////////////////////////////////////////////
399/// Return range of function
400
402{
403 xmin = fXmin;
404 xmax = fXmax;
405 ymin = fYmin;
406 ymax = fYmax;
407 zmin = fZmin;
408 zmax = fZmax;
409}
410
411
412////////////////////////////////////////////////////////////////////////////////
413/// Get value corresponding to X in array of fSave values
414
416{
417 if (fSave.size() < 9) return 0;
418 Int_t nsave = fSave.size() - 9;
419 Double_t xmin = fSave[nsave+0];
420 Double_t xmax = fSave[nsave+1];
421 Double_t ymin = fSave[nsave+2];
422 Double_t ymax = fSave[nsave+3];
423 Double_t zmin = fSave[nsave+4];
424 Double_t zmax = fSave[nsave+5];
425 Int_t npx = Int_t(fSave[nsave+6]);
426 Int_t npy = Int_t(fSave[nsave+7]);
427 Int_t npz = Int_t(fSave[nsave+8]);
428 Double_t x = xx[0];
429 Double_t dx = (xmax-xmin)/npx;
430 if (x < xmin || x > xmax) return 0;
431 if (dx <= 0) return 0;
432 Double_t y = xx[1];
433 Double_t dy = (ymax-ymin)/npy;
434 if (y < ymin || y > ymax) return 0;
435 if (dy <= 0) return 0;
436 Double_t z = xx[2];
437 Double_t dz = (zmax-zmin)/npz;
438 if (z < zmin || z > zmax) return 0;
439 if (dz <= 0) return 0;
440
441 //we make a trilinear interpolation using the 8 points surrounding x,y,z
442 Int_t ibin = TMath::Min(npx-1, Int_t((x-xmin)/dx));
443 Int_t jbin = TMath::Min(npy-1, Int_t((y-ymin)/dy));
444 Int_t kbin = TMath::Min(npz-1, Int_t((z-zmin)/dz));
445 Double_t xlow = xmin + ibin*dx;
446 Double_t ylow = ymin + jbin*dy;
447 Double_t zlow = zmin + kbin*dz;
448 Double_t t = (x-xlow)/dx;
449 Double_t u = (y-ylow)/dy;
450 Double_t v = (z-zlow)/dz;
451 Int_t k1 = (ibin ) + (npx+1)*((jbin ) + (npy+1)*(kbin ));
452 Int_t k2 = (ibin+1) + (npx+1)*((jbin ) + (npy+1)*(kbin ));
453 Int_t k3 = (ibin+1) + (npx+1)*((jbin+1) + (npy+1)*(kbin ));
454 Int_t k4 = (ibin ) + (npx+1)*((jbin+1) + (npy+1)*(kbin ));
455 Int_t k5 = (ibin ) + (npx+1)*((jbin ) + (npy+1)*(kbin+1));
456 Int_t k6 = (ibin+1) + (npx+1)*((jbin ) + (npy+1)*(kbin+1));
457 Int_t k7 = (ibin+1) + (npx+1)*((jbin+1) + (npy+1)*(kbin+1));
458 Int_t k8 = (ibin ) + (npx+1)*((jbin+1) + (npy+1)*(kbin+1));
459 Double_t r = (1-t)*(1-u)*(1-v)*fSave[k1] + t*(1-u)*(1-v)*fSave[k2] + t*u*(1-v)*fSave[k3] + (1-t)*u*(1-v)*fSave[k4] +
460 (1-t)*(1-u)*v*fSave[k5] + t*(1-u)*v*fSave[k6] + t*u*v*fSave[k7] + (1-t)*u*v*fSave[k8];
461 return r;
462}
463
464////////////////////////////////////////////////////////////////////////////////
465/// Return Integral of a 3d function in range [ax,bx],[ay,by],[az,bz]
466/// with a desired relative accuracy.
467
469{
470 Double_t a[3], b[3];
471 a[0] = ax;
472 b[0] = bx;
473 a[1] = ay;
474 b[1] = by;
475 a[2] = az;
476 b[2] = bz;
477 Double_t relerr = 0;
478 Int_t n = 3;
480 Int_t nfnevl, ifail;
481 Double_t result = IntegralMultiple(n,a,b,maxpts,epsrel,epsrel, relerr,nfnevl,ifail);
482 if (ifail > 0) {
483 Warning("Integral","failed for %s code=%d, maxpts=%d, epsrel=%g, nfnevl=%d, relerr=%g ",GetName(),ifail,maxpts,epsrel,nfnevl,relerr);
484 }
485 if (gDebug) {
486 Info("Integral","Integral of %s using %d and tol=%f is %f , relerr=%f nfcn=%d",GetName(),maxpts,epsrel,result,relerr,nfnevl);
487 }
488 return result;
489}
490
491////////////////////////////////////////////////////////////////////////////////
492/// Return kTRUE is the point is inside the function range
493
495{
496 if (x[0] < fXmin || x[0] > fXmax) return kFALSE;
497 if (x[1] < fYmin || x[1] > fYmax) return kFALSE;
498 if (x[2] < fZmin || x[2] > fZmax) return kFALSE;
499 return kTRUE;
500}
501
502////////////////////////////////////////////////////////////////////////////////
503/// Create a histogram for axis range.
504
506{
507 TH1* h = new TH3F("R__TF3",(char*)GetTitle(),fNpx,fXmin,fXmax
509 ,fNpz,fZmin,fZmax);
510 h->SetDirectory(nullptr);
511 return h;
512}
513
514////////////////////////////////////////////////////////////////////////////////
515/// Paint this 3-D function with its current attributes
516
518{
519
520 TString opt = option;
521 opt.ToLower();
522
523//- Create a temporary histogram and fill each channel with the function value
524 if (!fHistogram) {
525 fHistogram = new TH3F("R__TF3",(char*)GetTitle(),fNpx,fXmin,fXmax
527 ,fNpz,fZmin,fZmax);
528 fHistogram->SetDirectory(nullptr);
529 }
530
532
533 if (opt.Index("tf3") == kNPOS)
534 opt.Append("tf3");
535
536 fHistogram->Paint(opt.Data());
537}
538
539////////////////////////////////////////////////////////////////////////////////
540/// Set the function clipping box (for drawing) "off".
541
543{
545 fClipBox[0] = fClipBox[1] = fClipBox[2] = 0;
546}
547
548////////////////////////////////////////////////////////////////////////////////
549/// Save values of function in array fSave
550
552{
553 if (!fSave.empty()) fSave.clear();
554 Int_t npx = fNpx, npy = fNpy, npz = fNpz;
555 if ((npx < 2) || (npy < 2) || (npz < 2))
556 return;
557
558 Double_t dx = (xmax-xmin)/fNpx;
559 Double_t dy = (ymax-ymin)/fNpy;
560 Double_t dz = (zmax-zmin)/fNpz;
561 if (dx <= 0) {
562 dx = (fXmax-fXmin)/fNpx;
563 npx--;
564 xmin = fXmin + 0.5*dx;
565 xmax = fXmax - 0.5*dx;
566 }
567 if (dy <= 0) {
568 dy = (fYmax-fYmin)/fNpy;
569 npy--;
570 ymin = fYmin + 0.5*dy;
571 ymax = fYmax - 0.5*dy;
572 }
573 if (dz <= 0) {
574 dz = (fZmax-fZmin)/fNpz;
575 npz--;
576 zmin = fZmin + 0.5*dz;
577 zmax = fZmax - 0.5*dz;
578 }
579 Int_t nsave = (npx + 1)*(npy + 1)*(npz + 1);
580 fSave.resize(nsave + 9);
581 Double_t xv[3];
582 Double_t *pp = GetParameters();
583 InitArgs(xv,pp);
584 for (Int_t k = 0, l = 0; k <= npz; k++) {
585 xv[2] = zmin + dz*k;
586 for (Int_t j = 0; j <= npy; j++) {
587 xv[1] = ymin + dy*j;
588 for (Int_t i = 0; i <= npx; i++) {
589 xv[0] = xmin + dx*i;
590 fSave[l++] = EvalPar(xv, pp);
591 }
592 }
593 }
594 fSave[nsave+0] = xmin;
595 fSave[nsave+1] = xmax;
596 fSave[nsave+2] = ymin;
597 fSave[nsave+3] = ymax;
598 fSave[nsave+4] = zmin;
599 fSave[nsave+5] = zmax;
600 fSave[nsave+6] = npx;
601 fSave[nsave+7] = npy;
602 fSave[nsave+8] = npz;
603}
604
605////////////////////////////////////////////////////////////////////////////////
606/// Save primitive as a C++ statement(s) on output stream out
607
608void TF3::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/)
609{
610 char quote = '"';
611 TString f3Name(GetName());
612 out<<" "<<std::endl;
613 if (gROOT->ClassSaved(TF3::Class())) {
614 out<<" ";
615 } else {
616 out<<" TF3 *";
617 }
618 if (!fMethodCall) {
619 out<<f3Name.Data()<<" = new TF3("<<quote<<f3Name.Data()<<quote<<","<<quote<<GetTitle()<<quote<<","<<fXmin<<","<<fXmax<<","<<fYmin<<","<<fYmax<<","<<fZmin<<","<<fZmax<<");"<<std::endl;
620 } else {
621 out<<f3Name.Data()<<" = new TF3("<<quote<<f3Name.Data()<<quote<<","<<GetTitle()<<","<<fXmin<<","<<fXmax<<","<<fYmin<<","<<fYmax<<","<<fZmin<<","<<fZmax<<","<<GetNpar()<<");"<<std::endl;
622 }
623
624 if (GetFillColor() != 0) {
626 out<<" "<<f3Name.Data()<<"->SetFillColor(ci);" << std::endl;
627 else
628 out<<" "<<f3Name.Data()<<"->SetFillColor("<<GetFillColor()<<");"<<std::endl;
629 }
630 if (GetLineColor() != 1) {
632 out<<" "<<f3Name.Data()<<"->SetLineColor(ci);" << std::endl;
633 else
634 out<<" "<<f3Name.Data()<<"->SetLineColor("<<GetLineColor()<<");"<<std::endl;
635 }
636
637 if (GetNpx() != 30)
638 out<<" "<<f3Name.Data()<<"->SetNpx("<<GetNpx()<<");"<<std::endl;
639 if (GetNpy() != 30)
640 out<<" "<<f3Name.Data()<<"->SetNpy("<<GetNpy()<<");"<<std::endl;
641 if (GetNpz() != 30)
642 out<<" "<<f3Name.Data()<<"->SetNpz("<<GetNpz()<<");"<<std::endl;
643
644 if (GetChisquare() != 0)
645 out<<" "<<f3Name.Data()<<"->SetChisquare("<<GetChisquare()<<");"<<std::endl;
646
647 Double_t parmin, parmax;
648 for (Int_t i=0;i<GetNpar();i++) {
649 out<<" "<<f3Name.Data()<<"->SetParameter("<<i<<","<<GetParameter(i)<<");"<<std::endl;
650 out<<" "<<f3Name.Data()<<"->SetParError("<<i<<","<<GetParError(i)<<");"<<std::endl;
651 GetParLimits(i,parmin,parmax);
652 out<<" "<<f3Name.Data()<<"->SetParLimits("<<i<<","<<parmin<<","<<parmax<<");"<<std::endl;
653 }
654
655 if (GetXaxis()) GetXaxis()->SaveAttributes(out, f3Name.Data(), "->GetXaxis()");
656 if (GetYaxis()) GetYaxis()->SaveAttributes(out, f3Name.Data(), "->GetYaxis()");
657 if (GetZaxis()) GetZaxis()->SaveAttributes(out, f3Name.Data(), "->GetZaxis()");
658
659 out<<" "<<f3Name.Data()<<"->Draw("
660 <<quote<<option<<quote<<");"<<std::endl;
661}
662
663////////////////////////////////////////////////////////////////////////////////
664/// Set the function clipping box (for drawing) "on" and define the clipping box.
665/// xclip, yclip and zclip is a point within the function range. All the
666/// function values having x<=xclip and y<=yclip and z>=zclip are clipped.
667
669{
671 fClipBox[0] = xclip;
672 fClipBox[1] = yclip;
673 fClipBox[2] = zclip;
674}
675
676////////////////////////////////////////////////////////////////////////////////
677/// Set the number of points used to draw the function
678///
679/// The default number of points along x is 30 for 2-d/3-d functions.
680/// You can increase this value to get a better resolution when drawing
681/// pictures with sharp peaks or to get a better result when using TF3::GetRandom2
682/// the minimum number of points is 4, the maximum is 10000 for 2-d/3-d functions
683
685{
686 if (npz < 4) {
687 Warning("SetNpz","Number of points must be >=4 && <= 10000, fNpz set to 4");
688 fNpz = 4;
689 } else if(npz > 10000) {
690 Warning("SetNpz","Number of points must be >=4 && <= 10000, fNpz set to 10000");
691 fNpz = 10000;
692 } else {
693 fNpz = npz;
694 }
695 Update();
696}
697
698////////////////////////////////////////////////////////////////////////////////
699/// Initialize the upper and lower bounds to draw the function
700
702{
703 fXmin = xmin;
704 fXmax = xmax;
705 fYmin = ymin;
706 fYmax = ymax;
707 fZmin = zmin;
708 fZmax = zmax;
709 Update();
710}
711
712////////////////////////////////////////////////////////////////////////////////
713/// Stream an object of class TF3.
714
716{
717 if (R__b.IsReading()) {
718 UInt_t R__s, R__c;
719 Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
720 if (R__v > 0) {
721 R__b.ReadClassBuffer(TF3::Class(), this, R__v, R__s, R__c);
722 return;
723 }
724
725 } else {
726 Int_t saved = 0;
727 if (fType != EFType::kFormula && fSave.empty() ) { saved = 1; Save(fXmin,fXmax,fYmin,fYmax,fZmin,fZmax);}
728
729 R__b.WriteClassBuffer(TF3::Class(),this);
730
731 if (saved) { fSave.clear(); }
732 }
733}
734
735////////////////////////////////////////////////////////////////////////////////
736/// Return x^nx * y^ny * z^nz moment of a 3d function in range [ax,bx],[ay,by],[az,bz]
737/// \author Gene Van Buren <gene@bnl.gov>
738
740{
741 Double_t norm = Integral(ax,bx,ay,by,az,bz,epsilon);
742 if (norm == 0) {
743 Error("Moment3", "Integral zero over range");
744 return 0;
745 }
746
747 // define integrand function as a lambda : g(x,y,z)= x^(nx) * y^(ny) * z^(nz) * f(x,y,z)
748 auto integrand = [&](double *x, double *) {
749 return std::pow(x[0], nx) * std::pow(x[1], ny) * std::pow(x[2], nz) * this->EvalPar(x, nullptr);
750 };
751 // compute integral of g(x,y,z)
752 TF3 fnc("TF3_ExpValHelper", integrand, ax, bx, ay, by, az, bz, 0);
753 // set same points as current function to get correct max points when computing the integral
754 fnc.fNpx = fNpx;
755 fnc.fNpy = fNpy;
756 fnc.fNpz = fNpz;
757 return fnc.Integral(ax, bx, ay, by, az, bz, epsilon) / norm;
758}
759
760////////////////////////////////////////////////////////////////////////////////
761/// Return x^nx * y^ny * z^nz central moment of a 3d function in range [ax,bx],[ay,by],[az,bz]
762/// \author Gene Van Buren <gene@bnl.gov>
763
765{
766 Double_t norm = Integral(ax,bx,ay,by,az,bz,epsilon);
767 if (norm == 0) {
768 Error("CentralMoment3", "Integral zero over range");
769 return 0;
770 }
771
772 Double_t xbar = 0;
773 Double_t ybar = 0;
774 Double_t zbar = 0;
775 if (nx!=0) {
776 // compute first momentum in x
777 auto integrandX = [&](double *x, double *) { return x[0] * this->EvalPar(x, nullptr); };
778 TF3 fncx("TF3_ExpValHelperx", integrandX, ax, bx, ay, by, az, bz, 0);
779 fncx.fNpx = fNpx;
780 fncx.fNpy = fNpy;
781 fncx.fNpz = fNpz;
782 xbar = fncx.Integral(ax, bx, ay, by, az, bz, epsilon) / norm;
783 }
784 if (ny!=0) {
785 auto integrandY = [&](double *x, double *) { return x[1] * this->EvalPar(x, nullptr); };
786 TF3 fncy("TF3_ExpValHelpery", integrandY, ax, bx, ay, by, az, bz, 0);
787 fncy.fNpx = fNpx;
788 fncy.fNpy = fNpy;
789 fncy.fNpz = fNpz;
790 ybar = fncy.Integral(ax,bx,ay,by,az,bz,epsilon)/norm;
791 }
792 if (nz!=0) {
793 auto integrandZ = [&](double *x, double *) { return x[2] * this->EvalPar(x, nullptr); };
794 TF3 fncz("TF3_ExpValHelperz", integrandZ, ax, bx, ay, by, az, bz, 0);
795 fncz.fNpx = fNpx;
796 fncz.fNpy = fNpy;
797 fncz.fNpz = fNpz;
798 zbar = fncz.Integral(ax,bx,ay,by,az,bz,epsilon)/norm;
799 }
800 // define integrand function as a lambda : g(x,y)= (x-xbar)^(nx) * (y-ybar)^(ny) * f(x,y)
801 auto integrand = [&](double *x, double *) {
802 double xxx = (nx != 0) ? std::pow(x[0] - xbar, nx) : 1.;
803 double yyy = (ny != 0) ? std::pow(x[1] - ybar, ny) : 1.;
804 double zzz = (nz != 0) ? std::pow(x[2] - zbar, nz) : 1.;
805 return xxx * yyy * zzz * this->EvalPar(x, nullptr);
806 };
807 // compute integral of g(x,y, z)
808 TF3 fnc("TF3_ExpValHelper",integrand,ax,bx,ay,by,az,bz,0) ;
809 fnc.fNpx = fNpx;
810 fnc.fNpy = fNpy;
811 fnc.fNpz = fNpz;
812 return fnc.Integral(ax,bx,ay,by,az,bz,epsilon)/norm;
813}
#define b(i)
Definition RSha256.hxx:100
#define f(i)
Definition RSha256.hxx:104
#define a(i)
Definition RSha256.hxx:99
#define h(i)
Definition RSha256.hxx:106
RooAbsReal & function()
int Int_t
Definition RtypesCore.h:45
short Version_t
Definition RtypesCore.h:65
unsigned int UInt_t
Definition RtypesCore.h:46
constexpr Bool_t kFALSE
Definition RtypesCore.h:101
constexpr Ssiz_t kNPOS
Definition RtypesCore.h:124
constexpr Bool_t kTRUE
Definition RtypesCore.h:100
const char Option_t
Definition RtypesCore.h:66
#define ClassImp(name)
Definition Rtypes.h:377
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
char name[80]
Definition TGX11.cxx:110
float xmin
float ymin
float xmax
float ymax
Int_t gDebug
Definition TROOT.cxx:597
#define gROOT
Definition TROOT.h:407
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
#define gPad
Param Functor class for Multidimensional functions.
virtual Color_t GetFillColor() const
Return the fill area color.
Definition TAttFill.h:30
virtual Color_t GetLineColor() const
Return the line color.
Definition TAttLine.h:33
void SaveAttributes(std::ostream &out, const char *name, const char *subname) override
Save axis attributes as C++ statement(s) on output stream out.
Definition TAxis.cxx:710
Buffer base class used for serializing objects.
Definition TBuffer.h:43
virtual Version_t ReadVersion(UInt_t *start=nullptr, UInt_t *bcnt=nullptr, const TClass *cl=nullptr)=0
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=nullptr)=0
Bool_t IsReading() const
Definition TBuffer.h:86
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
static Bool_t SaveColor(std::ostream &out, Int_t ci)
Save a color with index > 228 as a C++ statement(s) on output stream out.
Definition TColor.cxx:2426
Int_t fNdim
Function dimension.
Definition TF1.h:266
virtual void GetParLimits(Int_t ipar, Double_t &parmin, Double_t &parmax) const
Return limits for parameter ipar.
Definition TF1.cxx:1942
TAxis * GetYaxis() const
Get y axis of the function.
Definition TF1.cxx:2411
virtual Double_t GetParError(Int_t ipar) const
Return value of parameter number ipar.
Definition TF1.cxx:1932
Double_t GetChisquare() const
Return the Chisquare after fitting. See ROOT::Fit::FitResult::Chi2()
Definition TF1.h:470
Double_t fXmin
Lower bounds for the range.
Definition TF1.h:263
std::unique_ptr< TMethodCall > fMethodCall
! Pointer to MethodCall in case of interpreted function
Definition TF1.h:284
virtual void Update()
Called by functions such as SetRange, SetNpx, SetParameters to force the deletion of the associated h...
Definition TF1.cxx:3619
TAxis * GetZaxis() const
Get z axis of the function. (In case this object is a TF2 or TF3)
Definition TF1.cxx:2422
virtual Int_t GetNpar() const
Definition TF1.h:507
TH1 * fHistogram
! Pointer to histogram used for visualisation
Definition TF1.h:283
virtual Double_t * GetParameters() const
Definition TF1.h:546
virtual void InitArgs(const Double_t *x, const Double_t *params)
Initialize parameters addresses.
Definition TF1.cxx:2482
virtual Double_t IntegralMultiple(Int_t n, const Double_t *a, const Double_t *b, Int_t maxpts, Double_t epsrel, Double_t epsabs, Double_t &relerr, Int_t &nfnevl, Int_t &ifail)
This function computes, to an attempted specified accuracy, the value of the integral.
Definition TF1.cxx:2851
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
Compute distance from point px,py to a function.
Definition TF1.cxx:1295
EFType fType
Definition TF1.h:268
virtual Double_t EvalPar(const Double_t *x, const Double_t *params=nullptr)
Evaluate function with given coordinates and parameters.
Definition TF1.cxx:1470
Int_t fNpx
Number of points used for the graphical representation.
Definition TF1.h:267
void ExecuteEvent(Int_t event, Int_t px, Int_t py) override
Execute action corresponding to one event.
Definition TF1.cxx:1538
std::vector< Double_t > fSave
Array of fNsave function values.
Definition TF1.h:277
virtual Double_t GetMinMaxNDim(Double_t *x, Bool_t findmax, Double_t epsilon=0, Int_t maxiter=0) const
Find the minimum of a function of whatever dimension.
Definition TF1.cxx:1725
std::vector< Double_t > fIntegral
! Integral of function binned on fNpx bins
Definition TF1.h:278
@ kFormula
Formula functions which can be stored,.
Definition TF1.h:255
virtual Int_t GetNpx() const
Definition TF1.h:516
Double_t fXmax
Upper bounds for the range.
Definition TF1.h:264
virtual Int_t GetNdim() const
Definition TF1.h:511
virtual Double_t GetParameter(Int_t ipar) const
Definition TF1.h:538
TAxis * GetXaxis() const
Get x axis of the function.
Definition TF1.cxx:2400
A 2-Dim function with parameters.
Definition TF2.h:29
void Copy(TObject &f2) const override
Copy this F2 to a new F2.
Definition TF2.cxx:210
Int_t fNpy
Number of points along y used for the graphical representation.
Definition TF2.h:34
Double_t fYmax
Upper bound for the range in y.
Definition TF2.h:33
Double_t fYmin
Lower bound for the range in y.
Definition TF2.h:32
Int_t GetNpy() const
Definition TF2.h:96
A 3-Dim function with parameters.
Definition TF3.h:28
void GetRange(Double_t &xmin, Double_t &xmax) const override
Return range of a 1-D function.
Definition TF3.h:139
static TClass * Class()
TF3()
F3 default constructor.
Definition TF3.cxx:37
void Paint(Option_t *option="") override
Paint this 3-D function with its current attributes.
Definition TF3.cxx:517
Int_t GetNpz() const
Definition TF3.h:91
virtual void GetRandom3(Double_t &xrandom, Double_t &yrandom, Double_t &zrandom, TRandom *rng=nullptr)
Return 3 random numbers following this function shape.
Definition TF3.cxx:342
Bool_t IsInside(const Double_t *x) const override
Return kTRUE is the point is inside the function range.
Definition TF3.cxx:494
void ExecuteEvent(Int_t event, Int_t px, Int_t py) override
Execute action corresponding to one event.
Definition TF3.cxx:209
Double_t fZmin
Lower bound for the range in z.
Definition TF3.h:31
TF3 & operator=(const TF3 &rhs)
Operator =.
Definition TF3.cxx:145
virtual void SetNpz(Int_t npz=30)
Set the number of points used to draw the function.
Definition TF3.cxx:684
virtual Double_t GetMaximumXYZ(Double_t &x, Double_t &y, Double_t &z)
Compute the X, Y and Z values corresponding to the maximum value of the function on its range.
Definition TF3.cxx:314
virtual Double_t CentralMoment3(Double_t nx, Double_t ax, Double_t bx, Double_t ny, Double_t ay, Double_t by, Double_t nz, Double_t az, Double_t bz, Double_t epsilon=0.000001)
Return x^nx * y^ny * z^nz central moment of a 3d function in range [ax,bx],[ay,by],...
Definition TF3.cxx:764
void Copy(TObject &f3) const override
Copy this F3 to a new F3.
Definition TF3.cxx:170
Double_t FindMinMax(Double_t *x, bool findmax) const override
Return minimum/maximum value of the function.
Definition TF3.cxx:230
~TF3() override
F3 default destructor.
Definition TF3.cxx:155
Bool_t fClipBoxOn
! Is clip box on
Definition TF3.h:34
virtual Double_t GetMinimumXYZ(Double_t &x, Double_t &y, Double_t &z)
Compute the X, Y and Z values corresponding to the minimum value of the function on its range.
Definition TF3.cxx:298
virtual Double_t Moment3(Double_t nx, Double_t ax, Double_t bx, Double_t ny, Double_t ay, Double_t by, Double_t nz, Double_t az, Double_t bz, Double_t epsilon=0.000001)
Return x^nx * y^ny * z^nz moment of a 3d function in range [ax,bx],[ay,by],[az,bz].
Definition TF3.cxx:739
Int_t fNpz
Number of points along z used for the graphical representation.
Definition TF3.h:33
void Draw(Option_t *option="") override
Draw this function with its current attributes.
Definition TF3.cxx:194
virtual void SetClippingBoxOn(Double_t xclip=0, Double_t yclip=0, Double_t zclip=0)
Set the function clipping box (for drawing) "on" and define the clipping box.
Definition TF3.cxx:668
virtual void SetClippingBoxOff()
Set the function clipping box (for drawing) "off".
Definition TF3.cxx:542
Double_t fZmax
Upper bound for the range in z.
Definition TF3.h:32
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
Compute distance from point px,py to a function.
Definition TF3.cxx:185
virtual Double_t Integral(Double_t ax, Double_t bx, Double_t ay, Double_t by, Double_t az, Double_t bz, Double_t epsrel=1.e-6)
Return Integral of a 3d function in range [ax,bx],[ay,by],[az,bz] with a desired relative accuracy.
Definition TF3.cxx:468
void Streamer(TBuffer &) override
Stream an object of class TF3.
Definition TF3.cxx:715
TH1 * CreateHistogram() override
Create a histogram for axis range.
Definition TF3.cxx:505
void Save(Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax) override
Save values of function in array fSave.
Definition TF3.cxx:551
Double_t GetSave(const Double_t *x) override
Get value corresponding to X in array of fSave values.
Definition TF3.cxx:415
Double_t fClipBox[3]
! Coordinates of clipbox
Definition TF3.h:35
void SetRange(Double_t xmin, Double_t xmax) override
Initialize the upper and lower bounds to draw the function.
Definition TF3.h:143
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save primitive as a C++ statement(s) on output stream out.
Definition TF3.cxx:608
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:59
virtual void SetDirectory(TDirectory *dir)
By default, when a histogram is created, it is added to the list of histogram objects in the current ...
Definition TH1.cxx:8905
TVirtualHistPainter * GetPainter(Option_t *option="")
Return pointer to painter.
Definition TH1.cxx:4488
void Paint(Option_t *option="") override
Control routine to paint any kind of histograms.
Definition TH1.cxx:6174
3-D histogram with a float per channel (see TH1 documentation)
Definition TH3.h:317
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:48
Mother of all ROOT objects.
Definition TObject.h:41
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:973
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition TObject.cxx:184
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:987
void MakeZombie()
Definition TObject.h:53
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition TObject.cxx:961
This is the base class for the ROOT Random number generators.
Definition TRandom.h:27
Double_t Rndm() override
Machine independent random number generator.
Definition TRandom.cxx:559
Basic string class.
Definition TString.h:139
void ToLower()
Change string to lower-case.
Definition TString.cxx:1182
const char * Data() const
Definition TString.h:376
TString & Append(const char *cs)
Definition TString.h:572
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:632
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition TString.h:651
virtual void ProcessMessage(const char *mess, const TObject *obj)=0
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Definition TMathBase.h:250
Double_t QuietNaN()
Returns a quiet NaN as defined by IEEE 754.
Definition TMath.h:902
Int_t Finite(Double_t x)
Check if it is finite with a mask in order to be consistent in presence of fast math.
Definition TMath.h:770
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Definition TMathBase.h:198
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Binary search in an array of n values to locate value.
Definition TMathBase.h:347
Double_t Infinity()
Returns an infinity as defined by the IEEE standard.
Definition TMath.h:917
TLine l
Definition textangle.C:4
auto * tt
Definition textangle.C:16