*
* KUMAC pour montrer l'estimation min Q2 d'un polynome
*
his/cre/1d 200 'Coefficient a0' 50  0.7 0.9
his/cre/1d 201 'Coefficient a1' 50 -2.0 2.0
his/cre/1d 202 'Coefficient a2' 50 -2.0 2.0
his/cre/1d 203 'Coefficient a3' 50 -2.0 2.0
his/cre/1d 300 'chi2 minimum  ' 50  0.0 25.
his/cre/1d 301 'chi2 attendu  ' 50  0.0 25.
*                             Generer echantillon de mesures polynome
vec/cre par(4)
vec/cre  x(10) r .0 .1 .2 .3 .4 .5 .6 .7 .8 .9
vec/cre  y(10) r .0 .0 .0 .0 .0 .0 .0 .0 .0 .0
vec/cre ey(10) r .015 .015 .015 .015 .015 .015 .015 .015 .015 .015
sigma s=0.015
sigma alph0=0.80 ; sigma alph1=0.60
sigma alph2=0.40 ; sigma alph3=0.20
sigma mu=alph0+alph1*x+alph2*x**2+alph3*x**3
ne = 1000
do i=1,[ne]
   sigma r1=rndm(array(10))
   sigma r2=rndm(array(10))
   sigma y=sin(2.*pi*r1)*sqrt(-2.*log(r2))*s+mu
*                             Pour chaque membre estimer les coefficients
   vec/fit x y ey p3 'qn' 4 par
   call hfill(200,par(1),0.,1.)
   call hfill(201,par(2),0.,1.)
   call hfill(202,par(3),0.,1.)
   call hfill(203,par(4),0.,1.)
   sigma r=y-(par(1)+par(2)*x+par(3)*x**2+par(4)*x**3)
   sigma chi2=vsum(r**2)/s**2
   call hfill(300,chi2,0.,1.)
enddo
*                             Calculer distribution de chi2
sigma c=array(50,0.25#24.25)
sigma fc=[ne]*.5*0.5*(c/2.)**2.*exp(-c/2.)/gamma(3.)
his/put/con 301 fc
*                             Distributions des coefficients
zone 2 2
opt stat
set stat 1111
his/plot 200
his/plot 201
his/plot 202
his/plot 203
wait
*                             Distributions de chi2
zone 1 2
his/plot 300
his/plot 301 'cs'