Subversion Repositories wimsdev

Rev

Rev 4161 | Go to most recent revision | Details | Last modification | View Log | RSS feed

Rev Author Line No. Line
20 reyssat 1
!if $wims_read_parm!=slib_header
2
 !goto proc
3
!endif
4
 
5
slib_title=Generation of Gaussian random data
6
slib_parms=3\
7
,m\
8
0,mu\
9
1,sigma
10
 
11
slib_author=Adeline Grelot and Bernadette PERRIN-RIOU
12
slib_out=random  data of m reals following a  normal law of parameters mu and  sigma.
13
slib_comment=Box Muller method
14
slib_example=10,0,2
15
!exit
16
 
17
:proc
18
!distribute item $wims_read_parm into  slib_m,slib_mu,slib_sigma
19
!default slib_mu=0
20
!default slib_sigma=1
21
 
22
slib_m2= $[floor($slib_m/2)]
23
slib_U1=!random 10^-20,1 repeat $[$slib_m2+1]
24
slib_U2=!random 0,1 repeat $[$slib_m2+1]
25
 
26
slib_w=!values sqrt(-2*log(x)) for x in $slib_U1
27
slib_w1=!values cos(2*Pi*x) for x in $slib_U2
28
slib_w2=!values sin(2*Pi*x) for x in $slib_U2
29
 
30
slib_V=$empty
31
!for slib_i=1 to $[$slib_m2 +1]
32
slib_v=!item  $slib_i of $slib_w
33
slib_v1=!item  $slib_i of $slib_w1
34
slib_v2=!item  $slib_i of $slib_w2
35
slib_v1=$[$slib_sigma*$slib_v1*$slib_v+$slib_mu] 
36
slib_v2=$[$slib_sigma*$slib_v2*$slib_v+$slib_mu] 
37
 
38
!if  $slib_i<=$slib_m2
39
slib_V=!append item $slib_v1 to $slib_V
40
slib_V=!append item $slib_v2 to $slib_V
41
!endif
42
!if $slib_i>$slib_m2 and   $[2*$slib_m2]!=$slib_m
43
slib_V=!append item $slib_v1 to $slib_V
44
!endif
45
!next slib_i
46
 
47
!if $slib_m >1 
48
slib_V=!shuffle $slib_V
49
!endif
50
 
51
slib_out=!trim $slib_V