How should I interpret the output of numpy.fft.rfft2?
Also note the ordering of the coefficients in the fft output:
According to the doc: by default the 1st element is the coefficient for 0 frequency component (effectively the sum or mean of the array), and starting from the 2nd we have coeffcients for the postive frequencies in increasing order, and starts from n/2+1 they are for negative frequencies in decreasing order. To have a view of the frequencies for a length-10 array:

np.fft.fftfreq(10)

the output is:

array([ 0. ,  0.1,  0.2,  0.3,  0.4, -0.5, -0.4, -0.3, -0.2, -0.1])

Use np.fft.fftshift(cf), where cf=np.fft.fft(array), the output is shifted so that it corresponds to this frequency ordering:

array([-0.5, -0.4, -0.3, -0.2, -0.1,  0. ,  0.1,  0.2,  0.3,  0.4])

which makes for sense for plotting.
In the 2D case it is the same. And the fft2 and rfft2 difference is as explained by others.
Obviously the rfft2 function simply computes the discrete fft of the input matrix. However how do I interpret a given index of the output? Given an index of the output, which Fourier coefficient am I looking at?
I am especially confused by the sizes of the output. For an n by n matrix, the output seems to be an n by (n/2)+1 matrix (for even n). Why does a square matrix ends up with a non-square fourier transform?