Orbit-4
Given x3, y3, z3, dx3, dy3, dz3, k, ie 
Find  Vel, a, th,e, p, ee, t, o, w, i
------------------
*******************************************************************************
x3 = x1 * Cos[O] - y1 * Sin[O]   =  -0.27098619163 AU
y3 = y2 * Cos[IE] - z2 * Sin[IE] =  -0.439148484293  AU
z3 = z2 * Cos[IE] + y2 * Sin[IE] =  -0.247105509699 AU

dx3 = dx1 * Cos[O] - dy1 * Sin[O]   =   1.4415181119E-4  AU/day
dy3 = dy2 * Cos[IE] - dz2 * Sin[IE] =  -6.38837169935E-3 AU/day
dz3 = dz2 * Cos[IE] + dy2 * Sin[IE] =  -2.09096438056E-3 AU/day

r   = SQRT[ x3^2 + y3^2 + z3^2 ]         = 0.572141626031 AU
Vel = SQRT[ dx3^2 + dy3^2 + dz3^2 ]      = 2.18642465407E-2 AU/day 

The sum of Kinetic Energy + Potential Energy is a constant.

The mass of the sun mm is one unit.

(1/2) * Vel^2 + (-k^2 * mm / r ) = (-k^2 * mm /( 2*a) )
a = 1 /[ 2/r-(Vel/k)^2 ] = .531875 AU


Set up vector s1 = [ x3 * I + y3 * J + z3* K ]
Set up vector s2 = [ dx3 I + dy3 * J + dz3* K ]
Cos[s4] = (s1 dot s2)/ ( Abs (s1) ]*Abs(s2) ) 
Abs (s1) = SQRT[ x3^2 + y3^2 + z3^2 ] = r
Abs (s2) = SQRT[ dx3^2 + dy3^2 + dz3^2 ] = Vel
s1 dot s2 = x3*dx3 +y3*dy3 + z3*dz3
s4= 79.3009561505 degrees

Vel * Sin[s4] =  2.14841561789E-2
Vel * Cos[s4] = 4.04910212669E-3 = drt

(drt* r^2/k)^2 = 0.5106 p  = semilatus rectum = vertical line from focus
to orbit

s5= p/r-1 = -0.1075563622747 = e*Cos[th] 

s6= drt*SQRT[p]/k] =0.168612179444 = e*Sin[th]

th= AtcTan[s5/56]   = 122.535231561 degrees true anomaly

p/(1+e)= q = perihelion AU = .4255

s6 = r*Cos[th]/a = -.37853467478 = Cos[ee] in radians

SQRT[s5^2 + s6^2] = 0.2 eccentricity = e

s7 = (r*Sin[th])/(a*SQRT[1-e^2]) = .925587110968  = Sin[ee] in radians

ee= ArcTan[s6/s7] = 1.95900897924 radians = eccentric anomaly

(ee-e*Sin[ee]*a^1.5/k= t = 40 = days past perihelion

Now find o,w, and i
*******************************************************************************
We will compute Px Qx Py Qy Pz Qz  which are combinations of the
o, w, and i elements

Sin[i] * Sin[w] = Pz * Cos[ie]-Py * Sin[ie]
Sin[i] * Cos[w] = Qz * Cos[ie]-Qy * Sin[ie]
Sin[o] = ( Py * Cos[w]-Qy*Sin[w] )/Cos[ie]
Cos[o] = Px * Cos[w]- Qx * Sin[w]
Cos[i] = -( Px * Sin[w] + Qx* Cos[w] )/Sin[o]

These equations lump the 4 rotations into w, i, o, and ie are combined into
the P and Q equations.

x3  = Px * r * Cos[th] + Qx * r* Sin[th]
xd3 = k/SQRT[p] *  [ -Px * Sin[th] + Qx * (Cos[th] + e) ]

y3  = Py * r * Cos[th] + Qy * r* Sin[th]
yd3 = k/SQRT[p] *  [ -Py * Sin[th] + Qy * (Cos[th] + e) ]

z3  = Pz * r * Cos[th] + Qz * r* Sin[th]
zd3 = k/SQRT[p] *  [ -Pz * Sin[th] + Qz * (Cos[th] + e) ]

s8 = r* Cos[th]
s9 = r * Sin[th]
s10 = k/SQRT[p] * Sin[th]
s11 = k/SQRT[p] * ( Cos[th]+e )
s12 = s8 * s11 - s9 * s10

Px = (x3 * s11 - xd3 * s9)/s12 = 0.173630530853
Py = (y3 * s11 - yd3 * s9)/s12 = -3.98586415568 E-2
Pz = (z3 * s11 - yd3 * s9)/s12 = 0.984003926544

Qx = ( dx3 * s8 - x3 * s10 )/s12 = -0.451038790789
Qy = ( dy3 * s8 - y3 * s10 )/s12 = 0.88500767159
Qz = ( dz3 * s8 - z3 * s10 )/s12 = 0.115435828233

s13 = Pz * Cos[ie] - Py * Sin[ie] = Sin[i] * Sin[w] = 0.918650051351
s14 = Qz * Cos[ie] - Qy * Sin[ie] = Sin[i] * Cos[w] = -.246151539338

====
Find IA, OA, WA from PQ,XYZ and IE
s1 = PZ * COS[ IE ] - PY * SIN[ IE ]  = 0.918650051364
S2 = QZ * COS[ IE ] - QY * SIN[ IE ]  = -0.246151539338
TAN[ WA ] = s1/s2 = (0.9510565516296 = S6 )(WA =  105 )

s3 = PX * COS[ WA ] - QX *  SIN[ WA ] (S3 = 0.390731128484 )
s4 = ( PY * COS[ WA ] - QY *  SIN[ WA ] )/COS[ IE ] = ( -0.920504853456 )
s4/s3= Tan[ OA ]  = 1 @ -67= QA
s5 = PX * -( SIN[ WA ] + QX * COS[ WA ] )/SIN [ OA ] = 0.30916994376
TAN[ IA ] = s5 / s6 -> 1 @ 72
WA = 105   OA = -67   IA = 72
Program E12 does this
======
This matches program E31
Here is another approach to this problem.
VEA = Velocity(AU/DAY) = SQRT[ DX3^2 + DY3^2 + DZ3^2 ]
R   =                  = SQRT[ X3^2 + Y3^2 + Z3^2 ]
A   = 2/R-(VEA/K)^2 = Semimajor Axis
PH = angle between radius and velocity vectors
COS[ PH ] = ([ X3 Y3 Z3 ] DOT [ DX3  DY3  DZ3 ])/( R * VEA )
d( TH )/ dT = DTHA = VEA* SIN[ PH ]/R^2
d( R )/ dT = DRTA = VEA * COS[ PH ]
PA = Semilatus Rectum = ( DTHA * R^2/K )^2
EA * COS[ THA ] = (PA/R) - 1    EA * SIN[ THA ] = DRTA * SQRT[ PA ]/K
This gives EA = Eccentricity and True Anomaly THA
QA = Perihelion = PA/(1+EA)  

COS[ EEA ] = ( R * COS[ THA ]/A) +EA
SIN[ EEA ] = ( R * SIN[ THA ]/(A * SQRT[1 - EA^2 ] )
This gives Eccentric Anomaly EEA

EEA - EA * SIN[ EEA ] = K * TA/ A^1.5    Kepler's Equation
TA is days since Perihelion passage
Find I, W, and O

Solve a simultaneous equation for P and Q, X Y Z
s1  = R * SIN[ THA ]
s2  = R * COS[ THA ]

X3A = PX * ( R * COS[ THA ] ) + QX * ( R * SIN[ THA ] )

Take first derivative

DX3 = PX * s3 - QX * (R * SIN[ THA ] * DTHA )

s3 = DRTA * COS[ THA ] - R * SIN[ THA ] *DTHA
s4 = DRTA * SIN[ THA ] + R * COS[ THA ] * DTHA


S1 = R * COS[ THA ]
S2 = R * SIN[ THA ]

[  S1 S2 ]   [ PX ]    [ X3A ]
[  S3 S4 ] * [ QX ]  = [ DX3 ]  solve simultaneous equation for PX,QX

Solve for PY,QY  and PZ,QZ
Program E12 derives the values of I, W, and O