How to generate Gaussian peak?

Particles, Track, Event, Run and Biasing
1

Please fill out the following information to help in answering your question, and also see tips for posting code snippets. If you don’t provide this information it will take more time to help with your problem!

_Geant4 Version: 11.1
_Operating System: Ubuntu 20.04
Compiler/Version:
CMake Version:

Dear all,
I want to simulate HPGe detector response by using Geant4. Does anyone have a way to generate Gauss peaks? I tried to use Root to read results and add a sigma factor, but it was incorrect. So is there any way to directly handle the results of each step or event?

Root file analysis:

void readNTuple()
{
    gROOT->Reset();
	gROOT->SetStyle("Plain");

	TFile* f = new TFile("defaultFilename.root", "READ");

	// list all Tree and Historgram
	f->ls();

    TTree* tree = (TTree*)f->Get("HPGe;1");

    tree->Print();

    double eDep;
    tree->SetBranchAddress("fEdep", &eDep);

    int entries;
    entries = tree->GetEntries();
    cout << "Number of entries: " << entries << endl;

    TH1D* eDepHist = new TH1D("HPGe Histogram", "Deposition Energies", 200, 0., 2.);
    TH1D* eDepFittedHist = new TH1D("HPGe Histogram", "Deposition Energies", 200, 0., 2.);

    double Ei;
    double sigma;

    for(int i=0; i<entries; ++i)
    {
        tree->GetEntry(i);

        Ei = (eDep + gRandom->Rndm());
        sigma = (0.5/100)*Ei; // 0.5%

        eDepHist->Fill(eDep);
        eDepFittedHist->Fill(gRandom->Gaus(Ei,sigma));
    }

    // Draw
    TCanvas* c1 = new TCanvas("c1", " ", 20, 20, 800, 1000);
    c1->Divide(1, 2);

    // Draw pad 1
	c1->cd(1);
	// gPad->SetLogy(1);
	eDepHist->SetTitle("HPGe Spectrum ");
	eDepHist->GetXaxis()->CenterTitle(true);
	eDepHist->GetYaxis()->CenterTitle(true);
	eDepHist->GetXaxis()->SetTitle("MeV");
	eDepHist->GetYaxis()->SetTitle("Counts");
	eDepHist->Draw("HIST");

    // Draw pad 2
	c1->cd(2);
    eDepFittedHist->SetTitle("HPGe Fitted Spectrum");
	eDepFittedHist->GetXaxis()->CenterTitle(true);
	eDepFittedHist->GetYaxis()->CenterTitle(true);
	eDepFittedHist->GetXaxis()->SetTitle("MeV");
	eDepFittedHist->GetYaxis()->SetTitle("Counts");
	eDepFittedHist->Draw("HIST");

}
2

Dear @hungbt1908 ,

I have not run your code but the logic of randomizing the bin and then adding an entry that has been “Gaussianized” seems correct. If you want to do this inside Geant4 one would basically do the same at the place where you are filling your Tree.

Note that you might want to save two different branches, one raw and one with the smearing added. This is so in case you realize you used the wrong detector resolution, you can still fix it in postanalysis as you are doing here that you claim it’s incorrect.

/Pico

3

If you aim to create your primary particle source with certain spatial, spectral or angular distribution specifications (e.g. Gaussian distribution, exponential distribution), you need to replace your G4ParticleGun with G4GeneralParticleSource.

You can then use macro-commands to define your source.

For example, a gaussian energy profile centered on 100 MeV with a certain deviation is defined as follows:

/gps/ene/type Gauss
/gps/ene/mono 100. MeV
/gps/ene/sigma 20. MeV

Cheers,
Younes

4

I use the chrono and random to generate the gauss position proton
I don’t know whether the way is correct.

image
image1107×530 51.4 KB

5

Thanks @pico, @Ziyang_He, @belmoussay