%% Code for Figure 5.7B
%% Carol Lucas

clear CT Cray

fi=0.01;kserca=0.2;vserca=100;vRyR=5;vleak=.2;sig=.02;
ka=0.4;kb=0.6;kc=.1;kcm=.1;

Cray=[.01:.01:3];

for k=1:length(Cray)
    C=Cray(k);
    winf=(1+(ka/C)^4+(C/kb)^3)/(1+(1/kc)+(ka/C)^4+(C/kb)^3);
    Po=winf*(1+(C/kb)^3)/(1+(ka/C)^4 + (C/kb)^3);
    CT(k)=(C+(C^2*vserca/(C^2+kserca^2))/(vRyR*Po+vleak))*sig+C;
end

figure(1);clf;plot(CT,Cray);xlabel('CT');ylabel('[Ca2+]');title('Figure 5.7 (B)');
axis([0 4 0 3]);