/*****************************************************************************
*                                                                            *
*   DEMODFM.SPL  Copyright (C) 2000-2013 DSP Development Corporation         *
*                               All Rights Reserved                          *
*                                                                            *
*   Author:      Randy Race                                                  *
*                                                                            *
*   Synopsis:    FM Demodulation using the Hilbert Transform                 *
*                                                                            *
*   Revisions:    8 Feb 2000  RRR  Creation                                  *
*                22 May 2013  RRR  optional fmin/fmax modulation range       *
*                                                                            *
*****************************************************************************/

#if @HELP_DEMODFM

    DEMODFM

    Purpose: Demodulates an FM series using the Hilbert Transform

    Syntax:  DEMODFM(s, fmin, fmax)

                 s - the input series

              fmin - Optional. A real, the minimum modulation frequency used
                     to scale the result. If not specified, no output scaling
                     is performed.

              fmax - Optional. A real, the maximum modulation frequency used
                     to scale the result. If not specified, no output scaling
                     is performed.

    Returns: A series or array

    Example:
             W1: gtriwave(1000,.001, 4)
             W2: W1 * 100 + 20
             W3: cos(integ(2*pi*w2))
             W4: demodfm(w3)

             W2 represents the scaled information signal and W3 is the
             resulting frequency modulated signal. The amplitude of W2
             determines the instantaneous frequency of W3. The instantaneous
             frequency of W3 ranges from:

                             min(w3) == 20 Hz

                                     to

                             max(w3) == 120 Hz


             W4 is the demodulated result.

    Example:
             W1: gtriwave(1000, .001, 4)
             W2: modfm(w1, 20, 100)
             W3: demodfm(w2)
             W4: demodfm(w2, 20, 100)
             
             Same as above except the minimum and maximum modulation
             frequencies are specified in W4 to scale the result.

    Remarks:
             DEMODFM uses HILB to calculate the Hilbert Transform.

             See MODFM to frequency modulate a series.

             See DEMODAM to demodulate an amplitude modulated series.

    See Also:
             Demodam
             Hilb
             Modfm

#endif


/* demodulate an FM series using the Hilbert transform */
ITERATE demodfm(s, fmin, fmax)
{
        local h, fm;

        if (argc < 3)
        {
                if (argc < 2)
                {
                        if (argc < 1) error("demodfm - input series required");
                        
                        /* no scaling */
                        fmin = -1;
                }

                fmax = fmin;
        }

        /* hilbert transform */
        h  = hilb(s) * exp(-i * 2 * pi * rate(s) * xvals(s));
        fm = deriv(unwrap(phase(h))) / (2 * pi);

        /* amplitude scaling */
        if (fmin > 0 && fmax > 0 && not(fmin == fmax))
        {
                fm = rescale(fm, min(s), max(s), fmin, fmax);
        }

        /* restore the offset */
        setxoffset(fm, xoffset(s));

        return(fm);
}