From 1554cab999d8cb42557c75173bd7641539a2ce29 Mon Sep 17 00:00:00 2001
From: Mauco03 <marco.gaibotti@skywarder.eu>
Date: Sat, 25 May 2024 09:03:29 +0200
Subject: [PATCH] [refactoring-ode][simulations] Added stdStability

---
 functions/simulations/stdStability.m | 53 ++++++++++++++++++++++++++--
 1 file changed, 51 insertions(+), 2 deletions(-)

diff --git a/functions/simulations/stdStability.m b/functions/simulations/stdStability.m
index cd315f4..be52786 100644
--- a/functions/simulations/stdStability.m
+++ b/functions/simulations/stdStability.m
@@ -1,4 +1,4 @@
-function stability = stdStability(rocket, env, wind, settings, wrapper)
+function stabilityOut = stdStability(rocket, env, wind, settings, wrapper)
 arguments
     rocket      Rocket
     env         Environment
@@ -43,10 +43,59 @@ W0 = [0; 0; 0];
 
 Y0 = [X0; V0; W0; Q0; 0];
 
-%% LAUNCHPAD
+%% FLIGHT
 t0 = 0;
 solution = ode113(@ballistic, [t0, tf], Y0, options, ...
     rocket, env, wind, [], 0, wrapper);
 
 stability = getOdeFcn(solution);
+
+%% GET EXIT CONDITION IN AERO FRAME
+Ypad = stability.state.Y(:, end);
+tPad = stability.state.time(end);
+
+alpha = stability.interp.alpha(end);
+beta = stability.interp.beta(end);
+alt = env.z0 + env.effectiveRampAltitude;
+mach = stability.interp.mach(end);
+xcg = interpLinear(rocket.motor.time, rocket.xCg, tPad);
+
+%%% create dissileMatcom input data struct
+dissileVars.alpha = sort([-alpha, alpha, 0, -1, 1]);
+dissileVars.beta = beta;
+dissileVars.alt = alt;
+dissileVars.mach = mach;
+dissileVars.xcg = xcg;
+dissileVars.hprot = [];
+
+input = createDissileInput(rocket, dissileVars);
+
+input.flagXCP = true;
+coeff = dissileMatcom(input);
+
+vars.useAlphaTot = false;
+vars.checkOpenRocket = false;
+vars.modHaack = true;
+
+if dissileVars.alpha(1) == alpha
+    iAlpha = 1;
+else
+    iAlpha = 5;
+end
+
+%%% compute the stability margin and coresponding quality indices
+[XCPvect, qualityIndex, stabilityCoeffs] = getXCP(coeff, iAlpha, 1, input, vars);
+stabilityOut.time = tPad;
+stabilityOut.state = Ypad;
+
+stabilityOut.wind.NED = stability.wind.NED(:, end);
+stabilityOut.wind.body = stability.wind.body(:, end);
+
+stabilityOut.alpha = alpha;
+stabilityOut.beta = beta;
+
+stabilityOut.XCPvect = XCPvect;
+stabilityOut.qualityIndex = qualityIndex;
+stabilityOut.stabilityCoeffs = stabilityCoeffs;
+stabilityOut.uncertanty = [];
 end
\ No newline at end of file
-- 
GitLab