diff options
Diffstat (limited to 'ride/elevation_outlier.go')
| -rw-r--r-- | ride/elevation_outlier.go | 336 |
1 files changed, 336 insertions, 0 deletions
diff --git a/ride/elevation_outlier.go b/ride/elevation_outlier.go new file mode 100644 index 0000000..aa4c7ab --- /dev/null +++ b/ride/elevation_outlier.go @@ -0,0 +1,336 @@ +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" +} |