3D matrix time 2D matrix

I need to do the following and I am wondering if they can be done in a much simple and efficient way.

s = size(relu);
M = s(1) * s(2);
N = s(3);

delta = single(zeros(s));
for i = 1 : s(1)
    for j = 1 : s(2)
        for k = 1 : s (3)
            if relu(i, j, k) > 0
                delta(i, j, k) = reshape(relu(i, j, :), [1, s(3)]) * gamma(:, k) / M / M / N / N;
            else
                delta(i, j, k) = 0;
            end
        end
    end
end

Answer 1

[Andrei Bobrov]

Let your gamma -> gamma1 (because gamma - name of the function from MATLAB)

s = size(relu);
M = s(1) * s(2);
N = s(3);

r1 = reshape(relu,[],s(3));
delta1 = reshape(r1*gamma1,s).*(relu > 0)/(M*N)^2;

[fireman]

It worked very well. Thanks,