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