49 int dim =
X.GetNrows();
51 for (
int i = 0; i<dim; ++i) {
X[i] = x[i] - p[k]; k++; }
52 for (
int i = 0; i<dim; ++i) {
53 CovMat(i,i) = p[k]*p[k];
56 for (
int i = 0; i<dim; ++i) {
57 for (
int j = i+1; j<dim; ++j) {
59 CovMat(i,j) = p[k]*
sqrt(CovMat(i,i)*CovMat(j,j));
60 CovMat(j,i) = CovMat(i,j);
69 double det = CovMat.Determinant();
71 Fatal(
"GausND",
"Determinant is <= 0 det = %f",det);
78 double fval =
std::exp( - 0.5 * CovMat.Similarity(
X) )/ norm;
81 std::cout <<
"det " << det << std::endl;
82 std::cout <<
"norm " << norm << std::endl;
83 std::cout <<
"fval " << fval << std::endl;
90 using namespace ROOT::Math;
92 void multidimSampling() {
95 std::cout <<
"DO NOT RUN WITH CINT:" << std::endl;
96 std::cout <<
"we are using a custom function which requires" << std::endl;
97 std::cout <<
"that this tutorial must be compiled with ACLIC" << std::endl;
101 const int N = 100000;
105 double xmin[] = {-10,-10,-10, -10};
106 double xmax[] = { 10, 10, 10, 10};
107 double par0[] = { 1., -1., 2, 0,
109 0.5,0.,0.,0.,0.,0.8 };
111 const int NPAR = DIM + DIM*(DIM+1)/2;
114 TF1 * f =
new TF1(
"functionND",gaus4d,0,1,14);
115 f->SetParameters(par0);
117 double x0[] = {0,0,0,0};
119 if (debug) f->EvalPar(x0,0);
123 for (
int i = 0; i <
NPAR; ++i ) {
124 if (i < DIM) f->SetParName(i, name.
Format(
"mu_%d",i+1) );
125 else if (i < 2*DIM) f->SetParName(i, name.
Format(
"sig_%d",i-DIM+1) );
126 else if (i < 2*DIM) f->SetParName(i, name.
Format(
"sig_%d",i-2*DIM+1) );
132 Info(
"multidimSampling",
"Default sampler %s is not available try with Foam ",
138 Error(
"multidimSampling",
"Foam sampler is not available - exit ");
144 bool ret = sampler->
Init();
146 std::vector<double> data1(DIM*N);
151 Error(
"Sampler::Init",
"Error initializing unuran sampler");
157 for (
int i = 0; i <
N; ++i) {
159 for (
int j = 0; j < DIM; ++j)
160 data1[N*j + i] = v[j];
166 TFile * file =
new TFile(
"multiDimSampling.root",
"RECREATE");
168 TTree * t1 =
new TTree(
"t1",
"Tree from Unuran");
169 t1->
Branch(
"x",x,
"x[4]/D");
170 for (
int i = 0; i <
N; ++i) {
171 for (
int j = 0; j < DIM; ++j) {
178 t1->
Draw(
"x[0]:x[1]:x[2]:x[3]",
"",
"candle");
179 TCanvas * c2 =
new TCanvas();
183 t1->
Draw(
"x[0]:x[1]");
185 t1->
Draw(
"x[0]:x[2]");
187 t1->
Draw(
"x[0]:x[3]");
189 t1->
Draw(
"x[1]:x[2]");
191 t1->
Draw(
"x[1]:x[3]");
193 t1->
Draw(
"x[2]:x[3]");
void Print(Option_t *option="") const
Print the real and cpu time passed between the start and stop events.
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
void Fatal(const char *location, const char *msgfmt,...)
static void SetDefaultSampler(const char *type)
virtual Int_t Fill()
Fill all branches.
const double * Sample()
sample one event and rerturning array x with coordinates
void Stop()
Stop the stopwatch.
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString...
void SetFunction(Function &func, unsigned int dim)
set the parent function distribution to use for sampling (generic case)
double pow(double, double)
void Info(const char *location, const char *msgfmt,...)
virtual Int_t Write(const char *name=0, Int_t option=0, Int_t bufsize=0)
Write this object to the current directory.
static ROOT::Math::DistSampler * CreateDistSampler(const std::string &samplerType="")
static method to create the distribution sampler class given a string specifying the type Supported s...
void SetRange(double xmin, double xmax, int icoord=0)
set range in a given dimension
TRObject operator()(const T1 &t1) const
Interface class for generic sampling of a distribution, i.e.
static const std::string & DefaultSampler()
virtual void Draw(Option_t *opt)
Default Draw method for all objects.
virtual Int_t Branch(TCollection *list, Int_t bufsize=32000, Int_t splitlevel=99, const char *name="")
Create one branch for each element in the collection.
virtual bool Init(const char *="")
initialize the generators with the given algorithm Implemented by derived classes who needs it (like ...
A TTree object has a header with a name and a title.
void Error(ErrorHandler_t func, int code, const char *va_(fmt),...)
Write error message and call a handler, if required.