/*
 *  dsp_rt.h
 *
 *   Only contains run-time functions suitable for inclusion in
 *   generated code.  No _S_i_m_S_t_r_u_c_t_ code is allowed in 
 *   this file. (The _'s are present because we have a test in 
 *   place to prohibit files containing that string from being 
 *   submitted to the src/rt directory!)
 *
 *  Copyright 2000 The MathWorks, Inc.
 *  $Revision: 1.5 $ $Date: 2000/05/18 21:05:05 $
 */
#ifndef dsp_rt_h
#define dsp_rt_h

#include <math.h>
#include <string.h>
#include "tmwtypes.h"
#include "simstruc_types.h"
#include "dsp_ic_rt.h" /* run-time IC handling utilities */

/*
 * ----------------------------------------------------------
 * Common constants
 * ----------------------------------------------------------
 */
#define DSP_PI        3.14159265358979323846
#define DSP_TWO_PI    6.283185307179586476925286766559005768394
/* xxx REMOVE THESE WHEN CONVENIENT */
#define PI_DOUBLE     (real64_T)3.14159265358979323846
#define TWO_PI_DOUBLE (real64_T)6.283185307179586476925286766559005768394

#define MIN_real64_T 2.225073858507201e-308
#define MAX_real64_T 1.7976931348623158e+308
#define MAX_real32_T 3.402823466e+38F
#define MAX_real_T   MAX_real64_T

#define EPS_real64_T 2.220446049250313e-016

/*
 * ----------------------------------------------------------
 * Common macros
 * ----------------------------------------------------------
 */
#ifndef MAX
#define MAX(a,b) ((a)>(b) ? (a) : (b))
#endif

#ifndef MIN
#define MIN(a,b) ((a)<(b) ? (a) : (b))
#endif


/*
 *  4 forms of complex multiply:
 */
#define CMULT_RE(X,Y)        ( (X).re * (Y).re - (X).im * (Y).im)
#define CMULT_IM(X,Y)        ( (X).re * (Y).im + (X).im * (Y).re)

#define CMULT_XCONJ_RE(X,Y)  ( (X).re * (Y).re + (X).im * (Y).im)
#define CMULT_XCONJ_IM(X,Y)  ( (X).re * (Y).im - (X).im * (Y).re)

#define CMULT_YCONJ_RE(X,Y)  ( (X).re * (Y).re + (X).im * (Y).im)
#define CMULT_YCONJ_IM(X,Y)  (-(X).re * (Y).im + (X).im * (Y).re)

#define CMULT_XYCONJ_RE(X,Y) ( (X).re * (Y).re - (X).im * (Y).im)
#define CMULT_XYCONJ_IM(X,Y) (-(X).re * (Y).im - (X).im * (Y).re)

/* Complex conjugate: */
#define CCONJ(X,Y)     \
{                      \
   (Y).re =  (X).re;   \
   (Y).im = -((X).im); \
}

/*
 *  Complex magnitude squared ( X * conj(X), or |X|^2 )
 *   CMAGSQ: arg is a complex struct
 *   CMAGSQp: arg is a pointer to a complex struct
 */
#define CMAGSQ(X)   ((X).re   * (X).re   + (X).im   * (X).im)
#define CMAGSQp(pX) ((pX)->re * (pX)->re + (pX)->im * (pX)->im)

/*
 * Quick-and-dirty (approximate) complex absolute value:
 */
#define CQABS(X) (fabs((X).re) + fabs((X).im))

/*
 *  Complex reciprocal: C = 1 / B  (A=1)
 */
#define CRECIP(B,C)                     \
{                                       \
    const real_T _s = 1.0 / CQABS(B);   \
    real_T  _d;                         \
    creal_T _bs;                        \
    _bs.re = (B).re * _s;               \
    _bs.im = (B).im * _s;               \
    _d = 1.0 / CMAGSQ(_bs);             \
    (C).re = ( _s * _bs.re) * _d;       \
    (C).im = (-_s * _bs.im) * _d;       \
}

/*
 * Complex division: C = A / B
 */
#define CDIV(A,B,C)                             \
{                                               \
    if ((B).im == 0.0) {                        \
	(C).re = (A).re / (B).re;               \
	(C).im = (A).im / (B).re;               \
    } else {                                    \
        const real_T _s = 1.0 / CQABS(B);       \
        real_T  _d;                             \
        creal_T _as, _bs;                       \
        _as.re = (A).re * _s;                   \
        _as.im = (A).im * _s;                   \
        _bs.re = (B).re * _s;                   \
        _bs.im = (B).im * _s;                   \
        _d = 1.0 / CMAGSQ(_bs);                 \
        (C).re = CMULT_YCONJ_RE(_as, _bs) * _d; \
        (C).im = CMULT_YCONJ_IM(_as, _bs) * _d; \
    }                                           \
}

/*
 *  Hypotenuse: c = sqrt(a^2 + b^2)
 */
#define CHYPOT(A,B,C)                          \
{                                              \
    if (fabs(A) > fabs(B)) {                   \
        real_T _tmp = (B)/(A);                 \
        (C) = (fabs(A)*sqrt(1+_tmp*_tmp));     \
    } else {                                   \
    if ((B) == 0.0) {                          \
            (C) = 0.0;                         \
        } else {                               \
            real_T _tmp = (A)/(B);             \
            (C) = (fabs(B)*sqrt(1+_tmp*_tmp)); \
        }                                      \
    }                                          \
}

/*
 *  Complex modulus: Y = abs(X)
 */
#define CABS(X,Y) CHYPOT((X).re, (X).im, (Y))


/*
 * ----------------------------------------------------------
 * Declare zero-crossing detector prototype:
 * ----------------------------------------------------------
 */
extern ZCEventType rt_ZCFcn(ZCDirection direction,
                            ZCSigState *prevSigState, 
                            real_T      zcSig);

#endif  /* dsp_rt_h */

/* [EOF] dsp_rt.h */
