diff --git a/functions/odeFunctions/ballistic.m b/functions/odeFunctions/ballistic.m index 0e16abd4f74355055c66f1839180c1d248e7e8fb..d08b3069d7274836e0f686c19963919facaec233 100644 --- a/functions/odeFunctions/ballistic.m +++ b/functions/odeFunctions/ballistic.m @@ -96,11 +96,13 @@ if theta < tb m = interpLinear(rocket.motor.time, rocket.mass, theta); T = interpLinear(rocket.motor.time, rocket.motor.thrust, theta); Pe = interpLinear(rocket.motor.time, rocket.motor.pe, theta); + xcg = interpLinear(rocket.motor.time, rocket.xcg, theta); T = T + rocket.motor.ae*(Pe - P); I = interpLinear(rocket.motor.time, rocket.inertia, theta); Idot = interpLinear(rocket.motor.time, rocket.inertiaDot, theta); else % for theta >= tb the flight condition is the empty one (no interpolation needed) + xcg = rocket.xcg(end); m = rocket.cutoffMass; T = 0; @@ -143,53 +145,20 @@ betaOut = beta; %% INTERPOLATE AERODYNAMIC COEFFICIENTS: -if abs(alpha)>25*pi/180 || abs(beta)>25*pi/180 - coeffsValues = interpN( rocket.coefficientsHighAOA.total,... - {rocket.coefficientsHighAOA.state.alphas, rocket.coefficientsHighAOA.state.machs, ... - rocket.coefficientsHighAOA.state.betas, rocket.coefficientsHighAOA.state.altitudes}, ... - [alpha, mach, beta, absoluteAltitude]); - angle0 = [alpha beta]; +if abs(alpha) > rocket.coefficients.state.alphas(end)*pi/180 || ... + abs(beta) > rocket.coefficients.state.betas(end)*pi/180 + coeffsValues = rocket.coefficientsHighAOA.get(alpha, mach, beta, absoluteAltitude, extension, xcg); else - %[coeffsValues, angle0] = rocket.interpCoeffs(theta, alpha, mach, beta, absoluteAltitude, opening); - [coeffsValues, angle0] = interpCoeffs(rocket.coefficients, theta, rocket.motor.cutoffTime, ... - alpha, mach, beta, absoluteAltitude, extension); + coeffsValues = rocket.coefficients.get(alpha, mach, beta, absoluteAltitude, extension, xcg); end % Retrieve Coefficients -CA = coeffsValues(1); CYB = coeffsValues(2); CY0 = coeffsValues(3); -CNA = coeffsValues(4); CN0 = coeffsValues(5); Cl = coeffsValues(6); -Clp = coeffsValues(7); Cma = coeffsValues(8); Cm0 = coeffsValues(9); -Cmad = coeffsValues(10); Cmq = coeffsValues(11); Cnb = coeffsValues(12); -Cn0 = coeffsValues(13); Cnr = coeffsValues(14); Cnp = coeffsValues(15); - -% Clb = coeffsValues(16); -% XCP_value = coeffsValues(16); - -% compute CN,CY,Cm,Cn (linearized with respect to alpha and beta): -alpha0 = angle0(1); beta0 = angle0(2); - -CN = (CN0 + CNA*(alpha - alpha0)); -CY = (CY0 + CYB*(beta - beta0)); -Cm = (Cm0 + Cma*(alpha - alpha0)); -Cn = (Cn0 + Cnb*(beta - beta0)); - -% XCPlon = Cm/CN; -if abs(alpha) <= pi/180 - XCPlon = Cma/CNA; -else - XCPlon = Cm/CN; -end - -% XCPlat = Cn/CY; -if abs(beta) <= pi/180 - XCPlat = Cnb/CYB; -else - XCPlat = Cn/CY; -end +CA = coeffsValues(1); CY = coeffsValues(2); CN = coeffsValues(3); +Cl = coeffsValues(4); Cm = coeffsValues(5); Cn = coeffsValues(6); +Clp = coeffsValues(7); Cmad = coeffsValues(8); Cmq = coeffsValues(9); +Cnr = coeffsValues(10); Cnp = coeffsValues(11); -% if Cn == 0 && CY == 0 -% XCPlat = -5; -% end +XCPtot = NaN; %% RAMP / FREE FLIGHT if altitude < environment.effectiveRampAltitude % No torque on the launchpad @@ -209,14 +178,10 @@ if altitude < environment.effectiveRampAltitude % No betaOut = NaN; fY = 0; fZ = 0; - XCPlon = NaN; - XCPlat = NaN; if T < Fg % No velocity untill T = Fg du = 0; end - - XCPtot = NaN; else %% FORCES % first computed in the body-frame reference system @@ -255,8 +220,6 @@ else CMaTot = cos(phi)*Cm - sin(phi)*Cn; if CFaTot ~= 0 XCPtot = CMaTot/CFaTot; - else - XCPtot = XCPlon; end dAngle = 0; @@ -331,19 +294,12 @@ if nargout == 2 || ~isempty(wrapper) parout.accelerometer.body_acc = ([-fX+T, fY, -fZ]')/m; parout.coeff.CA = CA; - parout.coeff.CYB = CYB; - parout.coeff.CNA = CNA; parout.coeff.Cl = Cl; parout.coeff.Clp = Clp; - %parout.coeff.Clb = Clb; - parout.coeff.Cma = Cma; parout.coeff.Cmad = Cmad; parout.coeff.Cmq = Cmq; - parout.coeff.Cnb = Cnb; parout.coeff.Cnr = Cnr; parout.coeff.Cnp = Cnp; - parout.coeff.XCPlon = XCPlon; - parout.coeff.XCPlat = XCPlat; parout.coeff.XCPtot = XCPtot; parout.uncertanty = [];