How to spline together Bezier curves to form a smoth closed curve? Given $k=m\cdot n$ points: $P_1,P_2,...,P_k$ (all points are two dimensional points), how can I spline together $m$ Bezier curves of $n$ degree to form a smooth closed curve?  
Denote $B_{i,j}(t)$ to be the Bezier curve defined by points $P_i$ through $P_j$ (meaning $P_i$ and $P_j$ are the edge points). I want to spline together the curves:  $B_{1,n+1}(t),B_{n+1,2n+1}(t),...,B_{(m-1)\cdot n+1,1}(t)$
For example, given the following: $P_1,P_2,...,P_8$

I want to find the Bezier curves $B_{1,3}(t),B_{3,5}(t),B_{5,7}(t),B_{7,1}(t)$ that will also create a nice smooth closed curve.
I know how to find Bezier curve of order $n$ given $n+1$ points, but the problem is how to spline them together, in a way that the edges (the points where two Bezier curves "meet") will be smooth.
In the example above, I want the curve to be smooth at points $P_1,P_3,P_5,P_7$.  
When interpolating (say, natural) cubic spline, we solve the above problem (smoothness in the knots) by adding constraints that states the derivative of every cubic spline at point $x_i$ must be equal to the derivative of the next cubic spline in $x_i$ (namely $S_i'(x_i)=S_{i+1}'(x)$, $\forall i=1,2,...,n-1$).
Can I do something similar here?

We can take a look at the problem for each dimension individually, since the cubic splines can be written as
$$B_i : [t_i, t_{i+1}] \to \mathbb R^2; ~~t\mapsto \pmatrix{a_3 t^3 + a_2 t^2 + a_1 t + a_0 \\ b_3 t^3 + b_2 t^2 + b_1 t + b_0}$$
For the $1D$-case we call the points $x_i$. Furthermore we assume that the number of points is even and at least four - this will ensure that we can form a closed curve and have "junctions" at the indices $2k-1$ and "internal points" at indices $2k$.
The individual cubic segments $B_i$ will then follow the constraints
$$\begin{align*}
B_i(t_i) &= x_{2i-1} \\
B_i(\frac12 (t_i + t_{i+1})) &= x_{2i} \\
B_i(t_{i+1}) &= x_{2i+1} \\
B_i'(t_{i+1}) - B_{i+1}'(t_{i+1}) & = 0
\end{align*}$$
Where the last two equations become
$$\begin{align*}
B_n(t_{n+1}) &= x_1 \\
B_n'(t_{n+1}) - B_1'(t_1) & = 0
\end{align*}$$
Where there are $2n$ knots and $n$ segments. These ensure the closedness.
All in all this gives us a system $4n$ linear equations and $4n$ unknowns ($4$ for each cubic). We can chose the points $\{t_i\}_{i=1}^{n+1}$ to our liking, they only affect the "speed" in wich the parameterization traverses the curve.
Patching the result together to a function
$$B : [t_1, t_{n+1}]\to\mathbb R$$
in a piecewise fashion ensures that every point $x_i$ is traversed at $t_i$ and $x_1$ is re-traversed at $t_{n+1}$, so $B$ is periodic. If we do this twice with the same choice of $\{t_i\}$ and the $y$-coordinates instead, we can join these two to a piecewise cubic $C^1$-curve in $\mathbb R^2$, each component given by the appropriate $B$.
