How to compute softmax and its gradient?

I am creating a simple two layer neural network where the activation function of the output layer will be softmax.

I have this for creating softmax in a numerically stable way

function g = softmax(z)

    dim = 1;

    s = ones(1, ndims(z));
    s(dim) = size(z, dim);

    maxz = max(z, [], dim);
    expz = exp(z-repmat(maxz, s));
    g = expz ./ repmat(sum(expz, dim), s);

z is a matrix that contains all of the data calculated by the previous layer one row at a time.

In order to compute the derivative of this though I will need to use the Kronecker delta but I am not sure how to do it.

Can someone provide me with a vectorized implementation for computing it in Matlab?

[usama pervaiz]

You can find here how to compute softmax of a matrix and its gradient http://peterroelants.github.io/posts/neural_network_implementation_intermezzo02/