An orbital speeds script

[Toti]

Hello,
this script puts the viewer above the selected body and computes tangential, angular and radial speeds for it, with respect to its "parent" body (e.g. Earth-Sol, ISS-Earth, Dactyl-Ida, etc.)

LUA precision limitations arise at high time speeds, or when viewer is really close to small bodies

Code: Select all

--*********************************************************************
--*********************************************************************
--            A script for computing speeds
--               coded by Toti
--*********************************************************************
--*********************************************************************

--*********************************************************************
--               constants
--*********************************************************************
DELAY         = 0.5
HORIG         = 0
VORIG         = -1

KM_MLY      = 9466411.842
SEC_DAY      = 86400

ANGLE         = 20
TAN_ANGLE      = math.tan(math.rad(ANGLE))   -- this one is used by viewBodyAbove(), but it's done here because of
                        -- performance
HOW_FAR      = 1.1               -- this one is used by viewBodyAbove(), it's a measure of how close
                        -- the viewer is from parent body
DELTA_S      = 1000            -- time increment in seconds (default: 100 sec. = 16.6 min. you can reduce
                        -- this in order to get smaller "differentials")
DELTA_D      = DELTA_S / SEC_DAY      -- time increment converted to days   (default: about 0.011 days)

obs         = celestia:getobserver()


--*********************************************************************
--      flags and status save and restoration functions
--*********************************************************************
function setRenderOn ()
-- sets the necessary renderflags on
   celestia:setrenderflags{constellations=false,eclipseshadows=true,galaxies=true,grid=false,lightdelay=false,markers=true,
         nightmaps=true,orbits=true,planets=true,ringshadows=true,smoothlines=false,stars=true}
end

function setRenderOff ()
-- sets the necessary renderflags off
   celestia:setrenderflags{markers=false,orbits=false}
   -- a 'labels' render flag would be handy:
   celestia:setlabelflags{asteroids=false,comets=false,moons=false,planets=false,spacecraft=false,constellations=false,
               galaxies=false,locations=false,stars=false}
end

function setOrbits (sType)
-- sets the orbit on for this particular body's type. (There is probably a direct way to do this):
   local t = {}
   t.Asteroid   = false
   t.Comet   = false
   t.Moon   = false
   t.Planet   = false
   t.Spacecraft= false

   if sType == "asteroid" then
      t.Asteroid   = true
   elseif sType == "comet" then
      t.Comet   = true
   elseif sType == "moon" then
      t.Moon   = true
   elseif sType == "planet" then
      t.Planet   = true
   elseif sType == "spacecraft" then
      t.Spacecraft= true
   end
   celestia:setorbitflags(t)
end

function setLabels (sType)
-- sets the label on for this particular body's type. (Again, there is probably a direct way to do this):
   local t = {}
   t.asteroids      = false
   t.comets      = false
   t.moons      = false
   t.planets      = false
   t.spacecraft   = false
   t.constellations   = false
   t.galaxies      = false
   t.locations      = false
   t.stars      = false

   if sType == "asteroid" then
      t.asteroids      = true
   elseif sType == "comet" then
      t.comets      = true
   elseif sType == "moon" then
      t.moons      = true
   elseif sType == "planet" then
      t.planets      = true
   elseif sType == "spacecraft" then
      t.spacecraft   = true
      t.planets      = true
   end
   celestia:setlabelflags(t)
end

function celestia_cleanup_callback()
-- restore user's settings
   celestia:setrenderflags(oldrf)
   celestia:setlabelflags(oldlf)
   celestia:setorbitflags(oldof)
end

function saveOldSettings()
--store user's settings
   oldrf      = celestia:getrenderflags()
   oldlf      = celestia:getlabelflags()
   oldof      = celestia:getorbitflags()
end

--*********************************************************************
--            some support functions
--*********************************************************************
function selectpair()
   local s   = celestia:getselection()
   local   p   

   if s:type() == "null" then         -- nothing selected
      s   = nil
   else
      if s:type() == "star" then       -- stars does not have parent, so select its first son
         p   = s                  
         s   = p:getchildren()[1]
         if s ~= nil then         -- this star does not have orbiting bodies
            celestia:select(s)
         end
      else                  -- the son is another class of body
         p   = s:getinfo()['parent']   
      end
   end
   
   if s   ~= nil then   
      setRenderOn()
      setOrbits(s:type())
      setLabels(s:type())
      celestia:unmarkall()
      s:mark ("#ff0000", "diamond", 30)
      p:mark ("#ffff00", "diamond", 50)
   else
      setRenderOff()
   end

   return s, p
end

function viewBodyAbove (norm, sPos, pPos)
   -- puts the observer perpendicular to the local orbit's plane, "above" the parent.
   -- because our line of sight will be locally orthogonal to the orbital plane, we can use simple trigonometrics to
   -- calculate the distance we need to be placed from parent, in order to see both son and parent with an angle = ANGLE
   local radiusVector = sPos-pPos
   local distance    = radiusVector:length() * HOW_FAR / TAN_ANGLE

   -- place the viewer
   obs:setposition ( pPos + (norm * distance) )

   -- keep the top of the screen oriented along radiusvector
   obs:lookat (pPos, radiusVector)

end


--*********************************************************************
--               speed functions
--*********************************************************************
function tspeed (p0, p1, dt)
   -- tangential speed
   return (p1-p0):length() * KM_MLY / dt
end

function aspeed (p0, p1, dt)
   -- angular speed
   return math.deg( math.acos( (p0*p1)/( p0:length() * p1:length() ) ) ) / dt
end

function rspeed (p0, p1, dt)
   -- radial speed
   return (p1:length() - p0:length()) * KM_MLY  / dt
end

--*********************************************************************
--               main script
--*********************************************************************
--save Celestia status
saveOldSettings()

repeat
   son, par = selectpair()
   
   if son ~= nil then
      parPos   = par:getposition()
      sonPos   = son:getposition()

      sonPosR1   = parPos - sonPos
      t1      = celestia:gettime()
      t0      = t1 - DELTA_D
      sonPosR0   = par:getposition(t0) - son:getposition(t0)   -- parent is moving

      ts      = tspeed (sonPosR0, sonPosR1, DELTA_S)
      as      = aspeed (sonPosR0, sonPosR1, DELTA_S)
      rs      = rspeed (sonPosR0, sonPosR1, DELTA_S)

      celestia:print(string.format("Speeds relative to " .. par:name() .. "\nTangential: %4.5f km/sec\nAngular: %14.15f deg/sec\nRadial: %4.5f km/sec", ts, as, rs), DELAY, HORIG, VORIG)

      -- compute a normal for the orbit. It must be done inside the loop because of local planar accuracy
      normal   = (sonPosR0^sonPosR1):normalize()

      -- if you want to move freely, just comment this line
      viewBodyAbove (normal, sonPos, parPos)

   end

   wait()

until 1==2


EDIT: I have modified the way speeds are computed (see below)

Please post bugs, suggestions and comments here.


Bye

[don]

Looks like a fun script Toti! Hope to try it out in the next couple of days.

Man, I need more time in a day so I can spend more time playing with scripts.

Hope to catch up with all the new scripts soon...

[Harry]

Nice...

Maybe you can improve the precision of the angular speed display by using positions which are 1 second apart (or ten, or ....). You can set t0 to current time - 10s, and use object:getposition(t0) instead of saving the values from the last loop. Does this help?

Harald

[don]

Hi Toti,

Can you explain how these velocity values are useful in Celestia?

Thanks!

[Toti]

...and use object:getposition(t0)...

I didn't know that you could pass parameters to this method. Thank you very much, this indeed helped: I have modified the script following your suggestions.
I also removed the abs() function from the angular speed.

Can you explain how these velocity values are useful in Celestia?

They show as precise numbers some things that you can intuitively perceive in Celestia, like how fast a body is moving in its orbit, when and how much does it accelerate, etc.

Example: for Sol-Earth
*] Tangential speed is measured along the orbital trajectory: how many kilometers of Earth's orbit path were done in the last second.

*] Earth sweeps 360 degrees in one year (an entire orbit around Sol). The angular speed tells you how many degrees has Earth swept in the last second, in its orbit around Sol.

*] The radial vector is a line starting at Sol's center and ending at Earth's center. If the orbit is circular and its origin is Sol's center, then this line has constant length.
But Earth's orbit is elliptical, so this line's length is constantly changing.
The radial speed measures how many kilometers the length of this line has varied in the last second.

Of course the above "definitions" can be applied to other bodies: you can compute angular speed of Dactyl's orbit around Ida, radial speed of ISS around Earth, etc.

You can try this, it is quite interesting:
Set time to 10-19-1989
Select Galileo
Speed up time (X 100,000 -- X 1,000,000)
Now you can see the speed behaviour during the mission: there are some spectacular accelerations during flybys: tangential speed boosts from 22 to 36 km/s in a few hours (only seconds for you); you can also watch the radial speed readout: it will tell you when the spacecraft is moving away from Sol (positive values) and when it is getting closer (negative values). It will show you also how fast Galileo is "climbing" or "falling" its orbit.
You can clearly see how tangential speed is increased when the body "falls" (i.e. the radial speed is negative), and how is decreased when it "climbs" (radial speed positive)
The values of this two speeds can help you appreciate how the spacecraft's kinetic energy is converted to potential energy and vice versa. (This is true for all bodies: there must be some internet pages where this effects are explained in an accessible way)
There is something more: the angular speed shows you a corollary of Kepler's Second Law: the spacecraft sweeps equal orbital areas in equal times. The radial vector gets shorter when the craft gets closer to Sol, so the swept angle must be larger to compensate. Thus, Galileo's angular speed increases as it gets closer to Sol.

Note: there is a point in the orbit when Galileo is captured by Jupiter (around Feb 1995), when the spacecraft's trajectory is abruptly distorted (this happens several times). You can see this because there are some local normal inversions (the viewer goes above and below the trajectory's plane in a fraction of a second).
This has two causes:
1) a limitation of the way that orbit's normals are calculated in the script (I quickly tested other methods but found too much precision reduction with them (there was a "trembling" everywhere))
2) actually, there is a local variation of the normal (i.e. the trajectory bends "up" or "down" when near Jupiter)

You can repeat the experiment with other bodies (try Halley, Voyager 1 and 2, EPS Eri b, Neptune and Pluto, etc.)

I hope this helps

Bye

[don]

Excellent explanations Toti, thank you! And some excellent experiments to try also!

Learning about our solar system and space sure is FUN in Celestia, and with the help of all the great folks here, like you.

[Toti]

Thank you very much for your compliment, Don.
I have modified the script again, and now it displays planet labels when a spacecraft is selected (so you can easily recognize flybys and other maneuvers)

And I have found a few pages with explanations on orbits and energy:
http://www.space.edu/projects/book/chapter5.html
http://www.go.ednet.ns.ca/~larry/orbits/orbits.html
http://kosmoi.com/Science/Physics/Mechanics/tpecp1.shtml#m3

Bye

[overex]

Thank you very much for your compliment, Don.
I have modified the script again, and now it displays planet labels when a spacecraft is selected (so you can easily recognize flybys and other maneuvers)

And I have found a few pages with explanations on orbits and energy:
http://www.space.edu/projects/book/chapter5.html
http://www.go.ednet.ns.ca/~larry/orbits/orbits.html
http://kosmoi.com/Science/Physics/Mechanics/tpecp1.shtml#m3

Bye

Toti great work !

[Toti]

Thanks, man.