OpenFin/backend/services/strategy_engine.py

"""
Generic N-leg (1-4) option strategy pricer: entry cost with real broker spread,
scenario repricing at a future horizon on a shocked vol surface, payoff curves,
greeks, and bounded-risk / non-directional checks.

A "leg" dict: {expiry_date, days_to_expiry, strike, option_type ("call"/"put"),
               position ("long"/"short"), quantity}
"""
import math
from typing import Any, Dict, List, Optional

import numpy as np
from scipy.optimize import minimize_scalar

from services.options_pricer import black_scholes
from services.option_chain import find_quote
from services.vol_surface import Surface, ScenarioSurface

DEFAULT_SPREAD_PCT = 0.05  # fallback relative bid/ask spread when no live quote is found
DEFAULT_CONTRACT_SIZE = 100_000  # notional per 1 contract/lot (e.g. a standard FX lot); "quantity" on a leg is the number of these


def to_native(obj: Any) -> Any:
    """Recursively converts numpy scalars (bool_, int64, float64, ...) to native Python
    types. Comparisons/aggregations over numpy-typed inputs (e.g. bid/ask sourced from a
    DB row that came back as a numpy type, or scipy's norm.cdf) can leave a stray
    numpy.bool_/numpy.float64 buried in a nested result — FastAPI's default JSON encoder
    doesn't know those types and fails ('X object is not iterable', then a secondary
    'vars() argument must have __dict__ attribute' from its own fallback). Applied once
    at the API boundary (price_combo/payoff_curves/optimizer results) rather than chasing
    the exact field through the whole pricing pipeline."""
    if isinstance(obj, dict):
        return {k: to_native(v) for k, v in obj.items()}
    if isinstance(obj, (list, tuple)):
        return [to_native(v) for v in obj]
    if isinstance(obj, np.generic):
        return obj.item()
    return obj


def _sign(leg: Dict[str, Any]) -> int:
    return 1 if leg.get("position", "long") == "long" else -1


def _intrinsic(S: float, K: float, option_type: str) -> float:
    return max(0.0, S - K) if option_type == "call" else max(0.0, K - S)


def _quote_spread_pct(quote: Optional[Dict[str, Any]]) -> float:
    if not quote or quote["bid"] <= 0 or quote["ask"] <= 0:
        return DEFAULT_SPREAD_PCT
    mid = (quote["bid"] + quote["ask"]) / 2
    if mid <= 0:
        return DEFAULT_SPREAD_PCT
    return (quote["ask"] - quote["bid"]) / mid


def entry_price(leg: Dict[str, Any], chain_slice: Dict[str, Any], surface_now: Surface, r: float) -> Dict[str, float]:
    """Real execution price (crossing the spread) + theoretical mid, for one leg today."""
    quote = find_quote(chain_slice, leg["expiry_date"], leg["strike"], leg["option_type"])
    T = max(leg["days_to_expiry"], 0.001) / 365
    sigma = surface_now.iv_at(leg["strike"], leg["days_to_expiry"])
    theo_mid = black_scholes(chain_slice["spot"], leg["strike"], T, r, sigma, leg["option_type"])["price"]

    if quote and quote["bid"] > 0 and quote["ask"] > 0:
        exec_price = quote["ask"] if leg.get("position", "long") == "long" else quote["bid"]
        mid = quote["mid"] or theo_mid
    else:
        spread = theo_mid * DEFAULT_SPREAD_PCT
        exec_price = theo_mid + spread / 2 if leg.get("position", "long") == "long" else max(0.0, theo_mid - spread / 2)
        mid = theo_mid

    return {"exec_price": exec_price, "mid": mid, "spread_pct": _quote_spread_pct(quote)}


def value_at(
    legs: List[Dict[str, Any]],
    S: float,
    eval_days_from_now: float,
    surface: Any,
    r: float,
    contract_size: float = DEFAULT_CONTRACT_SIZE,
) -> float:
    """Signed portfolio value (BS reprice for unexpired legs, intrinsic for expired ones)."""
    total = 0.0
    for leg in legs:
        remaining = leg["days_to_expiry"] - eval_days_from_now
        qty = leg.get("quantity", 1)
        sign = _sign(leg)
        if remaining <= 0:
            price = _intrinsic(S, leg["strike"], leg["option_type"])
        else:
            sigma = surface.iv_at(leg["strike"], remaining)
            price = black_scholes(S, leg["strike"], remaining / 365, r, sigma, leg["option_type"])["price"]
        total += sign * price * qty * contract_size
    return total


def greeks_at(legs: List[Dict[str, Any]], S: float, eval_days_from_now: float, surface: Any, r: float) -> Dict[str, float]:
    net = {"delta": 0.0, "gamma": 0.0, "theta": 0.0, "vega": 0.0}
    for leg in legs:
        remaining = max(leg["days_to_expiry"] - eval_days_from_now, 0.001)
        qty = leg.get("quantity", 1)
        sign = _sign(leg)
        sigma = surface.iv_at(leg["strike"], remaining)
        g = black_scholes(S, leg["strike"], remaining / 365, r, sigma, leg["option_type"])
        for k in net:
            net[k] += g[k] * qty * sign
    return {k: round(v, 4) for k, v in net.items()}


def price_combo(
    legs: List[Dict[str, Any]],
    chain_slice: Dict[str, Any],
    surface_now: Surface,
    surface_scenario: ScenarioSurface,
    horizon_days: int,
    r: float = 0.05,
    contract_size: float = DEFAULT_CONTRACT_SIZE,
    precise: bool = True,
) -> Dict[str, Any]:
    spot_now = chain_slice["spot"]
    spot_scenario = surface_scenario.spot

    entry_ref = 0.0
    entry_ref_mid = 0.0
    for leg in legs:
        ep = entry_price(leg, chain_slice, surface_now, r)
        sign = _sign(leg)
        qty = leg.get("quantity", 1)
        entry_ref += sign * ep["exec_price"] * qty * contract_size
        entry_ref_mid += sign * ep["mid"] * qty * contract_size

    # Scenario exit: apply each leg's own bid/ask spread (est. from entry quote) to the
    # theoretical scenario value, since we don't have a live quote for the future date.
    scenario_mid = 0.0
    scenario_exec = 0.0
    for leg in legs:
        remaining = max(leg["days_to_expiry"] - horizon_days, 0.001)
        qty = leg.get("quantity", 1)
        sign = _sign(leg)
        sigma = surface_scenario.iv_at(leg["strike"], remaining)
        theo = black_scholes(spot_scenario, leg["strike"], remaining / 365, r, sigma, leg["option_type"])["price"]
        quote = find_quote(chain_slice, leg["expiry_date"], leg["strike"], leg["option_type"])
        spread_pct = _quote_spread_pct(quote)
        exec_price = theo * (1 - spread_pct / 2) if leg.get("position", "long") == "long" else theo * (1 + spread_pct / 2)
        scenario_mid += sign * theo * qty * contract_size
        scenario_exec += sign * exec_price * qty * contract_size

    net_pnl = scenario_exec - entry_ref
    broker_cost = (entry_ref - entry_ref_mid) + (scenario_mid - scenario_exec)

    # Uses surface_scenario (not surface_now): for a single-expiry combo (condor,
    # butterfly, straddle...) all legs are simultaneously intrinsic at eval_days, so the
    # surface is never actually consulted there and this changes nothing. For a
    # calendar/diagonal spread, the far leg is still alive at the near leg's expiry and
    # DOES need a vol assumption to be priced — using surface_now there silently mixed
    # "today's vol" into a number sitting next to net_pnl (which uses the scenario's
    # shocked vol), producing a max_gain that could be below net_pnl. Pricing both with
    # the same scenario vol view keeps them consistent.
    bounded = check_bounded_risk(legs, entry_ref, surface_scenario, spot_now, r, contract_size, precise)
    delta_now = greeks_at(legs, spot_now, 0, surface_now, r)["delta"]
    delta_scenario = greeks_at(legs, spot_scenario, horizon_days, surface_scenario, r)["delta"]

    return to_native({
        "entry_cost": round(entry_ref, 2),
        "entry_cost_mid": round(entry_ref_mid, 2),
        "scenario_value": round(scenario_exec, 2),
        "scenario_value_mid": round(scenario_mid, 2),
        "net_pnl": round(net_pnl, 2),
        "broker_spread_cost": round(broker_cost, 2),
        "max_gain": bounded["max_gain"],
        "max_loss": bounded["max_loss"],
        "bounded_risk": bounded["bounded"],
        "greeks_now": greeks_at(legs, spot_now, 0, surface_now, r),
        "greeks_scenario": greeks_at(legs, spot_scenario, horizon_days, surface_scenario, r),
        "net_delta_now": delta_now,
        "net_delta_scenario": delta_scenario,
    })


def check_bounded_risk(
    legs: List[Dict[str, Any]], entry_ref: float, surface: Any, spot: float, r: float = 0.05,
    contract_size: float = DEFAULT_CONTRACT_SIZE, precise: bool = True,
) -> Dict[str, Any]:
    """
    Scan a wide log-spaced spot range at expiry and inspect both tails independently for LOSS
    vs GAIN direction. "Bounded risk" only requires the loss side to be capped — a long
    straddle (loss capped at the premium, gain uncapped) is a textbook risk-bounded,
    non-directional trade and must not be excluded just because its gain is open-ended.

    A tail is loss-bounded if moving further to that extreme does not make the P&L any worse
    than a point already well into that tail; symmetrically for gain-bounded.

    "At expiry" means the earliest expiry among the legs. For a single-expiry combo that's
    every leg simultaneously (pure intrinsic, no vol assumption involved at all). For a
    calendar/diagonal spread it's only the near leg — the far leg is still alive and needs
    `surface` to be priced, so max_gain/max_loss there is only as good as that vol input.

    `precise=False` skips only the 1-D refinement below (the expensive part — two
    minimize_scalar calls). The dense near-strike grid stays even in the fast path: a long
    straddle's max loss (or a butterfly's max gain) sits exactly AT a strike too, same as a
    calendar spread's peak, and the wide log-spaced tail grid alone is coarse enough near
    the money to miss it by an order of magnitude, not just a rounding difference (confirmed:
    on a real long straddle, tail-grid-only reported max_loss=-117 when the true V-bottom at
    the strike was -1140 — nearly 10x off, and wrong enough to silently let a candidate past
    a max_loss_cap filter it should have failed). The optimizer's bulk scan (hundreds of
    candidates) uses this fast path for speed; the single-candidate /price call refines
    further with minimize_scalar for pixel-perfect precision.
    """
    eval_days = min(l["days_to_expiry"] for l in legs)

    def f(s: float) -> float:
        return value_at(legs, s, eval_days, surface, r, contract_size) - entry_ref

    # Wide, log-spaced tail grid — used to detect whether the payoff flattens out (bounded)
    # toward either extreme, or keeps moving further away.
    tail_grid = np.geomspace(spot * 0.05, spot * 20, 300)
    tail_values = [f(float(p)) for p in tail_grid]

    tail_n = max(3, len(tail_values) // 30)
    tol = max(abs(entry_ref), 1.0) * 0.01
    lo_edge, lo_in = tail_values[0], tail_values[tail_n]
    hi_edge, hi_in = tail_values[-1], tail_values[-1 - tail_n]

    loss_bounded = (lo_edge >= lo_in - tol) and (hi_edge >= hi_in - tol)
    gain_bounded = (lo_edge <= lo_in + tol) and (hi_edge <= hi_in + tol)

    # Dense linear sweep across the legs' own strikes — needed even in the fast path (see
    # docstring above), only the final 1-D refinement is reserved for precise=True.
    strikes = [l["strike"] for l in legs]
    lo_k, hi_k = min(strikes) * 0.7, max(strikes) * 1.3
    near_grid = np.linspace(max(lo_k, spot * 0.05), min(hi_k, spot * 20), 400)
    near_values = [f(float(p)) for p in near_grid]
    combined_values = tail_values + near_values

    if not precise:
        return {
            "bounded": loss_bounded,
            "max_loss": round(min(combined_values), 2) if loss_bounded else None,
            "max_gain": round(max(combined_values), 2) if gain_bounded else None,
        }

    # Any FIXED grid — however dense — is a different finite sampling of the same
    # continuous curve than whatever grid a chart or another caller uses, so two "close but
    # not identical" readings of the same peak are pretty much guaranteed (confirmed: this
    # grid vs the payoff chart's own grid gave two different peak heights for the same
    # trade). Refine with a bounded 1-D optimizer around the grid's own best point instead.
    grid_all = np.concatenate([tail_grid, near_grid])
    values_all = np.concatenate([tail_values, near_values])
    order = np.argsort(grid_all)
    grid_sorted, values_sorted = grid_all[order], values_all[order]

    # ...then refine with a bounded 1-D optimizer in a NARROW bracket around that grid
    # point (the payoff at a fixed date/vol is smooth in spot, built from Black-Scholes).
    # A wide bracket is actively harmful here: tested directly on a real calendar spread,
    # minimize_scalar given a wide bound (±50% of the strike) converged on a flat false
    # optimum far from the true spike (-65 instead of +1170) because Brent's method isn't
    # guaranteed to find the global optimum over a non-unimodal interval. Anchoring the
    # bracket tightly around the grid's own best point keeps the search unimodal, and
    # taking max()/min() against the grid value means refinement can never do worse.
    def refine(is_max: bool) -> float:
        idx = int(np.argmax(values_sorted)) if is_max else int(np.argmin(values_sorted))
        grid_value = float(values_sorted[idx])
        step = grid_sorted[min(idx + 1, len(grid_sorted) - 1)] - grid_sorted[max(idx - 1, 0)]
        step = float(step) if step > 0 else spot * 0.001
        lo_b, hi_b = grid_sorted[idx] - 3 * step, grid_sorted[idx] + 3 * step
        if lo_b >= hi_b:
            return grid_value
        res = minimize_scalar((lambda s: -f(s)) if is_max else f, bounds=(float(lo_b), float(hi_b)), method="bounded")
        refined = -res.fun if is_max else res.fun
        return max(refined, grid_value) if is_max else min(refined, grid_value)

    return {
        "bounded": loss_bounded,
        "max_loss": round(refine(False), 2) if loss_bounded else None,
        "max_gain": round(refine(True), 2) if gain_bounded else None,
    }


def payoff_curve_expiry(legs: List[Dict[str, Any]], surface_now: Surface, spot: float, n: int = 100, contract_size: float = DEFAULT_CONTRACT_SIZE) -> List[Dict[str, float]]:
    eval_days = min(l["days_to_expiry"] for l in legs)
    prices = np.linspace(spot * 0.5, spot * 1.5, n)
    return [
        {"underlying": round(float(p), 2), "pnl": round(float(value_at(legs, float(p), eval_days, surface_now, 0.05, contract_size)), 2)}
        for p in prices
    ]


def expected_pnl_scenario(
    legs: List[Dict[str, Any]],
    surface_scenario: ScenarioSurface,
    horizon_days: int,
    r: float,
    entry_ref: float,
    n: int = 200,
    contract_size: float = DEFAULT_CONTRACT_SIZE,
) -> float:
    """
    Probability-weighted expected P&L at the scenario date: integrates the payoff over a
    risk-neutral lognormal density for the underlying, centered on the scenario spot with a
    variance driven by the scenario surface's ATM IV over the residual time to the nearest
    leg's expiry (the remaining uncertainty once the scenario date is reached).
    """
    spot_scenario = surface_scenario.spot
    remaining = max(min(l["days_to_expiry"] for l in legs) - horizon_days, 1)
    sigma = surface_scenario.iv_at(spot_scenario, remaining)
    T = remaining / 365
    sd = sigma * math.sqrt(T)
    if sd <= 1e-6:
        return value_at(legs, spot_scenario, horizon_days, surface_scenario, r, contract_size) - entry_ref

    mu = math.log(spot_scenario) + (r - 0.5 * sigma ** 2) * T
    grid = np.geomspace(spot_scenario * 0.15, spot_scenario * 4, n)
    log_grid = np.log(grid)
    density = np.exp(-0.5 * ((log_grid - mu) / sd) ** 2) / (grid * sd * math.sqrt(2 * math.pi))
    pnl = np.array([value_at(legs, float(s), horizon_days, surface_scenario, r, contract_size) - entry_ref for s in grid])

    numerator = np.trapz(pnl * density, grid)
    denominator = np.trapz(density, grid)
    return float(numerator / denominator) if denominator > 1e-12 else float(pnl.mean())


def payoff_curves(
    legs: List[Dict[str, Any]],
    chain_slice: Dict[str, Any],
    surface_now: Surface,
    surface_scenario: ScenarioSurface,
    horizon_days: int,
    r: float = 0.05,
    n: int = 100,
    contract_size: float = DEFAULT_CONTRACT_SIZE,
) -> Dict[str, Any]:
    spot = chain_slice["spot"]
    priced = price_combo(legs, chain_slice, surface_now, surface_scenario, horizon_days, r, contract_size)
    entry_ref = priced["entry_cost"]

    lo, hi = spot * 0.6, spot * 1.4
    # A calendar/ratio spread's payoff can spike sharply right at a strike — a plain
    # uniform sweep over the full 0.6x-1.4x range (n points) can straddle right over that
    # peak without ever sampling it (same issue fixed in check_bounded_risk). Blend a
    # coarse baseline (overall shape) with a dense window around the legs' own strikes.
    # The exact spot/scenario spot are forced in as sample points too — otherwise hovering
    # right on the "Spot"/"Scénario" reference line reads the nearest grid point, which can
    # sit meaningfully off the true value on a peak this steep, and disagree with the
    # entry_cost/net_pnl tiles (computed at the exact spot, not a grid sample).
    baseline = np.linspace(lo, hi, n)
    strikes = [l["strike"] for l in legs]
    lo_k, hi_k = max(min(strikes) * 0.9, lo), min(max(strikes) * 1.1, hi)
    near_strikes = np.linspace(lo_k, hi_k, n * 3)
    exact_points = np.array([spot, surface_scenario.spot] + strikes)
    prices = np.unique(np.concatenate([baseline, near_strikes, exact_points]))
    prices.sort()
    eval_days_expiry = min(l["days_to_expiry"] for l in legs)

    # Both curves share the scenario's volatility view (surface_scenario) — only the
    # evaluation DATE differs: "à échéance" prices the day the near leg expires (matching
    # the Max gain/Max perte tile, itself computed with surface_scenario for the same
    # reason — see check_bounded_risk), "à J+8" prices your chosen scenario horizon. Using
    # surface_now here instead would silently mix in today's un-shocked vol, producing a
    # curve whose peak doesn't match the Max gain tile right next to it.
    at_expiry = [
        {"underlying": round(float(p), 4), "pnl": round(float(value_at(legs, float(p), eval_days_expiry, surface_scenario, r, contract_size) - entry_ref), 2)}
        for p in prices
    ]
    at_scenario = [
        {"underlying": round(float(p), 4), "pnl": round(float(value_at(legs, float(p), horizon_days, surface_scenario, r, contract_size) - entry_ref), 2)}
        for p in prices
    ]
    return {"at_expiry": at_expiry, "at_scenario": at_scenario, **priced}