Why don't you rotate the points so that the axis are aligned with y and z axes and use parameterized equations to get the solution?
If you have a general ellipse at some custom CX,CY then an axis aligned parameterization would be:
EX = CX + RX*cos(t)
EY = CY + RY*sin(t)
From here you need to find the values of t that give CX and CY and the values that do this will be t = 0 and t = pi/2.
If the ellipse has been rotated, then you need to rotate it so that it is axis aligned. You can rotate something by the inverse of its angle quite easily with a standard rotation matrix (with -theta as your argument).
Overall though, once you get it so that the rotation axis are aligned to Y and Z, then given your points for EX and EY for two points, you obtain values of RX and RY by picking a value of t that eliminates the RX RY terms to get CX and CY (like mentioned above).