<p>The Gaussian and Laguerre orthogonal ensembles are fundamental to random matrix theory, and the marginal eigenvalue distributions are basic observable quantities often relevant to applications. Notwithstanding a long history, a formulation providing high precision numerical evaluations for <i>N</i> large enough to probe asymptotic regimes, has not been provided. An exception is for the largest eigenvalue, where there is a formalism due to Chiani which uses a combination of a formulation as a positive real valued Pfaffian, and a recursive computation of the matrix elements. We augment this strategy by introducing a generating function for the conditioned gap probabilities. A finite Fourier series approach is then used to extract the sequence of marginal eigenvalue distributions as a linear combination of complex valued Pfaffians, with the latter then evaluated using an efficient numerical procedure available in the literature due to Wimmer. Applications are given to illustrating various asymptotic formulas, local central limit theorems, and central limit theorems, as well as to probing finite size corrections. Further, our data indicates that the mean values of the marginal distributions interlace with the zeros of the Hermite polynomial (Gaussian ensemble) and a particular Laguerre polynomial (Laguerre ensemble).</p>