function logistic_orbit(ri,rf,Nsteps,N,N0)
% Orbit diagram for logistic map: parameter values from ri to rf, using
% Nsteps different values of r; eliminate the first N0 steps (transients)
% and plot the next N iterates.
%
step = (rf-ri)/(Nsteps-1);
x = 0.3;          % Initial value
LT = [];  RR = [];
r = ri-step;
for j = 1:Nsteps
    r = r + step;
    x = rand(1);
    Xn = [];    
    for i = 1:N0      % Eliminate transients - first N0 steps
        x1 = r*x*(1-x);
        x = x1;
    end
    for i = 1:N
        x1 = r*x*(1-x);
        Xn = [Xn x1];
        x = x1;
    end
    RR = [RR r];
    LT = [LT , Xn'];
end
figure(2)
clf;
plot(RR',LT','b.','MarkerSize',2)
axis([ri rf 0 1])