We present a method for scattered data approximation with subdivision surfaces which actually uses the true representation of the limit surface as a linear combination of smooth basis functions associated with the control vertices. This is unlike previous techniques which used only piecewise linear approximations of the limit surface. By this we can assign arbitrary parameterizations to the given sample points, including those generated by parameter correction. We present a robust and fast algorithm for exact closest point search on Loop surfaces by combining Newton iteration and non-linear minimization. Based on this we perform unconditionally convergent parameter correction to optimize the approximation with respect to the L^2 metric and thus we make a well-established scattered data fitting technique which has been available before only for B-spline surfaces, applicable to subdivision surfaces. Further we exploit the fact that the control mesh of a subdivision surface can have arbitrary connectivity to reduce the L^\infty error up to a certain user-defined tolerance by adaptively restructuring the control mesh. By employing iterative least squares solvers, we achieve acceptable running times even for large amounts of data and we obtain high quality approximations by surfaces with relatively low control mesh complexity compared to the number of sample points. Since we are using plain subdivision surfaces, there is no need for multiresolution detail coefficients and we do not have to deal with the additional overhead in data and computational complexity associated with them.