Re: [Orekit Users] CelestialBody Position extrapolation

[MAISONOBE Luc <luc.maisonobe@c-s.fr>]

Hi Quentin,


Quentin Nénon <q.nenon@gmail.com> a écrit :

> Hello Orekit users,
>
> I find this post interesting because Simon tried to extrapolate the
> position and velocity of celestial bodies in the solar system. We don't
> need to do it as JPL already did it for us for instance, but I tried to
> propagate the orbit of Venus (for instance) during 30 days and I have a
> difference in the position, compare to the translated JPL ephemerides, of
> more than 35,000 km !
>
> I was interested by this when I was comparing results of propagation of
> spacecrafts far from Earth to their "actual" positions, given by the
> Horizon JPL website. I then decided to do the same comparison on a simple
> body, like Venus. You will find enclosed a main class in which I am just
> using keplerian attractions from the Sun and the 8 other planets of the
> Solar System (assuming Pluto is a planet...). I know that I should add the
> solar radiation pressure to the propagator to be more realistic, but this
> doesn't work (or gives very strange results, such as Venus at a distance of
> 100AU from the Sun, using a sphere model with the Venus radius). Even if I
> do not take into account the solar radiation pressure, I do not expect to
> have more than 35,000km of difference after 30 days !
>
> I tried to do the same thing propagating around the Sun, but it almost
> gives the same results (not exactly, as the propagation around the Sun and
> the Earth are not equivalent for the moment, see the
> "PropagationFromDifferentPointOfViewPb" post).
>
> Does anyone have an idea of what I did wrong or what is not working ? This
> may also help the discussion about the problem of propagation from
> different points of view.

Here is a solution I propose to you for this specific problem. I have  
not yet checked if this also correct the former problem you already  
submitted last before.

What you are attempting here is not really propagating one body around  
another one with a few less important perturbing bodies. It is rather  
solving a full n-body problem in free space, where all bodies are  
consider at the same level.

For this to work, you should rather work in a *fixed* frame, not  
really bound to any specific body, and then consider the *absolute*  
attraction of all bodies.
For the fixed frame, I suggest to use ICRF (using  
FramesFactory.getICRF()), which is at the center of mass of the solar  
system. Of course, it is close to Sun since Sun is the more massive  
body around, but it is not exactly at the same position (due mainly to  
Jupiter and Saturn).

In order to compute the absolute gravity pull from each other body  
with respect to your propagated body, you need to build a custom  
BodyAttraction force model (you will find one in the new Junit test  
case I committed a few minutes ago,  
<https://www.orekit.org/forge/projects/orekit/repository/revisions/e59497c723485eebb9e4b3e15100a50fb2a4c107/diff/src/test/java/org/orekit/bodies/SolarBodyTest.java>, have a look at the org.orekit.bodies.SolarBodyTest class and more precisely at the testPropagationVsEphemeris method and BodyAttraction internal class. BodyAttraction was simply created from the original ThirdBodyAttraction, but removing the subtraction of the acceleration with respect to the frame, in order to just compute the absolute acceleration. BEWARE, this works only because we have selected the computation frame at the center of mass! We use this force model for all bodies, including Sun (but obviuosly not the body you propagate of course, which is Venus in your use  
case).

Now since all the gravity fields are considered in the new force model  
instances, and since our central point is an empty point, we don't  
want to have a "central" attraction. Unfortunately, Orekit does not  
allow you to propagate without a central attraction coefficient. So we  
use a quick and dirty hack and set the fictitious Venus orbit around  
IRF and the mu in the propagator to be a negligible value, so this  
impossible to suppress force will in fact be computed as very small.

All this has been set up in the testPropagationVsEphemeris test case.  
The end result is that after one month propagation, the propagator and  
the ephemeris are within 1400m to each other (much smaller than 35000  
km!).

Could you look if this would also work in your other case?

On a side note, you will see that the new test case which is based on  
your own test case has been heavily changed. Here are some random  
advices I can give you and which I have applied:

 - don't set tool stringent accuracy in the integrator, propagating
   orbits with 1e-5 m tolerance on position is too much, even for
   spacecraft around Earth, so its even worse for planets. In your
   case, I have put 1000m
 - don't set the tolerances as litteral arrays which cannot be read,
   use the static NumericalPropagator.tolerances method
 - don't clutter the code with very low level computation of array
   lengths, just rely on List instances and let them be allocated
   on the fly, it's simpler to read (of course, if you think you may
   run low in memory because you have huge collections to store,
   computing the size beforehand can still be useful, so my comment
   is not an absolute rule)
 - instead of building arrays in the step handler and process them
   later (in your case to compute the difference) to build yet another
   array, just build the final one on the fly (I have even completely
   removed it in my test case, because I do not plot it, I just do some
   asserts, but if you plot it just fill up the positionDifferenceNorm
   from within the step handler and don't use an intermediate
   pvResultTable
 - replace litteral constants like 86400 with Constants.JULIAN_DAY

Hope this helps,
Luc

> Thank you very much !
>
> Quentin
>
>
> 2013-03-29 19:48 GMT+01:00 Spörri Simon (s) <simon.spoerri@students.fhnw.ch>
> :
>
> >  Hi Luc,
> >
> >  thanks a lot for your help on even the most amateurish questions!
> > it works charmingly now.
> >
> >  cheers,
> > simon
> >
> >   From: MAISONOBE Luc <luc.maisonobe@c-s.fr>
> > Reply-To: <orekit-users@orekit.org>
> > Date: Fri, 29 Mar 2013 14:19:07 +0100
> > To: <orekit-users@orekit.org>
> > Subject: Re: [Orekit Users] CelestialBody Position extrapolation
> >
> >   simon.spoerri@students.fhnw.ch a écrit :
> >
> >  Dear all,
> >
> >
> >  Hi Simon,
> >
> >
> >  i'd like to create a discrete sampling of celestial body positions for a
> > complete orbit of all bodies. I tried both getting positions from the
> > celestialBody object directly and using a Propagator to get positions of
> > the
> > orbits. However i run into ephemerides constraints and cannot obtain
> > coordinates past around 2030 (Error: out of range date for EARTH_MOON
> > ephemerides: 2030-05-20T14:36:30.712).
> >
> >
> >  I'm not sure about your needs. Orekit uses JPL ephemerides (or INPOP
> > ephemerides) to compute position/velocity of celestial bodies. It does
> > not do any propagation by itself for these bodies.
> >
> >  The limitation you get is simply the limitation of the JPL ephemerides
> > you use. I guess you simply relied on the orekit-data.zip file
> > provided for convenience in the download section of the forge and did
> > not change the default set of data contained in this zip archive. If
> > you want to go past 2030, you can look inside the archive and you will
> > see it contains a file named unxp1962.406. This is a reduce JPL
> > ephemeris file. Remove this file and replace it by the original files
> > from JPL that you will find here:
> > <ftp://ssd.jpl.nasa.gov/pub/eph/planets/SunOS/de406/>. Each file
> > contains ephemerides for 300 years. If you retrieve all 20 files, you
> > will be able to propagate from year -3000 to year +3000. You can look
> > at the wiki page
> > <https://www.orekit.org/forge/projects/orekit/wiki/Configuration> for
> > further information on how Orekit does load such data.
> >
> >  Please note that JPL ephemerides are already sampled. The correspond
> > to large sets of Chebyshev polynomials. Adding another sampling layer
> > on top of that is probably not worth the trouble. Just using the
> > celestial body itself is probably sufficient for most purposes. The
> > current version of celestial bodies that use JPL ephemerides is also
> > thread safe (it was not thread safe in the 5.x series of the Orekit
> > library), and uses efficient cache mechanisms even when several
> > threads uses dates very wide apart from each other. This may be useful
> > in a web application as you cannot know beforehand if two different
> > clients will request positions around the same dates.
> >
> >  best regards,
> > Luc
> >
> >  I try to obtain the positions in the SolarSystemBarycenter Frame.
> > My question is if it is possible at all in Orekit to extrapolate
> > body positions
> > around 150 years into the future/past? I know that this is not the
> > main purpose
> > of Orekit and I do not require very accurate results, nevertheless i
> > wanted to
> > use Orekit for that for convenience...
> > And if it is possible, how should i do it?
> >
> >  thanks a lot in advance,
> > Simon
> >
> >
> >
> >
> >
> >  ----------------------------------------------------------------
> > This message was sent using IMP, the Internet Messaging Program.



----------------------------------------------------------------
This message was sent using IMP, the Internet Messaging Program.