I would like to use NBOs to probe d-orbital splitting of transition metal centers in various geometries, and I'm running into results I did not expect. I'm trying to find out if I'm asking NBO6 to perform the wrong subroutines or if there's a flaw in my understanding.

I've set up a calculation on Ni(CO)4 in its standard tetrahedral geometry, and run NBO6 with the keywords NPA NBO NLMO PLOT AONBO=W77. From the natural population analysis, nickel's valence natural atomic orbitals are as follows:

The d orbital energies don't work out nicely in a low energy e set (z^2 and x^2-y^2) and a higher energy t2 set (xz, yz, xy). Following the NBO subroutine (total non-Lewis electrons: 1.94857 or 2.320% of 94), the nickel's valence NBOs are as follows:

The NBOs do seem to appear in an e set (18 & 19 correspond to z^2 and x^2-y^2, respectively) and a t2 set (20, 21, & 22 correspond to yz, xz, and xy, respectively); however, the t2 set is lower in energy than the e set, which goes against expectation from crystal field theory. Am I missing something?