""" I didn't want to have to write this...I'm going to try to do this as quickly as possible from old, bad code Convenience function for solving lambert's problem """ function lamberts(planet1::Body,planet2::Body,leave::DateTime,arrive::DateTime) time_leave = utc2et(Dates.format(leave,"yyyy-mm-ddTHH:MM:SS")) time_arrive = utc2et(Dates.format(arrive,"yyyy-mm-ddTHH:MM:SS")) tof_req = time_arrive - time_leave state1 = [spkssb(planet1.id, time_leave, "ECLIPJ2000"); 0.0] state2 = [spkssb(planet2.id, time_arrive, "ECLIPJ2000"); 0.0] r1 = state1[1:3] ; r1mag = norm(r1) r2 = state2[1:3] ; r2mag = norm(r2) μ = Sun.μ cos_dθ = dot(r1,r2)/(r1mag*r2mag) dθ = atan(r2[2],r2[1]) - atan(r1[2],r1[1]) dθ = dθ > 2π ? dθ-2π : dθ dθ = dθ < 0.0 ? dθ+2π : dθ DM = abs(dθ) > π ? -1 : 1 A = DM * √(r1mag * r2mag * (1 + cos_dθ)) dθ == 0 || A == 0 && error("Can't solve Lambert's Problem") ψ, c2, c3 = 0, 1//2, 1//6 ψ_down = -4π ; ψ_up = 4π^2 y = r1mag + r2mag + (A*(ψ*c3 - 1)) / √(c2) ; χ = √(y/c2) tof = ( χ^3*c3 + A*√(y) ) / √(μ) i = 0 while abs(tof-tof_req) > 1e-2 y = r1mag + r2mag + (A*(ψ*c3 - 1)) / √(c2) while y/c2 <= 0 # println("You finally hit that weird issue... ") ψ += 0.1 if ψ > 1e-6 c2 = (1 - cos(√(ψ))) / ψ ; c3 = (√(ψ) - sin(√(ψ))) / √(ψ^3) elseif ψ < -1e-6 c2 = (1 - cosh(√(-ψ))) / ψ ; c3 = (-√(-ψ) + sinh(√(-ψ))) / √((-ψ)^3) else c2 = 1//2 ; c3 = 1//6 end y = r1mag + r2mag + (A*(ψ*c3 - 1)) / √(c2) end χ = √(y/c2) tof = ( c3*χ^3 + A*√(y) ) / √(μ) tof < tof_req ? ψ_down = ψ : ψ_up = ψ ψ = (ψ_up + ψ_down) / 2 if ψ > 1e-6 c2 = (1 - cos(√(ψ))) / ψ ; c3 = (√(ψ) - sin(√(ψ))) / √(ψ^3) elseif ψ < -1e-6 c2 = (1 - cosh(√(-ψ))) / ψ ; c3 = (-√(-ψ) + sinh(√(-ψ))) / √((-ψ)^3) else c2 = 1//2 ; c3 = 1//6 end i += 1 i > 500 && return [NaN,NaN,NaN],[NaN,NaN,NaN] end f = 1 - y/r1mag ; g_dot = 1 - y/r2mag ; g = A * √(y/μ) v0t = (r2 - f*r1)/g ; vft = (g_dot*r2 - r1)/g return v0t - state1[4:6], vft - state2[4:6], tof_req end