*
* KUMAC pour montrer l'estimation par observable optimisee
*                                      MP 3/95
*
for/file 20 asym.ps
meta     20 -111
*                             Histogrammes pour la presentation
opt stat
set stat 1111
his/cre/1d 100 'Distribution en cos(theta)' 50 -1. 1.
his/cre/1d 200 'Distribution en Obs Opt   ' 50 -.5 .5
*                             Generer distribution angulaire
ne = 20000
sigma asy=.50
sigma cth=rndm(array([ne]))*2. -1.
sigma pp =((1.+cth**2)+asy*cth)/(2.+asy)
sigma pm = rndm(array([ne]))
ns = 0
vec/cre ct([ne])
do i=1,[ne]
   if (pp([i])>pm([i])) then
      ns = [ns] + 1
      vec/copy cth([i])  ct([ns])
   endif
enddo
vec/cre cc([ns]) r
vec/copy ct(1:[ns]) cc
vec/hfill cc 100
*                             Observable optimisee
sigma g=cc/(1.+cc**2)
vec/hfill g 200
*                             Presentation
zone 1 2
his/plot 100
his/plot 200
*
for/close 20