Hello,

I have been trying to plot a shaded area (right y-axis) behind two lines (left y-axis), as shown in the code below. However, the shaded area always appears last, covering up the two lines. I’ve done this kind of graph before but have never had this issue. Best I can tell, it is because the right y-axis is plotted last. Is there a way to make the right y-axis plotted first? I can’t simply switch which data is on the axis, unfortunately: they need to be in this order. Thanks for any help!

clear; clc; close all;
Data = xlsread("FigureOfMeritResults.xlsx");
% Extract data
Date = Data(469:1908,1);
FoM = Data(469:1908,3);
GHI = Data(469:1908,4);
% Calculate the running average over 5 data points
FoM_avg = movmean(FoM, 5);
GHI_avg = movmean(GHI, 5);
dFoM = zeros(length(FoM),1) + mean(FoM);
ticks = {’12:00AM’,’4:00AM’,’8:00AM’,’12:00PM’,’4:00PM’,’8:00PM’,’12:00AM’};
fz = 15;
%%
fig = figure;
set(gcf,’Color’,’w’,’Units’,’inches’)
set(gcf,’Position’,[1,1,5,6])
GHI_color = [0.8 0.8 1];
yyaxis left;
p1 = scatter(Date, FoM, ‘Marker’, ‘.’, ‘SizeData’, 15, ‘MarkerFaceColor’, "#0072BD");
hold on;
p2 = scatter(Date, dFoM, ‘Marker’, ‘_’, ‘SizeData’, 15, ‘MarkerFaceColor’, ‘black’);
hold on;
ax = gca;
ax.YColor = ‘k’;
yyaxis right;
xData = [Date; flipud(Date)];
yData = [zeros(size(GHI)); flipud(GHI)];
fillHandle = fill(xData, yData, GHI_color,’EdgeColor’,’none’);
hold on;
ax = gca;
ax.YColor = ‘k’;
xlabel(‘Time of Day (hrs)’, ‘FontSize’, fz);
numTicks = length(ticks);
tickPositions = linspace(min(Date), max(Date), numTicks);
xticks(tickPositions);
xticklabels(ticks);

set(gca, ‘FontSize’, fz);
ax = gca;
box on;
set(ax,’BoxStyle’,’full’,’LineWidth’,2,’XGrid’,’off’,’XMinorTick’,’off’,’YGrid’,’off’,…
‘YMinorTick’,’on’);Hello,

I have been trying to plot a shaded area (right y-axis) behind two lines (left y-axis), as shown in the code below. However, the shaded area always appears last, covering up the two lines. I’ve done this kind of graph before but have never had this issue. Best I can tell, it is because the right y-axis is plotted last. Is there a way to make the right y-axis plotted first? I can’t simply switch which data is on the axis, unfortunately: they need to be in this order. Thanks for any help!

clear; clc; close all;
Data = xlsread("FigureOfMeritResults.xlsx");
% Extract data
Date = Data(469:1908,1);
FoM = Data(469:1908,3);
GHI = Data(469:1908,4);
% Calculate the running average over 5 data points
FoM_avg = movmean(FoM, 5);
GHI_avg = movmean(GHI, 5);
dFoM = zeros(length(FoM),1) + mean(FoM);
ticks = {’12:00AM’,’4:00AM’,’8:00AM’,’12:00PM’,’4:00PM’,’8:00PM’,’12:00AM’};
fz = 15;
%%
fig = figure;
set(gcf,’Color’,’w’,’Units’,’inches’)
set(gcf,’Position’,[1,1,5,6])
GHI_color = [0.8 0.8 1];
yyaxis left;
p1 = scatter(Date, FoM, ‘Marker’, ‘.’, ‘SizeData’, 15, ‘MarkerFaceColor’, "#0072BD");
hold on;
p2 = scatter(Date, dFoM, ‘Marker’, ‘_’, ‘SizeData’, 15, ‘MarkerFaceColor’, ‘black’);
hold on;
ax = gca;
ax.YColor = ‘k’;
yyaxis right;
xData = [Date; flipud(Date)];
yData = [zeros(size(GHI)); flipud(GHI)];
fillHandle = fill(xData, yData, GHI_color,’EdgeColor’,’none’);
hold on;
ax = gca;
ax.YColor = ‘k’;
xlabel(‘Time of Day (hrs)’, ‘FontSize’, fz);
numTicks = length(ticks);
tickPositions = linspace(min(Date), max(Date), numTicks);
xticks(tickPositions);
xticklabels(ticks);

set(gca, ‘FontSize’, fz);
ax = gca;
box on;
set(ax,’BoxStyle’,’full’,’LineWidth’,2,’XGrid’,’off’,’XMinorTick’,’off’,’YGrid’,’off’,…
‘YMinorTick’,’on’); Hello,

I have been trying to plot a shaded area (right y-axis) behind two lines (left y-axis), as shown in the code below. However, the shaded area always appears last, covering up the two lines. I’ve done this kind of graph before but have never had this issue. Best I can tell, it is because the right y-axis is plotted last. Is there a way to make the right y-axis plotted first? I can’t simply switch which data is on the axis, unfortunately: they need to be in this order. Thanks for any help!

clear; clc; close all;
Data = xlsread("FigureOfMeritResults.xlsx");
% Extract data
Date = Data(469:1908,1);
FoM = Data(469:1908,3);
GHI = Data(469:1908,4);
% Calculate the running average over 5 data points
FoM_avg = movmean(FoM, 5);
GHI_avg = movmean(GHI, 5);
dFoM = zeros(length(FoM),1) + mean(FoM);
ticks = {’12:00AM’,’4:00AM’,’8:00AM’,’12:00PM’,’4:00PM’,’8:00PM’,’12:00AM’};
fz = 15;
%%
fig = figure;
set(gcf,’Color’,’w’,’Units’,’inches’)
set(gcf,’Position’,[1,1,5,6])
GHI_color = [0.8 0.8 1];
yyaxis left;
p1 = scatter(Date, FoM, ‘Marker’, ‘.’, ‘SizeData’, 15, ‘MarkerFaceColor’, "#0072BD");
hold on;
p2 = scatter(Date, dFoM, ‘Marker’, ‘_’, ‘SizeData’, 15, ‘MarkerFaceColor’, ‘black’);
hold on;
ax = gca;
ax.YColor = ‘k’;
yyaxis right;
xData = [Date; flipud(Date)];
yData = [zeros(size(GHI)); flipud(GHI)];
fillHandle = fill(xData, yData, GHI_color,’EdgeColor’,’none’);
hold on;
ax = gca;
ax.YColor = ‘k’;
xlabel(‘Time of Day (hrs)’, ‘FontSize’, fz);
numTicks = length(ticks);
tickPositions = linspace(min(Date), max(Date), numTicks);
xticks(tickPositions);
xticklabels(ticks);

set(gca, ‘FontSize’, fz);
ax = gca;
box on;
set(ax,’BoxStyle’,’full’,’LineWidth’,2,’XGrid’,’off’,’XMinorTick’,’off’,’YGrid’,’off’,…
‘YMinorTick’,’on’); plot, area MATLAB Answers — New Questions

​

So I’ve been trying to simulate a parachute descent (very basic model, almost a mass body with a coefficient of drag), from a height. While simulating only in the axial direction(z axis), the results coincide with kinematics equations. However, when I give a non-zero initial condition, the outputs go crazy (definitely not correct) with periodic solutions to related terms. I seem to have checked my equations, they seem to be right. Is there anything with ‘odeset’ I can work this out around.
Zero non-axial initial condition:
x0 = [0; 0; -100; 0; 0; 1; 0; 0; 0; 0; 0; 0];

Non-zero IC:
x0 = [0; 0; -100; 0; 1; 1; 0; 0; 0; 0; 0; 0];

0

this is the code i’ve used for these plots. Please feel free to ask any clarifications
%% 6DOF fucntion definition
function [dxdt, aoa_series] = six_DOF_dynamics(t,x)

M_para = 5; %mass of the parachute in Kgs
g = 9.81; %gravitational acceleration in m/s^2
I_x = 1; I_y = 1; I_z = 3; %Inertia matrix diagonal principal of the parachute in kg-m^2
CD_para = 0.9; %coefficient of drag of parachute
CL_para = 0; %coefficient of lift of parachute
rho = 1; %density of the atmosphere in Kg/m^3
x_com = 0; y_com = 0; z_com = -5; %position of COM of the paarchute from the integrated system’s COM(body axis origin) in m
A = 10; %parachute area in m^2

%state variables defining
x_pos = x(1); y_pos = x(2); z_pos = x(3); %postion variables
u = x(4); v = x(5); w = x(6); %velocity variables
phi = x(7); theta = x(8); psi = x(9); %attitude/euler angles variables
p = x(10); q = x(11); r = x(12); %euler angular rate variables

%rate of change of postion -> velocity
dx_pos = u;
dy_pos = v;
dz_pos = w;

%Gravity moment of the parachute along x,y and z components
GravityMoment_para_x = M_para*g*cos(theta)*(y_com*cos(phi) – z_com*sin(phi));
GravityMoment_para_y = M_para*g*(z_com*sin(theta) – x_com*cos(theta)*cos(phi));
GravityMoment_para_z = M_para*g*(x_com*cos(theta)*sin(phi) + y_com*sin(theta));

V_total = sqrt(u^2 + v^2 + w^2); %net speed

persistent angle_of_attack %defining this variable to take it as a output outside ode45

alpha = acos(w/V_total); %angle of attack calculated in radians

angle_of_attack = [angle_of_attack; alpha]; %concatinating current iterated angle of attack in global variable

%Componetize CD and CL into tangential and normal force coefficients
CT = CD_para*cos(alpha) – CL_para*sin(alpha);
CN = CL_para*cos(alpha) + CD_para*sin(alpha);

% Aerodynamic force along x,y,z axes F = 0.5*rho*v^2*A
if (u == 0 & v == 0)
Aerodynamic_force_x = 0;
Aerodynamic_force_y = 0;
Aerodynamic_force_z = -0.5*rho*CT*V_total*V_total*A;
else
Aerodynamic_force_x = -0.5*rho*CN*V_total*V_total*A*(u/(u^2 + v^2)^0.5);
Aerodynamic_force_y = -0.5*rho*CN*V_total*V_total*A*(v/(u^2 + v^2)^0.5);
Aerodynamic_force_z = -0.5*rho*CT*V_total*V_total*A;
end

%Moment caused by aerodynamic force along x,y,x
Aerodynamic_moment_x = y_com*Aerodynamic_force_z – z_com*Aerodynamic_force_y;
Aerodynamic_moment_y = z_com*Aerodynamic_force_x – x_com*Aerodynamic_force_z;
Aerodynamic_moment_z = x_com*Aerodynamic_force_y – y_com*Aerodynamic_force_x;

%rate change of velocity definition
du = v*r – w*q + g*sin(theta) + Aerodynamic_force_x/M_para;
dv = -u*r + w*p + g*cos(theta)*sin(phi) + Aerodynamic_force_y/M_para;
dw = u*q – v*p + g*cos(theta)*cos(phi) + Aerodynamic_force_z/M_para;

%rate change of euler angles definition
dphi = p + (r*cos(phi) + q*sin(phi))*tan(theta);
dtheta = q*cos(phi) – r*sin(phi);
dpsi = (cos(theta)^-1)*(r*cos(phi) + q*sin(phi));

%rate change of angular rates definition
dp = (1/I_x)*(Aerodynamic_moment_x + GravityMoment_para_x);
dq = (1/I_y)*(Aerodynamic_moment_y + GravityMoment_para_y);
dr = (1/I_z)*(Aerodynamic_moment_z + GravityMoment_para_z);

dxdt = [dx_pos; dy_pos; dz_pos; du; dv; dw; dphi; dtheta; dpsi; dp; dq; dr];

%block to output angle of attack
if nargout>1
aoa_series = angle_of_attack;

end

end

tspan = [0 3];

x0 = [0; 0; -100; 0; 1; 1; 0; 0; 0; 0; 0; 0];

[t, x] = ode45(@six_DOF_dynamics, tspan, x0);

[~, aoa_series] = six_DOF_dynamics(t,x);

%adding alpha as a part of x

% Plot results
figure;
subplot(4,1,1);
plot(t, x(:,1:3));
title(‘Position’);
legend(‘x’, ‘y’, ‘z’);

subplot(4,1,2);
plot(t, x(:,4:6));
title(‘Velocity’);
legend(‘u’, ‘v’, ‘w’);

subplot(4,1,3);
plot(t, x(:,7:9));
title(‘Euler Angles’);
legend(‘phi’, ‘theta’, ‘psi’);

subplot(4,1,4);
plot(t, x(:,10:12));
title(‘Euler Rates’);
legend(‘p’, ‘q’, ‘r’);So I’ve been trying to simulate a parachute descent (very basic model, almost a mass body with a coefficient of drag), from a height. While simulating only in the axial direction(z axis), the results coincide with kinematics equations. However, when I give a non-zero initial condition, the outputs go crazy (definitely not correct) with periodic solutions to related terms. I seem to have checked my equations, they seem to be right. Is there anything with ‘odeset’ I can work this out around.
Zero non-axial initial condition:
x0 = [0; 0; -100; 0; 0; 1; 0; 0; 0; 0; 0; 0];

Non-zero IC:
x0 = [0; 0; -100; 0; 1; 1; 0; 0; 0; 0; 0; 0];

0

this is the code i’ve used for these plots. Please feel free to ask any clarifications
%% 6DOF fucntion definition
function [dxdt, aoa_series] = six_DOF_dynamics(t,x)

M_para = 5; %mass of the parachute in Kgs
g = 9.81; %gravitational acceleration in m/s^2
I_x = 1; I_y = 1; I_z = 3; %Inertia matrix diagonal principal of the parachute in kg-m^2
CD_para = 0.9; %coefficient of drag of parachute
CL_para = 0; %coefficient of lift of parachute
rho = 1; %density of the atmosphere in Kg/m^3
x_com = 0; y_com = 0; z_com = -5; %position of COM of the paarchute from the integrated system’s COM(body axis origin) in m
A = 10; %parachute area in m^2

%state variables defining
x_pos = x(1); y_pos = x(2); z_pos = x(3); %postion variables
u = x(4); v = x(5); w = x(6); %velocity variables
phi = x(7); theta = x(8); psi = x(9); %attitude/euler angles variables
p = x(10); q = x(11); r = x(12); %euler angular rate variables

%rate of change of postion -> velocity
dx_pos = u;
dy_pos = v;
dz_pos = w;

%Gravity moment of the parachute along x,y and z components
GravityMoment_para_x = M_para*g*cos(theta)*(y_com*cos(phi) – z_com*sin(phi));
GravityMoment_para_y = M_para*g*(z_com*sin(theta) – x_com*cos(theta)*cos(phi));
GravityMoment_para_z = M_para*g*(x_com*cos(theta)*sin(phi) + y_com*sin(theta));

V_total = sqrt(u^2 + v^2 + w^2); %net speed

persistent angle_of_attack %defining this variable to take it as a output outside ode45

alpha = acos(w/V_total); %angle of attack calculated in radians

angle_of_attack = [angle_of_attack; alpha]; %concatinating current iterated angle of attack in global variable

%Componetize CD and CL into tangential and normal force coefficients
CT = CD_para*cos(alpha) – CL_para*sin(alpha);
CN = CL_para*cos(alpha) + CD_para*sin(alpha);

% Aerodynamic force along x,y,z axes F = 0.5*rho*v^2*A
if (u == 0 & v == 0)
Aerodynamic_force_x = 0;
Aerodynamic_force_y = 0;
Aerodynamic_force_z = -0.5*rho*CT*V_total*V_total*A;
else
Aerodynamic_force_x = -0.5*rho*CN*V_total*V_total*A*(u/(u^2 + v^2)^0.5);
Aerodynamic_force_y = -0.5*rho*CN*V_total*V_total*A*(v/(u^2 + v^2)^0.5);
Aerodynamic_force_z = -0.5*rho*CT*V_total*V_total*A;
end

%Moment caused by aerodynamic force along x,y,x
Aerodynamic_moment_x = y_com*Aerodynamic_force_z – z_com*Aerodynamic_force_y;
Aerodynamic_moment_y = z_com*Aerodynamic_force_x – x_com*Aerodynamic_force_z;
Aerodynamic_moment_z = x_com*Aerodynamic_force_y – y_com*Aerodynamic_force_x;

%rate change of velocity definition
du = v*r – w*q + g*sin(theta) + Aerodynamic_force_x/M_para;
dv = -u*r + w*p + g*cos(theta)*sin(phi) + Aerodynamic_force_y/M_para;
dw = u*q – v*p + g*cos(theta)*cos(phi) + Aerodynamic_force_z/M_para;

%rate change of euler angles definition
dphi = p + (r*cos(phi) + q*sin(phi))*tan(theta);
dtheta = q*cos(phi) – r*sin(phi);
dpsi = (cos(theta)^-1)*(r*cos(phi) + q*sin(phi));

%rate change of angular rates definition
dp = (1/I_x)*(Aerodynamic_moment_x + GravityMoment_para_x);
dq = (1/I_y)*(Aerodynamic_moment_y + GravityMoment_para_y);
dr = (1/I_z)*(Aerodynamic_moment_z + GravityMoment_para_z);

dxdt = [dx_pos; dy_pos; dz_pos; du; dv; dw; dphi; dtheta; dpsi; dp; dq; dr];

%block to output angle of attack
if nargout>1
aoa_series = angle_of_attack;

end

end

tspan = [0 3];

x0 = [0; 0; -100; 0; 1; 1; 0; 0; 0; 0; 0; 0];

[t, x] = ode45(@six_DOF_dynamics, tspan, x0);

[~, aoa_series] = six_DOF_dynamics(t,x);

%adding alpha as a part of x

% Plot results
figure;
subplot(4,1,1);
plot(t, x(:,1:3));
title(‘Position’);
legend(‘x’, ‘y’, ‘z’);

subplot(4,1,2);
plot(t, x(:,4:6));
title(‘Velocity’);
legend(‘u’, ‘v’, ‘w’);

subplot(4,1,3);
plot(t, x(:,7:9));
title(‘Euler Angles’);
legend(‘phi’, ‘theta’, ‘psi’);

subplot(4,1,4);
plot(t, x(:,10:12));
title(‘Euler Rates’);
legend(‘p’, ‘q’, ‘r’); So I’ve been trying to simulate a parachute descent (very basic model, almost a mass body with a coefficient of drag), from a height. While simulating only in the axial direction(z axis), the results coincide with kinematics equations. However, when I give a non-zero initial condition, the outputs go crazy (definitely not correct) with periodic solutions to related terms. I seem to have checked my equations, they seem to be right. Is there anything with ‘odeset’ I can work this out around.
Zero non-axial initial condition:
x0 = [0; 0; -100; 0; 0; 1; 0; 0; 0; 0; 0; 0];

Non-zero IC:
x0 = [0; 0; -100; 0; 1; 1; 0; 0; 0; 0; 0; 0];

0

this is the code i’ve used for these plots. Please feel free to ask any clarifications
%% 6DOF fucntion definition
function [dxdt, aoa_series] = six_DOF_dynamics(t,x)

M_para = 5; %mass of the parachute in Kgs
g = 9.81; %gravitational acceleration in m/s^2
I_x = 1; I_y = 1; I_z = 3; %Inertia matrix diagonal principal of the parachute in kg-m^2
CD_para = 0.9; %coefficient of drag of parachute
CL_para = 0; %coefficient of lift of parachute
rho = 1; %density of the atmosphere in Kg/m^3
x_com = 0; y_com = 0; z_com = -5; %position of COM of the paarchute from the integrated system’s COM(body axis origin) in m
A = 10; %parachute area in m^2

%state variables defining
x_pos = x(1); y_pos = x(2); z_pos = x(3); %postion variables
u = x(4); v = x(5); w = x(6); %velocity variables
phi = x(7); theta = x(8); psi = x(9); %attitude/euler angles variables
p = x(10); q = x(11); r = x(12); %euler angular rate variables

%rate of change of postion -> velocity
dx_pos = u;
dy_pos = v;
dz_pos = w;

%Gravity moment of the parachute along x,y and z components
GravityMoment_para_x = M_para*g*cos(theta)*(y_com*cos(phi) – z_com*sin(phi));
GravityMoment_para_y = M_para*g*(z_com*sin(theta) – x_com*cos(theta)*cos(phi));
GravityMoment_para_z = M_para*g*(x_com*cos(theta)*sin(phi) + y_com*sin(theta));

V_total = sqrt(u^2 + v^2 + w^2); %net speed

persistent angle_of_attack %defining this variable to take it as a output outside ode45

alpha = acos(w/V_total); %angle of attack calculated in radians

angle_of_attack = [angle_of_attack; alpha]; %concatinating current iterated angle of attack in global variable

%Componetize CD and CL into tangential and normal force coefficients
CT = CD_para*cos(alpha) – CL_para*sin(alpha);
CN = CL_para*cos(alpha) + CD_para*sin(alpha);

% Aerodynamic force along x,y,z axes F = 0.5*rho*v^2*A
if (u == 0 & v == 0)
Aerodynamic_force_x = 0;
Aerodynamic_force_y = 0;
Aerodynamic_force_z = -0.5*rho*CT*V_total*V_total*A;
else
Aerodynamic_force_x = -0.5*rho*CN*V_total*V_total*A*(u/(u^2 + v^2)^0.5);
Aerodynamic_force_y = -0.5*rho*CN*V_total*V_total*A*(v/(u^2 + v^2)^0.5);
Aerodynamic_force_z = -0.5*rho*CT*V_total*V_total*A;
end

%Moment caused by aerodynamic force along x,y,x
Aerodynamic_moment_x = y_com*Aerodynamic_force_z – z_com*Aerodynamic_force_y;
Aerodynamic_moment_y = z_com*Aerodynamic_force_x – x_com*Aerodynamic_force_z;
Aerodynamic_moment_z = x_com*Aerodynamic_force_y – y_com*Aerodynamic_force_x;

%rate change of velocity definition
du = v*r – w*q + g*sin(theta) + Aerodynamic_force_x/M_para;
dv = -u*r + w*p + g*cos(theta)*sin(phi) + Aerodynamic_force_y/M_para;
dw = u*q – v*p + g*cos(theta)*cos(phi) + Aerodynamic_force_z/M_para;

%rate change of euler angles definition
dphi = p + (r*cos(phi) + q*sin(phi))*tan(theta);
dtheta = q*cos(phi) – r*sin(phi);
dpsi = (cos(theta)^-1)*(r*cos(phi) + q*sin(phi));

%rate change of angular rates definition
dp = (1/I_x)*(Aerodynamic_moment_x + GravityMoment_para_x);
dq = (1/I_y)*(Aerodynamic_moment_y + GravityMoment_para_y);
dr = (1/I_z)*(Aerodynamic_moment_z + GravityMoment_para_z);

dxdt = [dx_pos; dy_pos; dz_pos; du; dv; dw; dphi; dtheta; dpsi; dp; dq; dr];

%block to output angle of attack
if nargout>1
aoa_series = angle_of_attack;

end

end

tspan = [0 3];

x0 = [0; 0; -100; 0; 1; 1; 0; 0; 0; 0; 0; 0];

[t, x] = ode45(@six_DOF_dynamics, tspan, x0);

[~, aoa_series] = six_DOF_dynamics(t,x);

%adding alpha as a part of x

% Plot results
figure;
subplot(4,1,1);
plot(t, x(:,1:3));
title(‘Position’);
legend(‘x’, ‘y’, ‘z’);

subplot(4,1,2);
plot(t, x(:,4:6));
title(‘Velocity’);
legend(‘u’, ‘v’, ‘w’);

subplot(4,1,3);
plot(t, x(:,7:9));
title(‘Euler Angles’);
legend(‘phi’, ‘theta’, ‘psi’);

subplot(4,1,4);
plot(t, x(:,10:12));
title(‘Euler Rates’);
legend(‘p’, ‘q’, ‘r’); 6dof, parachute, periodic MATLAB Answers — New Questions

​