By integrating the truncated complex scalar gravitational motion equations for an anelastic, rotating, slightly elliptical Earth, the complex frequency dependent Earth transfer functions are computed directly. Unlike the conventional method, the effects of both oceanic loads and tidal currents are included via outer surface boundary conditions, all of which are expanded to second order in ellipticity. A modified ellipticity profile in second order accuracy for the non-hydrostatic Earth is obtained from Clairaut’s equation and the PREM Earth model by adjusting both the ellipticity of the core-mantle boundary and the global dynamical ellipticity to modern observations. The effects of different Earth models, anelastic models, and ocean models are computed and compared. The atmospheric contributions to prograde annual, retrograde annual and retrograde semiannual nutation are also included as oceanic effects. Finally, a complete new nutation series of more than 340 periods, including in-phase and out-of-phase parts of longitude and obliquity terms, for a more realistic Earth, is obtained and compared with other available nutation series and observations.