and Lambert scattering laws was used (see, e.g., Kaasalainen et al., 2002).

The quadratic drift in mean anomaly, $\Delta M_d$, was fitted as an independent parameter. It is the coefficient in the second term of the expansion of the time-variable mean anomaly: $$M (t) = M (t_0) + n (t-t_0) + \Delta M_d (t - t_0)^2, \label{dMd1}$$ where $$\Delta M_d = \frac{1}{2} \dot{n}, \label{dMd2}$$ where $n$ is the mean motion, $t$ is the time, and $t_0$ is the epoch. $\Delta M_d$ was stepped from $-15$ to $+15$ deg/yr$^2$ in the case of 1999 KW4 and from $-9$ to $+39$ deg/yr$^2$ in the case of 2001 SL9 and all other parameters were fitted at each step.[^1]

To reduce the complexity of the model, we estimated upper limits on the eccentricity of the mutual orbits by fitting the data from the best covered apparitions: the 2001 apparition for 1999 KW4 and the 2013 apparition for 2001 SL9. The model includes a precession of the line of apsides. The pericenter drift rate depends on the polar flattening of the primary (see Murray and Dermott, 1999, Eq. (6.249)), but as the polar flattenings are poorly constrained from the data (see Tables \[tablePropKW4\] and \[tablePropSL9\]), 