gRandom->BreitWigner trouble

Hi all:
I’m working on a simulation involving neutron-emitting states, whose energy I would like to populate with gRandom->BreitWigner
Unfortunately, while it seems to behave more or less as expected in CINT (using variously root versions 5.10 and 5.06), it breaks when I compile and run with g++ (version 3.3.5, the latest available on our campus cluster.)
I’m filling a branch in my Tree with gRandom->BreitWigner(energy,width), so I expect a BW distribution at the end. Instead I keep getting lots of bins filled to a constant value, clustered around the desired mean (like so):

    |            |     |  | | |||||  |     |        |                 |       

More of a density plot, and definately not useful to this simulation.
I’ve tried implementing the same thing by filling a ‘parent’ histogram and using h->GetRandom() instead, with the same result.

It’s quite simple, but not working right now!
Any suggestions for getting this to work with a slightly older compiler?
Thank you for the help!
Caleb

Here is the subroutine returning the energy:

float GetBe11State(struct Level Be11)
{
	int index = Be11.index;

	float Width = Be11.Lifetime;

	float E = Be11.E;
	
	float EState = (float)gRandom->BreitWigner(E,Width);

	#ifdef DEBUG

		cout << "GetBe11State w/index=" << index << "EBe11= " << E << " EState= " << EState << endl;

	#endif

	return ( EState );
}

I just found the TRandom3 class (it speeds up my simulation too, I’m pleased to say)

Apparantly the random generator (before, just using TRandom) had a periodicity depending on the seed… and sometimes amazingly small:

GetBe11State w/index=6EBe11= 8.816 EState= 8.845
GetBe11State w/index=6EBe11= 8.816 EState= 8.80491
GetBe11State w/index=6EBe11= 8.816 EState= 8.76096
GetBe11State w/index=6EBe11= 8.816 EState= 8.69105
GetBe11State w/index=6EBe11= 8.816 EState= 8.47117
GetBe11State w/index=6EBe11= 8.816 EState= 9.71773
GetBe11State w/index=6EBe11= 8.816 EState= 8.99769
GetBe11State w/index=6EBe11= 8.816 EState= 8.89603
GetBe11State w/index=6EBe11= 8.816 EState= 8.845
GetBe11State w/index=6EBe11= 8.816 EState= 8.80491
GetBe11State w/index=6EBe11= 8.816 EState= 8.76096
GetBe11State w/index=6EBe11= 8.816 EState= 8.69105

Cheers,
Caleb

The periodicity of TRandom is small, however this cannot explain the problem with all values in the same bin.
Could you post a concrete script showing this problem?
See also le class description of TRandom at:

root.cern.ch/root/htmldoc//TRandom.html

and also section “Improvements in the random Numbers” at

root.cern.ch/root/Version511.news.html

Rene

Rene,
I should have included a graphic earlier: one is now attached along with the script that generated it (this is a simplified version of my simulation).
The script (example.C) should be compiled with g++, and executed as

./example 19283612 (or any number to seed the random generator)

It creates filesave.root, with the TTree ‘mc’
In root, try:
f = new TFile(“filesave.root”)
h = new TH1F(“h”,"",20000,0,20)
mc->Draw(“EBe11>>h”)

Sometimes (apparently depending on the choice of random seed) I get a nice Breit Wigner distribution, other times (I assume because the random periodicity is shorter) I get something like the attached figure.

Fortunately, the easiest solution appears to be making the random generator rootrand a TRandom3 rather than a TRandom.

Cheers,
Caleb



example.C (1.74 KB)

I cannot reproduce your problem with your initial seed.
Are you sure about your seed? I see a starnge way (via atof)
to set your seed.

Rene

Here are some seeds that cause problems:

536870912
65536
1394606080
134217728

I’m choosing to pass the seed as an argument so that, when this runs as a distributed job, my server can give each client unique parameters for a simulation (including a new seed) each time (and the server can track what parameters went where).

Caleb

This is a known problem with TRandom,

see for example this note from J. Heinrich

www-cdf.fnal.gov/publications/cd … tfalls.pdf

seed which have the right bits sets to zero, reduce the period drammatically, and if you give a seed of 2**28 the period is only 2.

The new TRandom algorithm in ROOT 5.11.06 fixes this problem.
However, TRandom, has still a very short period (231 ~ 109) and
SHOULD NEVER be used in any statistical study.
Please use TRandom3,

Regards, Lorenzo

Lorenzo, Rene:
Thank you for the information and the informative link!
One concern: the link in the previous post mentions that seeds of 2**n are particularly bad choices of seed when using TRandom.
The attached script suggests that the problem still exists for TRandom3. (Although, if setting off the system clock, one could just check if the seed is even and then do TRandom3::SetSeed( TRandom3::GetSeed()+1 ) to avoid the danger).
Thank you again for your help!
Caleb
seed.C (332 Bytes)

Caleb,

you are right ! Also TRandom3 has the problem with these seeds,
in fact it is using the same multipliers as TRandom for initializing its state.
This is a known problem for the Mersenne and Twister generator and it has been fixed using a different initialization procedure.
The TRandom3 CVS version will be now fixed using this new procedure.

Thank you very much for reporting this problem !

Lorenzo