NavList:
A Community Devoted to the Preservation and Practice of Celestial Navigation and Other Methods of Traditional Wayfinding
Re: lunars with and without altitudes
From: Dave Walden
Date: 2006 Dec 02, 17:07 -0800
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 , send email to NavList-@fer3.com -~----------~----~----~----~------~----~------~--~---