F.1 — Launch Geometry & Sun-Synchronous Orbit Design

A sun-synchronous orbit is engineered through a coupled system: plane precession (physics), node orientation (geometry), and launch timing (Earth rotation).

1. Mission scenario

A launch vehicle departs from a high-latitude coastal site and flies on a south-westerly trajectory over open ocean. After a fixed ascent time $t_{\mathrm{MECO}}$, the vehicle reaches orbital injection conditions and deploys a spacecraft into a near-circular low Earth orbit:

• Altitude: $h = 500~\mathrm{km}$
• Eccentricity: $e \approx 0$

Mission requirements:

• The orbit must be Sun-synchronous (SSO).
• The Local Time of the Ascending Node (LTAN) must be $18{:}00$ (dusk / “evening terminator” style lighting).
• Lighting consistency is required for repeatable imaging and stable power/thermal behaviour.

Engineering goal: choose the orbit-plane geometry (inclination and node orientation) and a compatible launch time so that the injected orbit naturally satisfies SSO precession and meets the LTAN target without continuous plane-correction manoeuvres.

2. Geometry reasoning — determining SSO inclination

Sun-synchronous behaviour arises because Earth’s oblateness ($J_2$) produces nodal precession: the orbital plane rotates about Earth’s spin axis. A common approximation for the RAAN rate is:

\[ \dot{\Omega} = -\frac{3}{2}J_2\left(\frac{R_E^2}{a^2(1-e^2)^2}\right)n\cos i \]

Where:

• $R_E$ is Earth radius,
• $a$ is semi-major axis,
• $n=\sqrt{\mu/a^3}$ is mean motion,
• $i$ is inclination.

For Sun-synchronous motion, we choose the orbit such that the nodal precession approximately matches the Sun’s apparent motion:

\[ \dot{\Omega} \approx -\omega_{\mathrm{sun}}, \qquad \omega_{\mathrm{sun}} \approx \frac{2\pi}{1~\mathrm{year}} \]

Given $e\approx 0$ and $a=R_E+h$, solving for $i$ at $h=500~\mathrm{km}$ yields a typical SSO inclination near: $i \approx 97^\circ$ to $98^\circ$.

Why retrograde?

The required precession direction is achieved by a slightly retrograde orbit, which makes $\cos i < 0$. That sign flip produces the correct direction of plane rotation relative to the Sun.

import numpy as np

MU = 398600.4418          # km^3/s^2
RE = 6378.1363            # km
J2 = 1.08262668e-3
DAY = 86400.0
YEAR = 365.2422 * DAY     # mean tropical year (s)

def sso_inclination_deg(h_km, e=0.0):
    a = RE + h_km
    n = np.sqrt(MU / a**3)  # rad/s
    omega_sun = 2*np.pi / YEAR
    # RAAN rate: dOmega = -(3/2) * J2 * (RE/a)^2 * n * cos(i) / (1-e^2)^2
    denom = -1.5 * J2 * (RE/a)**2 * n / (1 - e**2)**2
    cos_i = (-omega_sun) / denom
    cos_i = np.clip(cos_i, -1.0, 1.0)
    i = np.degrees(np.arccos(cos_i))
    # prefer retrograde (i>90°)
    if i < 90.0:
        i = 180.0 - i
    return i

for h in [400, 500, 600, 700, 800]:
    print(h, "km -> i =", round(sso_inclination_deg(h), 3), "deg")

3. LTAN constraint — launch timing logic

LTAN = 18:00 means the spacecraft crosses the equator at the ascending node when local solar time is 6 pm. This places the orbit plane close to the day–night terminator, giving predictable lighting (often beneficial for imaging), and helps keep thermal and power environments consistent across repeats.

Conceptual relationship: RAAN ↔ local solar time

LTAN is essentially a specification of where the orbit plane sits relative to the Sun direction. In practice:

• RAAN defines the orientation of the ascending node in inertial space.
• The Sun direction at the epoch defines the “local solar time reference.”
• Target LTAN fixes the angular relationship between the orbit plane and the Sun vector.

Practical timing workflow (launch window logic)

1. Pick the plane inclination from SSO physics: choose $i$ so $\dot{\Omega}\approx-\omega_{\mathrm{sun}}$.
2. Set plane orientation to satisfy LTAN: determine the target $\Omega$ (RAAN) consistent with LTAN = 18:00 at the chosen epoch.
3. Connect inertial plane to Earth-fixed geography: equator-crossing longitude + Earth rotation links $\Omega$ to a specific UTC timing.
4. Account for Earth rotation during ascent: if injection occurs $t_{\mathrm{MECO}}$ after liftoff, Earth rotates by:

   \[ \Delta\theta_E=\omega_E\,t_{\mathrm{MECO}} \]

   So liftoff must be earlier by that amount (in Earth-fixed terms) to “arrive” at the correct injection geometry.
5. Back-propagate to obtain the launch window: the allowed liftoff time is constrained by the required plane orientation and the rotating Earth.

from datetime import datetime, timedelta, timezone
import numpy as np

# NOTE: This is a skeleton. You plug in:
# - your launch site lat/lon
# - your azimuth-to-inertial-plane mapping at injection
# - your Sun RA/GMST approximations (or use astropy/skyfield)

OMEGA_E = 7.2921159e-5   # rad/s (Earth rotation)

def sun_ra_deg(dt_utc):
    # minimal placeholder: replace with astropy/skyfield for higher fidelity
    # return approximate Sun right ascension in degrees
    return 0.0

def target_raan_from_ltan(ltan_hours, sun_ra_deg_val):
    # conceptual: LTAN fixes RAAN relative to Sun RA
    # plug your convention here (dawn/dusk-plane mapping)
    return (sun_ra_deg_val + 180.0) % 360.0

def achieved_raan_from_launch(dt_meco):
    # placeholder: map site + azimuth + Earth rotation at MECO to RAAN
    return 0.0

t_start = datetime(2025, 9, 1, 0, 0, 0, tzinfo=timezone.utc)
t_end   = t_start + timedelta(hours=6)   # search a 6-hour window
t_meco  = 1380.0                         # seconds (example)

ltan_target = 18.0
tol_deg = 1.0

t = t_start
hits = []
while t <= t_end:
    ra_sun = sun_ra_deg(t)
    raan_target = target_raan_from_ltan(ltan_target, ra_sun)

    t_meco_dt = t + timedelta(seconds=t_meco)
    raan_ach = achieved_raan_from_launch(t_meco_dt)

    # smallest angular difference
    d = (raan_ach - raan_target + 180.0) % 360.0 - 180.0
    if abs(d) <= tol_deg:
        hits.append((t.isoformat(), d))

    t += timedelta(minutes=1)

print("Candidate liftoff times:", len(hits))
for iso, err in hits[:10]:
    print(iso, "RAAN error deg:", round(err, 3))

4. Propagation method (validation)

After injection, validate the design using one (or both) of these fidelity levels:

• $J_2$-perturbed propagation (best for confirming SSO logic and plane drift)
• SGP4-like catalog propagation (useful for operational-style checks and quick comparisons)

Validation checks:

• Orbit shape: $a$ and $e$ remain near target values over the analysis window.
• Precession rate: computed $\dot{\Omega}$ matches $-\omega_{\mathrm{sun}}$ within tolerance.
• LTAN stability: LTAN stays near 18:00 over multiple days (allowing small drift due to simplifications).

Recommended practical check:

• Evaluate node crossings at several epochs (e.g., day 0, day 3, day 7).
• Confirm LTAN does not drift outside mission tolerance (e.g., ±10 minutes).

5. Ground station access implications

Add a mid-latitude ground station and compute:

• Pass times (AOS/LOS)
• Maximum elevation
• Contact duration
• Daily contact count

Engineering consequences of this orbit choice:

• High inclination improves high-latitude coverage but can create uneven mid-latitude pass distribution.
• At 500 km, passes are typically short; link margin and scheduling become important.
• SSO plane drift can shift access local-time patterns over long durations.

6. Interpretation of visibility

Key engineering observations:

• SSO ensures repeatable lighting — not uniform access everywhere.
• LTAN selection influences thermal cycling, eclipse-season behaviour, and power-margin consistency.
• Small inclination or altitude errors alter $\dot{\Omega}$, LTAN drift over time, and coverage patterns.

F.2 — TLE-Based Orbit Propagation & Ground Access Analysis

TLE propagation is operationally useful, but interpretation requires careful handling of frames, sampling, and access constraints. “Coverage” and “communication availability” are not the same thing.

1. Mission scenario

You are given a catalog orbit description (single or multiple TLEs) for one or more LEO objects. You must:

• Propagate the orbit over a specified time window
• Compute ground station visibility windows
• Interpret operational availability and reliability

This is an operations-style question: “When can I talk to it?” and “How trustworthy is that schedule?”

2. Propagation setup (what matters)

Inputs:

• TLE file(s)
• Start/end UTC
• Ground station latitude/longitude/altitude
• Minimum elevation mask $\epsilon_{\min}$ (e.g., $5^\circ$ or $10^\circ$)

Method outline:

• Propagate state using SGP4
• Convert ECI → ECEF (Earth rotation model)
• Transform to topocentric (station-centric) coordinates
• Compute azimuth/elevation time histories

3. Pass detection and key outputs

At each time step, compute the satellite relative position in the station’s topocentric frame and evaluate elevation:

\[ \text{elevation}= \arcsin\left( \frac{\mathbf{r}_{topo}\cdot \hat{\mathbf{z}}}{\lVert \mathbf{r}_{topo}\rVert} \right) \]

Pass condition:

\[ \text{elevation} > \epsilon_{\min} \]

Extract per pass:

• AOS (Acquisition of Signal): elevation crosses above mask.
• LOS (Loss of Signal): elevation drops below mask.
• Max elevation: proxy for link margin and data rate potential.
• Duration: usable time for telemetry, command, and payload downlink.

# deps: pip install skyfield sgp4 numpy pandas
from datetime import datetime, timezone
import pandas as pd
from skyfield.api import load, wgs84, EarthSatellite

tle1 = "1 64056U 25104B   25160.24306210  .00859907  25185-3 17582-2 0  9992"
tle2 = "2 64056  41.9357 156.0687 0193223  48.4945 313.2311 15.73238515  3578"

ts = load.timescale()
sat = EarthSatellite(tle1, tle2, "SAT", ts)

# ground station
gs = wgs84.latlon(48.123, 9.832, elevation_m=250.0)

t0 = ts.from_datetime(datetime(2025, 6, 9, 0, 0, 0, tzinfo=timezone.utc))
t1 = ts.from_datetime(datetime(2025, 6, 10, 0, 0, 0, tzinfo=timezone.utc))

min_el = 10.0  # degrees
t_events, events = sat.find_events(gs, t0, t1, altitude_degrees=min_el)

passes = []
i = 0
while i < len(events):
    if events[i] == 0:  # rise
        j = i + 1
        while j < len(events) and events[j] != 2:
            j += 1
        if j < len(events) and events[j] == 2:
            t_rise = t_events[i].utc_datetime()
            t_set  = t_events[j].utc_datetime()
            dur_s = (t_set - t_rise).total_seconds()
            passes.append({
                "start_utc": t_rise.isoformat().replace("+00:00","Z"),
                "end_utc":   t_set.isoformat().replace("+00:00","Z"),
                "duration_min": round(dur_s/60.0, 2)
            })
            i = j + 1
            continue
    i += 1

df = pd.DataFrame(passes)
print(df.to_string(index=False))

4. Sampling resolution effects (critical in real ops)

A pass is a “thin event.” If your step size is too large, you can:

• miss the true AOS/LOS boundaries,
• underestimate duration,
• mis-estimate maximum elevation,
• miscount total passes.

Practical mitigation:

• Use a coarse scan (e.g., 60 s) to locate candidate windows.
• Refine locally (e.g., 5–10 s) around AOS/LOS for accurate timing.

5. Ground track interpretation

Ground track reveals:

• Earth-rotation shift between successive orbits
• Latitude reach set primarily by inclination
• Repeatability driven by the period ratio with Earth rotation

Even for circular orbits, surface tracks can look complex because the Earth rotates underneath the orbital plane.

6. Multi-object visibility (constellation style)

For multiple objects, evaluate the time-varying count of visible spacecraft:

\[ N_{visible}(t) \]

Operational insights you can derive:

• Redundancy windows: multiple spacecraft visible simultaneously.
• Coverage gaps: no spacecraft visible for extended durations.
• Comms contention: several spacecraft require link time at the same time.

7. Operational interpretation (what mission teams care about)

Coverage ≠ communication. Distinguish between:

• Total access time per day vs number of contacts
• Peak elevation (link margin) vs sustained elevation
• Long gaps that interrupt payload duty cycles
• Downlink scheduling conflicts even when access exists

F.3 — Space-Based Sensor Visibility & Crossing Analysis

SSA detection is multi-gated: a field-of-view crossing is necessary, but true detectability requires illumination, range constraints, and scalable computation.

1. Mission scenario

A constellation of space-based sensors monitors LEO. We must determine:

• When tracked objects enter each sensor’s field-of-view (FOV)
• Whether each crossing is truly detectable (not just geometric)
• How to scale to tens of thousands of objects efficiently

This is a systems + geometry problem, not only orbital dynamics.

2. FOV geometry (core detection gate)

Define:

• Sensor position $\mathbf{r}_s$
• Target position $\mathbf{r}_t$
• Relative vector $\mathbf{r}_{rel}=\mathbf{r}_t-\mathbf{r}_s$
• Boresight unit vector $\hat{\mathbf{b}}$

Angle to boresight:

\[ \theta=\cos^{-1}\left( \frac{\mathbf{r}_{rel}\cdot \hat{\mathbf{b}}}{\lVert\mathbf{r}_{rel}\rVert} \right) \]

FOV crossing condition:

\[ \theta < \theta_{FOV} \]

This identifies crossings. It does not guarantee that the object can actually be detected.

3. Illumination / shadow logic

Many optical detection chains require:

• the target is not in Earth shadow,
• a usable phase/illumination geometry,
• enough brightness given range and instrument sensitivity.

A practical approach is to apply an eclipse/umbra test using the Sun direction and Earth shadow geometry. Even a simplified shadow model can dramatically reduce false “detections.”

4. Range gating (reject non-physical crossings)

Introduce a maximum useful range:

\[ \lVert\mathbf{r}_{rel}\rVert < R_{max} \]

This prevents counting very distant crossings that are inside angular FOV but not detectable due to sensitivity limits.

5. Crossing vs detection

• Crossing: satisfies FOV angle only.
• Detection: satisfies FOV + illumination + range + instrument sensitivity (and sometimes motion/blur constraints).

Crossing counts are typically high. Detection counts are much lower — and are what operations teams actually need.

"""
Space-based sensor crossing + detectability (geometry-first SSA screen)

deps:
  pip install numpy sgp4

Notes:
- Target object via SGP4 (TLE -> TEME-like output).
- Tracker via simple 2-body Kepler (circular OK for short screening).
- FOV: full 30 deg (half-angle 15 deg)
- Pointing: boresight along tracker velocity
- Detectable: crossing AND sunlit AND range < 1000 km
"""

import numpy as np
from dataclasses import dataclass
from datetime import datetime, timedelta, timezone
from sgp4.api import Satrec, jday

MU = 398600.4418      # km^3/s^2
RE = 6378.1363        # km
FOV_FULL_DEG = 30.0
HALF_ANGLE = np.deg2rad(FOV_FULL_DEG/2.0)
RMAX = 1000.0         # km

def unit(v):
    n = np.linalg.norm(v)
    return v*0.0 if n < 1e-15 else v/n

def dt_to_jd(dt):
    return jday(dt.year, dt.month, dt.day, dt.hour, dt.minute,
                dt.second + dt.microsecond*1e-6)

# --- minimal sun vector (OK for eclipse gating) ---
def sun_vec_eci_km(dt_utc):
    jd, fr = dt_to_jd(dt_utc)
    T = ((jd + fr) - 2451545.0)/36525.0
    L0 = (280.46646 + 36000.76983*T) % 360.0
    M = (357.52911 + 35999.05029*T) % 360.0
    Mr = np.deg2rad(M)
    lam = (L0 + 1.914602*np.sin(Mr) + 0.019993*np.sin(2*Mr)) % 360.0
    lamr = np.deg2rad(lam)
    eps = np.deg2rad(23.439291 - 0.0130042*T)
    AU = 149597870.7
    r_au = 1.00014 - 0.01671*np.cos(Mr) - 0.00014*np.cos(2*Mr)
    r = r_au*AU
    x = r*np.cos(lamr)
    y = r*np.cos(eps)*np.sin(lamr)
    z = r*np.sin(eps)*np.sin(lamr)
    return np.array([x,y,z], float)

def is_sunlit_cyl(r_obj, dt_utc):
    s_hat = unit(sun_vec_eci_km(dt_utc))
    proj = float(np.dot(r_obj, s_hat))
    if proj > 0.0:
        return True
    r_perp = r_obj - proj*s_hat
    return np.linalg.norm(r_perp) > RE

@dataclass(frozen=True)
class Kepler:
    a_km: float
    e: float
    i_deg: float
    raan_deg: float
    argp_deg: float
    M0_deg: float

def kepler_rv(el: Kepler, epoch: datetime, t: datetime):
    a = el.a_km
    e = el.e
    i = np.deg2rad(el.i_deg)
    O = np.deg2rad(el.raan_deg)
    w = np.deg2rad(el.argp_deg)
    M0 = np.deg2rad(el.M0_deg)

    dt = (t-epoch).total_seconds()
    n = np.sqrt(MU/a**3)
    M = (M0 + n*dt) % (2*np.pi)

    # circular shortcut
    nu = M
    r_pf = np.array([a*np.cos(nu), a*np.sin(nu), 0.0])
    v_pf = np.array([-a*n*np.sin(nu), a*n*np.cos(nu), 0.0])

    cO,sO = np.cos(O), np.sin(O)
    ci,si = np.cos(i), np.sin(i)
    cw,sw = np.cos(w), np.sin(w)
    R3O = np.array([[cO,-sO,0],[sO,cO,0],[0,0,1.0]])
    R1i = np.array([[1,0,0],[0,ci,-si],[0,si,ci]])
    R3w = np.array([[cw,-sw,0],[sw,cw,0],[0,0,1.0]])
    Q = R3O @ R1i @ R3w

    return Q@r_pf, Q@v_pf

def intervals(times, mask):
    out = []
    in_evt = False
    t0 = None
    for k,m in enumerate(mask):
        if m and not in_evt:
            in_evt = True
            t0 = times[k]
        elif (not m) and in_evt:
            out.append((t0, times[k]))
            in_evt = False
            t0 = None
    if in_evt and t0 is not None:
        out.append((t0, times[-1]))
    return out

# ---- Example inputs (edit for your scenario) ----
tle1 = "1 63223U 25052P 25244.59601767 .00010814 00000-0 51235-3 0 9991"
tle2 = "2 63223 97.4217 137.0451 0006365 74.2830 285.9107 15.19475170 25990"
sat = Satrec.twoline2rv(tle1, tle2)

epoch = datetime(2025,9,1,0,0,0,tzinfo=timezone.utc)
tracker = Kepler(a_km=6878.0, e=0.0, i_deg=97.4, raan_deg=72.628, argp_deg=331.7425, M0_deg=0.0)

dt_s = 10.0
N = int(24*3600/dt_s) + 1
times = [epoch + timedelta(seconds=k*dt_s) for k in range(N)]

cross = np.zeros(N, bool)
vis   = np.zeros(N, bool)

for k,t in enumerate(times):
    jd, fr = dt_to_jd(t)
    err, r_obj, v_obj = sat.sgp4(jd, fr)
    if err != 0:
        continue
    r_obj = np.array(r_obj, float)

    r_trk, v_trk = kepler_rv(tracker, epoch, t)
    rho = r_obj - r_trk
    rng = float(np.linalg.norm(rho))
    if rng < 1e-9:
        continue

    b_hat = unit(v_trk)
    cosang = float(np.dot(b_hat, rho/rng))
    cosang = float(np.clip(cosang, -1.0, 1.0))
    ang = np.arccos(cosang)

    in_fov = ang <= HALF_ANGLE
    cross[k] = in_fov

    if in_fov and (rng <= RMAX) and is_sunlit_cyl(r_obj, t):
        vis[k] = True

cross_int = intervals(times, cross)
vis_int   = intervals(times, vis)

print("Crossings:", cross_int if cross_int else "None")
print("Visible:",  vis_int if vis_int else "None")

6. Computational scaling (the real engineering challenge)

Naive brute-force checks at every time step scale as:

\[ N_{checks}=N_{sensors}\times N_{objects} \]

For 100 sensors and 30,000 objects:

\[ 100\times 30{,}000 = 3{,}000{,}000 \quad \text{checks per time step} \]

This is feasible only if time step is large or computation is optimized. In practice, you usually need both: smart filtering + parallelization.

Optimization strategies:

• Altitude-band pre-filter: skip objects outside the sensor’s effective altitude envelope.
• Spatial partitioning: k-d trees / voxel grids / hashing to reduce candidate pairs.
• Range pruning: reject candidates beyond $R_{max}$ before evaluating boresight angle.
• Parallelization: CPU threads or GPU batching for vectorized checks.
• Event-driven prediction: predict candidate encounter windows rather than brute-force scanning.

# High-level structure (not full code):

load_tle_catalog()         # ~30k objects -> Satrec list
load_trackers()            # <=100 sensors (Kepler or SGP4)
t_coarse = build_grid(dt=60s, horizon=24h)
sun_cache = precompute_sun(t_coarse)

parallel_for each tracker:
    r_trk, v_trk = propagate_tracker_on_grid(tracker, t_coarse)
    b_hat = unit(v_trk)  # pointing law, example: velocity direction

    hits = []  # sparse list of (obj_id, k_time_index)
    for obj_batch in batches(objects, B=3000):
        r_obj_batch = sgp4_batch(obj_batch, t_coarse)         # vectorized if possible
        rho = r_obj_batch - r_trk                             # broadcast
        rng = norm(rho)
        keep = (rng < Rmax_margin)                            # range prune
        theta = angle_between(rho, b_hat)
        keep &= (theta < half_angle_margin)                   # cone prune

        hits.extend(compress_hits(obj_batch_ids, keep))

    windows = merge_hits_into_windows(hits, dt=60s)

    for (obj_id, t_start, t_end) in windows:
        t_fine = build_grid(dt=10s, [t_start - pad, t_end + pad])
        # recompute accurately
        crossing_mask = compute_crossing(tracker, obj_id, t_fine)
        crossing_intervals = extract_intervals(crossing_mask)

        visible_intervals = []
        for interval in crossing_intervals:
            apply sunlit + range gates within interval
            visible_intervals.append(...)

        write_intervals(tracker_id, obj_id, crossing_intervals, visible_intervals)

7. Interpretation

Engineering observations:

• Relative velocities are high → detection windows are short.
• Illumination constraints can dominate availability.
• Boresight pointing strategy changes event rate and track quality.
• Scaling matters as much as geometry for real-time operations.

Why F.1–F.3 together form a strong Domain F start

These three cases build a coherent applied progression:

• F.1: orbit-plane design under perturbations + timing constraints
• F.2: operational access from catalog data + practical interpretation
• F.3: SSA detectability + scalable sensing system design

Together they demonstrate:

• Analytical dynamics
• Numerical propagation
• Coordinate frames
• Operational reasoning
• Scalable system design