# -*- coding: utf-8 -*-
"""Median Center Processing Algorithm."""
from __future__ import annotations

import logging
import numpy as np

from qgis.PyQt.QtCore import QVariant
from ._mixins import HelpUrlMixin
from qgis.core import (
    NULL,
    QgsFeature,
    QgsField,
    QgsFields,
    QgsPointXY,
    QgsGeometry,
    QgsWkbTypes,
    QgsProcessing,
    QgsProcessingAlgorithm,
    QgsProcessingException,
    QgsProcessingParameterFeatureSource,
    QgsProcessingParameterField,
    QgsProcessingParameterFeatureSink
)

from ..core.stats_engines import calculate_median_center
from ..core.layer_metadata import apply_output_metadata
from ..core.weights import geometry_centroid_point

from ._icons import algorithm_icon


logger = logging.getLogger("PlanX GeoStats Lab")


class MedianCenterAlgorithm(HelpUrlMixin, QgsProcessingAlgorithm):
    INPUT = "INPUT"
    WEIGHT_FIELD = "WEIGHT_FIELD"
    OUTPUT = "OUTPUT"

    def __init__(self):
        super().__init__()
        self.out_layer_id = None

    def name(self) -> str:
        return "median_center"

    def displayName(self) -> str:
        return "Median Center"

    def group(self) -> str:
        return "04 | Centers, Direction and Dispersion"

    def groupId(self) -> str:
        return "planx_center_direction_spread"

    def icon(self):
        return algorithm_icon("median_center")

    def createInstance(self):
        return MedianCenterAlgorithm()

    def shortHelpString(self) -> str:
        return (
            "Computes the geometric median (also called the 1-median, "
            "Fermat-Weber point, or spatial median): the point in continuous "
            "space - not necessarily any input feature - that minimizes the "
            "sum of (optionally weighted) Euclidean distances to all features. "
            "Solved via Weiszfeld's iterative algorithm, initialized from the "
            "weighted mean center, with a fixed convergence tolerance of 1e-6 "
            "and a maximum of 100 iterations; typical well-conditioned datasets "
            "converge in 15-40 iterations. If an iterate ever lands exactly on "
            "a data point (a division-by-zero risk), the implementation "
            "perturbs it slightly and continues.\n\n"
            "Output fields: median_x, median_y, and total_dist (the minimized "
            "objective value - the sum of weighted distances from the median "
            "center to every feature). Unlike the mean center, the geometric "
            "median has a breakdown point of 0.5: up to half the data can be "
            "arbitrarily displaced without destroying the estimate, versus a "
            "breakdown point of 0 for the mean.\n\n"
            "Use this whenever the dataset may contain outliers or extreme "
            "values that would distort the mean center, or when solving a "
            "genuine minisum facility-location problem (minimize total travel "
            "distance) with no constraint that the answer be an existing site. "
            "For very small datasets (fewer than 5 features), the median "
            "center rarely differs meaningfully from the mean center; when the "
            "result should correspond to an actual candidate site, compare "
            "against Central Feature instead."
        )

    def initAlgorithm(self, config=None):
        self.addParameter(
            QgsProcessingParameterFeatureSource(
                self.INPUT,
                "Input vector layer",
                [QgsProcessing.TypeVectorAnyGeometry]
            )
        )
        self.addParameter(
            QgsProcessingParameterField(
                self.WEIGHT_FIELD,
                "Weight field (optional)",
                parentLayerParameterName=self.INPUT,
                type=QgsProcessingParameterField.Numeric,
                optional=True
            )
        )
        self.addParameter(
            QgsProcessingParameterFeatureSink(
                self.OUTPUT,
                "Output median center layer"
            )
        )

    def processAlgorithm(self, parameters, context, feedback):
        source = self.parameterAsSource(parameters, self.INPUT, context)
        if source is None:
            raise QgsProcessingException("Invalid input layer source.")

        weight_field = self.parameterAsString(parameters, self.WEIGHT_FIELD, context)

        # Extract features
        x_coords = []
        y_coords = []
        weights = []

        has_weight = weight_field != ""
        if has_weight:
            field_idx = source.fields().lookupField(weight_field)
            if field_idx < 0:
                raise QgsProcessingException(f"Weight field '{weight_field}' not found.")

        feedback.pushInfo("Extracting coordinates...")
        total = source.featureCount() or 1
        for idx, f in enumerate(source.getFeatures()):
            if feedback.isCanceled():
                break

            geom = f.geometry()
            if geom is None or geom.isEmpty():
                continue

            centroid = geometry_centroid_point(geom)
            if centroid is None:
                continue
            x_coords.append(centroid.x())
            y_coords.append(centroid.y())

            w_val = 1.0
            if has_weight:
                val = f.attribute(weight_field)
                if val is not None and val != NULL and str(val) != 'NULL':
                    try:
                        w_val = float(val)
                    except (ValueError, TypeError):
                        pass
            weights.append(w_val)
            feedback.setProgress(int(40 * (idx / total)))

        if len(x_coords) == 0:
            raise QgsProcessingException("No features with valid geometries were found.")

        x_arr = np.array(x_coords)
        y_arr = np.array(y_coords)
        w_arr = np.array(weights)

        feedback.pushInfo("Iterating using Weiszfeld's algorithm...")
        med_x, med_y, total_dist = calculate_median_center(x_arr, y_arr, w_arr)

        # Output fields
        out_fields = QgsFields()
        out_fields.append(QgsField("median_x", QVariant.Double, len=15, prec=6))
        out_fields.append(QgsField("median_y", QVariant.Double, len=15, prec=6))
        out_fields.append(QgsField("total_dist", QVariant.Double, len=15, prec=6))

        # Setup sink
        (sink, dest_id) = self.parameterAsSink(
            parameters,
            self.OUTPUT,
            context,
            out_fields,
            QgsWkbTypes.Point,
            source.sourceCrs()
        )
        self.out_layer_id = dest_id

        out_feat = QgsFeature()
        out_feat.setGeometry(QgsGeometry.fromPointXY(QgsPointXY(med_x, med_y)))
        out_feat.setFields(out_fields)
        out_feat.setAttribute("median_x", med_x)
        out_feat.setAttribute("median_y", med_y)
        out_feat.setAttribute("total_dist", total_dist)

        sink.addFeature(out_feat)
        feedback.setProgress(100)

        return {self.OUTPUT: dest_id}

    def postProcessAlgorithm(self, context, feedback):
        if self.out_layer_id is None:
            return {}
        layer = context.getMapLayer(self.out_layer_id)
        if not layer:
            return {}
        apply_output_metadata(
            layer,
            "PlanX GeoStats median center output",
            {
                "median_x": "Median center X coordinate",
                "median_y": "Median center Y coordinate",
                "total_dist": "Total weighted distance minimized by the median center",
            },
            self.displayName(),
        )
        return {}
