defectpl.participation_ratio

High-level calculator

ParticipationRatioCalculator

ParticipationRatioCalculator(
    procar,
    defect_entry,
    defect_structure_info=None,
    poscar=None,
    cutoff_radius=3.5,
    use_pymatgen=True,
)

High-level calculator for electronic-state P-ratio and IPR.

Orchestrates file reading, neighbour resolution, the core computation, and serialisation to JSON / CSV.

Parameters:

Name	Type	Description	Default
procar	str or Path	Path to VASP PROCAR file.	required
defect_entry	str or Path	Path to defect_entry.json.	required
defect_structure_info	str or Path or None	Path to defect_structure_info.json. When supplied, neighbour indices are read directly from this file.	None
poscar	str or Path or None	Path to POSCAR/CONTCAR used as fallback for neighbour search.	None
cutoff_radius	float	Cutoff radius (Å) for the distance-based fallback neighbour search. Default: 3.5 Å.	3.5
use_pymatgen	bool	Use pymatgen's Procar parser when available. Default: True.	True

Examples:

>>> calc = ParticipationRatioCalculator(
...     procar="PROCAR",
...     defect_entry="defect_entry.json",
...     defect_structure_info="defect_structure_info.json",
... )
>>> result = calc.run()
>>> calc.to_json("participation_ratio.json")
>>> calc.to_csv("participation_ratio.csv")

Source code in defectpl\participation_ratio.py

def __init__(
    self,
    procar: str | Path,
    defect_entry: str | Path,
    defect_structure_info: Optional[str | Path] = None,
    poscar: Optional[str | Path] = None,
    cutoff_radius: float = 3.5,
    use_pymatgen: bool = True,
) -> None:
    self.procar_path = Path(procar)
    self.defect_entry_path = Path(defect_entry)
    self.dsi_path = Path(defect_structure_info) if defect_structure_info else None
    self.poscar_path = Path(poscar) if poscar else None
    self.cutoff_radius = cutoff_radius
    self.use_pymatgen = use_pymatgen
    self._result: Optional[dict] = None

result property

result

The computed result dictionary or None if :meth:run has not been called.

run

run()

Execute the full P-ratio / IPR calculation pipeline.

Returns:

Type	Description
dict	Nested result dictionary (spin → k-point → band → metrics). The same object is also stored as self.result.

Source code in defectpl\participation_ratio.py

def run(self) -> dict:
    """
    Execute the full P-ratio / IPR calculation pipeline.

    Returns
    -------
    dict
        Nested result dictionary (spin → k-point → band → metrics).
        The same object is also stored as ``self.result``.
    """
    # 1. Load defect metadata
    entry = self._load_defect_entry()
    defect_name = entry.get("name", entry.get("defect_name", "unknown"))
    defect_center = entry.get("defect_center", [0.5, 0.5, 0.5])

    # 2. Resolve neighbour atom indices
    poscar_path = self._find_poscar()
    neighbor_indices, neighbor_source = resolve_neighbors(
        dsi_path=self.dsi_path,
        poscar_path=poscar_path,
        defect_center_frac=defect_center,
        cutoff_radius=self.cutoff_radius,
    )
    logger.info("Neighbours resolved from: %s", neighbor_source)

    # 3. Parse PROCAR
    procar_data = read_procar(self.procar_path, use_pymatgen=self.use_pymatgen)
    logger.info(
        "PROCAR: %d spin(s), %d k-point(s), %d band(s), %d ions",
        procar_data["n_spins"],
        procar_data["n_kpoints"],
        procar_data["n_bands"],
        procar_data["n_ions"],
    )

    # 4. Compute metrics
    self._result = compute_participation_ratios(
        procar_data=procar_data,
        neighbor_indices=neighbor_indices,
        defect_name=defect_name,
        defect_center=defect_center,
    )
    self._result["neighbor_source"] = neighbor_source

    return self._result

to_json

to_json(path)

Serialise the result to a JSON file.

Parameters:

Name	Type	Description	Default
path	str or Path		required

Returns:

Type	Description
Path — the written file path.

Source code in defectpl\participation_ratio.py

def to_json(self, path: str | Path) -> Path:
    """
    Serialise the result to a JSON file.

    Parameters
    ----------
    path : str or Path

    Returns
    -------
    Path  — the written file path.
    """
    if self._result is None:
        raise RuntimeError("Call run() before to_json()")
    path = Path(path)
    path.parent.mkdir(parents=True, exist_ok=True)
    with open(path, "w") as fh:
        json.dump(self._result, fh, indent=2)
    logger.info("JSON written: %s", path)
    return path

to_csv

to_csv(path)

Write a flat CSV summary of all (spin, k-point, band) metrics.

Columns: spin, kpt, band, energy, occ, p_ratio, ipr, p_neighbors, p_total

Parameters:

Name	Type	Description	Default
path	str or Path		required

Returns:

Type	Description
Path

Source code in defectpl\participation_ratio.py

def to_csv(self, path: str | Path) -> Path:
    """
    Write a flat CSV summary of all (spin, k-point, band) metrics.

    Columns: spin, kpt, band, energy, occ, p_ratio, ipr,
             p_neighbors, p_total

    Parameters
    ----------
    path : str or Path

    Returns
    -------
    Path
    """
    import csv

    if self._result is None:
        raise RuntimeError("Call run() before to_csv()")

    path = Path(path)
    path.parent.mkdir(parents=True, exist_ok=True)

    fieldnames = [
        "spin",
        "kpt",
        "band",
        "energy",
        "occ",
        "p_ratio",
        "ipr",
        "p_neighbors",
        "p_total",
    ]
    rows = []
    for sp_label, kpt_dict in self._result["data"].items():
        spin_idx = int(sp_label.split("_")[1])
        for kpt_label, band_dict in kpt_dict.items():
            kpt_idx = int(kpt_label.split("_")[1])
            for band_label, vals in band_dict.items():
                band_idx = int(band_label.split("_")[1])
                rows.append(
                    {
                        "spin": spin_idx,
                        "kpt": kpt_idx,
                        "band": band_idx,
                        "energy": vals.get("energy"),
                        "occ": vals.get("occupation"),
                        "p_ratio": vals["p_ratio"],
                        "ipr": vals["ipr"],
                        "p_neighbors": vals["p_neighbors"],
                        "p_total": vals["p_total"],
                    }
                )

    with open(path, "w", newline="") as fh:
        writer = csv.DictWriter(fh, fieldnames=fieldnames)
        writer.writeheader()
        writer.writerows(rows)

    logger.info("CSV written: %s (%d rows)", path, len(rows))
    return path

top_localized

top_localized(n=10, metric='p_ratio')

Return the n most localised states sorted by metric.

Parameters:

Name	Type	Description	Default
n	int, optional (default 10)		10
metric	('p_ratio', 'ipr')		"p_ratio"

Returns:

Type	Description
list of dict, each with keys:	spin, kpt, band, energy, occ, p_ratio, ipr

Source code in defectpl\participation_ratio.py

def top_localized(self, n: int = 10, metric: str = "p_ratio") -> List[dict]:
    """
    Return the *n* most localised states sorted by ``metric``.

    Parameters
    ----------
    n : int, optional (default 10)
    metric : {"p_ratio", "ipr"}, optional

    Returns
    -------
    list of dict, each with keys:
        spin, kpt, band, energy, occ, p_ratio, ipr
    """
    if self._result is None:
        raise RuntimeError("Call run() first.")

    rows = []
    for sp_label, kpt_dict in self._result["data"].items():
        for kpt_label, band_dict in kpt_dict.items():
            kpt_idx = int(kpt_label.split("_")[1])
            for band_label, vals in band_dict.items():
                band_idx = int(band_label.split("_")[1])
                rows.append(
                    dict(
                        spin=sp_label,
                        kpt=kpt_idx,
                        band=band_idx,
                        energy=vals.get("energy"),
                        occ=vals.get("occupation"),
                        p_ratio=vals["p_ratio"],
                        ipr=vals["ipr"],
                    )
                )

    rows.sort(key=lambda x: -(x[metric] or 0))
    return rows[:n]

Plotting

plot_pr_vs_energy

plot_pr_vs_energy(
    result,
    metric="p_ratio",
    threshold=0.2,
    vbm=None,
    cbm=None,
    emin=None,
    emax=None,
    kpt_idx=0,
    title=None,
    out=None,
    ax=None,
    figsize=(8, 5),
    dpi=150,
)

Scatter plot of P-ratio or IPR versus Kohn-Sham energy.

Parameters:

Name	Type	Description	Default
result	dict	Nested PR result dict (from :func:compute_participation_ratios or loaded from participation_ratio.json).	required
metric	('p_ratio', 'ipr')	Quantity plotted on the Y-axis.	"p_ratio"
threshold	float	Horizontal dashed threshold reference line. Default 0.2.	0.2
vbm	float	Valence Band Maximum energy (eV). Shown as a vertical orange line.	None
cbm	float	Conduction Band Minimum energy (eV). Shown as a vertical green line.	None
emin	float	Lower energy boundary for filtering.	None
emax	float	Upper energy boundary for filtering.	None
kpt_idx	int	0-based k-point index to use. Default 0.	0
title	str	Plot title. Defaults to the defect name in result.	None
out	str or Path	If given, save the figure to this path.	None
ax	Axes	Inject into an existing Axes. When None a new figure is created.	None
figsize	tuple	Figure size (width, height) in inches.	(8, 5)
dpi	int	Resolution when saving.	150

Returns:

Type	Description
Axes

Source code in defectpl\participation_ratio.py

def plot_pr_vs_energy(
    result: dict,
    metric: str = "p_ratio",
    threshold: float = 0.2,
    vbm: Optional[float] = None,
    cbm: Optional[float] = None,
    emin: Optional[float] = None,
    emax: Optional[float] = None,
    kpt_idx: int = 0,
    title: Optional[str] = None,
    out: Optional[str | Path] = None,
    ax=None,
    figsize: tuple = (8, 5),
    dpi: int = 150,
):
    """
    Scatter plot of P-ratio or IPR versus Kohn-Sham energy.

    Parameters
    ----------
    result : dict
        Nested PR result dict (from :func:`compute_participation_ratios` or
        loaded from ``participation_ratio.json``).
    metric : {"p_ratio", "ipr"}
        Quantity plotted on the Y-axis.
    threshold : float
        Horizontal dashed threshold reference line.  Default 0.2.
    vbm : float, optional
        Valence Band Maximum energy (eV).  Shown as a vertical orange line.
    cbm : float, optional
        Conduction Band Minimum energy (eV).  Shown as a vertical green line.
    emin : float, optional
        Lower energy boundary for filtering.
    emax : float, optional
        Upper energy boundary for filtering.
    kpt_idx : int
        0-based k-point index to use.  Default 0.
    title : str, optional
        Plot title.  Defaults to the defect name in *result*.
    out : str or Path, optional
        If given, save the figure to this path.
    ax : matplotlib.axes.Axes, optional
        Inject into an existing Axes.  When None a new figure is created.
    figsize : tuple
        Figure size (width, height) in inches.
    dpi : int
        Resolution when saving.

    Returns
    -------
    matplotlib.axes.Axes
    """
    _, plt = _require_matplotlib()

    rows = flatten_pr_result(result, kpt_idx=kpt_idx)
    if emin is not None:
        rows = [r for r in rows if r["energy"] is not None and r["energy"] >= emin]
    if emax is not None:
        rows = [r for r in rows if r["energy"] is not None and r["energy"] <= emax]

    standalone = ax is None
    if standalone:
        fig, ax = plt.subplots(figsize=figsize)

    _scatter_spin_rows(ax, rows, metric=metric, x_key="energy")

    ax.axhline(
        threshold,
        color="gray",
        linestyle="--",
        linewidth=1,
        label=f"threshold = {threshold}",
    )

    if vbm is not None:
        ax.axvline(
            vbm,
            color="orange",
            linestyle=":",
            linewidth=1.2,
            label=f"VBM = {vbm:.3f} eV",
        )
    if cbm is not None:
        ax.axvline(
            cbm,
            color="green",
            linestyle=":",
            linewidth=1.2,
            label=f"CBM = {cbm:.3f} eV",
        )

    ax.set_xlabel("Energy (eV)")
    ax.set_ylabel(_METRIC_LABEL.get(metric, metric))
    ax.set_title(title or result.get("defect_name", "defect"))
    ax.legend(fontsize=8)
    ax.grid(True, alpha=0.3)

    if standalone and out is not None:
        ax.get_figure().tight_layout()
        ax.get_figure().savefig(str(out), dpi=dpi, bbox_inches="tight")
        plt.close(ax.get_figure())

    return ax

plot_pr_vs_band_index

plot_pr_vs_band_index(
    result,
    metric="p_ratio",
    threshold=0.2,
    emin=None,
    emax=None,
    kpt_idx=0,
    title=None,
    out=None,
    ax=None,
    figsize=(10, 5),
    dpi=150,
)

Scatter plot of P-ratio or IPR versus band index.

Band indices match the VASP band numbering (1-based). Separate spin channels are drawn in different colours; filled markers denote occupied states, open markers denote empty states.

Parameters:

Name	Type	Description	Default
result	dict	Nested PR result dict.	required
metric	('p_ratio', 'ipr')	Quantity plotted on the Y-axis.	"p_ratio"
threshold	float	Horizontal dashed threshold reference line. Default 0.2.	0.2
emin	float	Restrict displayed bands to those with energy ≥ emin (eV).	None
emax	float	Restrict displayed bands to those with energy ≤ emax (eV).	None
kpt_idx	int	0-based k-point index to use. Default 0.	0
title	str	Plot title. Defaults to the defect name in result.	None
out	str or Path	If given, save the figure to this path.	None
ax	Axes	Inject into an existing Axes.	None
figsize	tuple	Figure size in inches.	(10, 5)
dpi	int	Resolution when saving.	150

Returns:

Type	Description
Axes

Source code in defectpl\participation_ratio.py

def plot_pr_vs_band_index(
    result: dict,
    metric: str = "p_ratio",
    threshold: float = 0.2,
    emin: Optional[float] = None,
    emax: Optional[float] = None,
    kpt_idx: int = 0,
    title: Optional[str] = None,
    out: Optional[str | Path] = None,
    ax=None,
    figsize: tuple = (10, 5),
    dpi: int = 150,
):
    """
    Scatter plot of P-ratio or IPR versus band index.

    Band indices match the VASP band numbering (1-based).  Separate spin
    channels are drawn in different colours; filled markers denote occupied
    states, open markers denote empty states.

    Parameters
    ----------
    result : dict
        Nested PR result dict.
    metric : {"p_ratio", "ipr"}
        Quantity plotted on the Y-axis.
    threshold : float
        Horizontal dashed threshold reference line.  Default 0.2.
    emin : float, optional
        Restrict displayed bands to those with energy ≥ emin (eV).
    emax : float, optional
        Restrict displayed bands to those with energy ≤ emax (eV).
    kpt_idx : int
        0-based k-point index to use.  Default 0.
    title : str, optional
        Plot title.  Defaults to the defect name in *result*.
    out : str or Path, optional
        If given, save the figure to this path.
    ax : matplotlib.axes.Axes, optional
        Inject into an existing Axes.
    figsize : tuple
        Figure size in inches.
    dpi : int
        Resolution when saving.

    Returns
    -------
    matplotlib.axes.Axes
    """
    _, plt = _require_matplotlib()

    rows = flatten_pr_result(result, kpt_idx=kpt_idx)
    if emin is not None:
        rows = [r for r in rows if r["energy"] is not None and r["energy"] >= emin]
    if emax is not None:
        rows = [r for r in rows if r["energy"] is not None and r["energy"] <= emax]

    standalone = ax is None
    if standalone:
        fig, ax = plt.subplots(figsize=figsize)

    _scatter_spin_rows(ax, rows, metric=metric, x_key="band")

    ax.axhline(
        threshold,
        color="gray",
        linestyle="--",
        linewidth=1,
        label=f"threshold = {threshold}",
    )

    ax.set_xlabel("Band index")
    ax.set_ylabel(_METRIC_LABEL.get(metric, metric))
    ax.set_title(title or result.get("defect_name", "defect"))
    ax.legend(fontsize=8)
    ax.grid(True, alpha=0.3)

    if standalone and out is not None:
        ax.get_figure().tight_layout()
        ax.get_figure().savefig(str(out), dpi=dpi, bbox_inches="tight")
        plt.close(ax.get_figure())

    return ax

flatten_pr_result

flatten_pr_result(result, kpt_idx=0)

Flatten the nested PR result dictionary into a list of per-band row dicts.

Parameters:

Name	Type	Description	Default
result	dict	As returned by :func:compute_participation_ratios or loaded from participation_ratio.json.	required
kpt_idx	int	0-based k-point index to extract. Default 0 (Gamma point).	0

Returns:

Type	Description
list of dict	Each dict has keys: spin, band, energy, occ, p_ratio, ipr, p_neighbors, p_total.

Source code in defectpl\participation_ratio.py

def flatten_pr_result(result: dict, kpt_idx: int = 0) -> List[dict]:
    """
    Flatten the nested PR result dictionary into a list of per-band row dicts.

    Parameters
    ----------
    result : dict
        As returned by :func:`compute_participation_ratios` or loaded from
        ``participation_ratio.json``.
    kpt_idx : int, optional
        0-based k-point index to extract.  Default 0 (Gamma point).

    Returns
    -------
    list of dict
        Each dict has keys: ``spin``, ``band``, ``energy``, ``occ``,
        ``p_ratio``, ``ipr``, ``p_neighbors``, ``p_total``.
    """
    rows: List[dict] = []
    data = result.get("data", {})
    kpt_label = f"kpt_{kpt_idx + 1}"
    for sp_label, kpt_dict in data.items():
        if kpt_label not in kpt_dict:
            continue
        for band_label, vals in kpt_dict[kpt_label].items():
            band_idx = int(band_label.split("_")[1])
            rows.append(
                {
                    "spin": sp_label,
                    "band": band_idx,
                    "energy": vals.get("energy"),
                    "occ": vals.get("occupation", 0.0) or 0.0,
                    "p_ratio": vals.get("p_ratio", 0.0) or 0.0,
                    "ipr": vals.get("ipr", 0.0) or 0.0,
                    "p_neighbors": vals.get("p_neighbors", 0.0) or 0.0,
                    "p_total": vals.get("p_total", 0.0) or 0.0,
                }
            )
    return rows

Low-level functions

read_procar

read_procar(procar_path, use_pymatgen=True)

Parse a VASP PROCAR file and return projection arrays.

Tries pymatgen's Procar class first (if available and use_pymatgen=True); falls back to the built-in native parser.

Parameters:

Name	Type	Description	Default
procar_path	str or Path	Path to the PROCAR file.	required
use_pymatgen	bool	Attempt to use pymatgen's parser first. Set to False to force the native parser. Default: True.	True

Returns:

Type	Description
dict with keys:
``n_kpoints`` : int
``n_bands`` : int
``n_ions`` : int
``n_spins`` : int
``spins`` : list of spin sentinels
``site_proj`` : dict[spin, np.ndarray shape (nk, nb, ni)]	Site-projected weight summed over all orbital channels.
``eigenvalues`` : dict[spin, np.ndarray shape (nk, nb)]
``occupancies`` : dict[spin, np.ndarray shape (nk, nb)]
``kpoints`` : np.ndarray shape (nk, 3)
``weights`` : np.ndarray shape (nk,)

Raises:

Type	Description
FileNotFoundError	If the PROCAR file does not exist.

Source code in defectpl\participation_ratio.py

def read_procar(procar_path: str | Path, use_pymatgen: bool = True) -> dict:
    """
    Parse a VASP PROCAR file and return projection arrays.

    Tries pymatgen's ``Procar`` class first (if available and
    ``use_pymatgen=True``); falls back to the built-in native parser.

    Parameters
    ----------
    procar_path : str or Path
        Path to the PROCAR file.
    use_pymatgen : bool, optional
        Attempt to use pymatgen's parser first.  Set to ``False`` to force
        the native parser.  Default: ``True``.

    Returns
    -------
    dict with keys:

    ``n_kpoints`` : int
    ``n_bands``   : int
    ``n_ions``    : int
    ``n_spins``   : int
    ``spins``     : list of spin sentinels
    ``site_proj`` : dict[spin, np.ndarray shape (nk, nb, ni)]
        Site-projected weight summed over all orbital channels.
    ``eigenvalues`` : dict[spin, np.ndarray shape (nk, nb)]
    ``occupancies`` : dict[spin, np.ndarray shape (nk, nb)]
    ``kpoints``   : np.ndarray shape (nk, 3)
    ``weights``   : np.ndarray shape (nk,)

    Raises
    ------
    FileNotFoundError
        If the PROCAR file does not exist.
    """
    procar_path = Path(procar_path)
    if not procar_path.exists():
        raise FileNotFoundError(f"PROCAR not found: {procar_path}")

    if use_pymatgen:
        try:
            from pymatgen.io.vasp.outputs import Procar

            logger.info("Parsing PROCAR with pymatgen: %s", procar_path)
            p = Procar(str(procar_path))

            spins = list(p.data.keys())
            data_arr = p.data[spins[0]]  # (nk, nb, ni, n_orb)
            n_kpoints, n_bands, n_ions, _ = data_arr.shape

            site_proj = {spin: arr.sum(axis=-1) for spin, arr in p.data.items()}

            return {
                "n_kpoints": n_kpoints,
                "n_bands": n_bands,
                "n_ions": n_ions,
                "n_spins": len(spins),
                "spins": spins,
                "site_proj": site_proj,
                "eigenvalues": getattr(p, "eigenvalues", None),
                "occupancies": getattr(p, "occupancies", None),
                "kpoints": np.array(p.kpoints)
                if getattr(p, "kpoints", None) is not None
                else None,
                "weights": np.array(p.weights)
                if getattr(p, "weights", None) is not None
                else None,
            }
        except Exception as exc:
            logger.warning(
                "pymatgen Procar failed (%s); falling back to native parser.", exc
            )

    logger.info("Parsing PROCAR natively: %s", procar_path)
    return _parse_procar_native(procar_path)

compute_participation_ratios

compute_participation_ratios(
    procar_data,
    neighbor_indices,
    defect_name="unknown",
    defect_center=(0.5, 0.5, 0.5),
)

Compute P-ratio and IPR for every (spin, k-point, band) triplet.

Parameters:

Name	Type	Description	Default
procar_data	dict	As returned by :func:read_procar.	required
neighbor_indices	sequence of int	0-based atom indices that constitute the defect neighbourhood. Used in the P-ratio numerator.	required
defect_name	str		'unknown'
defect_center	sequence of 3 floats	Fractional coordinates of the defect centre (stored in output only).	(0.5, 0.5, 0.5)

Returns:

Type	Description
dict	Serialisable result dictionary (spin → k-point → band → metrics).

Source code in defectpl\participation_ratio.py

def compute_participation_ratios(
    procar_data: dict,
    neighbor_indices: Sequence[int],
    defect_name: str = "unknown",
    defect_center: Sequence[float] = (0.5, 0.5, 0.5),
) -> dict:
    """
    Compute P-ratio and IPR for every (spin, k-point, band) triplet.

    Parameters
    ----------
    procar_data : dict
        As returned by :func:`read_procar`.
    neighbor_indices : sequence of int
        0-based atom indices that constitute the defect neighbourhood.
        Used in the P-ratio numerator.
    defect_name : str, optional
    defect_center : sequence of 3 floats, optional
        Fractional coordinates of the defect centre (stored in output only).

    Returns
    -------
    dict
        Serialisable result dictionary (spin → k-point → band → metrics).
    """
    n_kpoints = procar_data["n_kpoints"]
    n_bands = procar_data["n_bands"]
    n_ions = procar_data["n_ions"]
    spins = procar_data["spins"]
    site_proj = procar_data["site_proj"]  # {spin: (nk, nb, ni)}
    eigenvalues = procar_data.get("eigenvalues")
    occupancies = procar_data.get("occupancies")

    neighbor_arr = np.array(sorted(set(int(i) for i in neighbor_indices)), dtype=int)

    # Build a human-readable spin label mapping
    def _spin_label(spin, idx: int) -> str:
        # Works with pymatgen Spin objects and our _Spin sentinels
        try:
            from pymatgen.electronic_structure.core import Spin

            if spin == Spin.up:
                return "spin_1"
            if spin == Spin.down:
                return "spin_2"
        except ImportError:
            pass
        # Fallback: use the index
        return f"spin_{idx + 1}"

    data_out: PRatioData = {}

    for spin_idx, spin in enumerate(spins):
        sp_label = _spin_label(spin, spin_idx)
        proj = site_proj[spin]  # (nk, nb, ni)
        data_out[sp_label] = {}

        for ik in range(n_kpoints):
            kpt_label = f"kpt_{ik + 1}"
            data_out[sp_label][kpt_label] = {}

            for ib in range(n_bands):
                band_label = f"band_{ib + 1}"

                p_i = proj[ik, ib, :]  # (n_ions,)
                p_total = float(p_i.sum())

                if p_total > 1e-12:
                    p_neighbors = (
                        float(p_i[neighbor_arr].sum()) if neighbor_arr.size else 0.0
                    )
                    p_ratio = p_neighbors / p_total
                    ipr = float((p_i**2).sum()) / (p_total**2)
                else:
                    p_neighbors = 0.0
                    p_ratio = 0.0
                    ipr = 0.0

                energy = occ = None
                if eigenvalues is not None and spin in eigenvalues:
                    energy = float(eigenvalues[spin][ik, ib])
                if occupancies is not None and spin in occupancies:
                    occ = float(occupancies[spin][ik, ib])

                data_out[sp_label][kpt_label][band_label] = {
                    "energy": energy,
                    "occupation": occ,
                    "p_ratio": round(p_ratio, 8),
                    "ipr": round(ipr, 8),
                    "p_neighbors": round(p_neighbors, 8),
                    "p_total": round(p_total, 8),
                }

    return {
        "defect_name": defect_name,
        "defect_center": list(defect_center),
        "neighbor_atom_indices": neighbor_arr.tolist(),
        "n_atoms": n_ions,
        "n_spins": len(spins),
        "n_kpoints": n_kpoints,
        "n_bands": n_bands,
        "data": data_out,
    }

resolve_neighbors

resolve_neighbors(
    dsi_path,
    poscar_path,
    defect_center_frac,
    cutoff_radius=3.5,
)

Resolve neighbour atom indices using the best available source.

Priority: 1. defect_structure_info.json (authoritative source) 2. Distance-based search on POSCAR/CONTCAR 3. Empty list (no neighbour info; P-ratio will be 0 for all states)

Parameters:

Name	Type	Description	Default
dsi_path	path or None		required
poscar_path	path or None		required
defect_center_frac	sequence of 3 floats		required
cutoff_radius	float(Å)		3.5

Returns:

Type	Description
(indices, source_description)

Source code in defectpl\participation_ratio.py

def resolve_neighbors(
    dsi_path: Optional[str | Path],
    poscar_path: Optional[str | Path],
    defect_center_frac: Sequence[float],
    cutoff_radius: float = 3.5,
) -> Tuple[List[int], str]:
    """
    Resolve neighbour atom indices using the best available source.

    Priority:
    1. ``defect_structure_info.json`` (authoritative source)
    2. Distance-based search on POSCAR/CONTCAR
    3. Empty list (no neighbour info; P-ratio will be 0 for all states)

    Parameters
    ----------
    dsi_path : path or None
    poscar_path : path or None
    defect_center_frac : sequence of 3 floats
    cutoff_radius : float (Å)

    Returns
    -------
    (indices, source_description)
    """
    if dsi_path is not None:
        indices = neighbors_from_defect_structure_info(dsi_path)
        if indices is not None:
            return indices, f"defect_structure_info.json ({len(indices)} atoms)"

    if poscar_path is not None and Path(poscar_path).exists():
        try:
            indices = neighbors_from_structure(
                poscar_path, defect_center_frac, cutoff_radius
            )
            return (
                indices,
                f"distance search r<{cutoff_radius:.1f}Å ({len(indices)} atoms)",
            )
        except Exception as exc:
            logger.warning("Distance-based neighbour search failed: %s", exc)

    logger.warning("No neighbour information found; P-ratio will be 0 for all states.")
    return [], "none (P-ratio will be 0)"

neighbors_from_defect_structure_info

neighbors_from_defect_structure_info(dsi_path)

Extract neighbour atom indices from a defect_structure_info.json file.

The function tries all key names used across known pydefect versions.

Parameters:

Name	Type	Description	Default
dsi_path	str or Path		required

Returns:

Type	Description
list of int (0-based POSCAR order) or None if not found / parse error.

Source code in defectpl\participation_ratio.py

def neighbors_from_defect_structure_info(dsi_path: str | Path) -> Optional[List[int]]:
    """
    Extract neighbour atom indices from a ``defect_structure_info.json`` file.

    The function tries all key names used across known pydefect versions.

    Parameters
    ----------
    dsi_path : str or Path

    Returns
    -------
    list of int (0-based POSCAR order) or None if not found / parse error.
    """
    dsi_path = Path(dsi_path)
    if not dsi_path.exists():
        return None
    try:
        with open(dsi_path) as fh:
            dsi = json.load(fh)
    except (json.JSONDecodeError, OSError) as exc:
        logger.warning("Could not read %s: %s", dsi_path, exc)
        return None

    for key in (
        "neighbor_atom_indices",
        "neighboring_atom_indices",
        "neighboring_sites",
        "neighbor_sites",
    ):
        if key in dsi:
            val = dsi[key]
            if val and isinstance(val[0], dict):
                return [int(v["index"]) for v in val]
            return [int(v) for v in val]

    # Fallback: look for a 'neighbors' list of dicts that carry an 'index'
    if "neighbors" in dsi:
        entries = dsi["neighbors"]
        if entries and isinstance(entries[0], dict) and "index" in entries[0]:
            return [int(e["index"]) for e in entries]

    return None