path: root/ride/elevation_outlier.go

package ride

import (
	"bytes"
	"fmt"
	"io"
	"math"
	"os"

	"github.com/tkrajina/gpxgo/gpx"
)

// ElevationOutlier describes a contiguous burst of unrealistic grades.
type ElevationOutlier struct {
	StartIdx  int     `json:"startIndex"`
	EndIdx    int     `json:"endIndex"`
	StartDist float64 `json:"startDistM"`
	EndDist   float64 `json:"endDistM"`
	MaxGrade  float64 `json:"maxGrade"`
	DeltaE    float64 `json:"deltaEM"`
}

// DefaultGradeThreshold is the grade beyond which a segment is considered
// unrealistic for road cycling (50% = 0.5). Real road walls rarely exceed 25-30%
// even on short ramps; 50% over <30m is a strong tunnel/GPS signal.
const DefaultGradeThreshold = 0.5

// MaxGapM is the maximum horizontal distance for a segment to be considered
// for outlier detection. Longer segments with high grade are kept (e.g. a
// sustained 10% climb over 500m is fine).
const MaxGapM = 30

// DetectElevationOutliers flags indices where |Δe|/Δd exceeds threshold within MaxGapM.
// Consecutive flagged indices are merged into intervals.
func DetectElevationOutliers(distances, elevations []float64, gradeThreshold float64) []ElevationOutlier {
	if gradeThreshold <= 0 {
		gradeThreshold = DefaultGradeThreshold
	}
	if len(distances) < 2 || len(distances) != len(elevations) {
		return nil
	}
	// Collect raw flagged indices (the end point of each bad segment).
	flagged := make([]int, 0)
	maxGradeByEnd := map[int]float64{}
	deltaByEnd := map[int]float64{}
	for i := 1; i < len(distances); i++ {
		dd := distances[i] - distances[i-1]
		if dd <= 0 || dd > MaxGapM {
			continue
		}
		de := elevations[i] - elevations[i-1]
		// Ignore tiny elevation jitter (e.g. 0.3m over 0.5m = 60% but noise).
		if math.Abs(de) < 1.5 {
			continue
		}
		grade := math.Abs(de) / dd
		if grade > gradeThreshold {
			flagged = append(flagged, i)
			maxGradeByEnd[i] = grade
			deltaByEnd[i] = de
		}
	}
	if len(flagged) == 0 {
		return nil
	}
	// Merge consecutive flagged indices into intervals.
	// A gap of 1 index means contiguous; allow a small gap of 2-3 non-flagged
	// points to still merge (tunnel burst often alternates).
	outliers := make([]ElevationOutlier, 0)
	start := flagged[0]
	prev := flagged[0]
	maxGrade := maxGradeByEnd[prev]
	sumDelta := deltaByEnd[prev]
	for _, idx := range flagged[1:] {
		if idx-prev <= 10 {
			// same burst (tunnel glitch spans ~125 points with small clean stretches)
			prev = idx
			if maxGradeByEnd[idx] > maxGrade {
				maxGrade = maxGradeByEnd[idx]
			}
			sumDelta += deltaByEnd[idx]
			continue
		}
		// close current burst: expand by one anchor point on each side when possible
		sIdx := max(start-1, 0)
		eIdx := min(prev+1, len(distances)-1)
		// For burst we report the flagged core, but anchors are for fix suggestion.
		// Keep outlier as core; fix will use anchors.
		outliers = append(outliers, ElevationOutlier{
			StartIdx:  sIdx,
			EndIdx:    eIdx,
			StartDist: distances[sIdx],
			EndDist:   distances[eIdx],
			MaxGrade:  maxGrade,
			DeltaE:    sumDelta,
		})
		start = idx
		prev = idx
		maxGrade = maxGradeByEnd[idx]
		sumDelta = deltaByEnd[idx]
	}
	sIdx := max(start-1, 0)
	eIdx := min(prev+1, len(distances)-1)
	outliers = append(outliers, ElevationOutlier{
		StartIdx:  sIdx,
		EndIdx:    eIdx,
		StartDist: distances[sIdx],
		EndDist:   distances[eIdx],
		MaxGrade:  maxGrade,
		DeltaE:    sumDelta,
	})
	return outliers
}

// FixZone is a suggested linear-interpolation interval for the CLI.
type FixZone struct {
	StartIdx int     `json:"startIndex"`
	EndIdx   int     `json:"endIndex"`
	StartKm  float64 `json:"startKm"`
	EndKm    float64 `json:"endKm"`
	Reason   string  `json:"reason"`
}

// SuggestFixZones converts outliers into fix zones with human-readable reason.
func SuggestFixZones(outliers []ElevationOutlier) []FixZone {
	zones := make([]FixZone, 0, len(outliers))
	for _, o := range outliers {
		zones = append(zones, FixZone{
			StartIdx: o.StartIdx,
			EndIdx:   o.EndIdx,
			StartKm:  o.StartDist / 1000,
			EndKm:    o.EndDist / 1000,
			Reason:   fmt.Sprintf("max grade %.0f%% Δe %.0fm", o.MaxGrade*100, o.DeltaE),
		})
	}
	return zones
}

// LinearFix interpolates elevations linearly between each zone's anchors by distance.
// It mutates elevations in place and returns the number of points modified.
func LinearFix(distances, elevations []float64, zones []FixZone) int {
	modified := 0
	for _, z := range zones {
		if z.StartIdx < 0 || z.EndIdx >= len(distances) || z.StartIdx >= z.EndIdx {
			continue
		}
		startDist := distances[z.StartIdx]
		endDist := distances[z.EndIdx]
		startEle := elevations[z.StartIdx]
		endEle := elevations[z.EndIdx]
		span := endDist - startDist
		if span <= 0 {
			continue
		}
		for i := z.StartIdx + 1; i < z.EndIdx; i++ {
			t := (distances[i] - startDist) / span
			elevations[i] = startEle + t*(endEle-startEle)
			modified++
		}
	}
	return modified
}

// FixRide mutates ride elevations by linearly fixing the given zones.
func FixRide(r *Ride, zones []FixZone) int {
	return LinearFix(r.Distances, r.Elevations, zones)
}

// ResolveIndex converts a user-provided km or index string to a ride index.
// Accepts "123" (index) or "4.5" / "4.5km" (distance in km).
func ResolveIndex(distances []float64, raw string) (int, error) {
	// Try index first if no "km" suffix and value is integer-like.
	var kmVal float64
	var isKm bool
	s := raw
	if len(s) > 2 && s[len(s)-2:] == "km" {
		s = s[:len(s)-2]
		isKm = true
	}
	// If isKm or contains dot, treat as km.
	if isKm || containsDot(s) {
		_, err := fmt.Sscan(s, &kmVal)
		if err != nil {
			return 0, fmt.Errorf("invalid distance %q", raw)
		}
		target := kmVal * 1000
		// nearest index by distance
		best := 0
		bestDiff := math.Abs(distances[0] - target)
		for i, d := range distances {
			diff := math.Abs(d - target)
			if diff < bestDiff {
				bestDiff = diff
				best = i
			}
		}
		return best, nil
	}
	var idx int
	_, err := fmt.Sscan(s, &idx)
	if err != nil {
		return 0, fmt.Errorf("invalid index %q", raw)
	}
	if idx < 0 || idx >= len(distances) {
		return 0, fmt.Errorf("index %d out of range [0,%d)", idx, len(distances))
	}
	return idx, nil
}

func containsDot(s string) bool {
	for _, c := range s {
		if c == '.' {
			return true
		}
	}
	return false
}

// EncodeGPX writes the ride as a GPX file. Used by elevation fix to persist
// corrected elevations while preserving coordinates and timestamps.
func EncodeGPX(r Ride, w io.Writer) error {
	gpxContent := &gpx.GPX{
		Version: "1.1",
		Creator: "biking_home",
	}
	track := gpx.GPXTrack{
		Type: "Ride",
	}
	segment := gpx.GPXTrackSegment{
		Points: make([]gpx.GPXPoint, r.Len()),
	}
	for i := 0; i < r.Len(); i++ {
		segment.Points[i] = gpx.GPXPoint{
			Point: gpx.Point{
				Latitude:  r.Coords[i].Lat,
				Longitude: r.Coords[i].Lon,
				Elevation: *gpx.NewNullableFloat64(r.Elevations[i]),
			},
			Timestamp: r.timestamps[i],
		}
	}
	track.Segments = append(track.Segments, segment)
	gpxContent.Tracks = append(gpxContent.Tracks, track)
	xmlBytes, err := gpxContent.ToXml(gpx.ToXmlParams{Version: "1.1", Indent: true})
	if err != nil {
		return err
	}
	_, err = w.Write(xmlBytes)
	return err
}

// FixRideFiles rewrites FIT and/or GPX files for ride with linearly fixed elevations.
// Zones are applied to the parsed ride, then both files are rewritten (if present)
// keeping a .bak copy. It does not touch the DB; caller should Save afterwards.
func FixRideFiles(ride Ride, zones []FixZone, gpxPath, fitPath string) error {
	// Apply fix to a copy so caller can inspect.
	fixed := ride
	// Copy elevations slice to avoid mutating original ride's backing array aliasing?
	elevCopy := make([]float64, len(ride.Elevations))
	copy(elevCopy, ride.Elevations)
	fixed.Elevations = elevCopy
	LinearFix(fixed.Distances, fixed.Elevations, zones)

	if gpxPath != "" {
		if err := backupAndWrite(gpxPath, func(w io.Writer) error { return EncodeGPX(fixed, w) }); err != nil {
			return fmt.Errorf("write gpx %q: %w", gpxPath, err)
		}
	}
	if fitPath != "" {
		if err := backupAndWrite(fitPath, func(w io.Writer) error { return EncodeFIT(fixed, w) }); err != nil {
			return fmt.Errorf("write fit %q: %w", fitPath, err)
		}
	}
	// If only one of the two exists, create the missing counterpart from fixed ride
	// to keep both in sync (optional, but helpful).
	if gpxPath != "" && fitPath == "" {
		// try to derive fit path
		derived := deriveFitPath(gpxPath)
		if err := writeNew(derived, func(w io.Writer) error { return EncodeFIT(fixed, w) }); err != nil {
			return err
		}
	}
	if fitPath != "" && gpxPath == "" {
		derived := deriveGpxPath(fitPath)
		if err := writeNew(derived, func(w io.Writer) error { return EncodeGPX(fixed, w) }); err != nil {
			return err
		}
	}
	return nil
}

func backupAndWrite(path string, write func(io.Writer) error) error {
	if _, err := os.Stat(path); err == nil {
		bak := path + ".bak"
		// Keep first backup only.
		if _, err := os.Stat(bak); os.IsNotExist(err) {
			data, err := os.ReadFile(path)
			if err != nil {
				return err
			}
			if err := os.WriteFile(bak, data, 0o600); err != nil {
				return err
			}
		}
	}
	var buf bytes.Buffer
	if err := write(&buf); err != nil {
		return err
	}
	return os.WriteFile(path, buf.Bytes(), 0o600)
}

func writeNew(path string, write func(io.Writer) error) error {
	if _, err := os.Stat(path); err == nil {
		return nil
	}
	var buf bytes.Buffer
	if err := write(&buf); err != nil {
		return err
	}
	return os.WriteFile(path, buf.Bytes(), 0o600)
}

func deriveFitPath(gpxPath string) string {
	if len(gpxPath) > 4 && gpxPath[len(gpxPath)-4:] == ".gpx" {
		return gpxPath[:len(gpxPath)-4] + ".fit"
	}
	return gpxPath + ".fit"
}

func deriveGpxPath(fitPath string) string {
	if len(fitPath) > 4 && fitPath[len(fitPath)-4:] == ".fit" {
		return fitPath[:len(fitPath)-4] + ".gpx"
	}
	return fitPath + ".gpx"
}