//
// File: AMJamming_H.cpp
//
// MATLAB Coder version            : 5.2
// C/C++ source code generated on  : 09-Mar-2023 16:57:20
//

// Include Files
#include "fun8_AMJamming_H.h"
#include "fun8_FFTImplementationCallback.h"
#include "fun8_AMJamming_H_data.h"
#include "fun8_randn.h"
#include "rt_nonfinite.h"
#include "fun8_std.h"
#include "coder_array.h"
#include <cmath>
#include <iostream>
#include "fun8_eml_rand_mt19937ar_stateful.h"

// Function Definitions
//
// 信号长度
//
// Arguments    : double N
//                double f_s
//                double B_n
//                double f_0
//                double m_A
//                coder::array<double, 2U> &s
// Return Type  : void
//

void noise_am_interference_initialize()
{
    eml_rand_mt19937ar_stateful_init_8();
    isInitialized_noise_am_interference = true;
}

void AMJamming_H(double N, double f_s, double B_n, double f_0,
                 double m_A, coder::array<double, 2U> &s)
{
  coder::array<creal_T, 2U> S_n;
  coder::array<creal_T, 1U> b_S_n;
  coder::array<creal_T, 1U> yCol;
  coder::array<double, 2U> b_b;
  coder::array<double, 2U> b_costab1q;
  coder::array<double, 2U> r;
  coder::array<double, 2U> sintabinv;
  double costab1q[2];
  double N_n;
  double b;
  double re;
  double y;
  int i;
  int k;
  int pmax;
  int pmin;
  if (!isInitialized_noise_am_interference) {
    noise_am_interference_initialize();
  }
  if (std::isnan(N)) {
    s.set_size(1, 1);
    s[0] = rtNaN;
  } else if (N < 1.0) {
    s.set_size(s.size(0), 0);
  } else if (std::isinf(N) && (1.0 == N)) {
    s.set_size(1, 1);
    s[0] = rtNaN;
  } else {
    pmax = static_cast<int>(std::floor(N - 1.0));
    s.set_size(1, pmax + 1);
    for (i = 0; i <= pmax; i++) {
      s[i] = static_cast<double>(i) + 1.0;
    }
  }
  b = 1.0 / f_s;
  //  调制噪声参数
  // 射频信号要产生的带宽B_n时，噪声要产生的带宽
  N_n = std::round(N * (B_n / (2.0 * f_s)));
  // 频谱上采样点数
  costab1q[0] = 1.0;
  costab1q[1] = static_cast<int>(N);
    coder::randn_8(costab1q, r);
  pmax = static_cast<int>(N);
  b_costab1q.set_size(1, static_cast<int>(N));
  for (i = 0; i < pmax; i++) {
    b_costab1q[i] = r[i];
  }
  costab1q[0] = 1.0;
  costab1q[1] = static_cast<int>(N);
    coder::randn_8(costab1q, r);
  pmax = static_cast<int>(N);
  b_b.set_size(1, static_cast<int>(N));
  for (i = 0; i < pmax; i++) {
    b_b[i] = r[i];
  }
  S_n.set_size(1, b_costab1q.size(1));
  pmax = b_costab1q.size(1);
  for (i = 0; i < pmax; i++) {
    S_n[i].re = b_costab1q[i] + 0.0 * b_b[i];
    S_n[i].im = b_b[i];
  }
  // 高斯白噪声频谱,中值0，标准差1，维度(1，N)
  i = static_cast<int>(((N - N_n) - 1.0) + (1.0 - N_n));
  for (pmax = 0; pmax < i; pmax++) {
    pmin = static_cast<int>(N_n + static_cast<double>(pmax)) - 1;
    S_n[pmin].re = 0.0;
    S_n[pmin].im = 0.0;
  }
  if (S_n.size(1) == 0) {
    S_n.set_size(1, 0);
  } else {
    int n;
    int pow2p;
    boolean_T useRadix2;
    useRadix2 = ((S_n.size(1) & (S_n.size(1) - 1)) == 0);
    pow2p = 1;
    if (useRadix2) {
      pmax = S_n.size(1);
    } else {
      n = (S_n.size(1) + S_n.size(1)) - 1;
      pmax = 31;
      if (n <= 1) {
        pmax = 0;
      } else {
        boolean_T exitg1;
        pmin = 0;
        exitg1 = false;
        while ((!exitg1) && (pmax - pmin > 1)) {
          k = (pmin + pmax) >> 1;
          pow2p = 1 << k;
          if (pow2p == n) {
            pmax = k;
            exitg1 = true;
          } else if (pow2p > n) {
            pmax = k;
          } else {
            pmin = k;
          }
        }
      }
      pow2p = 1 << pmax;
      pmax = pow2p;
    }
    N_n = 6.2831853071795862 / static_cast<double>(pmax);
    n = pmax / 2 / 2;
    b_costab1q.set_size(1, n + 1);
    b_costab1q[0] = 1.0;
    pmax = n / 2 - 1;
    for (k = 0; k <= pmax; k++) {
      b_costab1q[k + 1] = std::cos(N_n * (static_cast<double>(k) + 1.0));
    }
    i = pmax + 2;
    pmin = n - 1;
    for (k = i; k <= pmin; k++) {
      b_costab1q[k] = std::sin(N_n * static_cast<double>(n - k));
    }
    b_costab1q[n] = 0.0;
    if (!useRadix2) {
      n = b_costab1q.size(1) - 1;
      pmax = (b_costab1q.size(1) - 1) << 1;
      r.set_size(1, pmax + 1);
      b_b.set_size(1, pmax + 1);
      r[0] = 1.0;
      b_b[0] = 0.0;
      sintabinv.set_size(1, pmax + 1);
      for (k = 0; k < n; k++) {
        sintabinv[k + 1] = b_costab1q[(n - k) - 1];
      }
      i = b_costab1q.size(1);
      for (k = i; k <= pmax; k++) {
        sintabinv[k] = b_costab1q[k - n];
      }
      for (k = 0; k < n; k++) {
        r[k + 1] = b_costab1q[k + 1];
        b_b[k + 1] = -b_costab1q[(n - k) - 1];
      }
      i = b_costab1q.size(1);
      for (k = i; k <= pmax; k++) {
        r[k] = -b_costab1q[pmax - k];
        b_b[k] = -b_costab1q[k - n];
      }
    } else {
      n = b_costab1q.size(1) - 1;
      pmax = (b_costab1q.size(1) - 1) << 1;
      r.set_size(1, pmax + 1);
      b_b.set_size(1, pmax + 1);
      r[0] = 1.0;
      b_b[0] = 0.0;
      for (k = 0; k < n; k++) {
        r[k + 1] = b_costab1q[k + 1];
        b_b[k + 1] = b_costab1q[(n - k) - 1];
      }
      i = b_costab1q.size(1);
      for (k = i; k <= pmax; k++) {
        r[k] = -b_costab1q[pmax - k];
        b_b[k] = b_costab1q[k - n];
      }
      sintabinv.set_size(1, 0);
    }
    if (useRadix2) {
      pmax = S_n.size(1);
      b_S_n = S_n.reshape(pmax);
      coder::internal::fft::fun8_FFTImplementationCallback::r2br_r2dit_trig_impl(
          b_S_n, S_n.size(1), r, b_b, yCol);
      if (yCol.size(0) > 1) {
        N_n = 1.0 / static_cast<double>(yCol.size(0));
        pmax = yCol.size(0);
        for (i = 0; i < pmax; i++) {
          yCol[i].re = N_n * yCol[i].re;
          yCol[i].im = N_n * yCol[i].im;
        }
      }
    } else {
      pmax = S_n.size(1);
      b_S_n = S_n.reshape(pmax);
      coder::internal::fft::fun8_FFTImplementationCallback::dobluesteinfft(
          b_S_n, pow2p, S_n.size(1), r, b_b, sintabinv, yCol);
    }
    pmax = S_n.size(1);
    S_n.set_size(1, pmax);
    for (i = 0; i < pmax; i++) {
      S_n[i] = yCol[i];
    }
  }
  y = coder::b_std_8(S_n);
  S_n.set_size(1, S_n.size(1));
  pmax = S_n.size(1) - 1;
  for (i = 0; i <= pmax; i++) {
    double ai;
    N_n = S_n[i].re;
    ai = S_n[i].im;
    if (ai == 0.0) {
      re = N_n / y;
      N_n = 0.0;
    } else if (N_n == 0.0) {
      re = 0.0;
      N_n = ai / y;
    } else {
      re = N_n / y;
      N_n = ai / y;
    }
    S_n[i].re = re;
    S_n[i].im = N_n;
  }
  // 归一化
  // 产生需要调制的射频噪声干扰s_0
  N_n = coder::b_std_8(S_n) / m_A;
  re = 6.2831853071795862 * f_0;
  s.set_size(1, s.size(1));
  pmax = s.size(1) - 1;
  for (i = 0; i <= pmax; i++) {
    s[i] = re * s[i] * b;
  }
  pmax = s.size(1);
  for (k = 0; k < pmax; k++) {
    s[k] = std::cos(s[k]);
  }
  s.set_size(1, S_n.size(1));
  pmax = S_n.size(1) - 1;
  for (i = 0; i <= pmax; i++) {
    s[i] = (N_n + S_n[i].re) * s[i];
  }
}

//int main()
//{
//    coder::array<double, 2U> s;
//    double N=10000;
//    double f_s=100e6;
//    double B_n=2e6;
//    double f_0=10e6;
//    double m_A=1;
//    AMJamming_H(N, f_s, B_n, f_0, m_A, s);
//    // 输出结果
//    for (int i = 0; i < s.size(1); i++) {
//        std::cout << s[i] << " ";
//    }
//    std::cout << std::endl;
//    return 0;
//}