How do I remove vertical lines in graphs?

I have generated these graphs with this vertical lines, how do I remove them? When I generated with the same code in octave, there were no vertical lines.

Here follow the code I used and the graphs generted:

% Flutter Chapter B04 Appendix
% Sets up the aeroelastic matrices for binary aeroelastic model,
% performs eigenvalue solution at desired speeds and determines the frequencies
% and damping ratios
% plots V_omega and V_g trends
% and plots flutter conic solution
% Initialize variables
clear; close all; clc
% System parameters
EI = 4*10^7;
GJ = 8*10^6;
s = 7.5; % semi span
c = 2; % chord
m = 400; % unit mass / area of wing
mc = 400;
kappa_freq = 1.6602; % flapping freq in Hz
theta_freq = 3.0201; % pitch freq in Hz
betha_freq = 7.6763; % control freq in Hz 
xcm = 0.5*c; % position of centre of mass from nose
xf = 0.4*c; % position of flexural axis from nose
xh = 0.8*c;
e = xf/c - 0.25; % eccentricity between flexural axis and aero centre (1/4 chord)
velstart = 1; % lowest velocity
velend = 200; % maximum velocity
velinc =0.1; % velocity increment
rho = 1.225; % air density
Mthetadot = -1.2;% unsteady aero damping term
Mbethadot = -0.1;
%M = (m*c^2 - 2*m*c*xcm)/(2*xcm); % leading edge mass term
damping_Y_N = 0; % =1 if damping included =0 if not included
%if damping_Y_N == 1
% structural proportional damping inclusion C = alpha*M + beta * K
% then two freqs and damps must be defined
% set dampings to zero for no structural damping
% z1 = 0.0; % critical damping at first frequency
% z2 = 0.0; % critical damping at second frequency
% w1 = 2*2*pi; % first frequency
% w2 = 14*2*pi; % second frequency
% alpha = 2*w1*w2*(-z2*w1 + z1*w2)/ (w1*w1*w2*w2);
% beta = 2*(z2*w2-z1*w1) / (w2*w2 - w1*w1);
%end
% Set up system matrices
% Inertia matrix
a11= m*s*c/5; % I kappa
a12 = (m*s/4)*(((c^2)/2)-(xf*c)); %I kappa theta
a13 = (m*s/3)*((((c^2)-(xh^2))/2)-xh*(c-xh));
a21 = a12;
a22= (m*s/3)*(((c^3)/3)-(xf*c^2) + ((xf^2)*c));
a23 = (m*s/2)*((((c^3) - (xh^3))/3)-(xf+xh)*(((c^2) - (xh^2))/2)+(xf*xh*(c-xh)));
a31 = a13;
a32 = a23;
a33 = m*s*((((c^3) - (xh^3))/3)+((xh^2)*c) - (xh*(c^2)));
A=[a11,a12,a13;a21,a22,a23;a31,a32,a33];
% Structural stiffness matrix
% k1 = (kappa_freq*pi*2)^2*a11; % k kappa heave stiffness
% k2 = (theta_freq*pi*2)^2*a22; % k theta pitch stiffness
% k3 = (betha_freq*pi*2)^2*a33; % k betha control stiffness
k1 = (4*EI)/(s^3);
k2 = (GJ)/(s);
k3 = 1*(10^3)*s;
% %% Non-linearity with k2 = k_alpha*alpha + k2_alpha*alpha^3, where k2_alpha = 100*k_alpha.
% pitch = 0.02; % Pitch angle
% k_alpha = (theta_freq*pi*2)^2*a22;
% k2_alpha = 100*k_alpha;
% k2 = k2_alpha*pitch + k2_alpha*pitch^3;
%%
E = [k1 0 0; 0 k2 0;0 0 k3];
Er = 0.2; % Ratio of the control surface
d = (2*xh/c)-1;
T12 = sqrt(1-(d^2))*(2+d) + acos(d*(2*d+1));
T10 = sqrt(1-d^2) + acos(d);
aw = 2*pi;
bw = e*aw;
cw = -T12/2;
ac = (aw/pi)*(acos(1-2*Er)+2*sqrt(Er*(1-Er)));
bc = -(aw/pi)*(1-Er)*sqrt(Er*(1-Er));
cc = -(T12*T10)/(2*pi);
icount = 0;
Ce = zeros(3,3);
for V = velstart:velinc:velend % loop for different velocities
    icount = icount +1;
    B = rho*V*[c*s*aw/10,0,0;(-c^2*s)*bw/8,-c^3*s*Mthetadot/24,0; -(c^2*s)*cw/6, 0, -c^3*s*Mbethadot/8];
    if damping_Y_N == 0 % damping matrices
        C = B;
    else % =1 if damping included
        C = B + Ce;
        % Aero and structural damping
    end
    K = (rho*V^2*[0,c*s*aw/8,c*s*ac/6; 0,-c^2*s*e*bw/6, c^2*s*bc/4; 0, -c^2*s*cw/4, -c^2*s*cc/2]) + E; %aero / structural stiffness
    Mat = [[0,0,0; 0,0,0; 0,0,0],eye(3); -A\K,-A\C]; % set up 1st order eigenvalue solution matrix
    lambda = eig(Mat); % eigenvalue solution
    % Natural frequencies and damping ratios
    for jj = 1:length(lambda)
        im(jj) = imag(lambda(jj));
        re(jj) = real(lambda(jj));
        freq(jj,icount) = sqrt(re(jj)^2+im(jj)^2)/(2*pi);
        damp(jj,icount) = -100*re(jj)/freq(jj,icount);
        % freq(jj,icount) = freq(jj,icount)/(2*pi); % convert frequency to hertz
    end
    Vel(icount) = V;
end
% Plot frequencies and dampings vs speed
figure(1)
subplot(2,1,1); plot(Vel,freq,'k');
vaxis = axis; xlim = ([0 vaxis(2)]);
xlabel ('Air Speed (m/s) '); ylabel ('Freq (Hz)'); grid
subplot(2,1,2);
plot(Vel,damp/100,'k')
xlim = ([0 vaxis(2)]); axis([xlim ylim]);
xlabel ('Air Speed (m/s) '); ylabel ('Damping Ratio (%)'); grid

Answer 1

[Adam Danz]

If you zoom into one of those vertical lines you'll see that they aren't vertical and they are a part of your data.

Another way to see this is to only plot 2 of the 6 lines in the first subplot with different colors.

plot(Vel,freq([4,6],:)) % lines 4 and 6

If this pattern in your calculated data are not expected, I would start investigating here at this line and work your way up in debug mode.

freq(jj,icount) = sqrt(re(jj)^2+im(jj)^2)/(2*pi);

Answer 2

[Ronit]

The vertical lines appear to be gridlines on the plot. You can remove them by turning off the grid. To do this, you can add the following line of code before your plot commands:

grid off;

This turns off the grid display for all subsequent plots until you run grid on again. You can place this line of code anywhere before your plot commands in the script, but it's recommended to put it before the first plot command.