1      REAL FUNCTION gennch(df,xnonc)
2C**********************************************************************
3C
4C     REAL FUNCTION GENNCH( DF, XNONC )
5C           Generate random value of Noncentral CHIsquare variable
6C
7C
8C                              Function
9C
10C
11
12C     Generates random deviate  from the  distribution  of a  noncentral
13C     chisquare with DF degrees  of freedom and noncentrality  parameter
14C     XNONC.
15C
16C
17C                              Arguments
18C
19C
20C     DF --> Degrees of freedom of the chisquare
21C            (Must be >= 1.0)
22C                         REAL DF
23C
24C     XNONC --> Noncentrality parameter of the chisquare
25C               (Must be >= 0.0)
26C                         REAL XNONC
27C
28C
29C                              Method
30C
31C
32C     Uses fact that  noncentral chisquare  is  the  sum of a  chisquare
33C     deviate with DF-1  degrees of freedom plus the  square of a normal
34C     deviate with mean sqrt(XNONC) and standard deviation 1.
35C
36C**********************************************************************
37C     .. Scalar Arguments ..
38      REAL df,xnonc
39C     ..
40C     .. External Functions ..
41C     JJV changed these to call SGAMMA and SNORM directly
42C      REAL genchi,gennor
43C      EXTERNAL genchi,gennor
44      REAL sgamma,snorm
45      EXTERNAL sgamma,snorm
46C     ..
47C     .. Intrinsic Functions ..
48      INTRINSIC sqrt
49C     ..
50C     JJV changed abort to df < 1, and added case: df = 1
51C     .. 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
57C     JJV changed this to call SGAMMA and SNORM directly
58C      gennch = genchi(df-1.0) + gennor(sqrt(xnonc),1.0)**2
59
60 10   IF (df.GE.1.000001) GO TO 20
61C     JJV case DF = 1.0
62      gennch = (snorm() + sqrt(xnonc))**2
63      GO TO 30
64
65C     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
70