Please Note:
Our new Windows program Aw-Radw includes an advanced random number generator that includes the option to record elapsed time between hits i.e. the program will store the number of seconds between events in millionths of a second. Aw-Radw includes random number generator

Lacking Aw-Radw or a Windows PC, Aware Electronics customers can request a copy of our dos program r-time.exe which Professor Coisson uses in the following experiment. We will E-Mail it as an attached file.


Parma University, Physics Department

http://www.fis.unipr.it/

"Acqusition and treatment of data from a Poisson process"

R. Coisson (a), T.N. Ha (b)

(a) Parma University, (b) Institute of Materials Science, NCST Vietnam.

abstract

Poisson processes, i.e. processes with events happening at random times, have many applications in Physics. With the use of a simple Geiger counter with direct acquisition on a computer, we show how to test the validity of various theoretical relations on the Poisson statistics, with simple programs in the Matlab programming language.

(This experiment was performed in the framework of a collaboration on Physics laboratory teaching between Parma University, Italy, and Hanoi National University, Vietnam (Project co-financed by the Italian Ministry of Foreign Affairs (DGPC Uff. V), International Centre for Theoretical Physics (Trieste), Parma University and Hanoi National University)).

The Poisson process is a random function which for example can model a sequence of pulses arriving at random times. The statistics of this process is a subject that is taught in Statistics courses.

When statistical processes are taught at the blackboard, the students do not get an operational and intuitive understanding of their characteristics. For this we have elaborated a simple experiment using a small Geiger counter interfaced with a computer, where the students are invited to solve some problem and to compare the observed count distributions with the values expected from the theory.

The students are not given ready-made programs and formulas (they can use the textbooks on which they have studied Statistics), as only in this way they learn to be autonomous.

The results of theory that can be tested are for ex. the distribution of time intervals, which is an exponential with characteristic time equal to the average value of the intervals. And then the distribution of counts in a fixed time interval, which is (expressed as probability of counting n when the average number in that interval is m):

p(n,m)=m^n*exp(-m)/n!

But once the students have understood the method, they can for ex. calculate and measure the distribution of intervals between each count and the second one after it; as well as other exercises that the teacher may propose.

The Geiger counter is a small device (RM-60 of Aware Electronics, http://www.aw-el.com )

bought for 150 Euro, which is connected to a standard RS-232 bus (and powered by it). Counts are provided either by the dark counts, or by a net for gas lamps (containing a small amount of Th).

A program allows to create a file composed of a list of numbers, which are the times (in seconds) between each pulse and the following one.

This file can be open with the Matlab or Scilab programming language, so calculations on these data can be done, as well as plotting the results of data elaboration and of theoretical expectation.

Use can be made of the function "hist', which calculates the histogram of a set of data for a given number of channels. But the students might also be invited to write a program of their own to calculate histograms, an exercise which helps the understanding of the procedure and of programming.

When the histogram of the data is done, it should appear exponentially distributed.

Actually it does not always happen so, as sometimes there are some wrong, unreasonable numbers; and these should be cancelled, for example by a command to erase all numbers greater than a value 30 times the average.

For the comparison with theory, the two sets of data (experimental and theoretical) should be normalised to the same quantity: they can be either normalised to one (having then the meaning of probabilities) or to the total number N of counts (having then the meaning of expected or observed numbers).

The students should then write the theoretical exponential with the correct normalisation, and normalise the experimental data.

Then, if we call 't' (with average tm=mean(t);)the data and [a,b]=hist(t,k); is the command to calculate the histogram with k bins, a command to print the resulkt can be as follows:

plot(b,a,'o',b,norm*exp(-t/tm))

where 'norm' is the suitable calculated normalisation factor.

Analogous, but a bit more complicated, is the procedure to calculate the count distribution in fixed time intervals.

One first calculates for each count the time elapsed since the beginning:

T=cumsum(t);

Then [c,d]=hist(T,h); gives the numbers of counts in successive time intervals, with an average of m=N/h counts per interval.

then [f,g]=hist(d,0:max(c)) gives the Poisson distribution, to be compared, after suitable normalisation, with eq. 1.

In the following we list a possible program with the Matlab syntax, and the plots showing the results.


%
%
% EXAMPLE OF MATLAB PROGRAM TO TEST A POISSON PROCESS
%
%
% geiger.m
% takes time intervals from a random pulse sequence
% in file 'poisson.dat' created by r-time.exe (change for other filename!)
% compares time-interval statistics with theory (figure 1)
% compares distribution of counts in a fixed time interval with theory
% (Poisson statistics, figure 2)


load poisson.dat; % load the file "poisson.dat" (3 columns)
t=poisson(:,1); % we take the first column only
t(find(t>100))=[]; % erase some wrong data
lt=length(t);
mt=mean(t); % mean interval between successive counts
[a,b]=hist(t,sqrt(lt));
th=exp(-b/mt);
th=th/sum(th); % theoretical prediction for 'a'
figure(1)
plot(b,a/sum(a),'*',b,th); % compare time interval distrib. with theory
figure(2)
T=cumsum(t); % times of arrival of counts
m=3; % choose a time interval such that m=3
[c,d]=hist(T,lt/m); % counts in successive time intervals
ma=max(c);
[f,g]=hist(c,0:ma+1);
for k=1:4*m+1,
poi(k)=m^(k-1)*exp(-m)/factorial(k-1);
end % theoretical Poisson distribution
plot(g,f/sum(f),'o',g,poi); % test of Poisson distribution
interv.jpg
poisson.jpg