SimpleKSPOrbitalOperations0.002

-- SimpleKSPOrbitalOperations is a purely calculation based lua library, with some constants hardcoded for Kerbal Space Program

-- v0.002

-- v0.002 changes: Added measureSkooIspFromTmbiTmbtwnNsDriver(TimeBeforeInitial, TimeBetweenEach, NumberSlices), testSkooIspFromTmbiTmbtwnNsMpcsDriver(TimeBeforeInitial, TimeBetweenEach, NumberSlices, MaxPointingChangeStable), cotestSkooIspFromTmbiTmbtwnNs(TimeBeforeInitial, TimeBetweenEach, NumberSlices), cotestSkooIspFromTmbiTmbtwnNsMpcs(TimeBeforeInitial, TimeBetweenEach, NumberSlices, MaxPointingChangeStable)

-- to use, put

-- require "SimpleKSPOrbitalOperations"

-- in your lua script, or run it with

-- dofile("SimpleKSPOrbitalOperations")

-- Original Work Copyright 2012, Nadrek, insofar as any elements that aren't obvious, trivial, or mathematics are concerned (none of which are copyrightable)

-- licensed under your choice of GPLv2, GPLv3, or Creative Commons Attribution-ShareAlike 3.0 Unported (CC BY-SA 3.0)  http://www.gnu.org/licenses/gpl-2.0.html  or  http://www.gnu.org/licenses/gpl.html  or  http://creativecommons.org/licenses/by-sa/3.0/

-- KSP specific functions will have KSP in the name; other functions will not, and can be used for real life calculations as well.



-- This script requires SimpleKSPOrbitalMechanics to do the math, which in turn requires SimpleOrbitalMechanics to do the heavy mathematical lifting.

-- This script contains KSP specific functions and calls useful to and reliant on Kerbal Space Program, MechJeb, and MechJeb's Autom8 module.











--[[ NOTE: Undocumented commands:

mechjeb.stage() -- fires off a stage, like hitting space bar

mechjeb.controlRelease() -- deactivates the .attitudeTo() and perhaps other commands

mechjeb.thrustActivate(Percentage)

mechjeb.thrustDeactivate()

]]





require "SimpleKSPOrbitalMechanics"





local function launchSkooFromApPeInc(Ap, Pe, Inclination)

  -- Ap is the Apoapsis in meters

  -- Pe is the Periapsis in meters

  -- Inclination is the desired orbital inclination in degrees

  -- FUTURE include wait(time) or warp(time), but right now do print(vessel.time) wait(2) print(vessel.time) end    freezes the game for 2 seconds and then displays two identical times.



  -- REQUIRES globalSomMainBodyReset to have been called properly

  --   the simple way is to run globalSkomKSPMainBodyResetBasedOnGuessing() first!



  -- First, let's reset our main body; perhaps we've switch orbits recently!

  globalSkomKSPMainBodyResetBasedOnGuessing()



  -- Check to see if the Ap and Pe is obviously bad

  if Ap == (-1) or Pe == (-1) then

    print("ERROR: Ap or Pe is -1, which indicates an inability to get a prior calculation to come out right when orbiting " .. globalSomMainBodyName)

    return

  end





  print("Launching to a Periapsis of " .. Pe .. " meters with Inclination " .. Inclination .. " degrees at " .. vessel.time)

  mechjeb.launchTo(Pe, Inclination)

  wait(mechjeb.free)

  mechjeb.autoStageActivate()

  local StartTime = vessel.time -- in seconds, with a massive offset.  Maybe for use later.

  print("Followon launch circularization at " .. StartTime)

  mechjeb.circularize()

  wait(mechjeb.free)

  print("Changing Apoapsis to " .. Ap .. " meters starting at " .. vessel.time)

  mechjeb.changeAp(Ap)

  wait(mechjeb.free)

  print("Tweaking orbit Pe starting at " .. vessel.time)

  mechjeb.changePe(Pe)

  wait(mechjeb.free)

  print("Tweaking orbit Ap starting at " .. vessel.time)

  mechjeb.changeAp(Ap)

  wait(mechjeb.free)

  mechjeb.autoStageDeactivate()

  print("Launch complete at " .. vessel.time)

  printSkomKSPCurrentOrbitalElements()

end





local function launchSkooMolniyaFromIncMoe(Inclination, MarginOfError)

  -- Inclination is the desired orbital inclination in degrees, and should be 63.4, -63.4, 116.6, or -116.6 degrees

  --   to correctly counteract the gravitational effect of an oblate body, like Earth's fat equatorial bulge

  -- MarginOfError 0.001 means allow 0.1% for a "safety zone" above the minimum stable altitude, or below the Hill sphere of influence.

  --   Do this if you want to just barely scrape by!

  --   2 means allow a 200% "safety zone" above the minimum stable altitude, or below the Hill sphere of influence.

  --   Do this if you want to have a much less eccentric orbit - perhaps you want to keep certain orbits clear for other uses!

  --   If in doubt, try setting MarginOfError to 0.1

  -- REQUIRES globalSomMainBodyReset to have been called properly

  --   the simple way is to run globalSkomKSPMainBodyResetBasedOnGuessing() first!





    -- First, let's reset our main body; perhaps we've switch orbits recently!

  globalSkomKSPMainBodyResetBasedOnGuessing()



  -- Now, did we get a good Inclination?  If not, give a polite warning.

  if Inclination ~= (63.4) and Inclination ~= (-63.4) and Inclination ~= (116.6) and Inclination ~= (-116.6) then

    print "WARNING: Classic Molniya orbits must have an inclination of 63.4, -63.4, 116.6, or -116.6 degrees"

    print "   to correctly counteract the gravitational effect of an oblate body, like Earth's fat equatorial bulge"

    print "   Continuing to try anyway.  You have been politely warned!"

  end



    -- Now, calculate a Molniya orbit.

  local localAp, localPe = calcSkomMolniyaApPeFromMainBodyCurrentlyBeingOrbitedFromMoe(MarginOfError)



    -- Did we get one?  If we got -1 back, we didn't.

  if localAp == (-1) or localPe == (-1) then

    print("ERROR: Cannot calculate a stable Molniya orbit with this MarginOfError around " .. globalSomMainBodyName)

    return

  end





  launchSkooFromApPeInc(localAp, localPe, Inclination)

end



function colaunchSkooFromApPeInc(Ap, Pe, Inclination)

  local coToLaunchWith = coroutine.create(launchSkooFromApPeInc)

  coroutine.resume(coToLaunchWith, Ap, Pe, Inclination)

end





function colaunchSkooMolniyaFromIncMoe(Inclination, MarginOfError)

  local coToLaunchWith = coroutine.create(launchSkooMolniyaFromIncMoe)

  coroutine.resume(coToLaunchWith, Inclination, MarginOfError)

end





local function measureSkooIspFromTmbiTmbtwnNsDriver(TimeBeforeInitial, TimeBetweenEach, NumberSlices)

  -- added in v0.002

  -- NOTE: This is designed to be run asynchronously to vessel control.

  -- NOTE: This will only work anywhere close to properly with a small TimeBetweenEach,

  --   in a circular orbit, with the vessel held in prograde or retrograde directions.

  -- TimeBeforeInitial is the "spool up" delay in seconds we're going to use

  --   this is how much time we think it takes our engines to come to stable thrust

  --   i.e. 0.5 means set thrust to the desired value, then wait half a second, then start the measurement for this "slice"

  --   For the bell rocket engines in 0.16, 0.25 is a reasonable minimum for TimeBeforeInitial

  -- TimeBetweenEach is the time in seconds between initial and final measurements for each "slice"

  --   For the bell rocket engines in 0.16, 0.25 is a reasonable minimum for TimeBetweenEach

  -- NumberSlices is how many slices we're going to divide our thrust into, i.e.

  --   1 means 100 (100% max thrust) / 1 = one test at 100%.

  --   2 means 100/2 = 50% for the first test, and 100% for the second.

  --   4 means 100/4 = 25% for the first test, 50% for the second, 75% for the third, and 100% for the fourth.

  -- Nothing is returned per se, but the Specific Impulse, kg mass lost, and Delta-V per mass is printed at each slice for human consumption.

  -- While we could use the Skom reporting functions, because we need initial and final mass

  --   instead of just differences, this is almost easier and is certainly simpler to understand.

  if TimeBeforeInitial < 0.25 then

    print "WARNING: TimeBeforeInitial below 0.25 seconds may or may not return accurate results.  Test your particular implementation."

  end



  if TimeBetweenEach < 0.25 then

    print "WARNING: TimeBetweenEach below 0.25 seconds may or may not return accurate results.  Test your particular implementation."

  end



  -- Deactivate autostaging, which could result in inaccurate fuel estimates, since we

  --   assume that all mass lost is fuel.

  mechjeb.autoStageDeactivate()



  for i = 1, NumberSlices do

    local localThrustThisSlice = 100 * (i / NumberSlices)

    mechjeb.thrustActivate(localThrustThisSlice)



    wait(TimeBeforeInitial)



    local localInitialMs = vessel.mass

    local localInitialOv = vessel.speedOrbital



    wait(TimeBetweenEach)



    local localFinalMs = vessel.mass

    local localFinalOv = vessel.speedOrbital



    local localIsp = calcSomIspFromTsmTemDv(localInitialMs, localFinalMs, (localFinalOv - localInitialOv))



    print("At thrust: " .. localThrustThisSlice .. "% Isp is: " .. localIsp)

    print("   Used " .. (localInitialMs - localFinalMs)*1000 .. " kg fuel, for " .. (localFinalOv - localInitialOv) .. " DeltaV, or " .. (localFinalOv - localInitialOv)/((localInitialMs - localFinalMs)*1000) .. " DeltaV per kg fuel and " .. (localFinalOv - localInitialOv)/TimeBetweenEach .. " DeltaV per second")

  end -- for i = 1, NumberSlices do



  -- End test thrust.

  mechjeb.thrustDeactivate()



end





local function testSkooIspFromTmbiTmbtwnNsMpcsDriver(TimeBeforeInitial, TimeBetweenEach, NumberSlices, MaxPointingChangeStable)

  -- added in v0.002

  -- NOTE: This is designed to be run asynchronously to the main window.

  -- see measureSkooIspFromTmbiTmbtwnNsDriver for main usage

  -- MaxPointingChangeWeCallStable is the maximum degrees of change in heading + pitch

  --   we'll call "stable" for purposes of detecting when mechjeb.attitudeTo() is done

  --   i.e. 0.2 means 2/10ths of a degree change during TimeBetweenEach seconds

  -- This is a wrapper that circularizes our orbit, verifies

  --   a stable orbit, orients prograde, and then runs the Isp test



  -- REQUIRES globalSomMainBodyReset to have been called properly

  --   the simple way is to run globalSkomKSPMainBodyResetBasedOnGuessing() first!



  -- First, let's reset our main body; perhaps we've switch orbits recently!

  globalSkomKSPMainBodyResetBasedOnGuessing()



  -- Now, circularize, just in case.  Wait in TimeBetweenEach increments.

  print "Circularizing"

  mechjeb.circularize()

  while mechjeb.busy() == true do

    wait(TimeBetweenEach)

  end -- while mechjeb.busy() == true do



  print "Are We Stable"

  -- if we're in a stable orbit, continue

  if isSkomCurrentOrbitStable() ~= 0 then

    print "WARNING: Orbit not stable!  Test may be inaccurate and/or result in destruction of your vessel."

  end -- if isSkomCurrentOrbitStable() ~= 0 then



  -- WARNING:  The edge case of switching orbiting main bodies during the test run

  --   is NOT HANDLED, and will generate bogus results.

  --   It may also have the side effect of killing your crew.



  print "Prograde"

  mechjeb.attitudeTo("forward","ORBIT")

  -- Ok, we've asked mechjeb to point us prograde... but are we done pointing yet?

  --   mechjeb.busy() will stay true until mechjeb.controlRelease() is called.

  -- We're going to do a cheap hack here - if our heading or pitch is changing "too fast"

  --   keep waiting until it changes "slower".



  -- v0.002 NOTE: These _should_ be locals; however, in that case, they don't appear to update in the while loop.

  --   Therefore, we're going to use distinctive, function-specific names to avoid side effect with other functions.

  globalHdPrev_testSkooIspFromTmbiTmbtwnNsMpcsDriver = vessel.vesselHeading

  globalPiPrev_testSkooIspFromTmbiTmbtwnNsMpcsDriver = vessel.vesselPitch

  print "Before initial wait for prograde test"

  wait(TimeBetweenEach)

  print "After initial wait for prograde test, waiting on prograde completion"

  print "  If the thrusters don't fire on their own, click PRO GRAD on Smart A.S.S. and try a larger MaxPointingChangeStable next time"

  --print ("gt " .. (while1 > while2))

  -- If the sum of the changes in heading and pitch is over our threshold, let's keep waiting.

  while (((math.abs(vessel.vesselHeading - globalHdPrev_testSkooIspFromTmbiTmbtwnNsMpcsDriver) + math.abs(vessel.vesselPitch - globalPiPrev_testSkooIspFromTmbiTmbtwnNsMpcsDriver))) >= MaxPointingChangeStable) do

    globalHdPrev_testSkooIspFromTmbiTmbtwnNsMpcsDriver = vessel.vesselHeading

    globalPiPrev_testSkooIspFromTmbiTmbtwnNsMpcsDriver = vessel.vesselPitch

    wait(TimeBetweenEach)

  end

  -- one last wait, just to make sure we've settled out as best we can.

  wait(TimeBetweenEach)

  print "Prograde settled; initiating test measurements."

  -- Now, mechjeb should be holding us Prograde in a circular orbit, and we can start our timing measurements.

  --   These just print out to the console currently, as getting data out of coroutines is a pain.

  measureSkooIspFromTmbiTmbtwnNsDriver(TimeBeforeInitial, TimeBetweenEach, NumberSlices)



  -- Release the .attitudeTo hold we've put on the craft.

  mechjeb.controlRelease()

end







function cotestSkooIspFromTmbiTmbtwnNs(TimeBeforeInitial, TimeBetweenEach, NumberSlices)

  -- added in v0.002

  -- IF YOU USE THIS, START IN A CIRCULAR ORBIT, WITH SMART A.S.S. set to Prograde!!!

  -- Part of initial development.

  local coMeasure = coroutine.create(measureSkooIspFromTmbiTmbtwnNsDriver)

  coroutine.resume(coMeasure, TimeBeforeInitial, TimeBetweenEach, NumberSlices)

end



function cotestSkooIspFromTmbiTmbtwnNsMpcs(TimeBeforeInitial, TimeBetweenEach, NumberSlices, MaxPointingChangeStable)

  -- added in v0.002

  --[[ ex.

    cotestSkooIspFromTmbiTmbtwnNsMpcs(0.25, 0.25, 2, 0.02)

  ]]



  local coTest = coroutine.create(testSkooIspFromTmbiTmbtwnNsMpcsDriver)

  coroutine.resume(coTest, TimeBeforeInitial, TimeBetweenEach, NumberSlices, MaxPointingChangeStable)

end







function helpSkooSimpleKSPOrbitalOperations1()

  print "Usage: colaunchSkooMolniyaFromIncMoe(Inclination, MarginOfError)"

  print "   If you're currently landed in launch position, this will either launch you into a Molniya orbit, or"

  print "   give you a message stating that no stable Molniya orbit exists for this body and Margin of Error."

  print "   Inclination should be 63.4, -63.4, 116.6, or -116.6 degrees"

  print "   Inclination is the desired orbital inclination in degrees, and should be 63.4, -63.4, 116.6, or -116.6 degrees"

  print "     to correctly counteract the gravitational effect of an oblate body, like Earth's fat equatorial bulge"

  print "   MarginOfError 0.001 means allow 0.1% for a \"safety zone\" above the minimum stable altitude, or below the Hill sphere of influence."

  print "     Do this if you want to just barely scrape by!"

  print "     2 means allow a 200% \"safety zone\" above the minimum stable altitude, or below the Hill sphere of influence."

  print "     Do this if you want to have a much less eccentric orbit - perhaps you want to keep certain orbits clear for other uses!"

  print "     If in doubt, try setting MarginOfError to 0.1"

  print "   Ex: "

  print "     colaunchSkooMolniyaFromIncMoe(-63.4, 0.01)"

  print "   REQUIRES globalSomMainBodyReset to have been called properly"

  print "     the simple way is to run globalSkomKSPMainBodyResetBasedOnGuessing() first!"

  print "Usage: colaunchSkooFromApPeInc(Ap, Pe, Inclination)"

  print "   Ap is the Apoapsis in meters"

  print "   Pe is the Periapsis in meters"

  print "   Inclination is the desired orbital inclination in degrees"

  print "   Ex: 125km circular orbit, 0 degree inclination"

  print "     colaunchSkooFromApPeInc(125000, 125000, 0)"

  print "   REQUIRES globalSomMainBodyReset to have been called properly"

  print "     the simple way is to run globalSkomKSPMainBodyResetBasedOnGuessing() first!"

end



function helpSkooSimpleKSPOrbitalOperations2()

  print "Usage: cotestSkooIspFromTmbiTmbtwnNsMpcs(TimeBeforeInitial, TimeBetweenEach, NumberSlices, MaxPointingChangeStable)"

  print "   TimeBeforeInitial is the \"spool up\" delay in seconds we're going to use"

  print "     this is how much time we think it takes our engines to come to stable thrust"

  print "     i.e. 0.5 means set thrust to the desired value, then wait half a second, then start the measurement for this \"slice\""

  print "     For the bell rocket engines in 0.16, 0.25 is a reasonable minimum for TimeBeforeInitial"

  print "   TimeBetweenEach is the time in seconds between initial and final measurements for each \"slice\""

  print "     For the bell rocket engines in 0.16, 0.25 is a reasonable minimum for TimeBetweenEach"

  print "   NumberSlices is how many slices we're going to divide our thrust into, i.e."

  print "     1 means 100 (100% max thrust) / 1 = one test at 100%."

  print "     2 means 100/2 = 50% for the first test, and 100% for the second."

  print "     4 means 100/4 = 25% for the first test, 50% for the second, 75% for the third, and 100% for the fourth."

  print "   MaxPointingChangeWeCallStable is the maximum degrees of change in heading + pitch"

  print "     we'll call \"stable\" for purposes of detecting when mechjeb.attitudeTo() is done"

  print "     i.e. 0.2 means 2/10ths of a degree change during TimeBetweenEach seconds"

  print "   Ex: one half second spool up and measurements, at 50% and 100% thrust, with 0.04 degrees per second change = \"stable\""

  print "     cotestSkooIspFromTmbiTmbtwnNsMpcs(0.5, 0.5, 2, 0.02)"

end





print "Usage: helpSkooSimpleKSPOrbitalOperations1() through helpSkooSimpleKSPOrbitalOperations2()"

print "   gives detailed help on all Skoo functions."

print "   Skoo functions relate to Kerbal Space Program (KSP) operational functionality."