function jac = rbfjacob(net, x)

%RBFJACOB Evaluate derivatives of RBF network outputs with respect to inputs.

%

%	Description

%	G = RBFJACOB(NET, X) takes a network data structure NET and a matrix

%	of input vectors X and returns a three-index matrix G whose I, J, K

%	element contains the derivative of network output K with respect to

%	input parameter J for input pattern I.

%

%	See also

%	RBF, RBFGRAD, RBFBKP

%



%	Copyright (c) Ian T Nabney (1996-2001)



% Check arguments for consistency

errstring = consist(net, 'rbf', x);

if ~isempty(errstring);

  error(errstring);

end



if ~strcmp(net.outfn, 'linear')

  error('Function only implemented for linear outputs')

end



[y, z, n2] = rbffwd(net, x);



ndata = size(x, 1);

jac = zeros(ndata, net.nin, net.nout);

Psi = zeros(net.nin, net.nhidden);

% Calculate derivative of activations wrt n2

switch net.actfn

case 'gaussian'

  dz = -z./(ones(ndata, 1)*net.wi);

case 'tps'

  dz = 2*(1 + log(n2+(n2==0)));

case 'r4logr'

  dz = 2*(n2.*(1+2.*log(n2+(n2==0))));

otherwise

   error(['Unknown activation function ', net.actfn]);

end



% Ignore biases as they cannot affect Jacobian

for n = 1:ndata

  Psi = (ones(net.nin, 1)*dz(n, :)).* ...

    (x(n, :)'*ones(1, net.nhidden) - net.c');

  % Now compute the Jacobian

  jac(n, :, :) =  Psi * net.w2;

end



end