/*
 * sqrt.c
 *
 * Code generation for function 'sqrt'
 *
 * C source code generated on: Wed Aug 26 14:59:33 2015
 *
 */

/* Include files */
#include "rt_nonfinite.h"
#include "Select_Ground_Motions.h"
#include "sqrt.h"
#include "SP_1996.h"
#include "Select_Ground_Motions_data.h"

/* Function Definitions */
void b_sqrt(real_T *x)
{
  if (*x < 0.0) {
    emlrtPushRtStackR2012b(&md_emlrtRSI, emlrtRootTLSGlobal);
    e_eml_error();
    emlrtPopRtStackR2012b(&md_emlrtRSI, emlrtRootTLSGlobal);
  }

  *x = muDoubleScalarSqrt(*x);
}

void c_sqrt(emxArray_real_T *x)
{
  int32_T k;
  int32_T i27;
  for (k = 0; k < x->size[1]; k++) {
    if (x->data[(int32_T)(1.0 + (real_T)k) - 1] < 0.0) {
      emlrtPushRtStackR2012b(&md_emlrtRSI, emlrtRootTLSGlobal);
      e_eml_error();
      emlrtPopRtStackR2012b(&md_emlrtRSI, emlrtRootTLSGlobal);
    }
  }

  i27 = x->size[1];
  for (k = 0; k < i27; k++) {
    x->data[(int32_T)(1.0 + (real_T)k) - 1] = muDoubleScalarSqrt(x->data
      [(int32_T)(1.0 + (real_T)k) - 1]);
  }
}

void d_sqrt(emxArray_real_T *x)
{
  int32_T i28;
  int32_T k;
  i28 = x->size[0] * x->size[1];
  for (k = 0; k < i28; k++) {
    if (x->data[(int32_T)(1.0 + (real_T)k) - 1] < 0.0) {
      emlrtPushRtStackR2012b(&md_emlrtRSI, emlrtRootTLSGlobal);
      e_eml_error();
      emlrtPopRtStackR2012b(&md_emlrtRSI, emlrtRootTLSGlobal);
    }
  }

  i28 = x->size[0] * x->size[1];
  for (k = 0; k < i28; k++) {
    x->data[(int32_T)(1.0 + (real_T)k) - 1] = muDoubleScalarSqrt(x->data
      [(int32_T)(1.0 + (real_T)k) - 1]);
  }
}

void e_sqrt(creal_T *x)
{
  real_T absxi;
  real_T absxr;
  if (x->im == 0.0) {
    if (x->re < 0.0) {
      absxi = 0.0;
      absxr = muDoubleScalarSqrt(muDoubleScalarAbs(x->re));
    } else {
      absxi = muDoubleScalarSqrt(x->re);
      absxr = 0.0;
    }
  } else if (x->re == 0.0) {
    if (x->im < 0.0) {
      absxi = muDoubleScalarSqrt(-x->im / 2.0);
      absxr = -absxi;
    } else {
      absxi = muDoubleScalarSqrt(x->im / 2.0);
      absxr = absxi;
    }
  } else if (muDoubleScalarIsNaN(x->re) || muDoubleScalarIsNaN(x->im)) {
    absxi = rtNaN;
    absxr = rtNaN;
  } else if (muDoubleScalarIsInf(x->im)) {
    absxi = rtInf;
    absxr = x->im;
  } else if (muDoubleScalarIsInf(x->re)) {
    if (x->re < 0.0) {
      absxi = 0.0;
      absxr = rtInf;
    } else {
      absxi = rtInf;
      absxr = 0.0;
    }
  } else {
    absxr = muDoubleScalarAbs(x->re);
    absxi = muDoubleScalarAbs(x->im);
    if ((absxr > 4.4942328371557893E+307) || (absxi > 4.4942328371557893E+307))
    {
      absxr *= 0.5;
      absxi *= 0.5;
      absxi = muDoubleScalarHypot(absxr, absxi);
      if (absxi > absxr) {
        absxi = muDoubleScalarSqrt(absxi) * muDoubleScalarSqrt(1.0 + absxr /
          absxi);
      } else {
        absxi = muDoubleScalarSqrt(absxi) * 1.4142135623730951;
      }
    } else {
      absxi = muDoubleScalarSqrt((muDoubleScalarHypot(absxr, absxi) + absxr) *
        0.5);
    }

    if (x->re > 0.0) {
      absxr = 0.5 * (x->im / absxi);
    } else {
      if (x->im < 0.0) {
        absxr = -absxi;
      } else {
        absxr = absxi;
      }

      absxi = 0.5 * (x->im / absxr);
    }
  }

  x->re = absxi;
  x->im = absxr;
}

/* End of code generation (sqrt.c) */