In This Package:

DsG4GdCaptureGammas.cc File Reference

#include "DsG4GdCaptureGammas.h"
#include "Randomize.hh"
#include <TRandom.h>
#include <vector>
#include "G4Gamma.hh"
#include "TROOT.h"
#include "TFile.h"
#include "TH1F.h"
#include "TMath.h"

Include dependency graph for DsG4GdCaptureGammas.cc:

Go to the source code of this file.

Functions
static double	getRandomCLHEP (TH1 *hh)

---

Function Documentation

static double getRandomCLHEP

(

TH1 *

hh

)

[static]

Definition at line 105 of file DsG4GdCaptureGammas.cc.

{
    // return a random number distributed according the histogram bin contents.
    // This function checks if the bins integral exists. If not, the integral
    // is evaluated, normalized to one.
    // The integral is automatically recomputed if the number of entries
    // is not the same then when the integral was computed.
    // NB Only valid for 1-d histograms. Use GetRandom2 or 3 otherwise.
  
    if (hh->GetDimension()> 1) {
        G4cout << "Error in getRandomCLHEP: Function only valid for 1-d histograms"
               <<endl;
        return 0;
    }
    int nbinsx = hh->GetNbinsX();
    double integral;
  
    double *myIntegral = 0; //array
    myIntegral = hh->GetIntegral();
  
    if (myIntegral) {
        if (myIntegral[nbinsx+1] != hh->GetEntries()) {
            integral = hh->ComputeIntegral();
            //get the pointer to the integral array again
            myIntegral = hh->GetIntegral();
        } 
    } else {
        integral = hh->ComputeIntegral();
        //get the pointer to the integral array again
        myIntegral = hh->GetIntegral(); 
        if (integral == 0 || myIntegral == 0) return 0;
    }
  
    //Here enforce the CLHEP random number generator!
    //Double_t r1 = gRandom->Rndm();
    double r1 = HepRandom::getTheEngine()->flat();
  
    int ibin = TMath::BinarySearch(nbinsx,myIntegral,r1);
    double x = hh->GetBinLowEdge(ibin+1);
    if (r1 > myIntegral[ibin]) 
        x +=
            hh->GetBinWidth(ibin+1)*(r1-myIntegral[ibin])
            /(myIntegral[ibin+1] - myIntegral[ibin]);
    return x;
}

---