Over the last several years, a new theory of Nevanlinna-Pick interpolation with complexity constraint has been developed for scalar interpolants. In this paper we generalize this theory to the matrix-valued, case, also allowing for multiple interpolation points. We parameterize a class of interpolants consisting of "most interpotants" of no higher degree than the central solution in terms of spectral zeros. This is a complete parameterization, and for each choice of interpolant we provide a convex optimization problem for determining it. This is derived in the context of duality theory of mathematical programming. To solve the convex optimization problem, we employ a homotopy continuation technique previously developed for the scalar case. These results can be applied to many classes of engineering problems,and, to illustrate this, we provide some,examples. In particular, we apply our method to a benchmark problem in multivariate robust control. By constructing a controller satisfying all design specifications but having only half the McMillan degree of conventional H-infinity controllers, we demonstrate the advantage of the proposed method.