How to get the energy of the orbitals?

How to get number and energy of occupation and virtual orbitals in RHF method?

The following Psi4NumPy tutorial should help you: https://github.com/psi4/psi4numpy/blob/master/Moller-Plesset/MP2.py#L77

The MP2 method requires these quantities so they are naturally obtained in this python script.