Densely sampling orbit of parent planet if on its moon

src/core/modules/Planet.cpp

@@ -4876,9 +4876,12 @@ void Planet::drawOrbit(const StelCore* core)
 		if (!hasValidPositionalData(lastJDE, PositionQuality::OrbitPlotting))
 			return;
 	}
-
-	// Update the orbit positions to the current planet date.
-	computeOrbit();
+	bool fromMoonPerspective = false;
+	if (core->getCurrentPlanet()->pType == isMoon || core->getCurrentPlanet()->pType == isObserver)  // if I am a moon or observer of a planet
+	{
+		// only if we are also drawing my parent
+		fromMoonPerspective = core->getCurrentPlanet()->getParent() == this;
+	}
 
 	const StelProjectorP prj = core->getProjection(StelCore::FrameHeliocentricEclipticJ2000);
 	KeplerOrbit *keplerOrbit=static_cast<KeplerOrbit*>(orbitPtr);
@@ -4890,13 +4893,44 @@ void Planet::drawOrbit(const StelCore* core)
 
 	sPainter.setColor(getCurrentOrbitColor(), orbitFader.getInterstate());
 	Vec3d onscreen;
+
+	if (fromMoonPerspective) {
+		double dateJDE = lastJDE;
+		double calc_date;
+
+		Vec3d myPos = core->getCurrentPlanet()->getHeliocentricEclipticPos(dateJDE);
+		Vec3d myparentPos = getHeliocentricEclipticPos(dateJDE);
+		// orbital plane normal of parent planet, computed with the pos ~quarter of the orbit later
+		Vec3d orbitalNormal = myparentPos ^ getHeliocentricEclipticPos(dateJDE + deltaOrbitJDE*ORBIT_SEGMENTS/4);
+		orbitalNormal.normalize();
+		// project myPos into the orbital plane of parent planet using the normal of orbital plane
+		Vec3d myPosProjected = myPos - myPos.dot(orbitalNormal) * orbitalNormal;
+		// angle between myPos already projected onto the orbital plane of parent planet and the plane
+		// we should concentrate on this theta when sampling orbit
+		double theta = (atan2((myparentPos ^ myPosProjected).dot(orbitalNormal), myparentPos.dot(myPosProjected)))  / M_PI;
+
+		double f;
+		for(int d = 0; d < ORBIT_SEGMENTS; d++)
+		{
+			f = (d - ORBIT_SEGMENTS/2.) / (ORBIT_SEGMENTS/2.);
+			// make sure only sample half of the orbit forward and half backward
+			// the sampling spacing is trial and error, but power of 13 seems to be good (i.e., densely sample around the current date)
+			calc_date = dateJDE + pow(f, 13)*deltaOrbitJDE*ORBIT_SEGMENTS/2. + theta*deltaOrbitJDE*ORBIT_SEGMENTS/2;
+			orbit[d] = getEclipticPos(calc_date);
+		}
+	}
+	else {
+		// Update the orbit positions to the current planet date.
+		computeOrbit();
+		}
 	const Vec3d savePos = orbit[ORBIT_SEGMENTS/2];
 	if (closeOrbit)
 	{
-		if (!keplerOrbit || keplerOrbit->getEccentricity()<=0.3)
+		if ((!keplerOrbit || keplerOrbit->getEccentricity()<=0.3) && !fromMoonPerspective) {
 			// special case - use current Planet position as center vertex so that draws
 			// on its orbit all the time (since orbit is shown as segmented rather than smooth curve)
 			orbit[ORBIT_SEGMENTS/2]=getHeliocentricEclipticPos()+aberrationPush;
+			}
 		orbit[ORBIT_SEGMENTS]=orbit[0];
 	}
 	int nbIter = closeOrbit ? ORBIT_SEGMENTS : ORBIT_SEGMENTS-1;