POV-Ray : Newsgroups : povray.advanced-users : Calculating planet positions Server Time
28 Apr 2024 06:42:45 EDT (-0400)
  Calculating planet positions (Message 6 to 15 of 25)  
<<< Previous 5 Messages Goto Latest 10 Messages Next 10 Messages >>>
From: clipka
Subject: Re: Calculating planet positions
Date: 15 Sep 2018 19:03:19
Message: <5b9d8fb7$1@news.povray.org>
Am 15.09.2018 um 22:12 schrieb Bald Eagle:

>> "Modulos the mean anomaly so that -180 <= M <= +180 ..."
...

> I do believe that this is typically done with a macro, and occurs in the POV-Ray
> source code.
> 
> #macro Modulate (_X)
>      #while (_X < 180) #local _X = _X+360; #end
>      #While (_X > 180) #local _X = _X-360; #end
>      _X
> #end


That would be way too time-consuming for very large or very small values.


> Then of course, there's
> http://www.povray.org/documentation/view/3.6.2/458/
> 
> clamp(V, Min, Max). A function that limits a value to a specific range, if it
> goes outside that range it is "clamped" to this range, wrapping around. As the
> input increases or decreases outside the given range, the output will repeatedly
> sweep through that range, making a "sawtooth" waveform.
> Parameters:

That'll indeed do it:

    #include "math.inc"
    #local Foo = clamp(M,-180,180);


Alternatively:

    #local Foo = mod(M+180,360)-180;
    #if (Foo < -180)
      #local Foo = Foo + 260;
    #end

The post-processing `#if` branch is necessary because mod(X,Y) wraps
positive values into the range [0..Y), but negative values into the
range (-Y..0].

[At least on platforms where converting a floating-point number to an
integer rounds towards 0. That's the case for all contemporary platform
I'm aware of, and is an official prerequisite for compiling POV-Ray, but
the C++ standard would also allow for rounding towards negative infinity
instead.]

`clamp()` effectively does the same, but is implemented as a function,
which probably makes it faster than a macro or "in-line" code.


Post a reply to this message

From: Mike Horvath
Subject: Re: Calculating planet positions
Date: 16 Sep 2018 14:26:55
Message: <5b9ea06f$1@news.povray.org>
On 9/15/2018 4:12 PM, Bald Eagle wrote:
> clamp(V, Min, Max). A function that limits a value to a specific range, if it
> goes outside that range it is "clamped" to this range, wrapping around. As the
> input increases or decreases outside the given range, the output will repeatedly
> sweep through that range, making a "sawtooth" waveform.
> Parameters:
> 
> V = Input value.
> Min = Minimum of output range.
> Max = Maximum of output range.
> 
> 
> Bill
> 
> 

Clamp works, thanks!


Mike


Post a reply to this message

From: Mike Horvath
Subject: Re: Calculating planet positions
Date: 16 Sep 2018 15:00:25
Message: <5b9ea849$1@news.povray.org>
https://ssd.jpl.nasa.gov/txt/aprx_pos_planets.pdf

I am a little confused about the "Solution of Kepler's Equation" on Page 
2 of this document. Have I implemented it correctly do you think?

	#local Orrery_Temp_EStar						= 180/pi * Orrery_Temp_Eccentricity;
	#local Orrery_Temp_EccentricAnomaly				= Orrery_Temp_MeanAnomaly + 
Orrery_Temp_EStar * sind(Orrery_Temp_MeanAnomaly);
	#local Orrery_Temp_EccentricAnomalyTolerance	= 10e-6;
	#local Orrery_Temp_EccentricAnomalyDelta		= 0;
	#while (abs(Orrery_Temp_EccentricAnomalyDelta) > 
Orrery_Temp_EccentricAnomalyTolerance)
		#local Orrery_Temp_MeanAnomalyDelta				= Orrery_Temp_MeanAnomaly - 
(Orrery_Temp_EccentricAnomaly - Orrery_Temp_EStar * 
sind(Orrery_Temp_EccentricAnomaly));
		#local Orrery_Temp_EccentricAnomalyDelta		= 
Orrery_Temp_MeanAnomalyDelta/(1 - Orrery_Temp_Eccentricity * 
cosd(Orrery_Temp_EccentricAnomaly));
		#local Orrery_Temp_EccentricAnomaly				= Orrery_Temp_EccentricAnomaly 
+ Orrery_Temp_EccentricAnomalyDelta;
	#end



Mike


Post a reply to this message

From: Mike Horvath
Subject: Re: Calculating planet positions
Date: 16 Sep 2018 16:19:39
Message: <5b9ebadb@news.povray.org>
On 9/16/2018 3:00 PM, Mike Horvath wrote:
> https://ssd.jpl.nasa.gov/txt/aprx_pos_planets.pdf
> 
> I am a little confused about the "Solution of Kepler's Equation" on Page 
> 2 of this document. Have I implemented it correctly do you think?
> 
>      #local Orrery_Temp_EStar                        =
180/pi * 
> Orrery_Temp_Eccentricity;
>      #local Orrery_Temp_EccentricAnomaly                = 
> Orrery_Temp_MeanAnomaly + Orrery_Temp_EStar * 
> sind(Orrery_Temp_MeanAnomaly);
>      #local Orrery_Temp_EccentricAnomalyTolerance    = 10e-6;
>      #local Orrery_Temp_EccentricAnomalyDelta        = 0;
>      #while (abs(Orrery_Temp_EccentricAnomalyDelta) > 
> Orrery_Temp_EccentricAnomalyTolerance)
>          #local Orrery_Temp_MeanAnomalyDelta                = 
> Orrery_Temp_MeanAnomaly - (Orrery_Temp_EccentricAnomaly - 
> Orrery_Temp_EStar * sind(Orrery_Temp_EccentricAnomaly));
>          #local Orrery_Temp_EccentricAnomalyDelta        = 
> Orrery_Temp_MeanAnomalyDelta/(1 - Orrery_Temp_Eccentricity * 
> cosd(Orrery_Temp_EccentricAnomaly));
>          #local Orrery_Temp_EccentricAnomaly                = 
> Orrery_Temp_EccentricAnomaly + Orrery_Temp_EccentricAnomalyDelta;
>      #end
> 
> 
> 
> Mike


I did end up finding a bug in that code, but I also discovered that 
there's not a visible difference with the code working or not working. 
So the major problems I am having lie elsewhere.

:(

Mike


Post a reply to this message

From: Mike Horvath
Subject: Re: Calculating planet positions
Date: 16 Sep 2018 22:48:47
Message: <5b9f160f$1@news.povray.org>
I found the bug, and will be uploading a corrected version to the Object 
Collection shortly.


Mike


Post a reply to this message

From: Mike Horvath
Subject: Re: Calculating planet rotations
Date: 20 Sep 2018 20:16:56
Message: <5ba43878$1@news.povray.org>
On 9/15/2018 1:07 PM, Mike Horvath wrote:
> I'm trying to adapt the following PDF to POV-Ray code for my 
> planetarium/orrery.
> 
> https://ssd.jpl.nasa.gov/txt/aprx_pos_planets.pdf
> 

In Step 6 of that PDF, it gives an equation on how to convert between 
ecliptic coordinates and J2000 coordinates.

What is the inverse of that equation, and how do I change it into a 
series of `rotate` commands instead?

Thanks.



Mike


Post a reply to this message

From: Mike Horvath
Subject: Re: Calculating planet rotations
Date: 21 Sep 2018 00:10:31
Message: <5ba46f37$1@news.povray.org>
On 9/20/2018 8:17 PM, Mike Horvath wrote:
> On 9/15/2018 1:07 PM, Mike Horvath wrote:
>> I'm trying to adapt the following PDF to POV-Ray code for my 
>> planetarium/orrery.
>>
>> https://ssd.jpl.nasa.gov/txt/aprx_pos_planets.pdf
>>
> 
> In Step 6 of that PDF, it gives an equation on how to convert between 
> ecliptic coordinates and J2000 coordinates.
> 
> What is the inverse of that equation, and how do I change it into a 
> series of `rotate` commands instead?
> 
> Thanks.
> 
> 
> 
> Mike

I found the inverse on Wikipedia.

https://en.wikipedia.org/wiki/Ecliptic_coordinate_system#Conversion_between_celestial_coordinate_systems

Still not sure how to replace the equation with a series of `rotate` 
commands though.


Mike


Post a reply to this message

From: Mike Horvath
Subject: Re: Calculating planet rotations
Date: 25 Sep 2018 01:39:44
Message: <5ba9ca20$1@news.povray.org>
Found a document with a table of values:

https://astropedia.astrogeology.usgs.gov/download/Docs/WGCCRE/WGCCRE2009reprint.pdf

The pertinent variables are alpha, delta and W. What do I do with these 
variables? Dunno for sure. This document provides a formula:

https://depositonce.tu-berlin.de/bitstream/11303/6237/4/burmeister_steffi.pdf

I /think/ I need to rotate around the z-axis by W, around the x-axis by 
(90-δ), and around the z-axis again by (90+α), in that order. Lastly, an 
additional rotation (around the x-axis?) by 23.43928 degrees must be 
done to get the body out of the ICRF frame and into the ecliptic frame.

But I still don't know the starting conditions. I.e. before applying the 
transformations, should the globe's North Pole point upward? Should the 
intersection of the Prime Meridian and Equator lie along the x-axis? 
Also, am I wrong, and should I perform the inverse matrix calculations 
instead? Either way, I have not gotten results that match what I see in 
Celestia for the same Julian Date.


Mike


Post a reply to this message

From: Le Forgeron
Subject: Re: Calculating planet rotations
Date: 25 Sep 2018 03:06:49
Message: <5ba9de89$1@news.povray.org>
Le 25/09/2018 à 07:39, Mike Horvath a écrit :
> Found a document with a table of values:
> 
> https://astropedia.astrogeology.usgs.gov/download/Docs/WGCCRE/WGCCRE2009reprint.pdf 
> 
> 
> The pertinent variables are alpha, delta and W. What do I do with these 
> variables? Dunno for sure. This document provides a formula:
> 
> https://depositonce.tu-berlin.de/bitstream/11303/6237/4/burmeister_steffi.pdf 
> 

Read start of 1.3, alpha, delta & W are defined there.

> 
> I /think/ I need to rotate around the z-axis by W, around the x-axis by 
> (90-δ), and around the z-axis again by (90+α), in that order.

W is the "day" part of the planet, adjusting the prime meridian.
Right-handed system, so, yes, W is for z-axis when applied first.
Otherwise, it would be applied along the rotation axis of the planet 
once transformed from z by alpha & delta.

Figure 1.3 , page 11

 From the point gamma on ICRF equator, starting sphere, the axis of 
rotation of the planet is located at 90°+alpha, for an amount of 
90°-delta. (a single tilting of the planet axis, along a single 
perpendicular axis)

If I assert +x is gamma at start, +z north, you have to apply a rotation 
along v_rotate( y, alpha*z) of 90°-delta, then a rotation of W along the 
new pole axis.

You can of course transfer the W part before the other rotation, as it 
is simpler to use z.


> Lastly, an 
> additional rotation (around the x-axis?) by 23.43928 degrees must be 
> done to get the body out of the ICRF frame and into the ecliptic frame.

That's a change of referential, the matrix should be well-known. (i.e. I 
have no clue)
I do not know the x,y,z of ICRF(2) compared to the ecliptic plan and the 
vernal point.

> 
> But I still don't know the starting conditions. I.e. before applying the 
> transformations, should the globe's North Pole point upward? Should the 
> intersection of the Prime Meridian and Equator lie along the x-axis? 

Yes. x2.

> Also, am I wrong, and should I perform the inverse matrix calculations 
> instead? Either way, I have not gotten results that match what I see in 
> Celestia for the same Julian Date.
> 
> 
> Mike


Post a reply to this message

From: Mike Horvath
Subject: Re: Calculating planet rotations
Date: 26 Sep 2018 02:32:47
Message: <5bab280f@news.povray.org>
Thanks for taking a look!


> Figure 1.3 , page 11
> 
>  From the point gamma on ICRF equator, starting sphere, the axis of 
> rotation of the planet is located at 90°+alpha, for an amount of 
> 90°-delta. (a single tilting of the planet axis, along a single 
> perpendicular axis)
> 
> If I assert +x is gamma at start, +z north, you have to apply a rotation 
> along v_rotate( y, alpha*z) of 90°-delta, then a rotation of W along the 
> new pole axis.
> 

I can't quite understand your steps. Do you want me to do this?

   #declare sphere_coo = vrotate(y, (90+alpha_0) * z);
   #declare sphere_coo = vrotate(sphere_coo, (90-delta_0) * x);

OTOH, I have attached my scene so far if you would please look at it. 
The result is not similar to what I see in the program Celestia for the 
same Julian day. (I have attached a screenshot of that program as well.)


Mike


Post a reply to this message


Attachments:
Download 'do_not_delete.pov.txt' (5 KB) Download '2k_earth_daymap.jpg' (453 KB) Download 'celestia_screenshot_01.png' (130 KB)

Preview of image '2k_earth_daymap.jpg'
2k_earth_daymap.jpg

Preview of image 'celestia_screenshot_01.png'
celestia_screenshot_01.png


 

<<< Previous 5 Messages Goto Latest 10 Messages Next 10 Messages >>>

Copyright 2003-2023 Persistence of Vision Raytracer Pty. Ltd.