A class of n-dimensional Poisson systems reducible to an unperturbed harmonic oscillator shall be considered. In such case, perturbations leaving invariant a given symplectic leaf shall be investigated. Our purpose will be to analyze the bifurcation phenomena of periodic orbits as a result of these perturbations in the period annulus associated to the unperturbed harmonic oscillator. This is accomplished via the averaging theory up to an arbitrary order in the perturbation parameter ε. In that theory we shall also use both branching theory and singularity theory of smooth maps to analyze the bifurcation phenomena at points where the implicit function theorem is not applicable. When the perturbation is given by a polynomial family, the associated Melnikov functions are polynomial and tools of computational algebra based on Grobner basis are employed in order to ¨ reduce the generators of some polynomial ideals needed to analyze the bifurcation problem. When the most general perturbation of the harmonic oscillator by a quadratic perturbation field is considered, the complete bifurcation diagram (except at a high codimension subset) in the parameter space is obtained. Examples are given.