You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

1111 lines
24 KiB

  1. #include <math.h>
  2. #include "assemblage_includes.h"
  3. #ifdef USE_FLOAT
  4. #define MATH_COS cosf
  5. #define MATH_SIN sinf
  6. #else
  7. #define MATH_COS cos
  8. #define MATH_SIN sin
  9. #endif
  10. /*
  11. * The following include file is generated by the Prelude compiler
  12. * Theoretically we should include main node specific include
  13. * i.e. assemblage_vX.h but we know that every main node
  14. * assemblage, assemblage_v2, assemblage_v3, etc... share the
  15. * very same data type.
  16. */
  17. //#include "assemblage.h"
  18. /* Period/Frequency of the nodes */
  19. // const REAL_TYPE Ts_h = 1.0/50.0;
  20. // const REAL_TYPE Ts_K1 = 1.0/50.0;
  21. // const REAL_TYPE Ts_K2 = 1.0/50.0;
  22. // const REAL_TYPE Ts_f_Va = 1.0/100.0;
  23. // const REAL_TYPE Ts_f_Vz = 1.0/100.0;
  24. // const REAL_TYPE Ts_f_q = 1.0/100.0;
  25. // const REAL_TYPE Ts_f_az = 1.0/100.0;
  26. // const REAL_TYPE Ts_f_h = 1.0/100.0;
  27. const REAL_TYPE dt = 1.0f/200.0;
  28. const REAL_TYPE dt_de = 1.0/200.0;
  29. const REAL_TYPE dt_dx = 1.0/200.0;
  30. /* Controller parameters */
  31. /* Altitude hold */
  32. const REAL_TYPE Kp_h = 0.1014048;
  33. const REAL_TYPE Ki_h = 0.0048288;
  34. const REAL_TYPE h_switch = 50.0;
  35. // Setpoint commands
  36. REAL_TYPE Vz_c = -2.5;
  37. REAL_TYPE Va_c = 0.0;
  38. REAL_TYPE h_c = 10000;
  39. /* Va Speed controller */
  40. // const REAL_TYPE K1_intVa = 0.06018;
  41. // const REAL_TYPE K1_Va = -0.53115;
  42. // const REAL_TYPE K1_Vz = -0.08956;
  43. // const REAL_TYPE K1_q = 24.44890;
  44. const REAL_TYPE K1_intVa = 0.049802610664357;
  45. const REAL_TYPE K1_Va = -0.486813084356079;
  46. const REAL_TYPE K1_Vz = -0.077603095495388;
  47. const REAL_TYPE K1_q = 21.692383376322041;
  48. /* Vz Speed controller */
  49. // const REAL_TYPE K2_intVz = 0.0006545;
  50. // const REAL_TYPE K2_Vz = -0.0031107;
  51. // const REAL_TYPE K2_q = 0.4490749;
  52. // const REAL_TYPE K2_az = -0.0002038;
  53. const REAL_TYPE K2_intVz = 0.000627342822264;
  54. const REAL_TYPE K2_Vz = -0.003252836726554;
  55. const REAL_TYPE K2_q = 0.376071446897134;
  56. const REAL_TYPE K2_az = -0.001566907423747;
  57. /* Trimming parameters */
  58. const REAL_TYPE h_eq = 10000.0;
  59. const REAL_TYPE Va_eq = 230.0;
  60. const REAL_TYPE Vz_eq = 0.0;
  61. const REAL_TYPE alpha_eq = 0.026485847681737;
  62. const REAL_TYPE theta_eq = 0.026485847681737;
  63. /* Atmosphere parameters */
  64. const REAL_TYPE rho0 = 1.225;
  65. const REAL_TYPE g0 = 9.80665;
  66. const REAL_TYPE T0_0 = 288.15;
  67. const REAL_TYPE T0_h = -0.0065;
  68. const REAL_TYPE Rs = 287.05;
  69. /* Aircraft parameters */
  70. const REAL_TYPE masse = 57837.5;
  71. const REAL_TYPE I_y = 3781272.0;
  72. const REAL_TYPE S = 122.6;
  73. const REAL_TYPE cbar = 4.29;
  74. const REAL_TYPE CD_0 = 0.016;
  75. const REAL_TYPE CD_alpha = 2.5;
  76. const REAL_TYPE CD_deltae = 0.05;
  77. const REAL_TYPE CL_alpha = 5.5;
  78. const REAL_TYPE CL_deltae = 0.193;
  79. const REAL_TYPE alpha_0 = -0.05;
  80. const REAL_TYPE Cm_0 = 0.04;
  81. const REAL_TYPE Cm_alpha = -0.83;
  82. const REAL_TYPE Cm_deltae = -1.5;
  83. const REAL_TYPE Cm_q = -30;
  84. /* in-memory buffer for traces */
  85. //REAL_TYPE sample[SPL_SIZE][NBMAX_SAMPLE];
  86. //static unsigned long instant = 0;
  87. //static unsigned long sample_instant=0;
  88. #define FMTFLOAT "%5.15f"
  89. /* Va filter 100 Hz */
  90. REAL_TYPE
  91. Va_filter_100(REAL_TYPE Va) {
  92. static REAL_TYPE y = 0.0;
  93. static REAL_TYPE x1 = 0.0;
  94. static REAL_TYPE x2 = 0.0;
  95. static REAL_TYPE x1_tmp = 0.0;
  96. static REAL_TYPE x2_tmp = 0.0;
  97. static unsigned short debut = 1;
  98. /* 100 Hz coefficients */
  99. static REAL_TYPE a0 = 0.956543675476034;
  100. static REAL_TYPE a1 = -1.955578398054313;
  101. static REAL_TYPE b0 = 0.000479064865372430;
  102. static REAL_TYPE b1 = 0.000486212556348925;
  103. if (debut) {
  104. debut = 0;
  105. x1 = Va_eq * (1.0 + a1 - b1);
  106. x2 = Va_eq;
  107. }
  108. // Output
  109. y = x2;
  110. // State
  111. x1_tmp = - a0 * x2 + b0 * Va;
  112. x2_tmp = x1 - a1 * x2 + b1 * Va;
  113. // Update
  114. x1 = x1_tmp;
  115. x2 = x2_tmp;
  116. return y;
  117. } /* end of Va filter 100 Hz */
  118. /* Va filter 50 Hz */
  119. REAL_TYPE
  120. Va_filter_50(REAL_TYPE Va) {
  121. static REAL_TYPE y = 0.0;
  122. static REAL_TYPE x1 = 0.0;
  123. static REAL_TYPE x2 = 0.0;
  124. static REAL_TYPE x1_tmp = 0.0;
  125. static REAL_TYPE x2_tmp = 0.0;
  126. static unsigned short debut = 1;
  127. /* 50 Hz coefficients */
  128. static REAL_TYPE a0 = 0.914975803093201;
  129. static REAL_TYPE a1 = -1.911199519984605;
  130. static REAL_TYPE b0 = 0.001860178914816;
  131. static REAL_TYPE b1 = 0.001916104193780;
  132. if (debut) {
  133. debut = 0;
  134. x1 = Va_eq * (1.0 + a1 - b1);
  135. x2 = Va_eq;
  136. }
  137. // Output
  138. y = x2;
  139. // State
  140. x1_tmp = - a0 * x2 + b0 * Va;
  141. x2_tmp = x1 - a1 * x2 + b1 * Va;
  142. // Update
  143. x1 = x1_tmp;
  144. x2 = x2_tmp;
  145. return y;
  146. } /* end of Va filter 50 Hz */
  147. /* Va filter 33 Hz */
  148. REAL_TYPE
  149. Va_filter_33(REAL_TYPE Va) {
  150. static REAL_TYPE y = 0.0;
  151. static REAL_TYPE x1 = 0.0;
  152. static REAL_TYPE x2 = 0.0;
  153. static REAL_TYPE x1_tmp = 0.0;
  154. static REAL_TYPE x2_tmp = 0.0;
  155. static unsigned short debut = 1;
  156. /* 33 Hz coefficients */
  157. static REAL_TYPE a0 = 0.874036784828483;
  158. static REAL_TYPE a1 = -1.865563793814790;
  159. static REAL_TYPE b0 = 0.004141433623051;
  160. static REAL_TYPE b1 = 0.004331557390642;
  161. if (debut) {
  162. debut = 0;
  163. x1 = Va_eq * (1.0 + a1 - b1);
  164. x2 = Va_eq;
  165. }
  166. // Output
  167. y = x2;
  168. // State
  169. x1_tmp = - a0 * x2 + b0 * Va;
  170. x2_tmp = x1 - a1 * x2 + b1 * Va;
  171. // Update
  172. x1 = x1_tmp;
  173. x2 = x2_tmp;
  174. return y;
  175. } /* end of Va filter 33 Hz */
  176. /* Va filter 25 Hz */
  177. REAL_TYPE
  178. Va_filter_25(REAL_TYPE Va) {
  179. static REAL_TYPE y = 0.0;
  180. static REAL_TYPE x1 = 0.0;
  181. static REAL_TYPE x2 = 0.0;
  182. static REAL_TYPE x1_tmp = 0.0;
  183. static REAL_TYPE x2_tmp = 0.0;
  184. static unsigned short debut = 1;
  185. /* 25 Hz coefficients */
  186. static REAL_TYPE a0 = 0.837180720246048;
  187. static REAL_TYPE a1 = -1.822731999002980;
  188. static REAL_TYPE b0 = 0.007010380719078;
  189. static REAL_TYPE b1 = 0.007438340523990;
  190. if (debut) {
  191. debut = 0;
  192. x1 = Va_eq * (1.0 + a1 - b1);
  193. x2 = Va_eq;
  194. }
  195. // Output
  196. y = x2;
  197. // State
  198. x1_tmp = - a0 * x2 + b0 * Va;
  199. x2_tmp = x1 - a1 * x2 + b1 * Va;
  200. // Update
  201. x1 = x1_tmp;
  202. x2 = x2_tmp;
  203. return y;
  204. } /* end of Va filter 25 Hz */
  205. /* Vz filter 100 Hz */
  206. REAL_TYPE
  207. Vz_filter_100(REAL_TYPE Vz) {
  208. static REAL_TYPE y = 0.0;
  209. static REAL_TYPE x1 = 0.0;
  210. static REAL_TYPE x2 = 0.0;
  211. static REAL_TYPE x1_tmp = 0.0;
  212. static REAL_TYPE x2_tmp = 0.0;
  213. static unsigned short debut = 1;
  214. /* 100 Hz coefficients */
  215. static REAL_TYPE a0 = 0.956543675476034;
  216. static REAL_TYPE a1 = -1.955578398054313;
  217. static REAL_TYPE b0 = 0.000479064865372430;
  218. static REAL_TYPE b1 = 0.000486212556348925;
  219. if (debut) {
  220. debut = 0;
  221. x1 = 0.0;
  222. x2 = 0.0;
  223. }
  224. // Output
  225. y = x2;
  226. // State
  227. x1_tmp = - a0 * x2 + b0 * Vz;
  228. x2_tmp = x1 - a1 * x2 + b1 * Vz;
  229. // Update
  230. x1 = x1_tmp;
  231. x2 = x2_tmp;
  232. return y;
  233. } /* end of Vz filter 100 Hz */
  234. /* Vz filter 50 Hz */
  235. REAL_TYPE
  236. Vz_filter_50(REAL_TYPE Vz) {
  237. static REAL_TYPE y = 0.0;
  238. static REAL_TYPE x1 = 0.0, x2 = 0.0;
  239. static REAL_TYPE x1_tmp = 0.0;
  240. static REAL_TYPE x2_tmp = 0.0;
  241. static unsigned short debut = 1;
  242. /* 50 Hz coefficients */
  243. static REAL_TYPE a0 = 0.914975803093201;
  244. static REAL_TYPE a1 = -1.911199519984605;
  245. static REAL_TYPE b0 = 0.001860178914816;
  246. static REAL_TYPE b1 = 0.001916104193780;
  247. if (debut) {
  248. debut = 0;
  249. x1 = 0.0;
  250. x2 = 0.0;
  251. }
  252. // Output
  253. y = x2;
  254. // State
  255. x1_tmp = - a0 * x2 + b0 * Vz;
  256. x2_tmp = x1 - a1 * x2 + b1 * Vz;
  257. // Update
  258. x1 = x1_tmp;
  259. x2 = x2_tmp;
  260. return y;
  261. } /* end of Vz filter 50 Hz */
  262. /* Vz filter 33 Hz */
  263. REAL_TYPE
  264. Vz_filter_33(REAL_TYPE Vz) {
  265. static REAL_TYPE y = 0.0;
  266. static REAL_TYPE x1 = 0.0, x2 = 0.0;
  267. static REAL_TYPE x1_tmp = 0.0;
  268. static REAL_TYPE x2_tmp = 0.0;
  269. static unsigned short debut = 1;
  270. /* 33 Hz coefficients */
  271. static REAL_TYPE a0 = 0.874036784828483;
  272. static REAL_TYPE a1 = -1.865563793814790;
  273. static REAL_TYPE b0 = 0.004141433623051;
  274. static REAL_TYPE b1 = 0.004331557390642;
  275. if (debut) {
  276. debut = 0;
  277. x1 = 0.0;
  278. x2 = 0.0;
  279. }
  280. // Output
  281. y = x2;
  282. // State
  283. x1_tmp = - a0 * x2 + b0 * Vz;
  284. x2_tmp = x1 - a1 * x2 + b1 * Vz;
  285. // Update
  286. x1 = x1_tmp;
  287. x2 = x2_tmp;
  288. return y;
  289. } /* end of Vz filter 33 Hz */
  290. /* Vz filter 25 Hz */
  291. REAL_TYPE
  292. Vz_filter_25(REAL_TYPE Vz) {
  293. static REAL_TYPE y = 0.0;
  294. static REAL_TYPE x1 = 0.0;
  295. static REAL_TYPE x2 = 0.0;
  296. static REAL_TYPE x1_tmp = 0.0;
  297. static REAL_TYPE x2_tmp = 0.0;
  298. static unsigned short debut = 1;
  299. /* 25 Hz coefficients */
  300. static REAL_TYPE a0 = 0.837180720246048;
  301. static REAL_TYPE a1 = -1.822731999002980;
  302. static REAL_TYPE b0 = 0.007010380719078;
  303. static REAL_TYPE b1 = 0.007438340523990;
  304. if (debut) {
  305. debut = 0;
  306. x1 = 0.0;
  307. x2 = 0.0;
  308. }
  309. // Output
  310. y = x2;
  311. // State
  312. x1_tmp = - a0 * x2 + b0 * Vz;
  313. x2_tmp = x1 - a1 * x2 + b1 * Vz;
  314. // Update
  315. x1 = x1_tmp;
  316. x2 = x2_tmp;
  317. return y;
  318. } /* end of Vz filter 25 Hz */
  319. /* q filter 100 Hz */
  320. REAL_TYPE
  321. q_filter_100(REAL_TYPE q){
  322. static REAL_TYPE y = 0.0;
  323. static REAL_TYPE x1 = 0.0;
  324. static REAL_TYPE x2 = 0.0;
  325. static REAL_TYPE x1_tmp = 0.0;
  326. static REAL_TYPE x2_tmp = 0.0;
  327. static unsigned short debut = 1;
  328. /* 100 Hz coefficients */
  329. static REAL_TYPE a0 = 0.766000101841272;
  330. static REAL_TYPE a1 = -1.734903205885821;
  331. static REAL_TYPE b0 = 0.014857648981438;
  332. static REAL_TYPE b1 = 0.016239246974013;
  333. if (debut) {
  334. debut = 0;
  335. x1 = 0.0;
  336. x2 = 0.0;
  337. }
  338. // Output
  339. y = x2;
  340. // State
  341. x1_tmp = - a0 * x2 + b0 * q;
  342. x2_tmp = x1 - a1 * x2 + b1 * q;
  343. // Update
  344. x1 = x1_tmp;
  345. x2 = x2_tmp;
  346. return y;
  347. } /* end of q filter 100 Hz */
  348. /* q filter 50 Hz */
  349. REAL_TYPE
  350. q_filter_50(REAL_TYPE q){
  351. static REAL_TYPE y = 0.0;
  352. static REAL_TYPE x1 = 0.0;
  353. static REAL_TYPE x2 = 0.0;
  354. static REAL_TYPE x1_tmp = 0.0;
  355. static REAL_TYPE x2_tmp = 0.0;
  356. static unsigned short debut = 1;
  357. /* 50 Hz coefficients */
  358. static REAL_TYPE a0 = 0.586756156020839;
  359. static REAL_TYPE a1 = -1.477888930110354;
  360. static REAL_TYPE b0 = 0.049596808318647;
  361. static REAL_TYPE b1 = 0.059270417591839;
  362. if (debut) {
  363. debut = 0;
  364. x1 = 0.0;
  365. x2 = 0.0;
  366. }
  367. // Output
  368. y = x2;
  369. // State
  370. x1_tmp = - a0 * x2 + b0 * q;
  371. x2_tmp = x1 - a1 * x2 + b1 * q;
  372. // Update
  373. x1 = x1_tmp;
  374. x2 = x2_tmp;
  375. return y;
  376. } /* end of q filter 50 Hz */
  377. /* q filter 33 Hz */
  378. REAL_TYPE
  379. q_filter_33(REAL_TYPE q){
  380. static REAL_TYPE y = 0.0;
  381. static REAL_TYPE x1 = 0.0;
  382. static REAL_TYPE x2 = 0.0;
  383. static REAL_TYPE x1_tmp = 0.0;
  384. static REAL_TYPE x2_tmp = 0.0;
  385. static unsigned short debut = 1;
  386. /* 33 Hz coefficients */
  387. static REAL_TYPE a0 = 0.445839214374383;
  388. static REAL_TYPE a1 = -1.227970132817902;
  389. static REAL_TYPE b0 = 0.094268996251840;
  390. static REAL_TYPE b1 = 0.123600085304640;
  391. if (debut) {
  392. debut = 0;
  393. x1 = 0.0;
  394. x2 = 0.0;
  395. }
  396. // Output
  397. y = x2;
  398. // State
  399. x1_tmp = - a0 * x2 + b0 * q;
  400. x2_tmp = x1 - a1 * x2 + b1 * q;
  401. // Update
  402. x1 = x1_tmp;
  403. x2 = x2_tmp;
  404. return y;
  405. } /* end of q filter 33 Hz */
  406. /* q filter 25 Hz */
  407. REAL_TYPE
  408. q_filter_25(REAL_TYPE q){
  409. static REAL_TYPE y = 0.0;
  410. static REAL_TYPE x1 = 0.0;
  411. static REAL_TYPE x2 = 0.0;
  412. static REAL_TYPE x1_tmp = 0.0;
  413. static REAL_TYPE x2_tmp = 0.0;
  414. static unsigned short debut = 1;
  415. /* 25 Hz coefficients */
  416. static REAL_TYPE a0 = 0.344282786628352;
  417. static REAL_TYPE a1 = -1.010643377701049;
  418. static REAL_TYPE b0 = 0.137177088974822;
  419. static REAL_TYPE b1 = 0.196462319952482;
  420. if (debut) {
  421. debut = 0;
  422. x1 = 0.0;
  423. x2 = 0.0;
  424. }
  425. // Output
  426. y = x2;
  427. // State
  428. x1_tmp = - a0 * x2 + b0 * q;
  429. x2_tmp = x1 - a1 * x2 + b1 * q;
  430. // Update
  431. x1 = x1_tmp;
  432. x2 = x2_tmp;
  433. return y;
  434. } /* end of q filter 25 Hz */
  435. /* az filter 100 Hz */
  436. REAL_TYPE
  437. az_filter_100(REAL_TYPE az){
  438. static REAL_TYPE y = 0.0;
  439. static REAL_TYPE x1 = 0.0;
  440. static REAL_TYPE x2 = 0.0;
  441. static REAL_TYPE x1_tmp = 0.0;
  442. static REAL_TYPE x2_tmp = 0.0;
  443. static unsigned short debut = 1;
  444. /* 100 Hz coefficient */
  445. static REAL_TYPE a0 = 0.411240701442774;
  446. static REAL_TYPE a1 = -1.158045899830964;
  447. static REAL_TYPE b0 = 0.107849979167580;
  448. static REAL_TYPE b1 = 0.145344822444230;
  449. if (debut) {
  450. debut = 0;
  451. x1 = 0.0;
  452. x2 = 0.0;
  453. }
  454. // Output
  455. y = x2;
  456. // State
  457. x1_tmp = - a0 * x2 + b0 * az;
  458. x2_tmp = x1 - a1 * x2 + b1 * az;
  459. // Update
  460. x1 = x1_tmp;
  461. x2 = x2_tmp;
  462. return y;
  463. } /* end of az filter 100 Hz */
  464. /* az filter 50 Hz */
  465. REAL_TYPE
  466. az_filter_50(REAL_TYPE az){
  467. static REAL_TYPE y = 0.0;
  468. static REAL_TYPE x1 = 0.0;
  469. static REAL_TYPE x2 = 0.0;
  470. static REAL_TYPE x1_tmp = 0.0;
  471. static REAL_TYPE x2_tmp = 0.0;
  472. static unsigned short debut = 1;
  473. /* 50 Hz coefficients */
  474. static REAL_TYPE a0 = 0.169118914523145;
  475. static REAL_TYPE a1 = -0.518588903229759;
  476. static REAL_TYPE b0 = 0.229019233988375;
  477. static REAL_TYPE b1 = 0.421510777305010;
  478. if (debut) {
  479. debut = 0;
  480. x1 = 0.0;
  481. x2 = 0.0;
  482. }
  483. // Output
  484. y = x2;
  485. // State
  486. x1_tmp = - a0 * x2 + b0 * az;
  487. x2_tmp = x1 - a1 * x2 + b1 * az;
  488. // Update
  489. x1 = x1_tmp;
  490. x2 = x2_tmp;
  491. return y;
  492. } /* end of az filter 50 Hz */
  493. /* az filter 33 Hz */
  494. REAL_TYPE
  495. az_filter_33(REAL_TYPE az){
  496. static REAL_TYPE y = 0.0;
  497. static REAL_TYPE x1 = 0.0;
  498. static REAL_TYPE x2 = 0.0;
  499. static REAL_TYPE x1_tmp = 0.0;
  500. static REAL_TYPE x2_tmp = 0.0;
  501. static unsigned short debut = 1;
  502. /* 33 Hz coefficients */
  503. static REAL_TYPE a0 = 0.067700864731348;
  504. static REAL_TYPE a1 = -0.115832026705568;
  505. static REAL_TYPE b0 = 0.263451167882487;
  506. static REAL_TYPE b1 = 0.688417670143293;
  507. if (debut) {
  508. debut = 0;
  509. x1 = 0.0;
  510. x2 = 0.0;
  511. }
  512. // Output
  513. y = x2;
  514. // State
  515. x1_tmp = - a0 * x2 + b0 * az;
  516. x2_tmp = x1 - a1 * x2 + b1 * az;
  517. // Update
  518. x1 = x1_tmp;
  519. x2 = x2_tmp;
  520. return y;
  521. } /* end of az filter 33 Hz */
  522. /* az filter 25 Hz */
  523. REAL_TYPE
  524. az_filter_25(REAL_TYPE az){
  525. static REAL_TYPE y = 0.0;
  526. static REAL_TYPE x1 = 0.0;
  527. static REAL_TYPE x2 = 0.0;
  528. static REAL_TYPE x1_tmp = 0.0;
  529. static REAL_TYPE x2_tmp = 0.0;
  530. static unsigned short debut = 1;
  531. /* 25 Hz coefficients */
  532. static REAL_TYPE a0 = 0.028601207249487;
  533. static REAL_TYPE a1 = 0.069303378493245;
  534. static REAL_TYPE b0 = 0.228783762747218;
  535. static REAL_TYPE b1 = 0.869120822995514;
  536. if (debut) {
  537. debut = 0;
  538. x1 = 0.0;
  539. x2 = 0.0;
  540. }
  541. // Output
  542. y = x2;
  543. // State
  544. x1_tmp = - a0 * x2 + b0 * az;
  545. x2_tmp = x1 - a1 * x2 + b1 * az;
  546. // Update
  547. x1 = x1_tmp;
  548. x2 = x2_tmp;
  549. return y;
  550. } /* end of az filter 25 Hz */
  551. /* h filter 100 Hz*/
  552. REAL_TYPE
  553. h_filter_100(REAL_TYPE h){
  554. static REAL_TYPE y = 0.0;
  555. static REAL_TYPE x1 = 0.0;
  556. static REAL_TYPE x2 = 0.0;
  557. static REAL_TYPE x1_tmp = 0.0;
  558. static REAL_TYPE x2_tmp = 0.0;
  559. static unsigned short debut = 1;
  560. /* 100 Hz coefficients */
  561. static REAL_TYPE a0 = 0.766000101841272;
  562. static REAL_TYPE a1 = -1.734903205885821;
  563. static REAL_TYPE b0 = 0.014857648981438;
  564. static REAL_TYPE b1 = 0.016239246974013;
  565. if (debut) {
  566. debut = 0;
  567. x1 = h_eq * (1.0 + a1 - b1);
  568. x2 = h_eq;
  569. }
  570. // Output
  571. y = x2;
  572. // State
  573. x1_tmp = - a0 * x2 + b0 * h;
  574. x2_tmp = x1 - a1 * x2 + b1 * h;
  575. // Update
  576. x1 = x1_tmp;
  577. x2 = x2_tmp;
  578. return y;
  579. } /* end of h filter 100 Hz */
  580. /* h filter 50 Hz*/
  581. REAL_TYPE
  582. h_filter_50(REAL_TYPE h){
  583. static REAL_TYPE y = 0.0;
  584. static REAL_TYPE x1 = 0.0;
  585. static REAL_TYPE x2 = 0.0;
  586. static REAL_TYPE x1_tmp = 0.0;
  587. static REAL_TYPE x2_tmp = 0.0;
  588. static unsigned short debut = 1;
  589. /* 50 Hz coefficients */
  590. static REAL_TYPE a0 = 0.586756156020839;
  591. static REAL_TYPE a1 = -1.477888930110354;
  592. static REAL_TYPE b0 = 0.049596808318647;
  593. static REAL_TYPE b1 = 0.059270417591839;
  594. if (debut) {
  595. debut = 0;
  596. x1 = h_eq * (1.0 + a1 - b1);
  597. x2 = h_eq;
  598. }
  599. // Output
  600. y = x2;
  601. // State
  602. x1_tmp = - a0 * x2 + b0 * h;
  603. x2_tmp = x1 - a1 * x2 + b1 * h;
  604. // Update
  605. x1 = x1_tmp;
  606. x2 = x2_tmp;
  607. return y;
  608. } /* end of h filter 50 Hz */
  609. /* h filter 33 Hz*/
  610. REAL_TYPE
  611. h_filter_33(REAL_TYPE h){
  612. static REAL_TYPE y = 0.0;
  613. static REAL_TYPE x1 = 0.0;
  614. static REAL_TYPE x2 = 0.0;
  615. static REAL_TYPE x1_tmp = 0.0;
  616. static REAL_TYPE x2_tmp = 0.0;
  617. static unsigned short debut = 1;
  618. /* 33 Hz coefficients */
  619. static REAL_TYPE a0 = 0.445839214374383;
  620. static REAL_TYPE a1 = -1.227970132817902;
  621. static REAL_TYPE b0 = 0.094268996251840;
  622. static REAL_TYPE b1 = 0.123600085304640;
  623. if (debut) {
  624. debut = 0;
  625. x1 = h_eq * (1.0 + a1 - b1);
  626. x2 = h_eq;
  627. }
  628. // Output
  629. y = x2;
  630. // State
  631. x1_tmp = - a0 * x2 + b0 * h;
  632. x2_tmp = x1 - a1 * x2 + b1 * h;
  633. // Update
  634. x1 = x1_tmp;
  635. x2 = x2_tmp;
  636. return y;
  637. } /* end of h filter 33 Hz */
  638. /* h filter 25 Hz*/
  639. REAL_TYPE
  640. h_filter_25(REAL_TYPE h){
  641. static REAL_TYPE y = 0.0;
  642. static REAL_TYPE x1 = 0.0;
  643. static REAL_TYPE x2 = 0.0;
  644. static REAL_TYPE x1_tmp = 0.0;
  645. static REAL_TYPE x2_tmp = 0.0;
  646. static unsigned short debut = 1;
  647. /* 25 Hz coefficients */
  648. static REAL_TYPE a0 = 0.344282786628352;
  649. static REAL_TYPE a1 = -1.010643377701049;
  650. static REAL_TYPE b0 = 0.137177088974822;
  651. static REAL_TYPE b1 = 0.196462319952482; /**/
  652. if (debut) {
  653. debut = 0;
  654. x1 = h_eq * (1.0 + a1 - b1);
  655. x2 = h_eq;
  656. }
  657. // Output
  658. y = x2;
  659. // State
  660. x1_tmp = - a0 * x2 + b0 * h;
  661. x2_tmp = x1 - a1 * x2 + b1 * h;
  662. // Update
  663. x1 = x1_tmp;
  664. x2 = x2_tmp;
  665. return y;
  666. } /* end of h filter 25 Hz */
  667. /* Altitude hold controller 50 Hz */
  668. REAL_TYPE
  669. altitude_hold_50(REAL_TYPE h_f, REAL_TYPE h_c){
  670. static REAL_TYPE y = 0.0;
  671. static REAL_TYPE Ts_h = 1.0/50.0;
  672. static REAL_TYPE integrator = 532.2730285;
  673. if ((h_f - h_c) < -50) {
  674. // Output
  675. y = Vz_c;
  676. }
  677. else if ((h_f - h_c) > 50) {
  678. // Output
  679. y = -Vz_c;
  680. }
  681. else {
  682. // Output
  683. y = Kp_h * (h_f - h_c) + Ki_h * integrator;
  684. // State
  685. integrator += Ts_h * (h_f - h_c);
  686. }
  687. return y;
  688. }
  689. /* Altitude hold controller 33 Hz */
  690. REAL_TYPE
  691. altitude_hold_33(REAL_TYPE h_f, REAL_TYPE h_c){
  692. static REAL_TYPE y = 0.0;
  693. static REAL_TYPE Ts_h = 1.0/33.0;
  694. static REAL_TYPE integrator = 532.2730285;
  695. if ((h_f - h_c) < -50) {
  696. // Output
  697. y = Vz_c;
  698. }
  699. else if ((h_f - h_c) > 50) {
  700. // Output
  701. y = -Vz_c;
  702. }
  703. else {
  704. // Output
  705. y = Kp_h * (h_f - h_c) + Ki_h * integrator;
  706. // State
  707. integrator += Ts_h * (h_f - h_c);
  708. }
  709. return y;
  710. }
  711. /* Altitude hold controller 25 Hz */
  712. REAL_TYPE
  713. altitude_hold_25(REAL_TYPE h_f, REAL_TYPE h_c){
  714. static REAL_TYPE y = 0.0;
  715. static REAL_TYPE Ts_h = 1.0/25.0;
  716. static REAL_TYPE integrator = 532.2730285;
  717. if ((h_f - h_c) < -50) {
  718. // Output
  719. y = Vz_c;
  720. }
  721. else if ((h_f - h_c) > 50) {
  722. // Output
  723. y = -Vz_c;
  724. }
  725. else {
  726. // Output
  727. y = Kp_h * (h_f - h_c) + Ki_h * integrator;
  728. // State
  729. integrator += Ts_h * (h_f - h_c);
  730. }
  731. return y;
  732. }
  733. /* Altitude hold controller 10 Hz */
  734. REAL_TYPE
  735. altitude_hold_10(REAL_TYPE h_f, REAL_TYPE h_c){
  736. static REAL_TYPE y = 0.0;
  737. static REAL_TYPE Ts_h = 1.0/10.0;
  738. static REAL_TYPE integrator = 532.2730285;
  739. if ((h_f - h_c) < -50) {
  740. // Output
  741. y = Vz_c;
  742. }
  743. else if ((h_f - h_c) > 50) {
  744. // Output
  745. y = -Vz_c;
  746. }
  747. else {
  748. // Output
  749. y = Kp_h * (h_f - h_c) + Ki_h * integrator;
  750. // State
  751. integrator += Ts_h * (h_f - h_c);
  752. }
  753. return y;
  754. }
  755. /* Va Speed controller 50 Hz */
  756. REAL_TYPE
  757. Va_control_50(REAL_TYPE Va_f, REAL_TYPE Vz_f, REAL_TYPE q_f, REAL_TYPE Va_c){
  758. static REAL_TYPE y = 0.0;
  759. static REAL_TYPE Ts_K1 = 1.0/50.0;
  760. static REAL_TYPE integrator = 0.0;
  761. // Output
  762. y = K1_intVa * integrator + K1_Va * (Va_f - Va_eq) + K1_Vz * Vz_f + K1_q * q_f + delta_th_eq;
  763. // State
  764. integrator += Ts_K1 * (Va_c - Va_f + Va_eq);
  765. return y;
  766. }
  767. /* Va Speed controller 33 Hz */
  768. REAL_TYPE
  769. Va_control_33(REAL_TYPE Va_f, REAL_TYPE Vz_f, REAL_TYPE q_f, REAL_TYPE Va_c){
  770. static REAL_TYPE y = 0.0;
  771. static REAL_TYPE Ts_K1 = 1.0/33.0;
  772. static REAL_TYPE integrator = 0.0;
  773. // Output
  774. y = K1_intVa * integrator + K1_Va * (Va_f - Va_eq) + K1_Vz * Vz_f + K1_q * q_f + delta_th_eq;
  775. // State
  776. integrator += Ts_K1 * (Va_c - Va_f + Va_eq);
  777. return y;
  778. }
  779. /* Va Speed controller 25 Hz */
  780. REAL_TYPE
  781. Va_control_25(REAL_TYPE Va_f, REAL_TYPE Vz_f, REAL_TYPE q_f, REAL_TYPE Va_c){
  782. static REAL_TYPE y = 0.0;
  783. static REAL_TYPE Ts_K1 = 1.0/25.0;
  784. static REAL_TYPE integrator = 0.0;
  785. // Output
  786. y = K1_intVa * integrator + K1_Va * (Va_f - Va_eq) + K1_Vz * Vz_f + K1_q * q_f + delta_th_eq;
  787. // State
  788. integrator += Ts_K1 * (Va_c - Va_f + Va_eq);
  789. return y;
  790. }
  791. /* Va Speed controller 10 Hz */
  792. REAL_TYPE
  793. Va_control_10(REAL_TYPE Va_f, REAL_TYPE Vz_f, REAL_TYPE q_f, REAL_TYPE Va_c){
  794. static REAL_TYPE y = 0.0;
  795. static REAL_TYPE Ts_K1 = 1.0/10.0;
  796. static REAL_TYPE integrator = 0.0;
  797. // Output
  798. y = K1_intVa * integrator + K1_Va * (Va_f - Va_eq) + K1_Vz * Vz_f + K1_q * q_f + delta_th_eq;
  799. // State
  800. integrator += Ts_K1 * (Va_c - Va_f + Va_eq);
  801. return y;
  802. }
  803. /* Vz Speed controller 50 Hz */
  804. REAL_TYPE
  805. Vz_control_50(REAL_TYPE Vz_f, REAL_TYPE Vz_c, REAL_TYPE q_f, REAL_TYPE az_f){
  806. static REAL_TYPE y = 0.0;
  807. static REAL_TYPE Ts_K2 = 1.0/50.0;
  808. static REAL_TYPE integrator = 0.0;
  809. // Output
  810. y = K2_intVz * integrator + K2_Vz * Vz_f + K2_q * q_f + K2_az * az_f + delta_e_eq;
  811. // State
  812. integrator += Ts_K2 * (Vz_c - Vz_f);
  813. return y;
  814. }
  815. /* Vz Speed controller 33 Hz */
  816. REAL_TYPE
  817. Vz_control_33(REAL_TYPE Vz_f, REAL_TYPE Vz_c, REAL_TYPE q_f, REAL_TYPE az_f){
  818. static REAL_TYPE y = 0.0;
  819. static REAL_TYPE Ts_K2 = 1.0/33.0;
  820. static REAL_TYPE integrator = 0.0;
  821. // Output
  822. y = K2_intVz * integrator + K2_Vz * Vz_f + K2_q * q_f + K2_az * az_f + delta_e_eq;
  823. // State
  824. integrator += Ts_K2 * (Vz_c - Vz_f);
  825. return y;
  826. }
  827. /* Vz Speed controller 25 Hz */
  828. REAL_TYPE
  829. Vz_control_25(REAL_TYPE Vz_f, REAL_TYPE Vz_c, REAL_TYPE q_f, REAL_TYPE az_f){
  830. static REAL_TYPE y = 0.0;
  831. static REAL_TYPE Ts_K2 = 1.0/25.0;
  832. static REAL_TYPE integrator = 0.0;
  833. // Output
  834. y = K2_intVz * integrator + K2_Vz * Vz_f + K2_q * q_f + K2_az * az_f + delta_e_eq;
  835. // State
  836. integrator += Ts_K2 * (Vz_c - Vz_f);
  837. return y;
  838. }
  839. /* Vz Speed controller 10 Hz */
  840. REAL_TYPE
  841. Vz_control_10(REAL_TYPE Vz_f, REAL_TYPE Vz_c, REAL_TYPE q_f, REAL_TYPE az_f){
  842. static REAL_TYPE y = 0.0;
  843. static REAL_TYPE Ts_K2 = 1.0/10.0;
  844. static REAL_TYPE integrator = 0.0;
  845. // Output
  846. y = K2_intVz * integrator + K2_Vz * Vz_f + K2_q * q_f + K2_az * az_f + delta_e_eq;
  847. // State
  848. integrator += Ts_K2 * (Vz_c - Vz_f);
  849. return y;
  850. }
  851. /* Engine */
  852. REAL_TYPE
  853. engine(REAL_TYPE delta_th_c) {
  854. static REAL_TYPE y = delta_th_eq;
  855. static REAL_TYPE x1 = delta_th_eq;
  856. static REAL_TYPE x1_dot = 0.0;
  857. static REAL_TYPE tau = 0.75;
  858. // Output
  859. y = 26350.0 * x1;
  860. // State Equation
  861. x1_dot = -tau * x1 + tau * delta_th_c;
  862. // Update State
  863. x1 += dt_dx * x1_dot;
  864. return y;
  865. }
  866. /* Elevator */
  867. REAL_TYPE
  868. elevator(REAL_TYPE delta_e_c) {
  869. static REAL_TYPE y = delta_e_eq;
  870. static REAL_TYPE x1 = delta_e_eq;
  871. static REAL_TYPE x2 = 0.0;
  872. static REAL_TYPE x1_dot = 0.0;
  873. static REAL_TYPE x2_dot = 0.0;
  874. static REAL_TYPE omega = 25.0;
  875. static REAL_TYPE xi = 0.85;
  876. // Output
  877. y = x1;
  878. // State Equation
  879. x1_dot = x2;
  880. x2_dot = -omega * omega * x1 - 2.0 * xi * omega * x2 + omega * omega * delta_e_c;
  881. // Update State
  882. x1 += dt_de * x1_dot;
  883. x2 += dt_de * x2_dot;
  884. return y;
  885. }
  886. /* Flight dynamics */
  887. void aircraft_dynamics (REAL_TYPE delta_e, REAL_TYPE T, struct aircraft_dynamics_outs_t *outputs){
  888. static int debut = 1;
  889. static REAL_TYPE u = 0.0;
  890. static REAL_TYPE w = 0.0;
  891. static REAL_TYPE q = 0.0;
  892. static REAL_TYPE theta = 0.0;
  893. static REAL_TYPE h = 0.0;
  894. static REAL_TYPE u_dot = 0.0;
  895. static REAL_TYPE w_dot = 0.0;
  896. static REAL_TYPE q_dot = 0.0;
  897. static REAL_TYPE theta_dot = 0.0;
  898. static REAL_TYPE h_dot = 0.0;
  899. static REAL_TYPE CD = 0.0;
  900. static REAL_TYPE CL = 0.0;
  901. static REAL_TYPE Cm = 0.0;
  902. static REAL_TYPE Xa = 0.0;
  903. static REAL_TYPE Za = 0.0;
  904. static REAL_TYPE Ma = 0.0;
  905. static REAL_TYPE alpha = 0.0;
  906. static REAL_TYPE qbar = 0.0;
  907. static REAL_TYPE V = 0.0;
  908. static REAL_TYPE rho = 0.0;
  909. if (debut) {
  910. debut = 0;
  911. u = Va_eq * MATH_COS(theta_eq);
  912. w = Va_eq * MATH_SIN(theta_eq);
  913. q = 0.0;
  914. theta = theta_eq;
  915. h = h_eq;
  916. }
  917. rho = rho0 * pow(1.0 + T0_h / T0_0 * h,- g0 / (Rs * T0_h) - 1.0);
  918. alpha = atan(w/u);
  919. V = sqrt(u * u + w * w);
  920. qbar = 0.5 * rho * V * V;
  921. CL = CL_deltae * delta_e + CL_alpha * (alpha - alpha_0);
  922. CD = CD_0 + CD_deltae * delta_e + CD_alpha * (alpha - alpha_0) * (alpha - alpha_0);
  923. Cm = Cm_0 + Cm_deltae * delta_e + Cm_alpha * alpha + 0.5 * Cm_q * q * cbar / V;
  924. Xa = - qbar * S * (CD * MATH_COS(alpha) - CL * MATH_SIN(alpha));
  925. Za = - qbar * S * (CD * MATH_SIN(alpha) + CL * MATH_COS(alpha));
  926. Ma = qbar * cbar * S * Cm;
  927. // Output
  928. outputs -> Va = V;
  929. outputs -> Vz = w * MATH_COS(theta) - u * MATH_SIN(theta);
  930. outputs -> q = q;
  931. outputs -> az = g0 * MATH_COS(theta) + Za / masse;
  932. outputs -> h = h;
  933. // State Equation
  934. u_dot = - g0 * sin (theta) - q * w + (Xa + T) / masse;
  935. w_dot = g0 * MATH_COS (theta) + q * u + Za / masse;
  936. q_dot = Ma / I_y;
  937. theta_dot = q;
  938. h_dot = u * MATH_SIN(theta) - w * MATH_COS(theta);
  939. // Update State
  940. u += dt * u_dot;
  941. w += dt * w_dot;
  942. q += dt * q_dot;
  943. theta += dt * theta_dot;
  944. h += dt * h_dot;
  945. static REAL_TYPE Time = 0.0;
  946. // instant++;
  947. Time = Time + dt;
  948. }
  949. /*
  950. * The commanded altitude
  951. */
  952. REAL_TYPE input_h_c() {
  953. return h_c;
  954. }
  955. /*
  956. * The commanded airspeed
  957. */
  958. REAL_TYPE input_Va_c() {
  959. return Va_c;
  960. }
  961. void output_delta_th_c(REAL_TYPE delta_th_c){
  962. }
  963. void output_delta_e_c(REAL_TYPE delta_e_c){
  964. }