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

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

/* Variable Definitions */
static emlrtRSInfo nd_emlrtRSI = { 28, "baker_jayaram_correlation",
  "/home/sebenik/Projects/Ground_Motion_Selection/GMS_con_cms_1gmc_dynamicInput_501_jure/baker_jayaram_correlation.m"
};

static emlrtRSInfo od_emlrtRSI = { 37, "baker_jayaram_correlation",
  "/home/sebenik/Projects/Ground_Motion_Selection/GMS_con_cms_1gmc_dynamicInput_501_jure/baker_jayaram_correlation.m"
};

/* Function Definitions */
void b_eml_error(void)
{
  static char_T cv18[3][1] = { { 'l' }, { 'o' }, { 'g' } };

  emlrtPushRtStackR2012b(&jd_emlrtRSI, emlrtRootTLSGlobal);
  emlrtErrorWithMessageIdR2012b(emlrtRootTLSGlobal, &ee_emlrtRTEI,
    "Coder:toolbox:ElFunDomainError", 3, 4, 3, cv18);
  emlrtPopRtStackR2012b(&jd_emlrtRSI, emlrtRootTLSGlobal);
}

real_T baker_jayaram_correlation(real_T T1, real_T T2)
{
  real_T rho;
  real_T C2;
  real_T T_min;
  real_T T_max;
  real_T C3;
  real_T x;

  /*  Created by Jack Baker, 2/28/07 (updated 6/25/2007) */
  /*  Compute the correlation of epsilons for the NGA ground motion models */
  /*  */
  /*  The function is strictly emperical, fitted over the range the range 0.01s <= T1, T2 <= 10s */
  /*  */
  /*  Documentation is provided in the following document: */
  /*  Baker, J.W. and Jayaram, N., "Correlation of spectral acceleration values from NGA ground  */
  /*  motion models," Earthquake Spectra, (in review). */
  /*  INPUT */
  /*  */
  /*    T1, T2      = The two periods of interest. The periods may be equal, */
  /*                  and there is no restriction on which one is larger. */
  /*  */
  /*  INPUT */
  /*  */
  /*    rho         = The predicted correlation coefficient */
  C2 = 0.0;
  T_min = muDoubleScalarMin(T1, T2);
  T_max = muDoubleScalarMax(T1, T2);
  emlrtPushRtStackR2012b(&nd_emlrtRSI, emlrtRootTLSGlobal);
  C3 = T_max / muDoubleScalarMax(T_min, 0.109);
  if (C3 < 0.0) {
    emlrtPushRtStackR2012b(&pd_emlrtRSI, emlrtRootTLSGlobal);
    b_eml_error();
    emlrtPopRtStackR2012b(&pd_emlrtRSI, emlrtRootTLSGlobal);
  }

  x = muDoubleScalarCos(1.5707963267948966 - muDoubleScalarLog(C3) * 0.366);
  emlrtPopRtStackR2012b(&nd_emlrtRSI, emlrtRootTLSGlobal);
  if (T_max < 0.2) {
    C2 = 1.0 - 0.105 * (1.0 - 1.0 / (1.0 + muDoubleScalarExp(100.0 * T_max - 5.0)))
      * (T_max - T_min) / (T_max - 0.0099);
  }

  if (T_max < 0.109) {
    C3 = C2;
  } else {
    C3 = 1.0 - x;
  }

  emlrtPushRtStackR2012b(&od_emlrtRSI, emlrtRootTLSGlobal);
  if (C3 < 0.0) {
    emlrtPushRtStackR2012b(&md_emlrtRSI, emlrtRootTLSGlobal);
    e_eml_error();
    emlrtPopRtStackR2012b(&md_emlrtRSI, emlrtRootTLSGlobal);
  }

  C3 = (1.0 - x) + 0.5 * (muDoubleScalarSqrt(C3) - C3) * (1.0 +
    muDoubleScalarCos(3.1415926535897931 * T_min / 0.109));
  emlrtPopRtStackR2012b(&od_emlrtRSI, emlrtRootTLSGlobal);
  if (T_max <= 0.109) {
    rho = C2;
  } else if (T_min > 0.109) {
    rho = 1.0 - x;
  } else if (T_max < 0.2) {
    rho = muDoubleScalarMin(C2, C3);
  } else {
    rho = C3;
  }

  return rho;
}

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