The problem of scattering of harmonic plane acoustic waves by fluid spheroids (prolate and oblate) is addressed from an analytical approach. Mathematically, it consists in solving the Helmholtz equation in an unbounded domain with Sommerfeld radiation condition at infinity. The domain where propagation takes place is characterised by density and sound speed values $\rho_0$ and $c_0$, respectively, while $\rho_1$ and $c_1$ are the corresponding density and sound speed values of an immersed object that is responsible of the scattered field. Since Helmholtz equation is separable in prolate/oblate spheroidal coordinates, its exact solution for the scattered field can be expressed as an expansion on prolate/oblate spheroidal functions multiplied by coefficients whose values depend upon the boundary conditions verified at the medium-immersed fluid obstacle interface. The general case ($c_0 \neq c_1$) is cumbersome because it requires to solve successive matrix systems that are ill-conditioned when $c_1/c_0$ is far from unity. In this paper, a numerical implementation of the general exact solution that is valid for any range of eccentricity values and for $c_0 \neq c_1$, is provided. The high level solver code has been written in the Julia programming language while a software package recently released in the literature has been used to compute the spheroidal functions. Several limit cases (Dirichlet and Neumann boundary conditions, spheroid tending to sphere) have been satisfactorily verified using the implemented code. The corresponding example scripts can be downloaded from the authors' web (GitHub) site. Additionally, the new code has been used to extend results reported in the literature.