Multibody simulation software has entered the industrial user's market several years ago and is nowadays part of a CAE engineer's everyday life. Due to usability and time advantages in the modeling process, such a simulation software is typically based on modeling strategies using a redundant set of coordinates. Although such a redundant set simplifies the modeling process, modeling without a minimal set of coordinates forces the multibody dynamics into differential algebraic equations (DAEs). Unnecessary degrees of freedom (DOFs) are introduced easily through kinematic relations and especially when talking about planar problems, modeled in a 3D space, the number of unnecessary DOFs increases rapidly. Further, unnecessary constraint equations which address such kinematic chains or unstressed DOFs are typically introduced. The d-index three DAE systems of second order considered in this contribution are modeled in the freeware multibody simulation software FreeDyn using Euler parameters and solved in the open source software Scilab, applying a HHT solver algorithm. Due to the characteristic nonlinear system matrices, model reduction is not possible a-priori but needs at least one full system forward simulation to rest on. Inspired by repetitive simulations necessary for e.g. parameter identification, we focus on reducing such a highly redundant set of coordinates close to a minimal set of coordinates representation by eliminating both, physical and constraint DOFs. This contribution first gives a brief introduction to the underlying DAE system. In the second part model reduction of redundant physical DOFs using a modied linear Galerkin projection based on a velocity level Proper Orthogonal Decomposition (POD) is proposed. The subspace found by POD is further applied to identify unneeded constraint equations related to a.) - constraints addressing unneeded physical DOFs and b.) constraints which are orthogonal to the found subspace and, hence, satisfied by the projection matrix. The method is illustrated using a simple analytical example and further applied to a rigid V8 crankdrive consisting of 119 DOFs as well as 118 constraint equations. In the conclusion the proposed methods are discussed and an outlook to parameter identification, using a reduced adjoint system, and therewith arising special needs to the underlying subspace is given.