Let's make sure we understand splitting, before we talk about joining. To subdivide a Bézier curve into two, you use the deCasteljau algorithm, as illustrated in the figure below.
![enter image description here](https://i.stack.imgur.com/I9wKC.png)
Suppose we are given a curve defined by four control points $A$, $P$, $Q$, $G$, and a splitting parameter value $u \in [0,1]$ (which you called $t_{\text{cut}}$). To split the curve, we repeatedly divide edges in the ratio $u:v$, where $v=1-u$. So, we divide the edge $AP$ in the ratio $u:v$ to get point $B$, divide $QG$ to get $F$, and so on. The division is done using convex combinations of points, so, for example, $B = uP + vA$, and $F = uG+vQ$. This process generates the points $B$, $C$, $D$, $E$, $F$, as shown. The two curves resulting from the splitting are shown in pink and blue in the picture. The pink curve has control points $(A,B,C,D)$, and the blue one has control points $(D,E,F,G)$.
Now let's talk about how we might "undo" this splitting process. So, we are given the points $A$, $B$, $C$, $D$, $E$, $F$, $G$, and we want to get the points $P$ and $Q$. We simply have to run the deCasteljau "dividing" steps in reverse.
First we find the splitting parameter $u$ (your $t_{\text{cut}}$). We calculate $k = DE/CD$. Then $E-D = k(D-C)$, and a little arithmetic shows that
$$
D = \frac{1}{1+k}E + \frac{k}{1+k}C
$$
But if $D$ was produced by the deCasteljau algorithm, we also know that $D = uE + vC$, so we have $u = 1/(1+k)$ and $v=k/(1+k)$.
Next, we find $P$ so that $B$ divides $AP$ in the ratio $u:v$. This means that $B = uP + vA$, which implies that
$$
P = \frac{1}{u}B - \frac{v}{u}A = (1+k)B - kA
$$
Similarly, since $F = uG+vQ$, we have
$$
Q = \frac{1}{v}F - \frac{u}{v}G = \frac{1+k}{k}F - \frac{1}{k}G
$$
The $u$ and $v$ are just to relate all this to the deCasteljau algorithm, for purposes of understanding. In code, you can forget about $u$ and $v$, and just use $k$. The (pseudo) code is:
k = Length(E-D)/Length(D-C)
P = (1+k)*B - k*A
Q = ((1+k)*F - G)/k
Of course, you can execute this code with any two Bézier curves as input, regardless of whether they were produced by splitting. I think the output curve will then be a reasonable approximation of the two input ones, which gives another answer to a question you asked elsewhere.
If we use the input data
A = ( 94 , 242)
B = ( 121 , 183.5)
C = ( 272.5 , 183 )
D = ( 417.25 , 210.875)
E = ( 562 , 238.75 )
F = ( 700 , 295 )
G = ( 700 , 350 )
Then the output is
k = 1
P = 2*B - A = (242,367) - (94,242) = (148,125)
Q = 2*F - G = (1400, 590) - (700,350) = (700,240)
Another example, as given in the comments. The input data is:
A = ( 94 , 242 )
B = ( 107.5 , 212.75 )
C = ( 152.125 , 198 )
D = ( 211.46875 , 194.046875 )
E = ( 389.5 , 182.1875 )
F = ( 700 , 267.5 )
G = ( 700 , 350 )
This gives the following output, which is easy to check by hand calculations:
k = 3
P = (4*B - A) = (148 , 125)
Q = (4*F - G)/3 = (700 , 240)
So, actually, this is the same curve as in the previous example, but it is split at $u = 1/(1+k) = 0.25$ rather than at $u = 0.5$.
I don't know what is wrong with the code given in the question, but I suspect the calculation of "c2" in the last line is not working correctly. Maybe this is because division of a vector by a scalar is not properly supported. I suggest rewriting it as
r = 1.0/k
c2 = (1+r)*F - r*G
This is a common problem, of course, because people are always trying to convert Postscript fonts (which use cubic curves) into TrueType ones (which use quadratic curves).
To answer your question, I don't think there is any accepted "gold standard" approach.
Many of the algorithms you'll find on the net are heuristic. The developer thinks his results "look good", and that's the end of the analysis. If you want to be more rigorous, then you have to define "good" or "best" in order to find the best algorithm. Lets call the original cubic $C$ and the approximating piecewise quadratic $Q$. Then some possible criteria of goodness are:
- The maximum distance between $C$ and $Q$ is small
- The maximum value of $\|C(t) - Q(t)\|$ is small. This is not the same thing as #1.
- The maximum angle between $C'(t)$ and $Q'(t)$ is small
- Tangent discontinuities in $Q$ are small (perhaps even zero)
- $Q$ and $C$ have the same turning points (important for font hinting)
and so on.
No matter what criteria you decide to choose, the "test-and-reduce" meta-algorithm you mentioned is reasonable: take two points on the cubic, and use these points and their tangent vectors to construct a quadratic. If the quadratic is acceptable (according to whatever criteria you decided), you're done; if not, choose two points that are closer together. Rinse and repeat.
Regardless of what @Timo said, there is no way to exactly represent a cubic curve as a string of quadratic ones unless the cubic curve was, (by some miraculous stroke of luck) a disguised quadratic in the first place. The string of quadratics you construct will be merely be an approximation of the cubic; all you can do is ensure that it's a "good" approximation.
Best Answer
The calculations are roughly the same regardless of what line is used.
Suppose the line has equation $lx+my+nz = d$. In vector form, this is $$ \mathbf{A}\cdot \mathbf{X} = d, $$ where $\mathbf{A} = (l,m,n)$. Also, suppose the Bézier curve has equation $$ \mathbf{X}(t) = (1-t)^3\mathbf{P}_0 + 3t(1-t)^2\mathbf{P}_1 + 3t^2(1-t)\mathbf{P}_2 + t^3\mathbf{P}_3 $$ Substituting this into the previous equation, we get $$ (1-t)^3(\mathbf{A} \cdot \mathbf{P}_0) + 3t(1-t)^2(\mathbf{A} \cdot \mathbf{P}_1) + 3t^2(1-t)(\mathbf{A} \cdot \mathbf{P}_2) + t^3(\mathbf{A} \cdot \mathbf{P}_3) - d = 0 $$ This is a cubic equation that you need to solve for $t$. The $t$ values you obtain (at most 3 of them) are the parameter values on the Bézier curve at the intersection points. If any of these $t$ values is outside the interval $[0,1]$, you'll probably want to ignore it.
In your particular case, you know that the line passes through either the start-point or the end-point of the Bézier curve, so either $t=0$ or $t=1$ is a solution of the cubic equation, which makes things easier. Suppose $t=0$ is a solution. Then the cubic above can be written in the form $$ t(at^2 + bt + c) =0 $$ Now you only have to find $t$ values in $[0,1]$ that satisfy $at^2 + bt +c =0$, which is easy.