
figure(1)
clf

lw=2;

%doves.dat
%duckweed.dat
%houston.dat
%kenya.dat
%usapop.dat
%worldpop.dat

data=load('worldpop.dat');X=data(:,1);Y=data(:,2);
subplot(2,2,1)
X=X+1650;
t0=min(X); t1=max(X);
T=t0:(t1-t0)/1000:t1;
R1_3=log(Y(3)/Y(1))/(X(3)-X(1))
R1_5=log(Y(5)/Y(1))/(X(5)-X(1))
R1_9=log(Y(9)/Y(1))/(X(9)-X(1))
plot(X,Y,'o',T,Y(1)*exp(R1_3*(T-X(1))), ...
             T,Y(1)*exp(R1_5*(T-X(1))), ...
             T,Y(1)*exp(R1_9*(T-X(1))), ...
             'LineWidth',lw)
title('World pop')
legend('data','R1-3','R1-5','R1-9',2)
axis([t0 t1 0 max(Y)])
set(gca,'LineWidth',1.5);

data=load('usapop.dat');X=data(:,1);Y=data(:,2);
subplot(2,2,2)
X=X+1790;
t0=min(X); t1=max(X);
T=t0:(t1-t0)/1000:t1;
R1_3=log(Y(3)/Y(1))/(X(3)-X(1))
R1_5=log(Y(5)/Y(1))/(X(5)-X(1))
R1_9=log(Y(9)/Y(1))/(X(9)-X(1))
plot(X,Y,'o',T,Y(1)*exp(R1_3*(T-X(1))), ...
             T,Y(1)*exp(R1_5*(T-X(1))), ...
             T,Y(1)*exp(R1_9*(T-X(1))), ...
             'LineWidth',lw)
title('USA pop')
legend('data','R1-3','R1-5','R1-9',2)
axis([t0 t1 0 max(Y)])
set(gca,'LineWidth',1.5);

data=load('houston.dat');X=data(:,1);Y=data(:,2);
subplot(2,2,3)
t0=min(X); t1=max(X);
T=t0:(t1-t0)/1000:t1;
R1_3=log(Y(3)/Y(1))/(X(3)-X(1))
R1_5=log(Y(5)/Y(1))/(X(5)-X(1))
R1_9=log(Y(9)/Y(1))/(X(9)-X(1))
plot(X,Y,'o',T,Y(1)*exp(R1_3*(T-X(1))), ...
             T,Y(1)*exp(R1_5*(T-X(1))), ...
             T,Y(1)*exp(R1_9*(T-X(1))), ...
             'LineWidth',lw)
title('Houston pop')
legend('data','R1-3','R1-5','R1-9',2)
axis([t0 t1 0 max(Y)])
set(gca,'LineWidth',1.5);

data=load('doves.dat');X=data(:,1);Y=data(:,2);
subplot(2,2,4)
t0=min(X); t1=max(X);
T=t0:(t1-t0)/1000:t1;
R1_3=log(Y(3)/Y(1))/(X(3)-X(1))
R1_5=log(Y(5)/Y(1))/(X(5)-X(1))
R1_9=log(Y(9)/Y(1))/(X(9)-X(1))
plot(X,Y,'o',T,Y(1)*exp(R1_3*(T-X(1))), ...
             T,Y(1)*exp(R1_5*(T-X(1))), ...
             T,Y(1)*exp(R1_9*(T-X(1))), ...
             'LineWidth',lw)
title('Doves')
legend('data','R1-3','R1-5','R1-9',2)
axis([t0 t1 0 max(Y)])
set(gca,'LineWidth',1.5);

%print -djpeg  algebraic_fit_exp.jpg

% zoom:

%subplot(2,2,1);axis([1650 1900 0.5 1.7]);
%subplot(2,2,2);axis([1790 1900 0 80]);
%subplot(2,2,4);axis([1955 1961 0 60]);

%print -djpeg  algebraic_fit_exp_zoom.jpg
