For given functions $t\in [0,T] \mapsto (x_1(t),y_1(t))\in {\Bbb R}^2\setminus\{(0,0)\}$,
I don't think that there is a general local solution to get a continuous function $\Theta_1(t)$. Take e.g. $x_1(t)=\cos(t/2), y_1(t)=\sin(t/2)$, $t\in[0,20\pi]$.
Ideally, you might want that $\Theta_1(t)=2 \arctan(\tan(t/2))=t$ should be the solution. But your software will plot a saw tooth graph with 9 discontinuities. You would have to add more and more multiples of $2\pi$ as $t$ increases.
You need to have global information on the graph.
If, however, $x_1(t)$ and $y_1(t)$ comes from solving an ODE, then you may recover a continuous solution $\theta_1(t)$, simply by adding a differential equation for it to your ODE.
EDIT: Inspired by the post of @ingix, here is a suggestion using vector calculus (in a way it mimics an ODE). I use scilab (I don't know about Mathcad but assume it is close).
If $x,y$ are vectors of length $N$ describing the discretization of a smooth motion winding around the origin. Let $T=atan(x,y)$ be the arctan function taking values in $(-\pi,\pi]$. Let $DT=T(2:\$)-T(1:\$-1)$ be the vector of increments. Most of the values of $DT$ are close to zero but when the arctan jumps the value is close to $\pm 2\pi$. Calculating $R=round(DT/2 \; \pi)$ now gives values in $\{0,\pm 1\}$ and yields precisely the jumps in $T$. Using cumulative sums to recontruct the original path,
$$ ST= T(1) + [0; cumsum(DT) - 2\pi \; cumsum(R)] $$
then gives you the angular path, now with discontinuities removed. (You may then multiply by 2 to get your $\Theta$).
As an example I have taken some smooth curve in the plane that avoids the origin, calculated the atan along the curve ( taking values in $(-\pi,\pi]$ and giving the 'saw-tooth' blue curve). Then reconstructed the continuous angle (the green curve) from the blue curve data, using the above formulae.