% NOTE
% Before run this MALTAB scripts, you have to calclate 
% transient Ca2+ peak and averaged Ca2+ response with calc_CaRes.m

clear
Ca = load('AveCaRes_IP-7.5_0.1_-4.5Canum0_100_25000.txt');
%Ca = load('CaRes_all.txt');
Ca(:,1) = Ca(:,1) + 6;
Ca(:,4) = log10(Ca(:,4));

IP25 = log10(2.5)
example_points = [  -1,0.06;  -1,1;  -1,2;  -1,3;  -1,4;  -1,5;  ...
                              -1,6;  -1,7;  -1,8;  -1,9;  -1,10; ...
                  IP25,0.06;IP25,1;IP25,2;IP25,3;IP25,4;IP25,5;  ...
                            IP25,6;IP25,7;IP25,8;IP25,9;IP25,10; ...
                    +1,0.06;  +1,1;  +1,2;  +1,3;  +1,4;  +1,5;  ...
                              +1,6;  +1,7;  +1,8;  +1,9;  +1,10     ];

linecolor = [0,   0,   0;    %black
             0.5, 0,   1;    %violet
             0,   0,   1;    %blue
             0,   1,   1;    %cyan
             0,   1,   0;    %dark green
             0.5, 1,   0;    %light green
             1,   1,   0;    %yellow
             1,   0.5, 0;    %orange
             1,   0,   0;    %red
             1,   0,   1;    %magenta
             1,   0.5, 1];   %pink

xlin = linspace(-1.3, 2.3, 331);
ylin = linspace(0, 12, 241);
[X, Y] = meshgrid(xlin, ylin);
Z = griddata(Ca(:,1), Ca(:,3), Ca(:,4), X, Y, 'cubic');
hight = -1.3:0.1:+1.4;
[C, h, CF] = contourf(X, Y, Z, hight) %interpolatedaxis tight;
hold on
for i = 1:11
semilogx(example_points(i,1), example_points(i,2), '.','Color',linecolor(i,:),'MarkerSize',24)
semilogx(example_points(i+11,1), example_points(i+11,2), '.','Color',linecolor(i,:),'MarkerSize',24)
semilogx(example_points(i+22,1), example_points(i+22,2), '.','Color',linecolor(i,:),'MarkerSize',24)
end
axis([log10(0.05) log10(20) 0 12])
set(gca,'xtick',[log10(0.05) log10(0.06) log10(0.07) log10(0.08) log10(0.09) log10(0.1) log10(0.2) log10(0.3) log10(0.4) log10(0.5) log10(0.6) log10(0.7) log10(0.8) log10(0.9) log10(1) log10(2) log10(3) log10(4) log10(5) log10(6) log10(7) log10(8) log10(9) log10(10) log10(20)]);
set(gca,'ytick',[0 1 2 3 4 5 6 7 8 9 10 11 12]);
set(gca,'XTickLabel','|||||0.1|||||||||1|||||||||10|')
xlabel('[IP3] (\muM)')
ylabel('Transient Ca2+ Peak (\muM)')
title('Averaged Ca2+ Response (\muM)')
hold on
