let $\text{C} : \mathbb{R} \to \mathbb{R^2}, \quad t \mapsto (x_C(t),y_C(t))$ be the parametrization of our circle. Then:
$$\begin{pmatrix}
x_C(t) \\
y_C(t)
\end{pmatrix}=\begin{pmatrix}
a + R\cos t \\
b + R\sin t
\end{pmatrix}$$
Define also the function $f : \mathbb{R^2} \to \mathbb{R}, \quad (x,y) \mapsto (x(t) − a)^2 + (y(t) − b)^2 - R^2$, so that the Circle $C$ is also defined by equation $f(x,y)=0$.
following @N.Ciccoli : $C$ and $\gamma$ have $n^{th}$ order of contact at $t_0$ means that $(f\circ \gamma)(t)$ is $o((t − t_0 )^n )$ as $t \longrightarrow t_0$.
from my understanding, reason is: since $x(t)$ and $y(t)$ are both assumed twice differentiable and not vanishing together:
- if $C$ and $\gamma$ do have order of contact $0$ then $\gamma(t_0)=\text{C}(t_0)$, which means that $\gamma(t_0)=(x(t_0),y(t_0))=(x_0,y_0)$ satisfies the implicit equation of $C$, $$(x_0 − a)^2 + (y_0 − b)^2 - R^2=0 \tag{1}$$
- if they have order of contact $1$ then $\gamma'(t_0)=\text{C}'(t_0)=(\dot x_C(t_0),\dot y_C(t_0))$. But the implicit function theorem tells us that $2(x_C-a)\dot x_C+2(y_C-b)\dot y_C=0$, $\forall t \in \mathbb{R}$ and thus $2(x_0-a)\dot x_0+2(y_0-b)\dot y_0=0\implies$ $$(x_0-a)=-\frac{\dot y_0}{\dot x_0}(y_0-b) \tag{2}$$
- same reasoning for when they do have order of contact 2 yields: $2\ddot x_0(x_0-a)+2\dot x_0^2+2\ddot y_0(y_0-b)+2\dot y_0^2=0$ substituting $(x_0-a)$ using the previous equation (2): $$(y_0-b)=-\frac{\dot x_0^2+\dot y_0^2}{\ddot y_0-\frac{\ddot x_0}{\dot x_0} \dot y_0} \tag{3}$$
...ect, the same reasoning goes for any order of contact if we assume $\gamma$ to be smooth.
Now substituting both $(x_0-a)$ and $(y_0-b)$ in equation $(1)$ so that we get an equation with only R and derivatives of $x_0(t)$ and $y_0(t)$, first substituting $(x_0-a)$:
ps: I'll be omitting the variable here, so $\dot x$, $\dot y$, $\ddot x$, $\ddot y$ are all evaluated at $t_0$ in these calculations. Just abbreviation here.
$$\left(\frac{\dot y}{\dot x}\right)^2(y-b)^2+(y-b)^2=R^2$$ i.e $$(y-b)^2\left(\frac{\dot y^2}{\dot x^2}+1\right)=R^2\implies (y-b)^2 \times \frac{\dot y^2+\dot x^2}{\dot x^2}=R^2$$
Now substituting $(y-b)$: $$\frac{\left(\dot x^2+\dot y^2\right)^2}{\left(\ddot y-\frac{\ddot x}{\dot x} \dot y\right)^2} \times \frac{\left(\dot x^2+\dot y^2\right)}{\dot x^2}=R^2$$
i.e
$$R^2=\frac{\left(\dot x^2+\dot y^2\right)^3}{\left[\dot x\left(\ddot y-\frac{\ddot x}{\dot x} \dot y\right)\right]^2}=\frac{\left(\dot x^2+\dot y^2\right)^3}{\left(\dot x\ddot y-\ddot x \dot y\right)^2}$$
Again $\dot x$, $\dot y$, $\ddot x$, $\ddot y$ are all evaluated at $t_0$, they were only abbreviated.
Thus:
$$R_0=\frac{\left(\dot x_0^2+\dot y_0^2\right)^\left(\frac{3}{2}\right)}{\left|\dot x_0\ddot y_0-\ddot x_0 \dot y_0\right|}$$
Because R is dependent on $t_0$ too. From that we conclude $a$ and $b$:
$$\begin{pmatrix}
a \\
b
\end{pmatrix}=\begin{pmatrix}
x(t_0) - R(t_0)\cos t_0 \\
y(t_0) - R(t_0)\sin t_0
\end{pmatrix}$$
Since a, b and R are already determined there is no going further to order of contact $3$ I guess, the solution ends here, I think.