Seismic vulnerability modelling
================================

This example runs the full MDOF vulnerability workflow: backbone
fitting, EDP-IM prediction, storey drift and floor acceleration demand
estimation, and final loss computation via ``VMMDOF``.

Input files
-----------

* `slf.json <https://github.com/djura-risk-data-engineering/djura/blob/main/tests/vm/assets/vm/slf.json>`_
  — storey loss functions per storey and component group, one set per
  building. See the :doc:`SLF examples <slf>` for how to generate these.

All other inputs (backbone parameters, hazard curve coefficients) are
defined inline in the script below.

Running the example
-------------------

.. code-block:: python

    import asyncio
    import copy
    import json

    from djura.edp_im.predict import edp_im_batch
    from djura.vulnerability_modeller.backbone import Backbone, get_infill_spo
    from djura.vulnerability_modeller.demands import Demands
    from djura.vulnerability_modeller.vm_mdof import VMMDOF
    from djura.vulnerability_modeller.utilities import deep_merge

    # --- 1. EDP-IM model parameters (one entry per building) ----------------

    data_edp_im = [
        {
            "hysteresis": "bilin",
            "im_type": "sa",
            "method": "shahnazaryan-oreilly",
            "backboneMethod": "backbone",
            "backbone": {
                "damping": 0.02,
                "ductility": 3,
                "ductility_f": 6,
                "hardening_ratio": 0.05,
                "period": 1.0,
            },
        },
        {
            "hysteresis": "bilin",
            "im_type": "sa",
            "method": "shahnazaryan-oreilly",
            "backboneMethod": "backbone",
            "backbone": {
                "damping": 0.03,
                "ductility": 3,
                "ductility_f": 6,
                "hardening_ratio": 0.02,
                "period": 0.5,
            },
        },
        {
            "hysteresis": "infill",
            "im_type": "sa",
            "backboneMethod": "backbone",
            "backbone": {
                "period": 1.0,
                "c_y": 3.5,
                "c_rp": 2.0,
                "mu_h": 3.2,
                "mu_s": 3.2,
                "mu_rp": 3.5,
                "mu_ult": 10.0,
            },
        },
    ]

    # --- 2. Building geometry and backbone parameters -----------------------

    data_est = [
        {
            "backbone": {"heights": [3.5, 3.0, 3.0, 3.0, 3.0],
                         "dy": 0.03, "ductility_f": 6},
            "building-type": "frame",
        },
        {
            "backbone": {"heights": [4.5, 3.5, 3.5],
                         "dy": 0.035, "ductility_f": 6},
            "building-type": "frame",
        },
        {
            "backbone": {"heights": [3.5, 3.0, 3.0, 3.0, 3.0],
                         "dy": 0.02, "ductility_f": 6},
            "building-type": "infill",
        },
    ]

    # --- 3. Site hazard curve coefficients (power-law fit) ------------------

    hazard = {
        "pga": [0.0001144, 2.6348751, 0.2340722],
        "im":  [0.000298042, 2.551209397, 0.242236676],
    }

    # --- 4. Predict EDP-IM relationships ------------------------------------

    edp_ims = edp_im_batch(data_edp_im)

    # --- 5. Estimate storey drifts and floor accelerations ------------------

    psds = asyncio.run(Demands(data_est).estimate_drifts())

    data_acc = [deep_merge(d1, d2) for d1, d2 in zip(data_est, data_edp_im)]
    for i, d in enumerate(data_acc):
        d["pfa-model"] = "muho"
        d["hazard"] = hazard
        d["psd"] = (
            next((r for r in psds["successes"] if r["index"] == i), None
                 )["out"]["psd"]
            if i in psds["success_ids"] else []
        )
        d["sr"] = (
            next((r for r in edp_ims["successes"] if r["index"] == i), None
                 )["out"]["strength-ratios"]
            if i in edp_ims["success_ids"] else []
        )

    pfas = asyncio.run(Demands(data_acc).estimate_accelerations())

    # --- 6. Load storey loss functions and assign storeys -------------------

    with open("slf.json", encoding="utf-8") as f:
        base_slfs = json.load(f)

    rc = 400_000   # replacement cost in the same currency as the SLFs

    n = len(edp_ims["failure_ids"]) + len(edp_ims["success_ids"])
    slfs = [copy.deepcopy(base_slfs) for _ in range(n)]

    # Building 1: 5 storeys
    slfs[0]["1 - NS: PFA"]["Storey"] = [0, 1, 2, 3, 4, 5]
    slfs[0]["3 - NS: PSD"]["Storey"] = [1, 2, 3, 4, 5]
    slfs[0]["2 - S: PSD"]["Storey"]  = [1, 2, 3, 4, 5]

    # Building 2: 3 storeys
    slfs[1]["1 - NS: PFA"]["Storey"] = [0, 1, 2, 3]
    slfs[1]["3 - NS: PSD"]["Storey"] = [1, 2, 3]
    slfs[1]["2 - S: PSD"]["Storey"]  = [1, 2, 3]

    # Building 3 (infill): expected to fail — storey assignment not needed

    for i, slf in enumerate(slfs):
        slf["rc"] = rc
        slf["cases"] = [i]

    # --- 7. Run vulnerability assessment ------------------------------------

    data = {
        "slfs": slfs,
        "edp-ims": edp_ims,
        "psds": psds,
        "pfas": pfas,
    }

    result = asyncio.run(VMMDOF([data]).vulnerability())

    print("Status :", result["status"])
    print("Results:", result)

Output
------

``result`` is a dictionary with a ``"status"`` key (``"success"`` or
``"failure"``) and per-building loss metrics. Buildings for which the
EDP-IM or demand estimation step failed are reported in
``"failure_ids"`` and excluded from the loss calculation.

Workflow summary
----------------

The five steps mirror the djura vulnerability assessment pipeline:

1. **EDP-IM** — ``edp_im_batch`` fits strength-ratio curves to the
   backbone model for each building.
2. **Storey drifts** — ``Demands.estimate_drifts`` converts the backbone
   into peak storey drift profiles at each IM level.
3. **Floor accelerations** — ``Demands.estimate_accelerations`` derives
   PFA profiles using the drift results and the selected PFA model.
4. **SLFs** — pre-generated storey loss functions (see :doc:`slf`) map
   EDP demands to repair costs per storey.
5. **VMMDOF** — assembles the above into expected annual loss and
   vulnerability curves for each building.