 # NavList:

## A Community Devoted to the Preservation and Practice of Celestial Navigation and Other Methods of Traditional Wayfinding

Message:αβγ
Message:abc
 Add Images & Files Posting Code: Name: Email:
Re: lunars with and without altitudes
From: Dave Walden
Date: 2006 Dec 02, 17:07 -0800

```Again with apologies for the rough and ready presentation.  Same
explanation as in the past; if I wait 'til it's all cleaned up,
it'll probably never get posted.

Goal here is to produce an Assumed Position-like location for use in
the latitude by GMT and lunar distance method.  We might call it the
First order position (FOP), since it's based on the observed lunar
distances without refraction corrections.  (It can't include
refraction since we need altitude to calculate refraction, and we need
location to calculate altitude, and location is what we're trying to
find.)

We start by writing the equation of a cone whose vertex is at the
earth's center and whose axis points to the North Pole.  The cone
angle is the observed lunar distance to the first star center to
center.  We then rotate the equation of the cone so the axis points at
the first star of interest.  We next translate the cone to the center
of the moon.  Call this equation a(x,y,z).  Repeat for the second lunar
distance-star.  Call this equation b(x,y,z).  Finally, since we want
our FOP to be on the surface of the earth, we write
c(x,y,z)=x^2+y^2+z^2-1.  The point we seek lies on the surface of both
cones and on the surface of the earth.  We now have three equations
with three unknowns.  Since we have  higher order terms in x, y, and z,
we have a number of solutions, most imaginary.  We do have two real
solutions however, and these are the two of interest.  As has been
pointed out, one is on the side of the earth facing the moon and the
other on the opposite side.  Taking the one facing the moon, we now
have the starting point AP=FOP, for the calculation using
interpolation/extrapolation from four corners, or using the
azimuth/intercept method.

For anyone interested, below is a FORTRAN program to find the cone
equations, and a MAXIMA script to solve the system of equations.  The
example is that given by Frank.  (With thanks again!)  The answers at
the end are:             71.2W         38.8W

Be warned that errors/bugs may certainly still remain.

real rotm(3,3)
pi=4.*atan(1.0)
print*,pi
kang=0
c true ld1
theta=38.+(13.1+16.1)/60
aw=0.
np=1
tol=1.e-7
print*,'kang, theta, aw, np, tol'
print*,kang, theta, aw, np, tol
call aptcois (kang, theta, aw, np, tol, qc, qx, qy, qz,
&                    qxy, qyz, qzx, qxx, qyy, qzz, nerr)
print*,'qc, qx, qy, qz,qxy, qyz, qzx, qxx, qyy, qzz, nerr'
print999,qc, qx, qy, qz,qxy, qyz, qzx, qxx, qyy, qzz, nerr
999 format(10f8.5,i2)
n1=2
c star1 gha=3.035 dec=28.01167
th1=(90.-28.01167)
n2=3
th2=-3.035
n3=1
th3=0.
print*,'n1, th1, n2, th2, n3, th3, kang, tol'
print'(i5,f12.8)',n1, th1, n2, th2, n3, th3, ku, tol
call aptrots (n1, th1, n2, th2, n3, th3, kang, tol,
&                    rotm, nerr)
print*,'rotm,nerr'
print'(3f10.6)',rotm
print*,nerr
inv=0
ax=0.
ay=0.
az=0.
call aptrois (rotm, inv, ax, ay, az, qc, qx, qy, qz,
&                    qxy, qyz, qzx, qxx, qyy, qzz, tol, nerr)
print999,qc, qx, qy, qz,qxy, qyz, qzx, qxx, qyy, qzz, nerr
print*,'moon pos'
r=58.852

c moon1 dist=58.852 earth radii gha=47.1 dec=27.64833
ph=pi/180*(90-27.64833)
th=pi/180*(-47.1)
print*,'r,ph,th',r,ph,th
c x=rsin ph cos th  y=r sin phi sin th  z=rcos ph ph=90-dec th=360-gha
dx=r*sin(ph)*cos(th)
dy=r*sin(ph)*sin(th)
dz=r*cos(ph)
c      subroutine apttris (nopt, ax, ay, az, qc, qx, qy, qz,
c     &                    qxy, qyz, qzx, qxx, qyy, qzz, tol, nerr)
nopt=1
print*,'nopt= ',nopt
call apttris (nopt, dx, dy, dz, qc, qx, qy, qz,
&                    qxy, qyz, qzx, qxx, qyy, qzz, tol, nerr)
print*,'after apttris1, translate,  nopt= ',nopt
print*,qc,qx,qy,qz,qxy
print999,qc, qx, qy, qz,qxy, qyz, qzx, qxx, qyy, qzz, nerr

print*,'BEGIN STAR2'
c star2 gha 49.635  dec=16.525
c moon2 gha 46.38 dec=27.645
c pos= 69.6 dec=38.6
c true ld2 11.51992
theta=11.51992/2.
theta=11.51992
theta=11.665
theta=11.+(39.9-16.1)/60.
print*,'kang, theta, aw, np, tol'
print*,kang, theta, aw, np, tol
call aptcois (kang, theta, aw, np, tol, qc, qx, qy, qz,
&                    qxy, qyz, qzx, qxx, qyy, qzz, nerr)
print999,qc, qx, qy, qz,qxy, qyz, qzx, qxx, qyy, qzz, nerr
n1=2
th1=(90.-16.525)
n2=3
th2=-49.635
n3=1
th3=0.
print*,'n1, th1, n2, th2, n3, th3, kang, tol'
print'(i5,f12.8)',n1, th1, n2, th2, n3, th3, ku, tol
call aptrots (n1, th1, n2, th2, n3, th3, kang, tol,
&                    rotm, nerr)
print*,'rotm,nerr'
print'(3f10.6)',rotm
print*,nerr
call aptrois (rotm, inv, ax, ay, az, qc, qx, qy, qz,
&                    qxy, qyz, qzx, qxx, qyy, qzz, tol, nerr)
print999,qc, qx, qy, qz,qxy, qyz, qzx, qxx, qyy, qzz, nerr
ph=pi/180.*(90.-27.645)
th=pi/180.*(-46.38)
dx=r*sin(ph)*cos(th)
dy=r*sin(ph)*sin(th)
dz=r*cos(ph)
print*,'dx,dy,dz',dx,dy,dz
nopt=1
print*,'nopt= ',nopt
call apttris (nopt, dx, dy,dz, qc, qx, qy, qz,
&                    qxy, qyz, qzx, qxx, qyy, qzz, tol, nerr)
print*,'after apttris, translate,  nopt= ',nopt
print*,qc,qx,qy,qz,qxy
print999,qc, qx, qy, qz,qxy, qyz, qzx, qxx, qyy, qzz, nerr

end

C:\G77\bin>a.exe
3.14159274
kang, theta, aw, np, tol
0  38.4866676  0. 1  1.00000001E-007
qc, qx, qy, qz,qxy, qyz, qzx, qxx, qyy, qzz, nerr
0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.61270
0.61270-0.38730
0
n1, th1, n2, th2, n3, th3, kang, tol
2 61.98833084
3 -3.03500009
1  0.00000000
0  0.00000010
rotm,nerr
0.468993 -0.024866 -0.882852
0.052946  0.998597  0.000000
0.881614 -0.046743  0.469651
0
0.00000 0.00000 0.00000 0.00000 0.08242 0.04391-0.82810-0.16454
0.61052 0.39213
0
moon pos
r,ph,th  58.8520012  1.08824193 -0.822050095
nopt=  1
after apttris1, translate,  nopt=  1
15.5745239  37.4410362  42.5059509  9.64573669  0.0824193433
15.5745237.4410442.50595 9.64574 0.08242 0.04391-0.82810-0.16454
0.61052 0.39213
0
BEGIN STAR2
kang, theta, aw, np, tol
0  11.3966665  0. 1  1.00000001E-007
0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.96095
0.96095-0.03905
0
n1, th1, n2, th2, n3, th3, kang, tol
2 73.47499847
3-49.63499832
1  0.00000000
0  0.00000010
rotm,nerr
0.184215 -0.216720 -0.958696
0.761934  0.647655  0.000000
0.620904 -0.730463  0.284434
0
0.00000 0.00000 0.00000 0.00000 0.90709 0.41554-0.35321 0.57543
0.42738 0.88005
0
dx,dy,dz  35.9654045 -37.7410126  27.3068504
nopt=  1
after apttris, translate,  nopt=  1
2.90124512  2.48841476 -11.711647 -19.6767006  0.90709424
2.90125 2.48841**************** 0.90709 0.41554-0.35321 0.57543
0.42738 0.88005
0

%i248) batch("/mnt/floppy/simpline3");

batching #p/mnt/auto/floppy/simpline3.mac
(%i249) a(x, y, z) := 0.39213 z z + 0.61052 y y - 0.16454 x x - 0.8281
z x
+ 0.04391 y z + 0.082419 x y + 9.6457 z + 42.5059 y + 37.441 x +
15.5745
(%o249) a(x, y, z) := 0.39213 z z + 0.61052 y y - 0.16454 x x - 0.8281
z x
+ 0.04391 y z + 0.082419 x y + 9.6457 z + 42.5059 y + 37.441 x +
15.5745
(%i250) b(x, y, z) := 0.88005 z z + 0.42738 y y + 0.57543 x x - 0.35321
z x
+ 0.41554 y z + 0.90709 x y - 19.6767 z - 11.7116 y + 2.48841 x +
2.90125
(%o250) b(x, y, z) := 0.88005 z z + 0.42738 y y + 0.57543 x x - 0.35321
z x
+ 0.41554 y z + 0.90709 x y - 19.6767 z - 11.7116 y + 2.48841 x +
2.90125
(%i251)               c(x, y, z) := - 1 + z z + y y + x x
(%o251)               c(x, y, z) := - 1 + z z + y y + x x
(%i252)        xx : solve([a(x, y, z), b(x, y, z), c(x, y, z)])

`rat' replaced 15.5745 by 31149//2000 = 15.5745

`rat' replaced 37.441 by 37441//1000 = 37.441

`rat' replaced -0.16454 by -2579//15674 = -0.164540002552

`rat' replaced 42.5059 by 28819//678 = 42.50589970501475

`rat' replaced 0.082419 by 1277//15494 = 0.08241900090358

`rat' replaced 0.61052 by 4074//6673 = 0.61052000599431

`rat' replaced 9.6457 by 27143//2814 = 9.645700071073206

`rat' replaced -0.8281 by -6725//8121 = -0.8280999876862

`rat' replaced 0.04391 by 1813//41289 = 0.0439100002422

`rat' replaced 0.39213 by 2003//5108 = 0.39212999216915

`rat' replaced 2.90125 by 2321//800 = 2.90125

`rat' replaced 2.48841 by 13741//5522 = 2.488409996378124

`rat' replaced 0.57543 by 10806//18779 = 0.57543000159753

`rat' replaced -11.7116 by -25543//2181 = -11.7116001834021

`rat' replaced 0.90709 by 4657//5134 = 0.90708998831321

`rat' replaced 0.42738 by 4817//11271 = 0.42738000177447

`rat' replaced -19.6767 by -17650//897 = -19.6767001114827

`rat' replaced -0.35321 by -6415//18162 = -0.3532099988988

`rat' replaced 0.41554 by 5225//12574 = 0.41554000318117

`rat' replaced 0.88005 by 17601//20000 = 0.88005
(%o252) [[z = 186.9795605059265 - 71.03037545338685 %i,
y = 156.6552474191666 %i + 78.65861212645633,
x = 107.8904598096931 %i + 8.888126366487127],
[z = 71.03037545338685 %i + 186.9795605059265,
y = 78.65861212645621 - 156.6552474191664 %i,
x = 8.888126366487402 - 107.8904598096932 %i],
[z = - 62.46498747201671 %i - 14.20858464799212,
y = 20.0682027972187 - 60.980812743406 %i,
x = 4.009354556116049 %i + 83.86293856926758],
[z = 62.46498747201671 %i - 14.20858464799212,
y = 60.98081274340601 %i + 20.06820279721871,
x = 83.86293856926757 - 4.009354556116057 %i],
[z = 16.55536532251499 - 3.291714645530333 %i,
y = 11.28740475395483 %i - 2.988682285874108,
x = - 13.62692118616466 %i - 6.474683747745093],
[z = 3.291714645530333 %i + 16.55536532251499,
y = - 11.28740475395483 %i - 2.988682285874106,
x = 13.62692118616466 %i - 6.474683747745092],
[z = 0.62787950383934, y = - 0.73687926570958, x = 0.25055201974282],
[z = - 0.22476558190844, y = 0.45003531073446, x = - 0.86426198386331]]
(%i253)                        yy : seventh(xx)
(%o253) [z = 0.62787950383934, y = - 0.73687926570958, x =
0.25055201974282]
(%i254)                       zf : rhs(first(yy))
(%o254)                        0.62787950383934
(%i255)                      yf : rhs(second(yy))
(%o255)                       - 0.73687926570958
(%i256)                       xf : rhs(third(yy))
(%o256)                        0.25055201974282
yf
180 atan(--)
xf
(%i257)                     ev(------------, numer)
%pi
(%o257)                       - 71.22105587662567
180 acos(zf)
(%i258)                  ev(90 - ------------, numer)
%pi
(%o258)                        38.8938488718472
(%i259)                         yy : eighth(xx)
(%o259) [z = - 0.22476558190844, y = 0.45003531073446, x = -
0.86426198386331]
(%i260)                       zf : rhs(first(yy))
(%o260)                       - 0.22476558190844
(%i261)                      yf : rhs(second(yy))
(%o261)                        0.45003531073446
(%i262)                       xf : rhs(third(yy))
(%o262)                       - 0.86426198386331
yf
180 atan(--)
xf
(%i263)                     ev(------------, numer)
%pi
(%o263)                       - 27.50672813187197
180 acos(zf)
(%i264)                  ev(90 - ------------, numer)
%pi
(%o264)                       - 12.98909392363139
(%i265)

--~--~---------~--~----~------------~-------~--~----~
To post to this group, send email to NavList@fer3.com
To unsubscribe, send email to NavList-unsubscribe@fer3.com
-~----------~----~----~----~------~----~------~--~---

``` Browse Files

Drop Files ### NavList ### What is NavList? ### Join NavList

 Name: (please, no nicknames or handles) Email:
 Do you want to receive all group messages by email? Yes No
You can also join by posting. Your first on-topic post automatically makes you a member. ### Posting Code

Enter the email address associated with your NavList messages. Your posting code will be emailed to you immediately.
 Email: ### Email Settings

 Posting Code: ### Custom Index

 Subject: Author: Start date: (yyyymm dd) End date: (yyyymm dd)