jrkerns / pylinac

An image analysis library for medical physics
https://pylinac.readthedocs.io/en/latest/
MIT License
157 stars 99 forks source link

Protocol.ELEKTA symmetry calc (IEC point difference quotient) #490

Closed mitch-1211 closed 2 months ago

mitch-1211 commented 7 months ago

Describe the bug pylinac/field_analysis.py

elekta_protocol = {
    "symmetry": {"calc": symmetry_pdq_iec, "unit": "", "plot": plot_symmetry_pdq},
    "flatness": {"calc": flatness_dose_ratio, "unit": "", "plot": plot_flatness},
}
def plot_symmetry_pdq(instance, profile: SingleProfile, axis: plt.Axes) -> None:
    """Plotting of the symmetry point difference quotient."""

    def calc_sym(lt, rt, _) -> float:
        return max(abs(lt / rt), abs(rt / lt))

    _plot_sym_common(
        instance, calc_sym, profile, axis, label="Symmetry (AU)", padding=(0.05, 0.01)
    )
def symmetry_pdq_iec(profile: SingleProfile, in_field_ratio: float, **kwargs) -> float:
    """Symmetry calculation by way of PDQ IEC. See :ref:`elekta_protocol`.

    A negative value means the right side is higher. A positive value means the left side is higher.
    """
    field = profile.field_data(
        in_field_ratio=in_field_ratio,
        slope_exclusion_ratio=kwargs.get("slope_exclusion_ratio", 0.2),
    )
    field_values = field["field values"]

    def calc_sym(lt, rt) -> float:
        sym1 = lt / rt
        sym2 = rt / lt
        if abs(sym1) > abs(sym2):
            sign = np.sign(sym1)
        else:
            sign = np.sign(sym2)
        return max(abs(lt / rt), abs(rt / lt)) * sign

    sym_values = [calc_sym(lt, rt) for lt, rt in zip(field_values, field_values[::-1])]
    sym_idx = np.argmax(np.abs(sym_values))

    return sym_values[sym_idx]

The value returned is missing a factor of 100 i.e. should read return 100 * max(abs(lt / rt), abs(rt / lt)) * sign Also I don't think this value can ever be negative

To Reproduce Apply Protocol.ELEKTA to any field analysis. Values returned will be 1.009, 1.010 etc

Expected behavior Apply Protocol.ELEKTA to any field analysis. Values returned will be 100.9, 101.0 etc

jrkerns commented 7 months ago

For backwards compatibility I'll have to leave it the way it is and update the documentation. Thanks