mars_retrograde_by_o3.pde

// MarsRetrograde.pde — Processing(Java) float‑only version
// --------------------------------------------------------------
// Visualises heliocentric Earth/Mars orbits (left) and
// geocentric, unit‑circle view of Mars as seen from Earth (right).
// 全て float で計算・描画するので、Processing 標準の sin(), cos(), sqrt() が使えます。

// === Orbital parameters (astronomical units, days) ===
final float aE = 1.0;      // semi‑major axis Earth (AU)
final float eE = 0.0167;   // eccentricity Earth
final float TE = 365.256;  // sidereal period Earth (days)

final float aM = 1.523;    // semi‑major axis Mars (AU)
final float eM = 0.0934;   // eccentricity Mars
final float TM = 686.980;  // sidereal period Mars (days)

final float nE = TWO_PI / TE; // mean motion Earth (rad/day)
final float nM = TWO_PI / TM; // mean motion Mars

// === Simulation state ===
float tDays = 0;             // current epoch (days)
final float dt = 2;          // days per frame

// === Display constants ===
final float scaleAU = 120;   // pixels per AU (heliocentric view)
PVector centerHelio;         // set in setup()
PVector centerGeo;
final float geoRadius = 100; // pixel radius of unit circle

void settings() {
  size(800, 400);
  smooth(4);
}

void setup() {
  frameRate(30);
  centerHelio = new PVector(200, height / 2);
  centerGeo   = new PVector(600, height / 2);
}

void draw() {
  background(255);

  // --- Static guides ---
  drawHeliocentricGuides();
  drawGeocentricGuides();

  // --- Planet positions ---
  PVector earth = keplerPosition(aE, eE, nE, tDays, 0);
  PVector mars  = keplerPosition(aM, eM, nM, tDays, PI / 4); // give Mars a phase offset

  // --- Draw heliocentric planets ---
  pushMatrix();
  translate(centerHelio.x, centerHelio.y);
  drawPlanet(earth, color(0, 102, 204));  // Earth blue
  drawPlanet(mars,  color(204, 51, 51));  // Mars red
  popMatrix();

  // --- Geocentric representation ---
  PVector rel = PVector.sub(mars, earth); // Mars‑Earth vector (AU)
  rel.normalize();                        // unit vector
  PVector geoPos = PVector.mult(rel, geoRadius);
  geoPos.add(centerGeo);

  fill(255, 0, 204);
  noStroke();
  ellipse(geoPos.x, geoPos.y, 8, 8);

  // --- Advance simulation time ---
  tDays += dt;
}

// =====================================================
// Helper functions
// =====================================================
PVector keplerPosition(float a, float e, float n, float t, float M0) {
  // Returns heliocentric coordinates (AU) in orbital plane.
  float M = n * t + M0;    // mean anomaly (rad)
  float E = M;             // initial guess (eccentric anomaly)
  for (int i = 0; i < 5; i++) {
    E -= (E - e * sin(E) - M) / (1 - e * cos(E));
  }
  float x = a * (cos(E) - e);
  float y = a * sqrt(1 - e * e) * sin(E);
  return new PVector(x, y);
}

void drawHeliocentricGuides() {
  pushMatrix();
  translate(centerHelio.x, centerHelio.y);

  // Sun
  fill(255, 153, 0);
  noStroke();
  ellipse(0, 0, 10, 10);

  noFill();
  stroke(200);
  strokeWeight(1);
  // Earth orbit (circle guide)
  ellipse(0, 0, aE * 2 * scaleAU, aE * 2 * scaleAU);
  // Mars orbit (ellipse, perihelion at +x)
  float aPix = aM * scaleAU;
  float bPix = aM * sqrt(1 - eM * eM) * scaleAU;
  ellipse(-aM * eM * scaleAU, 0, aPix * 2, bPix * 2);

  popMatrix();
}

void drawGeocentricGuides() {
  noFill();
  stroke(180);
  ellipse(centerGeo.x, centerGeo.y, geoRadius * 2, geoRadius * 2);
}

void drawPlanet(PVector posAU, int c) {
  fill(c);
  noStroke();
  ellipse(posAU.x * scaleAU, posAU.y * scaleAU, 8, 8);
}