For what it's worth, I can post the code I ended up with (it's been a while, I can't really tell the whys of everything, now). [Note: judging by all this, in retrospect, you'd probably be better checking for intersect between two geodesics. On Earth at least, the line with constant longitude is a geodesic. Simply initialize it at equator for example, and with your desired longitude, towards azimuth 0.] Anyway. This is the function doing the main computation work: bool Ellipsoid::initPointByLonDiffFromRef(const...
For what it's worth, I can post the code I ended up with (it's been a while, I can't really tell the whys of everything, now). [Note: judging by all this, in retrospect, you'd probably be better checking for intersect between two geodesics. The line with constant longitude is a geodesic. Simply initialize it at equator for example, and with your desired longitude, towards azimuth 0.] Anyway. This is the function doing the main computation work: bool Ellipsoid::initPointByLonDiffFromRef(const LineEqRefParams&...
For what it's worth, I can post the code I ended up with (it's been a while, I can't really tell the whys of everything, now). [Note: judging by all this, in retrospect, you'd probably be better checking for intersect between two geodesics. The line with constant longitude is a geodesic. Simply initialize it at equator for example, and with your desired longitude, towards azimuth 0.] Anyway. This is the function doing the main computation work: bool Ellipsoid::initPointByLonDiffFromRef(const LineEqRefParams&...
In case anyone is interested, my attempt at derivating lambda above was quite wrong. (I miserad a w for omega in Charles Karney's paper). the omega derivative part is much more like : [edit] Omg'(sig) = sin(alp0) / (cos²(sig) + sin²(alp0) * sin²(sig))) Lambda as a function of sigma however is far less well-behaved than was the distance function : in case of geodesics which are quite close to meridians themselves, there are two cases where applying the Newton's method naively would fail : First, when...
In case anyone is interested, my attempt at derivating lambda above was quite wrong. (I miserad a w for omega in Charles Karney's paper). the omega derivative part is much more like : Omg'(sig) = sin(alp0) / (cos²(sig) * (1 + sin²(alp0) * tan²(sig))) That is, in case sig is not too close to +- pi/2 (or wraps of it). I'm still wondering why it works as is without considering a sign switch when sigma gets above pi/2 in absolute value, but... oh well. Lambda as a function of sigma however is far less...
Thanks a lot for pointing me the way. My lack of mathematical background makes those "obvious" methods hard to spot for me. so... omg(sig) = sqrt(1-e²cos²(bet(sig))) = sqrt(1-e²(1 - sin²(bet(sig)))) = sqrt(1-e²(1 - cos²(alp0) * sin²(sig))) omg'(sig) = e² * cos²(alp0) * sin(sig)*cos(sig) / omg(sig) lam(sig) = omg(sig) + f * sin(alp0) * I3(sig) lam'(sig) = omg'(sig) + f * sin(alp0) * I3'(sig) = omg'(sig) + f * sin(alp0) * (2-f) / (1 + (1-f) * sqrt(1 + k² + sin²(sig))) solving for (lam(sig) - lam0)...
Thanks a lot for pointing me the way. My lack of mathematical background makes those "obvious" methods hard to spot for me. so... omg(sig) = sqrt(1-e²cos²(bet(sig))) = sqrt(1-e²(1 - sin²(bet(sig)))) = sqrt(1-e²(1 - cos²(alp0) * sin²(sig))) omg'(sig) = e² * cos²(alp0) * sin(sig)*cos(sig) / omg(sig) lam(sig) = omg(sig) + f * sin(alp0) * I3(sig) lam'(sig) = omg'(sig) + f * sin(alp0) * I3'(sig) = omg'(sig) + f * sin(alp0) * (2-f) / (1 + (1-f) * sqrt(1 + k² + sin²(sig)) solving for (lam(sig) - lam0) ==...
Thanks a lot for pointing me the way. My lack of mathematical background makes those "obvious" methods hard to spot for me. so... omg(sig) = sqrt(1-e²cos²(bet(sig))) = sqrt(1-e²(1 - sin²(bet(sig)))) = sqrt(1-e²(1 - cos²(alp0) * sin²(sig))) omg'(sig) = e² * cos²(alp0) * sin(sig)*cos(sig) / omg(sig) lam(sig) = omg(sig) + f * sin(alp0) * I3(sig) lam'(sig) = omg'(sig) + f * sin(alp0) * I3'(sig) = omg'(sig) + f * sin(alp0) * (2-f) / (1 + (1-f) * sqrt(1 + k² + sin²(sig)) solving for (lam(sig) - lam0) ==...