GNU Octave  3.8.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
gennch.f
Go to the documentation of this file.
1  REAL FUNCTION gennch(df,xnonc)
2 C**********************************************************************
3 C
4 C REAL FUNCTION GENNCH( DF, XNONC )
5 C Generate random value of Noncentral CHIsquare variable
6 C
7 C
8 C Function
9 C
10 C
11 
12 C Generates random deviate from the distribution of a noncentral
13 C chisquare with DF degrees of freedom and noncentrality parameter
14 C XNONC.
15 C
16 C
17 C Arguments
18 C
19 C
20 C DF --> Degrees of freedom of the chisquare
21 C (Must be >= 1.0)
22 C REAL DF
23 C
24 C XNONC --> Noncentrality parameter of the chisquare
25 C (Must be >= 0.0)
26 C REAL XNONC
27 C
28 C
29 C Method
30 C
31 C
32 C Uses fact that noncentral chisquare is the sum of a chisquare
33 C deviate with DF-1 degrees of freedom plus the square of a normal
34 C deviate with mean sqrt(XNONC) and standard deviation 1.
35 C
36 C**********************************************************************
37 C .. Scalar Arguments ..
38  REAL df,xnonc
39 C ..
40 C .. External Functions ..
41 C JJV changed these to call SGAMMA and SNORM directly
42 C REAL genchi,gennor
43 C EXTERNAL genchi,gennor
44  REAL sgamma,snorm
45  EXTERNAL sgamma,snorm
46 C ..
47 C .. Intrinsic Functions ..
48  INTRINSIC sqrt
49 C ..
50 C JJV changed abort to df < 1, and added case: df = 1
51 C .. Executable Statements ..
52  IF (.NOT. (df.LT.1.0.OR.xnonc.LT.0.0)) go to 10
53  WRITE (*,*) 'DF < 1 or XNONC < 0 in GENNCH - ABORT'
54  WRITE (*,*) 'Value of DF: ',df,' Value of XNONC',xnonc
55  CALL xstopx('DF < 1 or XNONC < 0 in GENNCH - ABORT')
56 
57 C JJV changed this to call SGAMMA and SNORM directly
58 C gennch = genchi(df-1.0) + gennor(sqrt(xnonc),1.0)**2
59 
60  10 IF (df.GE.1.000001) go to 20
61 C JJV case DF = 1.0
62  gennch = (snorm() + sqrt(xnonc))**2
63  go to 30
64 
65 C JJV case DF > 1.0
66  20 gennch = 2.0*sgamma((df-1.0)/2.0) + (snorm() + sqrt(xnonc))**2
67  30 RETURN
68 
69  END