function [ssL,sL,d] = chol_band5(ssB,sB,dB)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Compute the Cholesky factorization of a symmetric positive 
% definite banded matrix B with bandwidth 5.  
% The factorization is B = L*D*L' where L is a lower triangular
% matrix with 1s on the diagonal and two non-zero subdiagonals.
% D is a diagonal matrix.
%
% INPUTS:
%      ssB    Nx1 vector of sub-subdiagonal elements of B (uses 1st N-2)
%      sB     Nx1 vector of subdiagonal elements of B (uses 1st N-1)
%      dB     Nx1 vector of diagonal elements of B
%
% OUTPUT:
%      ssL    Nx1 vector containing the sub-sub-diagonal elements of L (uses 1st N-2)
%      sL     Nx1 vector containing the sub-diagonal elements of L (uses 1st N-1)
%      d      Nx1 vector containing the diagonals of D
%
% Written by Karen Chiswell on March 19, 2006
% Based on description of algorithm in Green & Silverman (1994), Sec 2.6
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

N = length(dB);

%       initialize output vectors
d = dB;
sL = sB;
ssL = ssB;

%       equate first three elements of B (B11, B21 and B22)
%       D1 = B11, L21 = B21/D1,  D2 = B22 - L21*L21*D1
d(1) = dB(1);
sL(1) = sB(1)/d(1);
d(2) = dB(2) - sL(1)*sL(1)*d(1);

%       now for the remaining elements
%       Li,i-2 = Bi,i-2/Di-2
%       Li,i-1 = (Bi,i-1 - Li-1,i-2*Li,i-2*Di-2)/Di-1
%       Di = Bii - Li,i-1*Li,i-1*Di-1 - Li,i-2*Li,i-2*Di-2
for i = 3:N
   ssL(i-2) = ssB(i-2)/d(i-2);
   sL(i-1)  = ( sB(i-1) - sL(i-2)*ssL(i-2)*d(i-2) ) / d(i-1);
   d(i)     = dB(i) - sL(i-1)*sL(i-1)*d(i-1) - ssL(i-2)*ssL(i-2)*d(i-2);
end   % loop on i

% end of chol_band5.m