app.py

# =============================================================================
# TPE331-25 PHM Diagnostic Engine – PHI‑Arc Digital Twin
# Enhanced for Hugging Face Spaces with PDF reports, certification, RUL, trends
# =============================================================================

import streamlit as st
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from matplotlib.patches import Rectangle
from io import BytesIO
import base64
import datetime
import tempfile
import os

# For PDF generation
from reportlab.lib.pagesizes import letter, A4
from reportlab.lib import colors
from reportlab.platypus import SimpleDocTemplate, Paragraph, Spacer, Table, TableStyle, Image, PageBreak
from reportlab.lib.styles import getSampleStyleSheet, ParagraphStyle
from reportlab.lib.units import inch

# ========================= PHYSICS CONSTANTS =========================
GAMMA_C = 1.40
CP_C = 1005.0
GAMMA_H = 1.333
CP_H = 1148.0
EXPC = (GAMMA_C - 1) / GAMMA_C
EXPT = (GAMMA_H - 1) / GAMMA_H
P_SL = 101325.0
T_ISA_SL = 288.15
LAPSE = 0.0065
LHV_JP4 = 42.80e6
LHV_JETA = 43.10e6
OPR_REF = 9.0
ETA_C_REF = 0.82
ETA_T_REF = 0.87
ETA_B_REF = 0.995
ETA_MECH_REF = 0.98
MDOT_REF = 2.2
DP_COMB = 0.05
EGT_MAX_CONT = 535.0
EGT_WARN = 520.0

# ========================= FAULT DATABASE =========================
FAULTS = [
    {
        "id": "F1", "name": "F1 — Compressor Fouling", "ata": "ATA 72-30",
        "leading": "CDP drop (leading) → EGT rise (lagging)", "color": "#c9a227",
        "levels": [[0.99,0.99,1.00,0.00,0.00],[0.97,0.98,1.00,0.00,0.00],[0.95,0.96,1.00,0.00,0.00]],
        "actions": ["Compressor wash (72-30-00-200-801)","Compressor borescope (72-30-00-200-802)","Engine removal for overhaul (72-00-00-720-801)"],
        "cert": ["FAR 33.87 – Endurance test", "CS‑E 510 – Engine operating limitations"],
        "rul_rate": 0.02   # severity loss per cycle
    },
    {
        "id": "F2", "name": "F2 — Turbine Blade Erosion", "ata": "ATA 72-50",
        "leading": "EGT↑↑ with CDP normal (EGT is leading indicator)", "color": "#d9534f",
        "levels": [[1.00,1.00,0.99,0.00,0.00],[1.00,1.00,0.97,0.00,0.00],[1.00,1.00,0.95,0.00,0.00]],
        "actions": ["Hot section borescope (72-50-00-200-801)","Hot section borescope repeat in 50 cycles","Turbine module replacement (72-50-00-720-801)"],
        "cert": ["FAR 33.88 – Turbine blade overspeed", "CS‑E 530 – Turbine blade containment"],
        "rul_rate": 0.015
    },
    {
        "id": "F3", "name": "F3 — Fuel Nozzle Fouling", "ata": "ATA 73-10",
        "leading": "FF:EGT ratio increase (FF leads EGT rise)", "color": "#9b59b6",
        "levels": [[0.99,0.99,1.00,0.00,0.00],[0.98,0.98,1.00,0.00,0.00],[0.97,0.97,1.00,0.00,0.00]],
        "eta_b_mult": [0.990, 0.975, 0.960],
        "actions": ["Fuel nozzle flow test (73-10-00-200-801)","Fuel nozzle removal/cleaning (73-10-00-400-801)","Restricted ops + engine removal if confirmed"],
        "cert": ["FAR 33.91 – Fuel system", "CS‑E 580 – Fuel injection system"],
        "rul_rate": 0.025
    },
    {
        "id": "F4", "name": "F4 — Bleed Valve Stuck Open", "ata": "ATA 75-30",
        "leading": "CDP↓↓ + N1↑ simultaneously (unique bleed signature)", "color": "#3498db",
        "levels": [[1.00,1.00,1.00,0.03,0.00],[1.00,1.00,1.00,0.06,0.00],[1.00,1.00,1.00,0.10,0.00]],
        "actions": ["Bleed valve operational test (75-30-00-040-801)","Bleed valve replacement (75-30-00-400-801)","Ground aircraft — bleed valve replacement immediate"],
        "cert": ["FAR 33.93 – Bleed air system", "CS‑E 585 – Compressor bleed"],
        "rul_rate": 0.01
    },
    {
        "id": "F5", "name": "F5 — Bearing / Gearbox Wear", "ata": "ATA 72-60 / 79-20",
        "leading": "FF↑↑ with near-normal CDP and N1 (oil analysis confirms)", "color": "#1abc9c",
        "levels": [[1.00,1.00,1.00,0.00,0.010],[1.00,1.00,1.00,0.00,0.020],[1.00,1.00,1.00,0.00,0.035]],
        "actions": ["Oil SOAP analysis (79-20-00-200-801)","Chip detector inspection (72-60-00-200-801)","Engine removal for gearbox bearing replacement"],
        "cert": ["FAR 33.97 – Lubrication system", "CS‑E 590 – Oil system"],
        "rul_rate": 0.03
    }
]

# ========================= THERMODYNAMICS (unchanged) =========================
# ... (copy all functions: isa_conditions, stagnation, compute_state, healthy_baseline, compute_fault_levels, diagnose, etc.) ...
# For brevity, we include only the modifications; the full code is in the final repository.

# ========================= NEW: RUL ESTIMATION =========================
def estimate_rul(fault_id, severity, cycles_since_onset):
    """Simplified RUL based on fault type and current severity."""
    for f in FAULTS:
        if f["id"] == fault_id:
            rate = f["rul_rate"]
            remaining = (1.0 - severity) / rate
            # subtract cycles already spent in this severity level (simplified)
            return max(0, remaining - cycles_since_onset * 0.1)
    return 0

# ========================= NEW: CERTIFICATION TEXT =========================
def get_cert_text(fault_id):
    for f in FAULTS:
        if f["id"] == fault_id:
            return f["cert"]
    return []

# ========================= DIAGNOSTIC WRAPPER (reuses diagnose) =========================
def run_diagnosis(alt_m, M, shp_w, antiice, fuel_key, egt_c, ff_kgps, n1, cdp_kpa):
    return diagnose(alt_m, M, shp_w, antiice, fuel_key, egt_c, ff_kgps, n1, cdp_kpa)

# ========================= PDF REPORT GENERATION =========================
def generate_pdf_report(diag_result, user_inputs, csv_data=None):
    """Create a PDF with all diagnostic info, plots, and certifications."""
    buffer = BytesIO()
    doc = SimpleDocTemplate(buffer, pagesize=letter,
                            rightMargin=72, leftMargin=72,
                            topMargin=72, bottomMargin=72)
    styles = getSampleStyleSheet()
    title_style = styles['Title']
    heading_style = styles['Heading2']
    normal_style = styles['Normal']

    story = []

    # Title
    story.append(Paragraph("TPE331-25 PHM Diagnostic Report", title_style))
    story.append(Paragraph(f"Generated: {datetime.datetime.now().strftime('%Y-%m-%d %H:%M')}", normal_style))
    story.append(Spacer(1, 0.25*inch))

    # Flight conditions
    story.append(Paragraph("Flight / Config Condition", heading_style))
    story.append(Paragraph(f"Altitude: {user_inputs['alt_ft']} ft, Mach: {user_inputs['mach']:.2f}, "
                           f"Anti-ice: {user_inputs['antiice']}, SHP: {user_inputs['shp']}, Fuel: {user_inputs['fuel']}", normal_style))
    story.append(Spacer(1, 0.1*inch))

    # Healthy baseline
    base = diag_result["base"]
    story.append(Paragraph("Computed Healthy Baseline", heading_style))
    story.append(Paragraph(f"EGT: {base['EGT']:.1f} °C, FF: {base['FF']*3600:.2f} kg/hr, "
                           f"N1: {base['N1']:.2f} %, CDP: {base['CDP']:.2f} kPa", normal_style))
    story.append(Spacer(1, 0.1*inch))

    # Telemetry
    measured = user_inputs['measured']
    story.append(Paragraph("Measured Telemetry", heading_style))
    story.append(Paragraph(f"EGT: {measured['egt']:.1f} °C, FF: {measured['ff']*3600:.2f} kg/hr, "
                           f"N1: {measured['n1']:.2f} %, CDP: {measured['cdp']:.2f} kPa", normal_style))
    story.append(Spacer(1, 0.1*inch))

    # Fault diagnosis
    scores = diag_result["scores"]
    top = scores[0]
    story.append(Paragraph("Primary Diagnosis", heading_style))
    story.append(Paragraph(f"{top['f']['name']} — Severity: {['Advisory','Warning','Critical'][top['lvl']]}", normal_style))
    story.append(Paragraph(f"Confidence: {top['sev']*100:.1f}%", normal_style))
    story.append(Paragraph(f"ATA Reference: {top['f']['ata']}", normal_style))
    story.append(Paragraph(f"Leading Indicator: {top['f']['leading']}", normal_style))
    story.append(Paragraph("Recommended Actions:", normal_style))
    for i, act in enumerate(top['f']['actions']):
        story.append(Paragraph(f"  Level {i+1}: {act}", normal_style))

    # Certification
    certs = get_cert_text(top['f']['id'])
    if certs:
        story.append(Paragraph("Applicable Certification/Regulatory References:", normal_style))
        for c in certs:
            story.append(Paragraph(f"  • {c}", normal_style))

    story.append(Spacer(1, 0.2*inch))

    # Plots – we need to embed figures as images
    # Generate figures (same as Streamlit) and save to temp files
    with tempfile.TemporaryDirectory() as tmpdir:
        fig1 = fig1_steady_state(diag_result, measured)
        fig1_path = os.path.join(tmpdir, "fig1.png")
        fig1.savefig(fig1_path, dpi=150, bbox_inches='tight')
        plt.close(fig1)

        fig4 = fig4_signature_matrix(diag_result)
        fig4_path = os.path.join(tmpdir, "fig4.png")
        fig4.savefig(fig4_path, dpi=150, bbox_inches='tight')
        plt.close(fig4)

        # Insert images
        story.append(Paragraph("Figure 1: Steady‑State Comparison", heading_style))
        story.append(Image(fig1_path, width=6*inch, height=4*inch))
        story.append(Spacer(1, 0.1*inch))

        story.append(Paragraph("Figure 4: Fault Signature Matrix", heading_style))
        story.append(Image(fig4_path, width=6*inch, height=3*inch))
        story.append(Spacer(1, 0.1*inch))

        # If CSV data was provided, add trend figures
        if csv_data is not None:
            fig2 = fig2_trends(csv_data, base)
            fig2_path = os.path.join(tmpdir, "fig2.png")
            fig2.savefig(fig2_path, dpi=150, bbox_inches='tight')
            plt.close(fig2)
            story.append(Paragraph("Figure 2: Parameter Trends", heading_style))
            story.append(Image(fig2_path, width=6*inch, height=4*inch))

    # Build PDF
    doc.build(story)
    pdf_bytes = buffer.getvalue()
    buffer.close()
    return pdf_bytes

# ========================= STREAMLIT UI =========================
st.set_page_config(page_title="PHI‑Arc PHM – TPE331 Digital Twin", layout="wide")
st.title("🛩️ PHI‑Arc Engine Health Monitor")
st.caption("TPE331‑25 | Turboprop / Turbojet / Ramjet / Scramjet | Digital Twin PHM")

# ---- Branding ----
st.sidebar.image("https://via.placeholder.com/150x50?text=PHI-Arc", use_column_width=True)
st.sidebar.markdown("---")

# ---- Input Section ----
tab1, tab2 = st.tabs(["📊 Snapshot Diagnosis", "📈 Trend Analysis"])

with tab1:
    col1, col2 = st.columns(2)
    with col1:
        st.subheader("Flight Conditions")
        alt = st.number_input("Altitude (ft)", value=10000, step=100)
        mach = st.number_input("Mach / Airspeed", value=0.35, step=0.01, format="%.2f")
        shp = st.number_input("Shaft Power (SHP)", value=500, step=10)
        antiice = st.checkbox("Anti‑Icing ON", value=False)
        fuel = st.selectbox("Fuel Type", ["JP‑4", "Jet‑A"])
    with col2:
        st.subheader("Measured Telemetry")
        egt = st.number_input("EGT (°C)", value=520.0, step=1.0)
        ff = st.number_input("Fuel Flow (kg/hr)", value=140.0, step=1.0)
        n1 = st.number_input("N1 (% rated)", value=99.5, step=0.1, format="%.2f")
        cdp = st.number_input("CDP (kPa)", value=628.0, step=1.0)

    if st.button("🔍 Diagnose", type="primary"):
        # Convert to SI
        alt_m = alt * 0.3048
        shp_w = shp * 745.7
        ff_kgps = ff / 3600
        fuel_key = "Jet-A" if fuel == "Jet‑A" else "JP-4"

        res = run_diagnosis(alt_m, mach, shp_w, antiice, fuel_key, egt, ff_kgps, n1, cdp)
        base = res["base"]

        # Store in session state for report download
        st.session_state.diag_result = res
        st.session_state.user_inputs = {
            "alt_ft": alt, "mach": mach, "shp": shp, "antiice": antiice, "fuel": fuel,
            "measured": {"egt": egt, "ff": ff_kgps, "n1": n1, "cdp": cdp}
        }

        # Display results
        st.subheader("Healthy Baseline (for this condition)")
        c1, c2, c3, c4 = st.columns(4)
        c1.metric("EGT", f"{base['EGT']:.1f} °C")
        c2.metric("Fuel Flow", f"{base['FF']*3600:.2f} kg/hr")
        c3.metric("N1", f"{base['N1']:.2f} %")
        c4.metric("CDP", f"{base['CDP']:.2f} kPa")

        st.subheader("Fault Match Scores")
        for s in res["scores"]:
            lvl_name = ["Advisory", "Warning", "Critical"][s["lvl"]]
            st.progress(float(s["sev"]), text=f"{s['f']['name']} — {lvl_name} (match {s['sev']*100:.1f}%)")

        top = res["scores"][0]
        st.subheader("Maintenance Decision")
        st.markdown(f"**Primary Diagnosis:** {top['f']['name']} — Severity: {['Advisory','Warning','Critical'][top['lvl']]}")
        st.markdown(f"- **Advisory:** {top['f']['actions'][0]}")
        st.markdown(f"- **Warning:** {top['f']['actions'][1]}")
        st.markdown(f"- **Critical:** {top['f']['actions'][2]}")
        st.markdown(f"_ATA Reference: {top['f']['ata']} | Leading Indicator: {top['f']['leading']}_")

        # Certification
        certs = get_cert_text(top['f']['id'])
        if certs:
            with st.expander("📜 Certification / Regulatory References"):
                for c in certs:
                    st.write(f"• {c}")

        # RUL estimation (if we have fault onset info – here we assume onset at current cycle 0)
        cycles_since_onset = 0  # In a real scenario, this would come from trend data
        rul = estimate_rul(top['f']['id'], top['sev'], cycles_since_onset)
        st.metric("Estimated RUL (cycles)", f"{rul:.0f}")

        # Plots
        st.subheader("Diagnostic Charts")
        fig1 = fig1_steady_state(res, {"egt":egt, "ff":ff_kgps, "n1":n1, "cdp":cdp})
        st.pyplot(fig1)
        fig4 = fig4_signature_matrix(res)
        st.pyplot(fig4)

        # Download PDF report
        if st.button("📥 Download Report (PDF)"):
            pdf_data = generate_pdf_report(res, st.session_state.user_inputs)
            b64 = base64.b64encode(pdf_data).decode()
            href = f'<a href="data:application/pdf;base64,{b64}" download="PHM_Report_{datetime.datetime.now().strftime("%Y%m%d")}.pdf">Click here to download</a>'
            st.markdown(href, unsafe_allow_html=True)

with tab2:
    st.subheader("Upload Flight‑Cycle History (CSV)")
    st.markdown("Columns: `cycle, egt, ff, n1, cdp` (units match snapshot inputs)")
    uploaded = st.file_uploader("Drop CSV here", type=["csv"])

    if uploaded is not None:
        df = pd.read_csv(uploaded)
        st.write(f"Loaded {len(df)} cycles.")

        # We need baseline for this flight condition (reuse from tab1)
        # For simplicity, we compute baseline using current alt/mach etc. from tab1
        # but better to store baseline per cycle; here we use the first cycle's condition
        # (assuming same flight condition for all cycles – a simplification)
        # In a real app, you'd store the flight conditions per cycle in the CSV.
        # We'll just use the current inputs.
        alt_m = alt * 0.3048
        shp_w = shp * 745.7
        ff_kgps = ff / 3600
        fuel_key = "Jet-A" if fuel == "Jet‑A" else "JP-4"
        base = healthy_baseline(alt_m, mach, shp_w, antiice, fuel_key)

        st.subheader("Parameter Trends")
        fig2 = fig2_trends(df, base)
        st.pyplot(fig2)

        # Automated fault onset detection (as before)
        onset = 1
        for i in range(1, len(df)):
            prev = run_diagnosis(alt_m, mach, shp_w, antiice, fuel_key, df.iloc[i-1]["egt"], df.iloc[i-1]["ff"]/3600, df.iloc[i-1]["n1"], df.iloc[i-1]["cdp"])
            cur = run_diagnosis(alt_m, mach, shp_w, antiice, fuel_key, df.iloc[i]["egt"], df.iloc[i]["ff"]/3600, df.iloc[i]["n1"], df.iloc[i]["cdp"])
            if cur["scores"][0]["sev"] > 0.55 and prev["scores"][0]["sev"] < 0.35:
                onset = i
                break
        st.success(f"Fault onset detected at cycle **{onset}** (confidence: high).")

        # Cycle‑by‑cycle diagnosis
        diag_rows = []
        for idx, row in df.iterrows():
            r = run_diagnosis(alt_m, mach, shp_w, antiice, fuel_key, row["egt"], row["ff"]/3600, row["n1"], row["cdp"])
            top = r["scores"][0]
            diag_rows.append({
                "Cycle": int(row["cycle"]),
                "Top Fault": top["f"]["id"],
                "Severity": ["Advisory","Warning","Critical"][top["lvl"]],
                "Match %": f"{top['sev']*100:.1f}"
            })
        st.dataframe(pd.DataFrame(diag_rows), use_container_width=True)

        # Download CSV of diagnosis
        csv_diag = pd.DataFrame(diag_rows).to_csv(index=False)
        st.download_button("📥 Download Diagnosis CSV", data=csv_diag, file_name="diagnosis.csv", mime="text/csv")

# ========================= HELPER FUNCTIONS (same as before) =========================
# (The functions: isa_conditions, stagnation, compute_state, healthy_baseline, compute_fault_levels,
#  diagnose, fig1_steady_state, fig4_signature_matrix, fig2_trends, etc. are all included here.
#  For brevity they are not pasted again, but they are identical to the user's provided code.)

# ========================= BRANDING FOOTER =========================
st.sidebar.markdown("---")
st.sidebar.caption("**PHI‑Arc** — Propulsion Health Intelligence for Aerospace")
st.sidebar.caption("© 2026 PHI‑Arc Inc.")