fun8_AMJamming_H.cpp 7.2 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284
  1. //
  2. // File: AMJamming_H.cpp
  3. //
  4. // MATLAB Coder version : 5.2
  5. // C/C++ source code generated on : 09-Mar-2023 16:57:20
  6. //
  7. // Include Files
  8. #include "fun8_AMJamming_H.h"
  9. #include "fun8_FFTImplementationCallback.h"
  10. #include "fun8_AMJamming_H_data.h"
  11. #include "fun8_randn.h"
  12. #include "rt_nonfinite.h"
  13. #include "fun8_std.h"
  14. #include "coder_array.h"
  15. #include <cmath>
  16. #include <iostream>
  17. #include "fun8_eml_rand_mt19937ar_stateful.h"
  18. // Function Definitions
  19. //
  20. // 信号长度
  21. //
  22. // Arguments : double N
  23. // double f_s
  24. // double B_n
  25. // double f_0
  26. // double m_A
  27. // coder::array<double, 2U> &s
  28. // Return Type : void
  29. //
  30. void noise_am_interference_initialize()
  31. {
  32. eml_rand_mt19937ar_stateful_init_8();
  33. isInitialized_noise_am_interference = true;
  34. }
  35. void AMJamming_H(double N, double f_s, double B_n, double f_0,
  36. double m_A, coder::array<double, 2U> &s)
  37. {
  38. coder::array<creal_T, 2U> S_n;
  39. coder::array<creal_T, 1U> b_S_n;
  40. coder::array<creal_T, 1U> yCol;
  41. coder::array<double, 2U> b_b;
  42. coder::array<double, 2U> b_costab1q;
  43. coder::array<double, 2U> r;
  44. coder::array<double, 2U> sintabinv;
  45. double costab1q[2];
  46. double N_n;
  47. double b;
  48. double re;
  49. double y;
  50. int i;
  51. int k;
  52. int pmax;
  53. int pmin;
  54. if (!isInitialized_noise_am_interference) {
  55. noise_am_interference_initialize();
  56. }
  57. if (std::isnan(N)) {
  58. s.set_size(1, 1);
  59. s[0] = rtNaN;
  60. } else if (N < 1.0) {
  61. s.set_size(s.size(0), 0);
  62. } else if (std::isinf(N) && (1.0 == N)) {
  63. s.set_size(1, 1);
  64. s[0] = rtNaN;
  65. } else {
  66. pmax = static_cast<int>(std::floor(N - 1.0));
  67. s.set_size(1, pmax + 1);
  68. for (i = 0; i <= pmax; i++) {
  69. s[i] = static_cast<double>(i) + 1.0;
  70. }
  71. }
  72. b = 1.0 / f_s;
  73. // 调制噪声参数
  74. // 射频信号要产生的带宽B_n时,噪声要产生的带宽
  75. N_n = std::round(N * (B_n / (2.0 * f_s)));
  76. // 频谱上采样点数
  77. costab1q[0] = 1.0;
  78. costab1q[1] = static_cast<int>(N);
  79. coder::randn_8(costab1q, r);
  80. pmax = static_cast<int>(N);
  81. b_costab1q.set_size(1, static_cast<int>(N));
  82. for (i = 0; i < pmax; i++) {
  83. b_costab1q[i] = r[i];
  84. }
  85. costab1q[0] = 1.0;
  86. costab1q[1] = static_cast<int>(N);
  87. coder::randn_8(costab1q, r);
  88. pmax = static_cast<int>(N);
  89. b_b.set_size(1, static_cast<int>(N));
  90. for (i = 0; i < pmax; i++) {
  91. b_b[i] = r[i];
  92. }
  93. S_n.set_size(1, b_costab1q.size(1));
  94. pmax = b_costab1q.size(1);
  95. for (i = 0; i < pmax; i++) {
  96. S_n[i].re = b_costab1q[i] + 0.0 * b_b[i];
  97. S_n[i].im = b_b[i];
  98. }
  99. // 高斯白噪声频谱,中值0,标准差1,维度(1,N)
  100. i = static_cast<int>(((N - N_n) - 1.0) + (1.0 - N_n));
  101. for (pmax = 0; pmax < i; pmax++) {
  102. pmin = static_cast<int>(N_n + static_cast<double>(pmax)) - 1;
  103. S_n[pmin].re = 0.0;
  104. S_n[pmin].im = 0.0;
  105. }
  106. if (S_n.size(1) == 0) {
  107. S_n.set_size(1, 0);
  108. } else {
  109. int n;
  110. int pow2p;
  111. boolean_T useRadix2;
  112. useRadix2 = ((S_n.size(1) & (S_n.size(1) - 1)) == 0);
  113. pow2p = 1;
  114. if (useRadix2) {
  115. pmax = S_n.size(1);
  116. } else {
  117. n = (S_n.size(1) + S_n.size(1)) - 1;
  118. pmax = 31;
  119. if (n <= 1) {
  120. pmax = 0;
  121. } else {
  122. boolean_T exitg1;
  123. pmin = 0;
  124. exitg1 = false;
  125. while ((!exitg1) && (pmax - pmin > 1)) {
  126. k = (pmin + pmax) >> 1;
  127. pow2p = 1 << k;
  128. if (pow2p == n) {
  129. pmax = k;
  130. exitg1 = true;
  131. } else if (pow2p > n) {
  132. pmax = k;
  133. } else {
  134. pmin = k;
  135. }
  136. }
  137. }
  138. pow2p = 1 << pmax;
  139. pmax = pow2p;
  140. }
  141. N_n = 6.2831853071795862 / static_cast<double>(pmax);
  142. n = pmax / 2 / 2;
  143. b_costab1q.set_size(1, n + 1);
  144. b_costab1q[0] = 1.0;
  145. pmax = n / 2 - 1;
  146. for (k = 0; k <= pmax; k++) {
  147. b_costab1q[k + 1] = std::cos(N_n * (static_cast<double>(k) + 1.0));
  148. }
  149. i = pmax + 2;
  150. pmin = n - 1;
  151. for (k = i; k <= pmin; k++) {
  152. b_costab1q[k] = std::sin(N_n * static_cast<double>(n - k));
  153. }
  154. b_costab1q[n] = 0.0;
  155. if (!useRadix2) {
  156. n = b_costab1q.size(1) - 1;
  157. pmax = (b_costab1q.size(1) - 1) << 1;
  158. r.set_size(1, pmax + 1);
  159. b_b.set_size(1, pmax + 1);
  160. r[0] = 1.0;
  161. b_b[0] = 0.0;
  162. sintabinv.set_size(1, pmax + 1);
  163. for (k = 0; k < n; k++) {
  164. sintabinv[k + 1] = b_costab1q[(n - k) - 1];
  165. }
  166. i = b_costab1q.size(1);
  167. for (k = i; k <= pmax; k++) {
  168. sintabinv[k] = b_costab1q[k - n];
  169. }
  170. for (k = 0; k < n; k++) {
  171. r[k + 1] = b_costab1q[k + 1];
  172. b_b[k + 1] = -b_costab1q[(n - k) - 1];
  173. }
  174. i = b_costab1q.size(1);
  175. for (k = i; k <= pmax; k++) {
  176. r[k] = -b_costab1q[pmax - k];
  177. b_b[k] = -b_costab1q[k - n];
  178. }
  179. } else {
  180. n = b_costab1q.size(1) - 1;
  181. pmax = (b_costab1q.size(1) - 1) << 1;
  182. r.set_size(1, pmax + 1);
  183. b_b.set_size(1, pmax + 1);
  184. r[0] = 1.0;
  185. b_b[0] = 0.0;
  186. for (k = 0; k < n; k++) {
  187. r[k + 1] = b_costab1q[k + 1];
  188. b_b[k + 1] = b_costab1q[(n - k) - 1];
  189. }
  190. i = b_costab1q.size(1);
  191. for (k = i; k <= pmax; k++) {
  192. r[k] = -b_costab1q[pmax - k];
  193. b_b[k] = b_costab1q[k - n];
  194. }
  195. sintabinv.set_size(1, 0);
  196. }
  197. if (useRadix2) {
  198. pmax = S_n.size(1);
  199. b_S_n = S_n.reshape(pmax);
  200. coder::internal::fft::fun8_FFTImplementationCallback::r2br_r2dit_trig_impl(
  201. b_S_n, S_n.size(1), r, b_b, yCol);
  202. if (yCol.size(0) > 1) {
  203. N_n = 1.0 / static_cast<double>(yCol.size(0));
  204. pmax = yCol.size(0);
  205. for (i = 0; i < pmax; i++) {
  206. yCol[i].re = N_n * yCol[i].re;
  207. yCol[i].im = N_n * yCol[i].im;
  208. }
  209. }
  210. } else {
  211. pmax = S_n.size(1);
  212. b_S_n = S_n.reshape(pmax);
  213. coder::internal::fft::fun8_FFTImplementationCallback::dobluesteinfft(
  214. b_S_n, pow2p, S_n.size(1), r, b_b, sintabinv, yCol);
  215. }
  216. pmax = S_n.size(1);
  217. S_n.set_size(1, pmax);
  218. for (i = 0; i < pmax; i++) {
  219. S_n[i] = yCol[i];
  220. }
  221. }
  222. y = coder::b_std_8(S_n);
  223. S_n.set_size(1, S_n.size(1));
  224. pmax = S_n.size(1) - 1;
  225. for (i = 0; i <= pmax; i++) {
  226. double ai;
  227. N_n = S_n[i].re;
  228. ai = S_n[i].im;
  229. if (ai == 0.0) {
  230. re = N_n / y;
  231. N_n = 0.0;
  232. } else if (N_n == 0.0) {
  233. re = 0.0;
  234. N_n = ai / y;
  235. } else {
  236. re = N_n / y;
  237. N_n = ai / y;
  238. }
  239. S_n[i].re = re;
  240. S_n[i].im = N_n;
  241. }
  242. // 归一化
  243. // 产生需要调制的射频噪声干扰s_0
  244. N_n = coder::b_std_8(S_n) / m_A;
  245. re = 6.2831853071795862 * f_0;
  246. s.set_size(1, s.size(1));
  247. pmax = s.size(1) - 1;
  248. for (i = 0; i <= pmax; i++) {
  249. s[i] = re * s[i] * b;
  250. }
  251. pmax = s.size(1);
  252. for (k = 0; k < pmax; k++) {
  253. s[k] = std::cos(s[k]);
  254. }
  255. s.set_size(1, S_n.size(1));
  256. pmax = S_n.size(1) - 1;
  257. for (i = 0; i <= pmax; i++) {
  258. s[i] = (N_n + S_n[i].re) * s[i];
  259. }
  260. }
  261. //int main()
  262. //{
  263. // coder::array<double, 2U> s;
  264. // double N=10000;
  265. // double f_s=100e6;
  266. // double B_n=2e6;
  267. // double f_0=10e6;
  268. // double m_A=1;
  269. // AMJamming_H(N, f_s, B_n, f_0, m_A, s);
  270. // // 输出结果
  271. // for (int i = 0; i < s.size(1); i++) {
  272. // std::cout << s[i] << " ";
  273. // }
  274. // std::cout << std::endl;
  275. // return 0;
  276. //}