Multiplying a list of matrices

I am attempting to create a 'Abeles matrix formalism' model to analyse some experimental data - I have attached a wiki link to this for reference so that you can see what I an attempting to achieve. http://en.wikipedia.org/wiki/Transfer-matrix_method_(optics)#Abeles_matrix_formalism

The crux of my issue is that I am unable to multiply four sets of matrices against each other, as in: A[1]*B[1]*C[1]*D[1], A[2]*B[2]*C[2]*D[2], ..., A[n]*B[n]*C[n]*D[n]. I then need to store the results as individual matrices of their own - each matrix represents the corresponding momentum transfer value from Qmin:Qstep:Qmax.

I believe I have sorted this out to an extent by using the BSX function, i.e AB = bsxfun(@times,A,B), ABC = bsxfun(@times,AB,C), ABCD = bsxfun(@times,ABC,D), although when I attempt to carry out the final step; R = abs((ABCD(2,1)./ABCD(1,1)).^2) I end up with a single value rather than a value of R for each Q value.

I have heard some suggestions about creating cell arrays ect. although I am not sure how to approach this. Any help or recommended reading will be appreciated.

my 'test' code is: %import data fid = fopen('run_22208_09.dat'); A = textscan(fid,'%f%f%f',270,'headerlines',0,'delimiter',','); %A = textscan(fid,'%f%f%f',270);

NQ = size(A{1,1}); NQ = NQ(1); Qmin = A{1,1}(1); Qmax = A{1,1}(NQ); Qstep = A{1,1}(2) - A{1,1}(1); fclose('all');

s0 = 2e-6; s1 = 10e-6; s2 = 6e-6; s3 = 4e-6; s4 = 8e-6; sn = 12e-6;

r1 = 2; r2 = 10; r3 = 3; r4 = 7; t1 = 10; t2 = 45; t3 = 5; t4 = 20;

Q=Qmin:Qstep:Qmax;

k = 2.*Q; k1 = (((k).^2) - 4.*pi.*(s1 - s0)).^0.5; k2 = (((k).^2) - 4.*pi.*(s2 - s0)).^0.5; k3 = (((k).^2) - 4.*pi.*(s3 - s0)).^0.5; k4 = (((k).^2) - 4.*pi.*(s4 - s0)).^0.5; kn = (((k).^2) - 4.*pi.*(sn - s0)).^0.5;

% %layer1 layer1 = ((k1 - k2)./(k1 + k2)).*(exp(-2.*k1.*k2.*(r1.^2))); beta1 = (sqrt(-1)).*k1.*t1;

%layer2 layer2 = ((k2 - k3)./(k2 + k3)).*(exp(-2.*k2.*k3.*(r2.^2))); beta2 = (sqrt(-1)).*k2.*t2;

%layer3 layer3 = ((k3 - k4)./(k3 + k4)).*(exp(-2.*k3.*k4.*(r3.^2))); beta3 = (sqrt(-1)).*k3.*t3;

%layer4 layer4 = ((k4 - kn)./(k4 + kn)).*(exp(-2.*k4.*kn.*(r4.^2))); beta4 = (sqrt(-1)).*k4.*t4;

%general matrix C1 = [exp(beta1),layer1.*(exp(beta1));layer1.*exp(-beta1),exp(-beta1)] C2 = [exp(beta2),layer2.*(exp(beta2));layer2.*exp(-beta2),exp(-beta2)]; C3 = [exp(beta3),layer3.*(exp(beta3));layer3.*exp(-beta3),exp(-beta3)]; C4 = [exp(beta4),layer4.*(exp(beta4));layer4.*exp(-beta4),exp(-beta4)];

% CA = bsxfun(@times,C1,C2) % CB = bsxfun(@times,CA,C3); % C = bsxfun(@times,CB,C4)

% R = abs((C(2,1)./C(1,1)).^2) % scalingfactor = 1 % R = R.*scalingfactor.*1e8

[EManStudent]

%import data

fid = fopen('run_22208_09.dat');

A = textscan(fid,'%f%f%f',270,'headerlines',0,'delimiter',',');

NQ = size(A{1,1});

NQ = NQ(1);

Qmin = A{1,1}(1);

Qmax = A{1,1}(NQ);

Qstep = A{1,1}(2) - A{1,1}(1);

fclose('all');

s0 = 2e-6; s1 = 10e-6;

s2 = 6e-6;

s3 = 4e-6;

s4 = 8e-6;

sn = 12e-6;

r1 = 2;

r2 = 10;

r3 = 3;

r4 = 7;

t1 = 10;

t2 = 45;

t3 = 5;

t4 = 20;

Q=Qmin:Qstep:Qmax;

k = 2.*Q;

k1 = (((k).^2) - 4.*pi.*(s1 - s0)).^0.5;

k2 = (((k).^2) - 4.*pi.*(s2 - s0)).^0.5;

k3 = (((k).^2) - 4.*pi.*(s3 - s0)).^0.5;

k4 = (((k).^2) - 4.*pi.*(s4 - s0)).^0.5;

kn = (((k).^2) - 4.*pi.*(sn - s0)).^0.5;

% %layer1

layer1 = ((k1 - k2)./(k1 + k2)).*(exp(-2.*k1.*k2.*(r1.^2)));

beta1 = (sqrt(-1)).*k1.*t1;

%layer2

layer2 = ((k2 - k3)./(k2 + k3)).*(exp(-2.*k2.*k3.*(r2.^2)));

beta2 = (sqrt(-1)).*k2.*t2;

%layer3

layer3 = ((k3 - k4)./(k3 + k4)).*(exp(-2.*k3.*k4.*(r3.^2)));

beta3 = (sqrt(-1)).*k3.*t3;

%layer4

layer4 = ((k4 - kn)./(k4 + kn)).*(exp(-2.*k4.*kn.*(r4.^2)));

beta4 = (sqrt(-1)).*k4.*t4;

%general matrix

C1 = [exp(beta1),layer1.*(exp(beta1));layer1.*exp(-beta1),exp(-beta1)]

C2 = [exp(beta2),layer2.*(exp(beta2));layer2.*exp(-beta2),exp(-beta2)];

C3 = [exp(beta3),layer3.*(exp(beta3));layer3.*exp(-beta3),exp(-beta3)];

C4 = [exp(beta4),layer4.*(exp(beta4));layer4.*exp(-beta4),exp(-beta4)];

% CA = bsxfun(@times,C1,C2)

% CB = bsxfun(@times,CA,C3);

% C = bsxfun(@times,CB,C4)

% R = abs((C(2,1)./C(1,1)).^2)