P015 — Adaptive Multi-Rate GPS/IMU Signal Cleaning

import sys
from pathlib import Path
import yaml
import duckdb
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

NB_DIR = Path.cwd()
TELEFORGE_ROOT = NB_DIR.parent.parent.parent.parent
NOSTOS_ROOT = TELEFORGE_ROOT.parent / "nostos"
sys.path.insert(0, str(NOSTOS_ROOT / "src"))

from nostos.context import TelematicsContext
from nostos.stages.d0_trip_detector import TripDetectorStage
from nostos.stages.d1_gps_cleaner import GPSCleanerStage

print("Stages loaded")

clermont_dir = TELEFORGE_ROOT / "datasets" / "fr_clermont_proto_2025-09"
with open(clermont_dir / "manifest.yaml") as f:
    mf_clermont = yaml.safe_load(f)
parquet = (clermont_dir / mf_clermont['data_files'][0]['path']).resolve()
df_clermont = pd.read_parquet(parquet)
print(f"Clermont : {len(df_clermont):,} samples, {df_clermont['trip_idx'].nunique()} trips")

ctx = TelematicsContext(cfg={}, df=df_clermont.copy(),
                        meta={'device_id': '?', 'country': 'FR', 'hz': 10})
TripDetectorStage().run(ctx)
GPSCleanerStage().run(ctx)

clermont_zupt_pct = 100 * ctx.df['zupt'].mean() if 'zupt' in ctx.df.columns else None
clermont_valid = ctx.df.copy() if 'zupt' in ctx.df.columns else None
print(f"Clermont ZUPT fraction: {clermont_zupt_pct:.1f}%" if clermont_zupt_pct else 'no zupt column')

conn = duckdb.connect(str(NOSTOS_ROOT / "data" / "flespi" / "storage" / "telemetry.duckdb"), read_only=True)
us_meta = conn.execute("""
    SELECT trip_id, device_id, distance_km, n_samples, parquet_path
    FROM trips
    WHERE carrier_state = 'mounted_driving'
      AND distance_km > 1.0
      AND ts_start >= '2026-04-10 11:06:00'
    ORDER BY distance_km DESC
    LIMIT 5
""").fetchdf()
conn.close()

us_zupt_results = []
for _, row in us_meta.iterrows():
    pq = Path(row['parquet_path'])
    if not pq.exists():
        alt = NOSTOS_ROOT / "data" / "flespi" / "trips" / str(row['device_id']) / f"{row['trip_id']}.parquet"
        if alt.exists(): pq = alt
    if not pq.exists(): continue
    trip_df = pd.read_parquet(pq)
    ctx_us = TelematicsContext(cfg={}, df=trip_df.copy(),
                               meta={'device_id': '?', 'country': 'US', 'hz': 10})
    TripDetectorStage().run(ctx_us)
    GPSCleanerStage().run(ctx_us)
    zupt_pct = 100 * ctx_us.df['zupt'].mean() if 'zupt' in ctx_us.df.columns else 0.0
    us_zupt_results.append({
        'trip': row['trip_id'][-20:],
        'distance_km': row['distance_km'],
        'n_samples': len(trip_df),
        'zupt_pct': zupt_pct,
    })

us_zupt_df = pd.DataFrame(us_zupt_results)
print(f"us_greensboro : {len(us_zupt_df)} trips analysés")
us_zupt_df

# Clermont ZUPT par trip
clermont_zupt_per_trip = None
if clermont_valid is not None and 'trip_id' in clermont_valid.columns:
    clermont_zupt_per_trip = (
        clermont_valid.groupby('trip_id')['zupt'].mean() * 100
    ).reset_index(drop=True)

summary = pd.DataFrame([
    {
        'dataset': 'Clermont prototype',
        'n_trips': df_clermont['trip_idx'].nunique() if 'trip_idx' in df_clermont else 0,
        'total_samples': len(df_clermont),
        'zupt_mean_pct': clermont_zupt_pct,
        'zupt_median_per_trip': (
            float(clermont_zupt_per_trip.median()) if clermont_zupt_per_trip is not None else None
        ),
    },
    {
        'dataset': 'us_greensboro Teltonika',
        'n_trips': len(us_zupt_df),
        'total_samples': int(us_zupt_df['n_samples'].sum()) if len(us_zupt_df) else 0,
        'zupt_mean_pct': float(us_zupt_df['zupt_pct'].mean()) if len(us_zupt_df) else None,
        'zupt_median_per_trip': float(us_zupt_df['zupt_pct'].median()) if len(us_zupt_df) else None,
    },
])
summary.round(2)

fig, axes = plt.subplots(1, 2, figsize=(12, 5))

# Left: ZUPT per trip, both datasets
ax = axes[0]
if clermont_zupt_per_trip is not None:
    ax.bar(range(len(clermont_zupt_per_trip)), clermont_zupt_per_trip.values,
           color='#0066CC', alpha=0.75, label=f'Clermont (n={len(clermont_zupt_per_trip)})')
if len(us_zupt_df):
    offset = len(clermont_zupt_per_trip) if clermont_zupt_per_trip is not None else 0
    ax.bar(range(offset, offset + len(us_zupt_df)), us_zupt_df['zupt_pct'].values,
           color='#E87700', alpha=0.75, label=f'us_greensboro (n={len(us_zupt_df)})')
ax.set_xlabel('trip index (all datasets)')
ax.set_ylabel('ZUPT fraction (%)')
ax.set_title('P015 ZUPT fraction per trip')
ax.legend()
ax.grid(True, alpha=0.3, axis='y')

# Right: summary bars
ax = axes[1]
x_labels = ['Clermont\nproto', 'us_greensboro\nTeltonika']
means = [clermont_zupt_pct or 0,
         float(us_zupt_df['zupt_pct'].mean()) if len(us_zupt_df) else 0]
bars = ax.bar(x_labels, means, color=['#0066CC', '#E87700'], alpha=0.85)
ax.set_ylabel('ZUPT fraction (%, mean)')
ax.set_title('ZUPT : signature métier')
ax.grid(True, alpha=0.3, axis='y')
for bar, v in zip(bars, means):
    ax.text(bar.get_x() + bar.get_width()/2, v + 1, f'{v:.1f}%',
            ha='center', va='bottom', fontsize=12, fontweight='bold')

plt.tight_layout()
plt.savefig('p015_cross_dataset_zupt.png', dpi=120, bbox_inches='tight')
plt.show()