function [BufferCa0, BufferCa2, Calcium]  = Get_ApproximantByName_AnyRestingCa(Method, lambda2, nu2, eps, gamma, CaInfty, lineType)

% Compute an approximation for concentrations of
%        Ca2+ and 2:1 buffer near a point channel / source
%        Case of nonzero background / resting [Ca2+] at infinity, CaInfty
%
% Input variables: 
%  (1)   Method string: {'RBA', 'ExpPade', 'ExpExp'}
%  (2-6) Dimensional parameters: lambda2, nu2, eps, gamma, CaInfty
%  (7)   lineType: MATLAB plot linetype string, e.g. 'b--'
%
% Return arrays: BufferCa0 = [Free buffer]        vs. distance from channel 
%                BufferCa2 = [Fully bound buffer] vs. distance from channel 
%                Calcium   = [Ca2+]               vs. distance from channel
%                R         = Array of distances from channel see Eq. 28)
%
%           Victor Matveev and Yinbo Chen (NJIT) March 2021
%  (GPL) Free to use, share, modify, with appropriate attribution
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

R    = logspace(-3, 2, 100);

bT = 1 + eps * CaInfty*(2 + CaInfty);
cT = CaInfty * ( 1 + nu2*eps*(1 + CaInfty) );

lambda1 = lambda2*gamma;
nu1     = eps*nu2;
q       = 1 / (bT + nu1*(1 + CaInfty*(2 + eps*CaInfty)));

uInf = CaInfty*(2 + CaInfty);
vInf = CaInfty^2;

u1x  = 2*q*(1 + CaInfty);
v1x  = 2*q*CaInfty*(1 + eps*CaInfty);

a1        = sqrt( (eps/lambda1)^2 + 1/(q*lambda1*(1 + CaInfty)) ) - eps/lambda1;
u         = uInf + u1x * (1 - exp(-a1*R)) ./ R;
BufferCa0 = bT - eps*u;

u0   = uInf + a1*u1x;
F    = (u0 - vInf) / v1x;

[BCa0_RBA, BCa2_RBA, Ca_RBA] = Get_RBA_AnyRestingCa(nu2, eps, gamma, CaInfty, '');

switch Method

    case {'RBA', 1}
	
		BufferCa0 = BCa0_RBA;
		BufferCa2 = BCa2_RBA;
		
    case {'ExpPade', 2}
        
        B = ( v1x + sqrt(1 + 8*lambda2*F) ) / (2*F);
        v = vInf + v1x ./ (B + R);
		BufferCa2 = eps*v;
		
    case {'ExpExp', 3} 
        
        a2 = (-1 + sqrt(1 + 4*lambda2*F) )/ (2*lambda2);
        v = vInf + v1x*(1 - exp(-a2*R)) ./ R;
		BufferCa2 = eps*v;
        
    otherwise
            disp([Method, 'Is an unknown method.']);
            return;
end

BufferCa1 = bT - BufferCa2 - BufferCa0;
Calcium   = 1./R + cT - nu2*(BufferCa1 + 2*BufferCa2)/2;

IND = (Calcium < Ca_RBA);
if numel(IND)
    Calcium(IND) = Ca_RBA(IND);
end

if ~isempty(lineType)
    
    subplot(411); hold on;
    set(gca, 'fontsize', 14);
    title(['[Free buffer] b (\lambda_2=', num2str(lambda2), '\nu_2=', num2str(nu2), ' \epsilon=', num2str(eps), ' \gamma=', num2str(gamma),' c_\infty=', num2str(CaInfty), ')']);
    plot(R, BufferCa0, lineType, 'linewidth', 2);
    set(gca, 'xscale', 'log');
    xlabel('r'); ylabel('b');
    axis tight;

    subplot(412); hold on;
    set(gca, 'fontsize', 14);
    title(['[Single-bound buffer] b*']);
    plot(R, BufferCa1, lineType, 'linewidth', 2);
    set(gca, 'xscale', 'log');
    xlabel('r'); ylabel('b*');
    axis tight;
    
    subplot(413); hold on;
    set(gca, 'fontsize', 14);
    title(['[Fully-bound buffer] b**']);
    plot(R, BufferCa2, lineType, 'linewidth', 2);
    set(gca, 'xscale', 'log');
    xlabel('r'); ylabel('b**');
    axis tight;
    
    subplot(414); hold on;
    set(gca, 'fontsize', 14);
    title(['[Ca2+] Calcium']);
    plot(R, Calcium, lineType, 'linewidth', 2);
    set(gca, 'xscale', 'log', 'yscale', 'log');
    xlabel('r'); ylabel('Calcium');
    axis tight;
end
