We study general parabolic equations of the form ut = divA(x, t, u, Du) + div (|F|Fp-2) + f whose principal part depends on the solution itself. The vector field A is assumed to have small mean oscillation in x, measurable in t, Lipschitz continuous in u, and its growth in Du is like the p-Laplace operator. We establish interior Calderon-Zygmund estimates for locally bounded weak solutions to the equations when p > 2n/(n + 2). This is achieved by employing a perturbation method together with developing a twoparameter technique and a new compactness argument. We also make crucial use of the intrinsic geometry method by DiBenedetto (Degenerate parabolic equations, Springer, New York, 1993) and the maximal function free approach by Acerbi and Mingione