'Fit data to this model.
Parameters
model : `~astropy.modeling.FittableModel`
model to fit to x, y, z
x : array
input coordinates
y : array
input coordinates
z : array (optional)
input coordinates
weights : array (optional)
weights
maxiter : int
maximum number of iterations
acc : float
Relative error desired in the approximate solution
epsilon : float
A suitable step length for the forward-difference
approximation of the Jacobian (if model.fjac=None). If
epsfcn is less than the machine precision, it is
assumed that the relative errors in the functions are
of the order of the machine precision.
estimate_jacobian : bool
If False (default) and if the model has a fit_deriv method,
it will be used. Otherwise the Jacobian will be estimated.
If True, the Jacobian will be estimated in any case.
Returns
model_copy : `~astropy.modeling.FittableModel`
a copy of the input model with parameters set by the fitter'

@fitter_unit_support
def __call__(self, model, x, y, z=None, weights=None, maxiter=DEFAULT_MAXITER, acc=DEFAULT_ACC, epsilon=DEFAULT_EPS, estimate_jacobian=False):

    from scipy import optimize
    model_copy = _validate_model(model, self.supported_constraints)
    farg = ((model_copy, weights) + _convert_input(x, y, z))
    if ((model_copy.fit_deriv is None) or estimate_jacobian):
        dfunc = None
    else:
        dfunc = self._wrap_deriv
    (init_values, _) = _model_to_fit_params(model_copy)
    (fitparams, cov_x, dinfo, mess, ierr) = optimize.leastsq(self.objective_function, init_values, args=farg, Dfun=dfunc, col_deriv=model_copy.col_fit_deriv, maxfev=maxiter, epsfcn=epsilon, xtol=acc, full_output=True)
    _fitter_to_model_params(model_copy, fitparams)
    self.fit_info.update(dinfo)
    self.fit_info[u'cov_x'] = cov_x
    self.fit_info[u'message'] = mess
    self.fit_info[u'ierr'] = ierr
    if (ierr not in [1, 2, 3, 4]):
        warnings.warn(u"The    fit    may    be    unsuccessful;    check    fit_info['message']    for    more    information.", AstropyUserWarning)
    if ((len(y) > len(init_values)) and (cov_x is not None)):
        sum_sqrs = np.sum((self.objective_function(fitparams, *farg) ** 2))
        dof = (len(y) - len(init_values))
        self.fit_info[u'param_cov'] = ((cov_x * sum_sqrs) / dof)
    else:
        self.fit_info[u'param_cov'] = None
    return model_copy
