/*****************************************************************************
*                                                                            *
*   BESSELH.SPL  Copyright (C) 2009 DSP Development Corporation              *
*                               All Rights Reserved                          *
*                                                                            *
*   Author:      Randy Race                                                  *
*                                                                            *
*   Synopsis:    Hv(z) bessel function of fractional order                   *
*                                                                            *
*   Revisions:    8 Nov 2009  RRR  Creation                                  *
*                                                                            *
*****************************************************************************/

#include 

#if @HELP_BESSELH

    BESSELH

    Purpose: Evaluates the Hankel function Hv(z) for any real order.

    Syntax:  BESSELH(v, K, z, opt)

                 v - A real or real series, the order. The order must
                     be real but need not be an integer.

                 K - Optional, an integer. Either 1 or 2, the type of
                     Hankel function. Defaults to 1.

                 z - Any scalar or series. The input value, can be complex.

               opt - Optional. An integer, the scaling method.

                       0: no scaling (default)
                       1: scale by exp(-i*z) for K == 1
                          scale by exp( i*z) for K == 2

    Returns:
                                                (K)
             A scalar or series, the value of Hv   (z), where v is the order,
             z is the input and K is either 1 or 2.

    Example:
             besselh(0, 3)
        
                                                           (1)
             returns -0.260052 + 0.376850i, the value of H0   (3).

    Example:
             W1: 0..0.2..1
             W2: besselh(1, w1)

                         (1)
             Evaluates H1   (z) for z between 0 and 1. W2 contains the 
             complex series:

             {NA, 
              0.099501 - 3.323825i,
              0.196027 - 1.780872i,
              0.286701 - 1.260391i,
              0.368842 - 0.978144i,
              0.440051 - 0.781213i}

    Example:
             besselh(3..9, 0..0.2..10)

             Evaluates the Hankel function (Bessel function of the
             third kind for orders 3 through 9 with inputs from 0 to
             10. Each column of the result contains the output for the
             specified order.

    Remarks:
             Bessel functions are solutions to the differential equation:

                 2
              2 d y     dy     2   2
             z  --- + z -- + (z - v ) y = 0
                  2
                dz      dz

             Jv(z) as a solution of the first kind. Yv(z) is a linearly 
             independent solution of the second kind.

                                   (K)
             The Hankel functions H   v(z) are related to the Bessel functions
             with:

               (1)
             Hv   (z) = Jv(z) + i Yz(v)

               (2)
             Hv   (z) = Jv(z) - i Yz(v)

             Hankel functions are also known as Bessel functions of the
             third kind.


             BESSELH is based on a FORTRAN library written by Donald E. Amos.

   See Also:
             Besselh
             Besseli
             Besselk
             Bessely
             Jn
             Yn

 References:
             [1] Abramowitz and Stegun
                 Handbook of Mathematical Functions (9th printing 1970)
                 US Gov. Printing Office
                 Section 9.1.1, 9.1.89, 9.12

             [2] Amos, D.E.
                 A Subroutine Package for Bessel Functions of a Complex 
                 Argument and Nonnegative Order
                 Sandia National Laboratory Report
                 SAND85-1018, May, 1985. 

             [3] Amos, D.E. 
                 A Portable Package for Bessel Functions of a Complex
                 Argument and Nonnegative Order
                 Trans. Math. Software, 1986

#endif


/* bessel function of the first kind */
ITERATE besselh(v, K, z, opt)
{
        local y, scalarval;

        (v, K, z, opt) = besselh_parse_args(v, K, z, opt);

        /* is input a single value */
        scalarval = not(isarray(v)) && not(isarray(z));

        /* use built-in */
        y = besselx("h", {z}, {v}, opt);

        if (scalarval) y = y[1];

        return(y);
}


besselh_parse_args(v, K, z, opt)
{
        if (argc < 4)
        {
                opt = 0;
                
                if (argc < 3)
                {
                        z = K;
                        K = 1;
                        
                        if (argc < 2)
                        {
                                if (argc < 1) error("besselh - input required");
                                
                                z = v;
                                v = 0;
                        }

                }
        }
        
        return(v, K, z, opt);
}