/*
 * File: CalcAtmosCOESA.c
 *
 * Code generated for Simulink model 'SIL'.
 *
 * Model version                  : 1.470
 * Simulink Coder version         : 9.0 (R2018b) 24-May-2018
 * C/C++ source code generated on : Thu Oct 22 10:25:15 2020
 */

#include "rtwtypes.h"
#include "COESA76DATA_MACRO.h"
#include "COESAATMOS_MACRO.h"
#include "rt_powd_snf.h"
#include <math.h>
#include "CalcAtmosCOESA.h"

/*
 *     Using cached pressure and temperature tables, find the
 *     working interval and perform logarithmic interpolation.
 */
void CalcAtmosCOESA(const real_T *altitude, real_T *temp, real_T *pressure,
                    real_T *density, real_T *speedofsound, real_T *temperature76,
                    real_T *pressureRatio76, int_T numPoints)
{
  int_T i;
  for (i=0; i < numPoints; i++) {
    int_T bottom = 0;
    int_T top = NUM1976PTS-1;
    int_T idx;

    /* Find altitude interval using binary search
     *
     * Deal with the extreme cases first:
     *   if altitude <= altitude76[bottom] then return idx = bottom
     *   if altitude >= altitude76[top]    then return idx = top
     */
    if (altitude[i] <= altitude76[bottom]) {
      idx = bottom;
    } else if (altitude[i] >= altitude76[top]) {
      idx = NUM1976PTS-2;
    } else {
      for (;;) {
        idx = (bottom + top)/2;
        if (altitude[i] < altitude76[idx]) {
          top = idx - 1;
        } else if (altitude[i] >= altitude76[idx+1]) {
          bottom = idx + 1;
        } else {
          /* we have altitude76[idx] <= altitude[i] < altitude76[idx+1],
           * so break and just use idx
           */
          break;
        }
      }
    }

    /* Interval has been obtained, now do linear temperature
     * interpolation and log pressure interpolation.
     */
    if (tempGradient76[idx] != 0.0 ) {
      temp[i] = temperature76[idx] +
        tempGradient76[idx] * (altitude[i] - altitude76[idx]);
      pressure[i] = PRESSURE0 * pressureRatio76[idx] *
        (rt_powd_snf(temperature76[idx]/temp[i], GMR/tempGradient76[idx]));
    } else {
      temp[i] = temperature76[idx];
      pressure[i] = PRESSURE0 * pressureRatio76[idx] *
        exp((-GMR)*(altitude[i] - altitude76[idx]) / temperature76[idx]);
    }

    density[i] = pressure[i] / ((R_HAT/MOL_WT)*temp[i]);
    speedofsound[i] = sqrt(GAMMA*temp[i]*(R_HAT/MOL_WT));
  }
}

/*
 * File trailer for generated code.
 *
 * [EOF]
 */