global N beta gamma
% Get Data
M=79;
Crange=['C1:C',num2str(M)];
C = xlsread('COVID19data.xlsx', 1, Crange);

% Constants
%t0=0;
%tf=300;
N=3.3E8;
NI=1;

Recovery = 14.0;
gamma = 1.0/Recovery;
R0 = 3.45;
beta = R0*gamma;

%b=[gamma, beta]
%x=1:M;
%R=getsir(b,x)
% need [t,y]=ode15s('sircovf',[t0,tf],[N-NI,NI,0]);

% Nonlinear fit
x=1:M;
b=[gamma, beta];
opts = statset('MaxIter',100);
bnew = nlinfit(x',C,'getsir',b,opts);

gamma=bnew(1);
beta=bnew(2);
N=3.3E8;
NI=1;

[t,y]=ode15s('sircovf',[0,300],[N-NI,NI,0]);
plot(t,y(:,2),1:M,C,'bo')
axis([0,100,0,1E6])
xlabel('Time (days) since Jan 22')
ylabel('Population')
title('Function Fit to COVID-19 Infected Data')