#!/usr/bin/env python3
"""
Unified Glycogen Depletion Model + Gamification Score
"""

import csv
import re
import statistics
from datetime import datetime, timedelta
from collections import defaultdict

# Load CGM data - flexible column parsing
cgm_data = []
with open('/home/ben/downloads/Blood Glucose-2026-02-19-2026-03-17.csv', 'r') as f:
    reader = csv.DictReader(f)
    cols = reader.fieldnames
    print(f"CSV columns: {cols}")
    
    # Find the right columns
    time_col = None
    bg_col = None
    for c in cols:
        if c.lower() in ('date/time', 'timestamp (yyyy-mm-ddthh:mm:ss)') or (('date' in c.lower() and 'time' in c.lower()) and 'meal' not in c.lower()):
            time_col = c
        if 'glucose' in c.lower() or 'blood' in c.lower():
            if 'mg' in c.lower() or bg_col is None:
                bg_col = c
    
    print(f"Using: time='{time_col}', bg='{bg_col}'")
    
    for row in reader:
        try:
            ts_str = row[time_col].strip()
            # Try multiple formats
            for fmt in ['%Y-%m-%dT%H:%M:%S', '%Y-%m-%d %H:%M:%S', '%m/%d/%Y %H:%M']:
                try:
                    ts = datetime.strptime(ts_str, fmt)
                    break
                except ValueError:
                    continue
            else:
                continue
            
            bg_str = row[bg_col].strip()
            if not bg_str:
                continue
            bg = float(bg_str)
            cgm_data.append((ts, bg))
        except Exception as e:
            pass

cgm_data.sort(key=lambda x: x[0])
print(f"Loaded {len(cgm_data)} CGM readings from {cgm_data[0][0]} to {cgm_data[-1][0]}")

# Parse meal log
with open('/home/ben/ava/ben/health/meal-log.md', 'r') as f:
    meal_text = f.read()

MONTH_MAP = {
    'jan': '01', 'january': '01', 'feb': '02', 'february': '02',
    'mar': '03', 'march': '03', 'apr': '04', 'april': '04',
    'may': '05', 'jun': '06', 'june': '06', 'jul': '07', 'july': '07',
    'aug': '08', 'august': '08', 'sep': '09', 'september': '09',
    'oct': '10', 'october': '10', 'nov': '11', 'november': '11',
    'dec': '12', 'december': '12'
}

meals = []
current_date = None

for line in meal_text.split('\n'):
    line = line.strip()
    
    date_match = re.match(r'^#+\s+(\w+)\s+(\d{1,2})\s*$', line)
    iso_match = re.match(r'^#+\s+(\d{4}-\d{2}-\d{2})', line)
    
    if iso_match:
        current_date = iso_match.group(1)
        continue
    elif date_match:
        month_str = date_match.group(1).lower()
        day = date_match.group(2)
        month = MONTH_MAP.get(month_str)
        if month:
            current_date = f"2026-{month}-{int(day):02d}"
            continue
    
    meal_match = re.match(r'^(\d{1,2}:\d{2})\s*(am|pm)?\s+(.+)', line, re.I)
    if meal_match and current_date:
        time_str = meal_match.group(1)
        food = meal_match.group(3).strip()
        hour, minute = map(int, time_str.split(':'))
        try:
            dt = datetime(2026, int(current_date[5:7]), int(current_date[8:10]), hour, minute)
            meals.append((dt, food))
        except:
            pass

meals.sort(key=lambda x: x[0])
print(f"Parsed {len(meals)} meals from {len(set(m[0].strftime('%Y-%m-%d') for m in meals))} days")

# Exercise data
exercises = [
    {"date": "2026-02-20", "type": "VirtualRide", "duration_min": 34, "kj": 227, "start": "06:00"},
    {"date": "2026-03-01", "type": "VirtualRide", "duration_min": 29, "kj": 227, "start": "06:00"},
    {"date": "2026-03-02", "type": "VirtualRide", "duration_min": 30, "kj": 227, "start": "06:00"},
    {"date": "2026-03-07", "type": "VirtualRide", "duration_min": 35, "kj": 277, "start": "06:00"},
    {"date": "2026-03-08", "type": "VirtualRide", "duration_min": 33, "kj": 228, "start": "06:00"},
    {"date": "2026-03-09", "type": "VirtualRide", "duration_min": 40, "kj": 315, "start": "06:00"},
    {"date": "2026-03-09", "type": "WeightTraining", "duration_min": 21, "kj": 166, "start": "13:00"},
    {"date": "2026-03-14", "type": "VirtualRide", "duration_min": 60, "kj": 449, "start": "10:00"},
    {"date": "2026-03-16", "type": "VirtualRide", "duration_min": 38, "kj": 270, "start": "06:00"},
]

# Group by date
meals_by_date = defaultdict(list)
for t, food in meals:
    meals_by_date[t.strftime('%Y-%m-%d')].append((t, food))

exercises_by_date = defaultdict(list)
for ex in exercises:
    exercises_by_date[ex['date']].append(ex)

cgm_by_date = defaultdict(list)
for t, bg in cgm_data:
    cgm_by_date[t.strftime('%Y-%m-%d')].append((t, bg))

def overnight_bg(date_str):
    readings = cgm_by_date.get(date_str, [])
    overnight = [bg for t, bg in readings if 4 <= t.hour < 6]
    return sum(overnight) / len(overnight) if overnight else None

def daily_mean_bg(date_str):
    readings = cgm_by_date.get(date_str, [])
    return sum(bg for _, bg in readings) / len(readings) if readings else None

def daily_tir(date_str, low=70, high=140):
    readings = cgm_by_date.get(date_str, [])
    if not readings:
        return None
    in_range = sum(1 for _, bg in readings if low <= bg <= high)
    return in_range / len(readings) * 100

def estimate_carbs(food):
    food_lower = food.lower()
    zero = ['water', 'coffee', 'tea', 'olipop', 'coke zero', 'diet', 'zero-sugar',
            'zero sugar', 'poppi', 'chamomile', 'creatine', 'beet root']
    for z in zero:
        if z in food_lower:
            return 2
    
    carbs = 0
    high_carb = {'rice': 50, 'pasta': 45, 'spaghetti': 45, 'noodle': 40, 'pancake': 50,
                 'mamaliga': 35, 'quinoa': 30, 'corn ': 25, 'cornbread': 20,
                 'cookie': 25, 'ice cream': 30, 'cake': 30, 'cobbler': 25,
                 'syrup': 30, 'pudding': 25, 'mashed potato': 25, 'hash brown': 20,
                 'sweet potato': 25, 'dumpling': 20, 'lassi': 20}
    med_carb = {'bread': 20, 'toast': 20, 'croissant': 25, 'pastry': 20, 'roll': 15,
                'pita': 20, 'tortilla': 15, 'wrap': 15, 'burrito': 35,
                'taco': 15, 'sandwich': 25, 'colac': 25, 'pizza': 25, 'pizzolli': 20,
                'fig': 20, 'banana': 25, 'apple': 15, 'fruit': 15, 'mango': 20,
                'berry': 8, 'blueberry': 8, 'strawberry': 5, 'cherry': 10, 'peach': 10,
                'yogurt': 12, 'milk': 12, 'honey': 15, 'jam': 10, 'applesauce': 15,
                'juice': 25, 'orange juice': 25, 'sugar': 10, 'chocolate': 15,
                'popcorn': 10, 'protein bar': 15, 'rice cake': 10, 'granola': 12,
                'oat': 15, 'quiche': 15, 'casserole': 15, 'mac': 30, 'macaroni': 30,
                'bean': 15, 'lentil': 15, 'pinto': 12, 'stew': 10, 'soup': 8}
    low_carb = {'steak': 0, 'sirloin': 0, 'chicken': 0, 'beef': 0, 'pork': 0,
                'fish': 0, 'shrimp': 0, 'egg': 2, 'cheese': 1, 'bacon': 0,
                'salad': 3, 'protein powder': 3, 'salmon': 0, 'sardine': 0,
                'tilapia': 0, 'gyro': 5, 'sausage': 1, 'ham': 1,
                'avocado': 3, 'peanut butter': 4, 'macadamia': 2, 'nut': 3,
                'broccoli': 3, 'mushroom': 2, 'green bean': 4, 'vegetable': 5,
                'carrot': 5, 'cucumber': 2, 'sauerkraut': 2, 'hummus': 5,
                'cream cheese': 1, 'heavy cream': 0, 'butter': 0}
    
    if 'nicnac' in food_lower:
        carbs += 35
    for item, g in high_carb.items():
        if item in food_lower:
            carbs += g
    for item, g in med_carb.items():
        if item in food_lower:
            carbs += g
    
    if carbs == 0:
        for item, g in low_carb.items():
            if item in food_lower:
                carbs += g
        if carbs == 0:
            carbs = 15
    
    return min(carbs, 150)

# ============================================================
# UNIFIED GLYCOGEN MODEL
# ============================================================

def model_day(date_str):
    day_meals = meals_by_date.get(date_str, [])
    day_exercises = exercises_by_date.get(date_str, [])
    
    liver_glyc = 50.0
    
    prev_date = (datetime.strptime(date_str, '%Y-%m-%d') - timedelta(days=1)).strftime('%Y-%m-%d')
    prev_meals = meals_by_date.get(prev_date, [])
    if prev_meals:
        last_meal_time = max(t for t, _ in prev_meals)
        midnight = datetime.strptime(f"{date_str} 00:00", "%Y-%m-%d %H:%M")
        hours_since = (midnight - last_meal_time).total_seconds() / 3600
        last_food = [food for t, food in prev_meals if t == last_meal_time][0]
        last_carbs = estimate_carbs(last_food)
        post_meal_liver = min(80, 30 + last_carbs * 0.4)
        liver_glyc = max(0, post_meal_liver - hours_since * 5)
    
    hourly = []
    hours_therapeutic = 0
    hours_deep = 0
    depletion_integral = 0.0
    
    for hour in range(24):
        h_start = datetime(int(date_str[:4]), int(date_str[5:7]), int(date_str[8:10]), hour)
        h_end = h_start + timedelta(hours=1)
        
        # Exercise depletion
        for ex in day_exercises:
            ex_start = datetime.strptime(f"{ex['date']} {ex['start']}", "%Y-%m-%d %H:%M")
            ex_end = ex_start + timedelta(minutes=ex['duration_min'])
            overlap_start = max(h_start, ex_start)
            overlap_end = min(h_end, ex_end)
            if overlap_end > overlap_start:
                overlap_min = (overlap_end - overlap_start).total_seconds() / 60
                kj_per_min = ex['kj'] / ex['duration_min']
                carb_burn = kj_per_min * overlap_min * 0.24 * 0.6 / 4
                liver_burn = carb_burn * 0.2
                liver_glyc = max(0, liver_glyc - liver_burn)
        
        meals_this_hour = [(t, f) for t, f in day_meals if t.hour == hour]
        if not meals_this_hour:
            liver_glyc = max(0, liver_glyc - 5)
        
        for t, food in meals_this_hour:
            carbs = estimate_carbs(food)
            liver_refill = carbs * 0.45
            for ex in day_exercises:
                ex_end_t = datetime.strptime(f"{ex['date']} {ex['start']}", "%Y-%m-%d %H:%M") + timedelta(minutes=ex['duration_min'])
                gap = (t - ex_end_t).total_seconds() / 3600
                if 0 <= gap <= 1:
                    liver_refill *= 0.6
            liver_glyc = min(90, liver_glyc + liver_refill)
        
        hourly.append((hour, round(liver_glyc, 1)))
        
        if liver_glyc < 30:
            hours_therapeutic += 1
            depletion_integral += (30 - liver_glyc)
        if liver_glyc < 15:
            hours_deep += 1
    
    eating_window = 0
    if day_meals:
        first = min(t for t, _ in day_meals)
        last = max(t for t, _ in day_meals)
        eating_window = (last - first).total_seconds() / 3600
    
    total_carbs = sum(estimate_carbs(f) for _, f in day_meals)
    
    return {
        'date': date_str,
        'hours_therapeutic': hours_therapeutic,
        'hours_deep': hours_deep,
        'depletion_integral': depletion_integral,
        'glycogen_trace': hourly,
        'n_meals': len(day_meals),
        'n_exercises': len(day_exercises),
        'eating_window': eating_window,
        'total_carbs_est': total_carbs,
        'overnight_bg': overnight_bg(date_str),
        'daily_mean': daily_mean_bg(date_str),
        'tir': daily_tir(date_str),
    }

# Run model
all_dates = sorted(set(t.strftime('%Y-%m-%d') for t, _ in cgm_data))
results = []
for d in all_dates:
    results.append(model_day(d))

for i in range(len(results) - 1):
    results[i]['next_overnight_bg'] = results[i+1]['overnight_bg']
for r in results:
    r.setdefault('next_overnight_bg', None)

# ============================================================
# PRINT RESULTS
# ============================================================
print()
print("=" * 85)
print("UNIFIED GLYCOGEN DEPLETION MODEL: Daily Summary")
print("=" * 85)
print()
print(f"{'Date':<12} {'TxHrs':>5} {'Deep':>4} {'DeplInt':>7} {'Meals':>5} {'Ex':>2} {'EatWin':>6} {'Carbs':>5} {'OvnBG':>5} {'Mean':>5} {'TIR':>5}")
print("-" * 72)
for r in results:
    ovn = f"{r['overnight_bg']:.0f}" if r['overnight_bg'] else "  -"
    mean = f"{r['daily_mean']:.0f}" if r['daily_mean'] else "  -"
    tir = f"{r['tir']:.0f}%" if r['tir'] is not None else "  -"
    ew = f"{r['eating_window']:.1f}h" if r['eating_window'] > 0 else "FAST"
    print(f"{r['date']:<12} {r['hours_therapeutic']:>5} {r['hours_deep']:>4} {r['depletion_integral']:>7.0f} {r['n_meals']:>5} {r['n_exercises']:>2} {ew:>6} {r['total_carbs_est']:>5} {ovn:>5} {mean:>5} {tir:>5}")

# ============================================================
# CORRELATIONS
# ============================================================
def corr(pairs):
    if len(pairs) < 3:
        return float('nan')
    x = [p[0] for p in pairs]
    y = [p[1] for p in pairs]
    mx, my = statistics.mean(x), statistics.mean(y)
    sx, sy = statistics.stdev(x), statistics.stdev(y)
    if sx == 0 or sy == 0:
        return 0
    cov = sum((xi - mx) * (yi - my) for xi, yi in zip(x, y)) / (len(x) - 1)
    return cov / (sx * sy)

print()
print("=" * 85)
print("CORRELATIONS")
print("=" * 85)

print("\nvs SAME-DAY MEAN BG:")
p = [(r['hours_therapeutic'], r['daily_mean']) for r in results if r['daily_mean'] is not None]
print(f"  Therapeutic hours:     r = {corr(p):>7.3f}  (n={len(p)})")
p = [(r['depletion_integral'], r['daily_mean']) for r in results if r['daily_mean'] is not None]
print(f"  Depletion integral:    r = {corr(p):>7.3f}  (n={len(p)})")
p = [(r['total_carbs_est'], r['daily_mean']) for r in results if r['daily_mean'] is not None]
print(f"  Estimated carbs:       r = {corr(p):>7.3f}  (n={len(p)})")

print("\nvs NEXT-DAY OVERNIGHT BG:")
p = [(r['hours_therapeutic'], r['next_overnight_bg']) for r in results if r['next_overnight_bg'] is not None]
print(f"  Therapeutic hours:     r = {corr(p):>7.3f}  (n={len(p)})")
p = [(r['depletion_integral'], r['next_overnight_bg']) for r in results if r['next_overnight_bg'] is not None]
print(f"  Depletion integral:    r = {corr(p):>7.3f}  (n={len(p)})")

# Time trend
p = [(i, r['overnight_bg']) for i, r in enumerate(results) if r['overnight_bg'] is not None]
print(f"\n  Time trend vs overnight: r = {corr(p):>7.3f}  (n={len(p)})")

# Therapeutic ratio
print("\nvs THERAPEUTIC RATIO (depletion / carbs):")
p = [(r['depletion_integral'] / max(r['total_carbs_est'], 1), r['daily_mean']) 
     for r in results if r['daily_mean'] is not None and r['n_meals'] > 0]
print(f"  Ratio vs same-day mean:    r = {corr(p):>7.3f}  (n={len(p)})")

# ============================================================
# EXERCISE vs FASTING COMPARISON
# ============================================================
print()
print("=" * 85)
print("UNIFIED MODEL: Exercise vs Fasting as Depletion Methods")
print("=" * 85)
print()

ex_days = [r for r in results if r['n_exercises'] > 0 and r['daily_mean'] is not None]
rest_days = [r for r in results if r['n_exercises'] == 0 and r['daily_mean'] is not None and r['n_meals'] > 0]

if ex_days:
    print(f"Exercise days (n={len(ex_days)}):")
    print(f"  avg therapeutic hrs: {statistics.mean([r['hours_therapeutic'] for r in ex_days]):.1f}")
    print(f"  avg daily mean: {statistics.mean([r['daily_mean'] for r in ex_days]):.1f}")
    print(f"  avg depletion integral: {statistics.mean([r['depletion_integral'] for r in ex_days]):.0f}")

if rest_days:
    print(f"\nRest days with meals (n={len(rest_days)}):")
    print(f"  avg therapeutic hrs: {statistics.mean([r['hours_therapeutic'] for r in rest_days]):.1f}")
    print(f"  avg daily mean: {statistics.mean([r['daily_mean'] for r in rest_days]):.1f}")
    print(f"  avg depletion integral: {statistics.mean([r['depletion_integral'] for r in rest_days]):.0f}")

# ============================================================
# GAMIFICATION SCORE
# ============================================================
print()
print("=" * 85)
print("GAMIFICATION: DAILY DEPLETION SCORE")
print("=" * 85)
print()

streak = 0
scores = []

for r in results:
    base = r['hours_therapeutic']
    depth_bonus = r['hours_deep'] * 0.5
    
    ex_hours = sum(ex['duration_min'] / 60 for ex in exercises_by_date.get(r['date'], []))
    ex_bonus = ex_hours * 2
    
    carb_penalty = 0
    protein_bonus = 0
    for t, food in meals_by_date.get(r['date'], []):
        est = estimate_carbs(food)
        if est > 40:
            carb_penalty += 3
        elif est > 25:
            carb_penalty += 1.5
        elif est > 15:
            carb_penalty += 0.5
        
        food_lower = food.lower()
        protein_words = ['steak', 'chicken', 'beef', 'pork', 'fish', 'shrimp', 'egg',
                        'salmon', 'protein', 'gyro', 'sardine', 'tilapia', 'turkey',
                        'sausage', 'bacon', 'ham']
        if any(pw in food_lower for pw in protein_words):
            protein_bonus += 1
    
    clean_close = 0
    day_meals_list = meals_by_date.get(r['date'], [])
    if day_meals_list:
        last_meal = max(t for t, _ in day_meals_list)
        if last_meal.hour < 19:
            clean_close += 2
        elif last_meal.hour < 20:
            clean_close += 1
    else:
        clean_close += 3
    
    score_raw = base + depth_bonus + ex_bonus + protein_bonus - carb_penalty + clean_close
    if score_raw >= 15:
        streak += 1
    else:
        streak = 0
    streak_bonus = min(streak * 1, 5)
    
    total_score = max(0, score_raw + streak_bonus)
    
    scores.append({
        'date': r['date'],
        'score': total_score,
        'base': base,
        'depth': depth_bonus,
        'exercise': ex_bonus,
        'protein': protein_bonus,
        'carb_pen': carb_penalty,
        'clean': clean_close,
        'streak_b': streak_bonus,
        'streak': streak,
        'overnight_bg': r['overnight_bg'],
        'daily_mean': r['daily_mean'],
        'tir': r['tir'],
    })

print(f"{'Date':<12} {'SCORE':>5} {'Base':>4} {'Dep+':>4} {'Ex+':>4} {'Pro+':>4} {'Crb-':>4} {'Cln+':>4} {'Strk':>4} {'OvnBG':>5} {'Mean':>5}")
print("-" * 68)
for s in scores:
    ovn = f"{s['overnight_bg']:.0f}" if s['overnight_bg'] else "  -"
    mean = f"{s['daily_mean']:.0f}" if s['daily_mean'] else "  -"
    print(f"{s['date']:<12} {s['score']:>5.1f} {s['base']:>4} {s['depth']:>4.1f} {s['exercise']:>4.1f} {s['protein']:>4} {s['carb_pen']:>4.1f} {s['clean']:>4} {s['streak_b']:>4.1f} {ovn:>5} {mean:>5}")

# Score correlations
print()
print("SCORE CORRELATIONS:")
p = [(s['score'], s['overnight_bg']) for s in scores if s['overnight_bg'] is not None]
print(f"  Score vs same-day overnight BG:  r = {corr(p):>7.3f} (n={len(p)})")
p = [(s['score'], s['daily_mean']) for s in scores if s['daily_mean'] is not None]
print(f"  Score vs same-day mean BG:       r = {corr(p):>7.3f} (n={len(p)})")

# Weekly
print()
print("WEEKLY SCORES:")
weekly = defaultdict(list)
for s in scores:
    wk = datetime.strptime(s['date'], '%Y-%m-%d').isocalendar()[1]
    weekly[wk].append(s)

for wk in sorted(weekly):
    wk_scores = weekly[wk]
    avg_score = statistics.mean([s['score'] for s in wk_scores])
    total_score = sum(s['score'] for s in wk_scores)
    ovns = [s['overnight_bg'] for s in wk_scores if s['overnight_bg']]
    means = [s['daily_mean'] for s in wk_scores if s['daily_mean']]
    avg_ovn = statistics.mean(ovns) if ovns else 0
    avg_mean = statistics.mean(means) if means else 0
    print(f"  Week {wk}: avg score {avg_score:>5.1f}, total {total_score:>6.1f}, overnight {avg_ovn:.0f}, mean {avg_mean:.0f} (n={len(wk_scores)})")

# Top/bottom
print()
print("TOP 5 DAYS:")
by_score = sorted(scores, key=lambda s: s['score'], reverse=True)
for s in by_score[:5]:
    ovn = f"{s['overnight_bg']:.0f}" if s['overnight_bg'] else "-"
    mean = f"{s['daily_mean']:.0f}" if s['daily_mean'] else "-"
    tir = f"{s['tir']:.0f}%" if s['tir'] else "-"
    print(f"  {s['date']}: score {s['score']:.1f} — overnight {ovn}, mean {mean}, TIR {tir}")

print()
print("BOTTOM 5 DAYS:")
for s in by_score[-5:]:
    ovn = f"{s['overnight_bg']:.0f}" if s['overnight_bg'] else "-"
    mean = f"{s['daily_mean']:.0f}" if s['daily_mean'] else "-"
    tir = f"{s['tir']:.0f}%" if s['tir'] else "-"
    print(f"  {s['date']}: score {s['score']:.1f} — overnight {ovn}, mean {mean}, TIR {tir}")

# Prescriptive
print()
print("=" * 85)
print("PRESCRIPTIVE TARGETS")
print("=" * 85)
print()

good_days = [s for s in scores if s['daily_mean'] and s['daily_mean'] < 115]
ok_days = [s for s in scores if s['daily_mean'] and 115 <= s['daily_mean'] < 120]
bad_days = [s for s in scores if s['daily_mean'] and s['daily_mean'] >= 120]

if good_days:
    print(f"Days with mean <115 (GOOD):  avg score {statistics.mean([s['score'] for s in good_days]):.1f}, n={len(good_days)}")
if ok_days:
    print(f"Days with mean 115-120 (OK): avg score {statistics.mean([s['score'] for s in ok_days]):.1f}, n={len(ok_days)}")
if bad_days:
    print(f"Days with mean >=120 (BAD):  avg score {statistics.mean([s['score'] for s in bad_days]):.1f}, n={len(bad_days)}")

# Glycogen traces for key days
print()
print("GLYCOGEN TRACES (selected days):")
for d in ['2026-03-14', '2026-03-07', '2026-03-09', '2026-02-24', '2026-03-03']:
    r = next((x for x in results if x['date'] == d), None)
    if r:
        trace = r['glycogen_trace']
        bar = ""
        for h, g in trace:
            if g < 15:
                bar += "█"  # deep therapeutic
            elif g < 30:
                bar += "▓"  # therapeutic  
            elif g < 50:
                bar += "░"  # moderate
            else:
                bar += "·"  # high
        print(f"  {d}: [{bar}] txh={r['hours_therapeutic']} deep={r['hours_deep']} meals={r['n_meals']} ex={r['n_exercises']}")
        print(f"           00      06      12      18      24")
print()
print("Legend: █=deep(<15g) ▓=therapeutic(<30g) ░=moderate(<50g) ·=high(≥50g)")