From $ROOTSYS/tutorials/foam/foam_demo.C

// Demonstrate the TFoam class.
//
//  To run this macro type from CINT command line
//
//  root [0] gSystem->Load("libFoam.so")
//  root [1] .x foam_demo.C+
//Author: Stascek Jadach
//____________________________________________________________________________

#if defined(__CINT__) && !defined(__MAKECINT__)
{
   std::cout << "Using ACliC to run this macro since it uses custom classes" << std::endl;
   TString macroFileName = gSystem->UnixPathName(__FILE__);
   gSystem->CompileMacro(macroFileName, "k");
   foam_demo();
}
#else

#include "Riostream.h"
#include "TFile.h"
#include "TFoam.h"
#include "TH1.h"
#include "TMath.h"
#include "TFoamIntegrand.h"
#include "TRandom3.h"

class TFDISTR: public TFoamIntegrand {
public:
  TFDISTR(){};
  Double_t Density(int nDim, Double_t *Xarg){
  // Integrand for mFOAM
  Double_t Fun1,Fun2,R1,R2;
  Double_t pos1=1e0/3e0;
  Double_t pos2=2e0/3e0;
  Double_t Gam1= 0.100e0;  // as in JPC
  Double_t Gam2= 0.100e0;  // as in JPC
  Double_t sPi = sqrt(TMath::Pi());
  Double_t xn1=1e0;
  Double_t xn2=1e0;
  int i;
  R1=0;
  R2=0;
  for(i = 0 ; i<nDim ; i++){
    R1=R1+(Xarg[i] -pos1)*(Xarg[i] -pos1);
    R2=R2+(Xarg[i] -pos2)*(Xarg[i] -pos2);
    xn1=xn1*Gam1*sPi;
    xn2=xn2*Gam2*sPi;
  }
  R1   = sqrt(R1);
  R2   = sqrt(R2);
  Fun1 = exp(-(R1*R1)/(Gam1*Gam1))/xn1;  // Gaussian delta-like profile
  Fun2 = exp(-(R2*R2)/(Gam2*Gam2))/xn2;  // Gaussian delta-like profile
  return 0.5e0*(Fun1+ Fun2);
}
  ClassDef(TFDISTR,1) //Class of testing functions for FOAM
};
ClassImp(TFDISTR)

Int_t foam_demo()
{
  TFile RootFile("foam_demo.root","RECREATE","histograms");
  long   loop;
  Double_t MCresult,MCerror,MCwt;
//=========================================================
  long NevTot   =     50000;   // Total MC statistics
  Int_t  kDim   =         2;   // total dimension
  Int_t  nCells   =     500;   // Number of Cells
  Int_t  nSampl   =     200;   // Number of MC events per cell in build-up
  Int_t  nBin     =       8;   // Number of bins in build-up
  Int_t  OptRej   =       1;   // Wted events for OptRej=0; wt=1 for OptRej=1 (default)
  Int_t  OptDrive =       2;   // (D=2) Option, type of Drive =0,1,2 for TrueVol,Sigma,WtMax
  Int_t  EvPerBin =      25;   // Maximum events (equiv.) per bin in buid-up
  Int_t  Chat     =       1;   // Chat level
//=========================================================
  TRandom *PseRan   = new TRandom3();  // Create random number generator
  TFoam   *FoamX    = new TFoam("FoamX");   // Create Simulator
  TFoamIntegrand    *rho= new TFDISTR();
  PseRan->SetSeed(4357);
//=========================================================
  cout<<"*****   Demonstration Program for Foam version "<<FoamX->GetVersion()<<"    *****"<<endl;
  FoamX->SetkDim(        kDim);      // Mandatory!!!
  FoamX->SetnCells(      nCells);    // optional
  FoamX->SetnSampl(      nSampl);    // optional
  FoamX->SetnBin(        nBin);      // optional
  FoamX->SetOptRej(      OptRej);    // optional
  FoamX->SetOptDrive(    OptDrive);  // optional
  FoamX->SetEvPerBin(    EvPerBin);  // optional
  FoamX->SetChat(        Chat);      // optional
//===============================
  FoamX->SetRho(rho);
  FoamX->SetPseRan(PseRan);
  FoamX->Initialize(); // Initialize simulator
  FoamX->Write("FoamX");     // Writing Foam on the disk, TESTING PERSISTENCY!!!
//===============================
  long nCalls=FoamX->GetnCalls();
  cout << "====== Initialization done, entering MC loop" << endl;
//======================================================================
  //cout<<" About to start MC loop: ";  cin.getline(question,20);
  Double_t *MCvect =new Double_t[kDim]; // vector generated in the MC run
//======================================================================
  TH1D  *hst_Wt = new TH1D("hst_Wt" ,  "Main weight of Foam",25,0,1.25);
  hst_Wt->Sumw2();
//======================================================================
  for(loop=0; loop<NevTot; loop++){
//===============================
    FoamX->MakeEvent();           // generate MC event
//===============================
    FoamX->GetMCvect( MCvect);
    MCwt=FoamX->GetMCwt();
    hst_Wt->Fill(MCwt,1.0);
    if(loop<15){
      cout<<"MCwt= "<<MCwt<<",  ";
      cout<<"MCvect= ";
      for ( Int_t k=0 ; k<kDim ; k++) cout<<MCvect[k]<<" "; cout<< endl;
    }
    if( ((loop)%100000)==0 ){
      cout<<"   loop= "<<loop<<endl;
    }
  }
//======================================================================
//======================================================================
  cout << "====== Events generated, entering Finalize" << endl;

  hst_Wt->Print("all");
  Double_t eps = 0.0005;
  Double_t Effic, WtMax, AveWt, Sigma;
  Double_t IntNorm, Errel;
  FoamX->Finalize(   IntNorm, Errel);     // final printout
  FoamX->GetIntegMC( MCresult, MCerror);  // get MC intnegral
  FoamX->GetWtParams(eps, AveWt, WtMax, Sigma); // get MC wt parameters
  Effic=0; if(WtMax>0) Effic=AveWt/WtMax;
  cout << "================================================================" << endl;
  cout << " MCresult= " << MCresult << " +- " << MCerror << " RelErr= "<< MCerror/MCresult << endl;
  cout << " Dispersion/<wt>= " << Sigma/AveWt << endl;
  cout << "      <wt>/WtMax= " << Effic <<",    for epsilon = "<<eps << endl;
  cout << " nCalls (initialization only) =   " << nCalls << endl;
  cout << "================================================================" << endl;

  delete [] MCvect;
  //
  RootFile.ls();
  RootFile.Write();
  RootFile.Close();
  cout << "***** End of Demonstration Program  *****" << endl;

  return 0;
} // end of Demo

#endif

 foam_demo.C:1
 foam_demo.C:2
 foam_demo.C:3
 foam_demo.C:4
 foam_demo.C:5
 foam_demo.C:6
 foam_demo.C:7
 foam_demo.C:8
 foam_demo.C:9
 foam_demo.C:10
 foam_demo.C:11
 foam_demo.C:12
 foam_demo.C:13
 foam_demo.C:14
 foam_demo.C:15
 foam_demo.C:16
 foam_demo.C:17
 foam_demo.C:18
 foam_demo.C:19
 foam_demo.C:20
 foam_demo.C:21
 foam_demo.C:22
 foam_demo.C:23
 foam_demo.C:24
 foam_demo.C:25
 foam_demo.C:26
 foam_demo.C:27
 foam_demo.C:28
 foam_demo.C:29
 foam_demo.C:30
 foam_demo.C:31
 foam_demo.C:32
 foam_demo.C:33
 foam_demo.C:34
 foam_demo.C:35
 foam_demo.C:36
 foam_demo.C:37
 foam_demo.C:38
 foam_demo.C:39
 foam_demo.C:40
 foam_demo.C:41
 foam_demo.C:42
 foam_demo.C:43
 foam_demo.C:44
 foam_demo.C:45
 foam_demo.C:46
 foam_demo.C:47
 foam_demo.C:48
 foam_demo.C:49
 foam_demo.C:50
 foam_demo.C:51
 foam_demo.C:52
 foam_demo.C:53
 foam_demo.C:54
 foam_demo.C:55
 foam_demo.C:56
 foam_demo.C:57
 foam_demo.C:58
 foam_demo.C:59
 foam_demo.C:60
 foam_demo.C:61
 foam_demo.C:62
 foam_demo.C:63
 foam_demo.C:64
 foam_demo.C:65
 foam_demo.C:66
 foam_demo.C:67
 foam_demo.C:68
 foam_demo.C:69
 foam_demo.C:70
 foam_demo.C:71
 foam_demo.C:72
 foam_demo.C:73
 foam_demo.C:74
 foam_demo.C:75
 foam_demo.C:76
 foam_demo.C:77
 foam_demo.C:78
 foam_demo.C:79
 foam_demo.C:80
 foam_demo.C:81
 foam_demo.C:82
 foam_demo.C:83
 foam_demo.C:84
 foam_demo.C:85
 foam_demo.C:86
 foam_demo.C:87
 foam_demo.C:88
 foam_demo.C:89
 foam_demo.C:90
 foam_demo.C:91
 foam_demo.C:92
 foam_demo.C:93
 foam_demo.C:94
 foam_demo.C:95
 foam_demo.C:96
 foam_demo.C:97
 foam_demo.C:98
 foam_demo.C:99
 foam_demo.C:100
 foam_demo.C:101
 foam_demo.C:102
 foam_demo.C:103
 foam_demo.C:104
 foam_demo.C:105
 foam_demo.C:106
 foam_demo.C:107
 foam_demo.C:108
 foam_demo.C:109
 foam_demo.C:110
 foam_demo.C:111
 foam_demo.C:112
 foam_demo.C:113
 foam_demo.C:114
 foam_demo.C:115
 foam_demo.C:116
 foam_demo.C:117
 foam_demo.C:118
 foam_demo.C:119
 foam_demo.C:120
 foam_demo.C:121
 foam_demo.C:122
 foam_demo.C:123
 foam_demo.C:124
 foam_demo.C:125
 foam_demo.C:126
 foam_demo.C:127
 foam_demo.C:128
 foam_demo.C:129
 foam_demo.C:130
 foam_demo.C:131
 foam_demo.C:132
 foam_demo.C:133
 foam_demo.C:134
 foam_demo.C:135
 foam_demo.C:136
 foam_demo.C:137
 foam_demo.C:138
 foam_demo.C:139
 foam_demo.C:140
 foam_demo.C:141
 foam_demo.C:142
 foam_demo.C:143
 foam_demo.C:144
 foam_demo.C:145
 foam_demo.C:146
 foam_demo.C:147
 foam_demo.C:148
 foam_demo.C:149
 foam_demo.C:150
 foam_demo.C:151