#! /usr/bin/env python
#
def gennch ( df, xnonc ):

#*****************************************************************************80
#
## GENNCH generates a noncentral Chi-Square random deviate.
#
#  Discussion:
#
#    This procedure generates a random deviate from the  distribution of a
#    noncentral chisquare with DF degrees of freedom and noncentrality parameter
#    XNONC.
#
#    It uses the fact that the noncentral chisquare is the sum of a chisquare
#    deviate with DF-1 degrees of freedom plus the square of a normal
#    deviate with mean XNONC and standard deviation 1.
#
#    A subtle ambiguity arises in the original formulation:
#
#      value = genchi ( arg1 ) + ( gennor ( arg2, arg3 ) ) ^ 2
#
#    because the compiler is free to invoke either genchi or gennor
#    first, both of which alter the random number generator state,
#    resulting in two distinct possible results.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    03 September 2018
#
#  Author:
#
#    Original FORTRAN77 version by Barry Brown, James Lovato.
#    Python version by John Burkardt.
#
#  Parameters:
#
#    Input, real DF, the degrees of freedom.
#    1.0 < DF.
#
#    Input, real XNONC, the noncentrality parameter.
#    0.0 <= XNONC.
#
#    Output, real VALUE, a random deviate from the distribution.
#
  import numpy as np
  from genchi import genchi
  from gennor import gennor
  from sys import exit

  if ( df <= 1.0 ):
    print ( ' ' )
    print ( 'GENNCH - Fatal error!' )
    print ( '  DF <= 1.' )
    exit ( 'GENNCH - Fatal error!' )

  if ( xnonc < 0.0 ):
    print ( ' ' )
    print ( 'GENNCH - Fatal error!' )
    print ( '  XNONC < 0.0.' )
    exit ( 'GENNCH - Fatal error!' )

  arg1 = df - 1.0
  arg2 = np.sqrt ( xnonc )
  arg3 = 1.0

  t1 = genchi ( arg1 )
  t2 = gennor ( arg2, arg3 )
  value = t1 + t2 * t2

  return value

def gennch_test ( phrase ):

#*****************************************************************************80
#
## GENNCH_TEST tests GENNCH, which generates noncentral Chi-Square deviates.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    03 September 2018
#
#  Author:
#
#    John Burkardt
#
  import numpy as np
  from genunf import genunf
  from initialize import initialize
  from phrtsd import phrtsd
  from set_initial_seed import set_initial_seed
  from stats import stats
  from trstat import trstat

  n = 1000

  print ( '' )
  print ( 'GENNCH_TEST' )
  print ( '  GENNCH generates noncentral Chi-square deviates.' )
#
#  Initialize the generators.
#
  initialize ( )
#
#  Set the seeds based on the phrase.
#
  seed1, seed2 = phrtsd ( phrase )
#
#  Initialize all generators.
#
  set_initial_seed ( seed1, seed2 )
#
#  Select the parameters at random within a given range.
#
  low = 2.0
  high = 10.0
  df = genunf ( low, high )

  low = 0.0
  high = 2.0
  xnonc = genunf ( low, high )

  print ( '' )
  print ( '  N = %d' % ( n ) )
  print ( '' )
  print ( '  Parameters:' )
  print ( '' )
  print ( '  DF =    %g' % ( df ) )
  print ( '  XNONC = %g' % ( xnonc ))
#
#  Generate N samples.
#
  array = np.zeros ( n )
  for i in range ( 0, n ):
    array[i] = gennch ( df, xnonc )
#
#  Compute statistics on the samples.
#
  av, var, xmin, xmax = stats ( array, n )
#
#  Request expected value of statistics for this distribution.
#
  pdf = 'nch'
  param = np.array ( [ df, xnonc ] )
  avtr, vartr = trstat ( pdf, param )

  print ( '' )
  print ( '  Sample data range:          %14.6g  %14.6g' % ( xmin, xmax ) )
  print ( '  Sample mean, variance:      %14.6g  %14.6g' % ( av,   var ) )
  print ( '  Distribution mean, variance %14.6g  %14.6g' % ( avtr, vartr ) )
#
#  Terminate.
#
  print ( '' )
  print ( 'GENNCH_TEST' )
  print ( '  Normal end of execution.' )
  return

if ( __name__ == '__main__' ):
  from timestamp import timestamp
  timestamp ( )
  phrase = 'Randomizer'
  gennch_test ( phrase )
  timestamp ( )