using Microsoft.Extensions.Logging;
using Microsoft.ML;
using Microsoft.ML.Data;

using OpenCvSharp;

using GSS2.Core.Extensions;
using GSS2.Core.Analysis.AmineContent.Database;

using MathNet.Numerics.LinearAlgebra;
using MathNet.Numerics.LinearAlgebra.Double;

namespace GSS2.Core.Analysis.AmineContent;

public partial class AnalyzerService
{
    //TODO: Перенести все гиперпараметры в appsettings.json
    public const int CHANNELS_COUNT = 3;
    public const int IMAGES_COUNT = 12;
    public const int FEATURES_LENGTH = (CHANNELS_COUNT + 3) * IMAGES_COUNT;
    public const int EFFECTIVE_IMAGE_SIZE = 200;
    public const int ANALYSIS_IMAGE_SIZE = 1000;
    public const int SEPARATOR_ROI_THRESHOLD = 5;

    public const float CONTOUR_MAX_SIZE = 100;
    public const float CONTOUR_MAX_AREA = (float)Math.PI * (CONTOUR_MAX_SIZE * CONTOUR_MAX_SIZE) / 5.0f;
    public const float CONTOUR_MIN_SIZE = 10;
    public const float CONTOUR_MIN_AREA = (float)Math.PI * (CONTOUR_MIN_SIZE * CONTOUR_MIN_SIZE) / 5.0f / 5.0f;

    public const int SEPARATOR_X = 519;
    public const int SEPARATOR_Y = 17;
    public const int SEPARATOR_WIDTH = 2969;
    public const int SEPARATOR_HEIGHT = 2965;

    public const double FILTERING_THRESHOLD = 5.0;

    private class ClusterizationData
    {
        [ColumnName("Label")]
        public required string Label { get; set; }
        [ColumnName("Features")]
        [VectorType(FEATURES_LENGTH)]
        public required float[] Features { get; set; }
    }
    private class ClusterizationPrediction
    {
        [ColumnName("PredictedLabel")]
        public string Label { get; set; }
        public float Probability { get; set; }
        public float[] Score { get; set; }

        public ClusterizationPrediction()
        {
            Label = "";
            Probability = 0;
            Score = new float[2];
        }
    }

    private class RegressionData
    {
        [ColumnName("Label")]
        public required float Value { get; set; }
        [ColumnName("Features")]
        [VectorType(FEATURES_LENGTH)]
        public required float[] Features { get; set; }
    }
    private class RegressionPrediction
    {
        [ColumnName("Score")]
        public float Value { get; set; }

        public RegressionPrediction()
        {
            Value = -1;
        }
    }

    private readonly ILogger<AnalyzerService> _logger;
    private readonly Context _context;

    private MLContext? _ml = null;
    private PredictionEngine<ClusterizationData, ClusterizationPrediction>? _clusterizationEngine = null;
    private PredictionEngine<RegressionData, RegressionPrediction>? _regressionEngine = null;

    public bool Initialized { get; private set; } = false;
    public int InitializedBrandId { get; private set; } = 0;
    public int InitializedSeparatorId { get; private set; } = 0;

    public AnalyzerService(ILogger<AnalyzerService> logger, Context Context)
    {
        _logger = logger;
        _context = Context;
    }

    public async Task Initialize(BrandRecord brand, SeparatorRecord separator, CancellationToken cancellationToken = default)
    {
        // TODO: Повторная инициализация если изменились калибровочные данные

        if (Initialized &&
            InitializedBrandId == brand.Id &&
            InitializedSeparatorId == separator.Id
        )
            return;

        Initialized = false;

        _logger.LogInformation("Инициализация");

        _ml = new MLContext();
        _ml.Log += (_, ea) =>
        {
            if (ea.Kind >= Microsoft.ML.Runtime.ChannelMessageKind.Info)
                _logger.LogInformation("ML: [{}] {}", ea.Kind, ea.Message);
        };

        try
        {
            await InitializeMlClusterizationEngine(brand, separator, cancellationToken);
            await InitializeMlRegressionEngine(brand, cancellationToken);
        }
        catch (Exception ex)
        {
            _logger.LogError(ex, "Исключение во время инициализации моделей");
            return;
        }

        Initialized = true;
        InitializedBrandId = brand.Id;
        InitializedSeparatorId = separator.Id;
    }
    private async Task InitializeMlClusterizationEngine(BrandRecord brand, SeparatorRecord separator, CancellationToken cancellationToken)
    {
        _logger.LogInformation("Инициализация модели кластеризации для марки: \"{} - {}\" и сепаратора: \"{}\"", brand.BrandName, brand.MixtureName, separator.Type);

        if (_ml is null)
            throw new InvalidOperationException("Анализатор не инициализирован");
        cancellationToken.ThrowIfCancellationRequested();

        _logger.LogInformation("Подготовка данных");
        var sampleRecords = _context.SampleRecords
            .Where(r => r.BrandId == brand.Id);

        _logger.LogInformation("Найдено калибровочных записей проб для марки: {}", sampleRecords.Count());
        if (sampleRecords.Count() == 0)
            throw new InvalidDataException($"Не найдено калибровочных записей для марки: \"{brand.BrandName} - {brand.MixtureName}\"");
        cancellationToken.ThrowIfCancellationRequested();

        if (separator.Features.Count() == 0)
            throw new InvalidDataException($"Не найдено калибровочных признаков для сепаратора: \"{separator.Type}\"");
        if (separator.Features.Count() % FEATURES_LENGTH != 0)
            throw new InvalidDataException($"Неожидаемое количество калибровочных признаков для сепаратора: \"{separator.Type}\"");
        var separatorDataSet = separator.Features.Chunk(FEATURES_LENGTH)
            .Select(v => new ClusterizationData
            {
                Label = "separator",
                Features = v.ToArray()
            })
            .Shuffle();

        _logger.LogInformation("Загружено векторов для сепаратора: {}", separatorDataSet.Count());
        if (separatorDataSet.Count() == 0)
            throw new InvalidDataException($"Не найдено калибровочных данных для сепаратора: \"{separator.Type}\"");
        cancellationToken.ThrowIfCancellationRequested();

        foreach (var sampleRecord in sampleRecords)
        {
            if (sampleRecord.Features.Count() == 0)
                _logger.LogError("Не найдено калибровочных признаков для пробы: \"{}\"", sampleRecord.SampleName);
            if (sampleRecord.Features.Count() % FEATURES_LENGTH != 0)
                _logger.LogError("Неожидаемое количество калибровочных признаков для пробы: \"{}\"", sampleRecord.SampleName);
        }

        var sampleDataSet = sampleRecords
            .Where(r => r.Features.Count() > 0 && r.Features.Count() % FEATURES_LENGTH == 0)
            .SelectMany(r => r.Features.Chunk(FEATURES_LENGTH))
            .Select(f => new ClusterizationData
            {
                Label = "sample",
                Features = f.ToArray()
            })
            .Shuffle();

        _logger.LogInformation("Загружено векторов для пробы: {}", sampleDataSet.Count());
        if (separatorDataSet.Count() == 0)
            throw new InvalidDataException($"Не найдено калибровочных данных для марки: \"{brand.BrandName} - {brand.MixtureName}\"");
        cancellationToken.ThrowIfCancellationRequested();

        var trainDataSet = Enumerable.Concat(separatorDataSet, sampleDataSet);
        _logger.LogInformation("Загружено векторов для обучения: {}", trainDataSet.Count());
        var trainData = _ml.Data.LoadFromEnumerable(trainDataSet);
        cancellationToken.ThrowIfCancellationRequested();

        _logger.LogInformation("Создание конвейера");
        var pipeline = _ml.Transforms.Conversion.MapValueToKey("Label")
            .Append(_ml.MulticlassClassification.Trainers.LbfgsMaximumEntropy())
            .Append(_ml.Transforms.Conversion.MapKeyToValue("PredictedLabel"));
        cancellationToken.ThrowIfCancellationRequested();

        _logger.LogInformation("Обучение модели");
        var model = pipeline.Fit(trainData);
        cancellationToken.ThrowIfCancellationRequested();

        _clusterizationEngine = _ml.Model.CreatePredictionEngine<ClusterizationData, ClusterizationPrediction>(model);
        model.Dispose();
    }
    private async Task InitializeMlRegressionEngine(BrandRecord brand, CancellationToken cancellationToken)
    {
        _logger.LogInformation("Инициализация модели регрессии для марки \"{} - {}\"", brand.BrandName, brand.MixtureName);

        if (_ml is null)
            throw new InvalidOperationException("Анализатор не инициализирован");
        cancellationToken.ThrowIfCancellationRequested();

        _logger.LogInformation("Подготовка данных");
        var sampleRecords = _context.SampleRecords
            .Where(r => r.BrandId == brand.Id);

        _logger.LogInformation("Найдено калибровочных записей проб для марки: {}", sampleRecords.Count());
        if (sampleRecords.Count() == 0)
            throw new InvalidDataException($"Не найдено калибровочных записей для марки: \"{brand.BrandName} - {brand.MixtureName}\"");
        cancellationToken.ThrowIfCancellationRequested();

        foreach (var sampleRecord in sampleRecords)
        {
            if (sampleRecord.Features.Count() == 0)
                _logger.LogError("Не найдено калибровочных признаков для пробы: \"{}\"", sampleRecord.SampleName);
            if (sampleRecord.Features.Count() % FEATURES_LENGTH != 0)
                _logger.LogError("Неожидаемое количество калибровочных признаков для пробы: \"{}\"", sampleRecord.SampleName);
        }

        _logger.LogInformation("Подготовка данных");
        var trainDataSet = sampleRecords
            .Where(r => r.Features.Count() > 0 && r.Features.Count() % FEATURES_LENGTH == 0)
            .Where(r => r.MeasuredContent >= 0 || r.MixtureActualRate >= 0)
            .SelectMany(r =>
                r.Features
                    .Chunk(FEATURES_LENGTH)
                    .Select(f =>
                        new RegressionData
                        {
                            Value = (float)(r.MeasuredContent >= 0 ? r.MeasuredContent : r.MixtureActualRate),
                            Features = f.ToArray()
                        }
                )
            )
            .Shuffle()
            .ToList();

        _logger.LogInformation("Загружено векторов для обучения: {}", trainDataSet.Count());
        if (trainDataSet.Count() == 0)
            throw new InvalidDataException($"Не найдено калибровочных данных для марки: \"{brand.BrandName} - {brand.MixtureName}\"");
        cancellationToken.ThrowIfCancellationRequested();

        int rowCount = trainDataSet.Count();
        double[][] featuresArray = trainDataSet.Select(d => d.Features.Select(f => (double)f).ToArray()).ToArray();
        Matrix<double> matrix = DenseMatrix.OfRowArrays(featuresArray);
        for (int i = 0; i < IMAGES_COUNT; i++)
        {
            matrix = matrix.RemoveColumn(i * 3);
            matrix = matrix.RemoveColumn(i * 3);
            matrix = matrix.RemoveColumn(i * 3);
        }
        cancellationToken.ThrowIfCancellationRequested();

        Vector<double> meanVector = matrix.ColumnSums() / matrix.RowCount;
        Matrix<double> covarianceMatrix = MatrixCovariance(matrix, meanVector);
        Matrix<double> invCovarianceMatrix = covarianceMatrix.Inverse();
        trainDataSet = trainDataSet
            .Where((data, index) =>
            {
                var row = matrix.Row(index);
                double distance = CalculateMahalanobis(row, meanVector, invCovarianceMatrix);

                return distance <= FILTERING_THRESHOLD;
            })
            .ToList();
        _logger.LogInformation("Векторов для обучения после фильтрации: {}", trainDataSet.Count());
        var trainData = _ml.Data.LoadFromEnumerable(trainDataSet);
        cancellationToken.ThrowIfCancellationRequested();

        _logger.LogInformation("Создание конвейера");
        var pipeline = _ml.Regression.Trainers.FastTree();
        cancellationToken.ThrowIfCancellationRequested();

        _logger.LogInformation("Обучение модели");
        var model = pipeline.Fit(trainData);
        _regressionEngine = _ml.Model.CreatePredictionEngine<RegressionData, RegressionPrediction>(model);

        model.Dispose();
    }

    public async Task<(Mat avgImage, Rect roi, Mat mask, IEnumerable<Point[]> contours, Mat result, string separatorType, string sampleBrand)?> Analyze(IEnumerable<ImageData> imagesData)
    {
        if (!Initialized)
            return null;

        if (imagesData.Count() != IMAGES_COUNT)
            return null;

        try
        {
            var segmentationResult = await Segment(imagesData);
            if (segmentationResult is null)
                return null;

            var (avgImage, roi, mask, contours, separatorType, sampleBrand) = segmentationResult.Value;

            var images = new List<(Mat mat, double visible, double uv365, double uv254)>();
            foreach (var imageData in imagesData
                .OrderBy(i => i.VisibleIntensity)
                .ThenBy(i => i.Uv365Intensity)
                .ThenBy(i => i.Uv254Intensity))
                images.Add((
                    LoadAndPrepare(imageData.Path, roi, ANALYSIS_IMAGE_SIZE),
                    imageData.VisibleIntensity,
                    imageData.Uv365Intensity,
                    imageData.Uv254Intensity
                ));

            var result = new Mat(new Size(ANALYSIS_IMAGE_SIZE, ANALYSIS_IMAGE_SIZE), MatType.CV_32FC1, new Scalar(0));
            var values = new float[ANALYSIS_IMAGE_SIZE * ANALYSIS_IMAGE_SIZE];

            for (var row = 0; row < ANALYSIS_IMAGE_SIZE; row++)
                for (var col = 0; col < ANALYSIS_IMAGE_SIZE; col++)
                {
                    if (mask.At<byte>(row, col) == 0)
                        continue;
                    var vector = new float[FEATURES_LENGTH];
                    foreach (var (i, image) in images.Enumerate())
                    {
                        vector[i * (CHANNELS_COUNT + 3) + 0] = (float)image.visible;
                        vector[i * (CHANNELS_COUNT + 3) + 1] = (float)image.uv365;
                        vector[i * (CHANNELS_COUNT + 3) + 2] = (float)image.uv254;
                        var color = image.mat.At<Vec3b>(row, col);
                        vector[i * (CHANNELS_COUNT + 3) + 3] = color.Item0 / 255.0f;
                        vector[i * (CHANNELS_COUNT + 3) + 4] = color.Item1 / 255.0f;
                        vector[i * (CHANNELS_COUNT + 3) + 5] = color.Item2 / 255.0f;
                    }
                    var prediction = _regressionEngine?.Predict(new RegressionData
                    {
                        Value = -1,
                        Features = vector
                    });
                    if (prediction is null)
                        continue;
                    if (prediction.Value >= 0)
                    {
                        result.Set(row, col, prediction.Value);
                        values[row * ANALYSIS_IMAGE_SIZE + col] = prediction.Value;
                    }
                }

            // foreach (var contour in contours)
            // {
            //     var contourMask = new Mat(new Size(ANALYSIS_IMAGE_SIZE, ANALYSIS_IMAGE_SIZE), MatType.CV_8UC1, new Scalar(0));
            //     var contourDistance = new Mat(new Size(ANALYSIS_IMAGE_SIZE, ANALYSIS_IMAGE_SIZE), MatType.CV_32FC1, new Scalar(0));
            //     Cv2.DrawContours(contourMask, [contour], -1, new Scalar(255), -1);
            //     Cv2.DistanceTransform(contourMask, contourDistance, DistanceTypes.L2, DistanceTransformMasks.Mask5);
            //     Cv2.Pow(contourDistance, 2, contourDistance);
            //     Cv2.MinMaxLoc(contourDistance, out double minVal, out double maxVal);
            //     contourDistance.ConvertTo(contourDistance, MatType.CV_32FC1, 1 / (maxVal - minVal), 1 * minVal / (maxVal - minVal));
            //     // Cv2.ConvertScaleAbs(contourDistance, contourDistance, -1, 1);
            //     Cv2.BitwiseNot(contourMask, contourMask);
            //     contourDistance.SetTo(new Scalar(1), contourMask);
            //     Cv2.Multiply(result, contourDistance, result);
            //     contourMask.Dispose();
            //     contourDistance.Dispose();
            // }

            // var maskDistance = new Mat();
            // var resultDistance = new Mat();
            // var temp = new Mat();
            // Cv2.DistanceTransform(mask, maskDistance, DistanceTypes.L2, DistanceTransformMasks.Mask5);

            // Cv2.MinMaxLoc(result, out double minValue, out double maxValue, out _, out _, mask);
            // result.ConvertTo(temp, MatType.CV_8UC1, 255 / (maxValue - 0), 255 * 0 / (maxValue - 0));
            // Cv2.DistanceTransform(temp, resultDistance, DistanceTypes.L2, DistanceTransformMasks.Mask5);

            // Cv2.Multiply(maskDistance, resultDistance, temp);
            // Cv2.Threshold(temp, temp, 0, 255, ThresholdTypes.Binary);

            // for (var row = 0; row < ANALYSIS_IMAGE_SIZE; row++)
            //     for (var col = 0; col < ANALYSIS_IMAGE_SIZE; col++)
            //         if (mask.At<byte>(row, col) == 0)
            //         {
            //             result.Set(row, col, 0);
            //             mask.Set(row, col, 0);
            //         }

            // maskDistance.Dispose();
            // resultDistance.Dispose();
            // temp.Dispose();

            // var oldMean = Cv2.Mean(result, mask).Val0;
            // float p001 = values[values.Count(v => v == 0) + 1];
            // float p999 = values[(int)(values.Length * 0.999)];
            // Cv2.Threshold(result, result, p999, p999, ThresholdTypes.Trunc);
            // var zeroMask = result.LessThanOrEqual(0);
            // result.SetTo(p001, zeroMask);
            // zeroMask.Dispose();
            // var newMean = Cv2.Mean(result, mask).Val0;
            // var scale = oldMean / newMean;
            // Cv2.ConvertScaleAbs(result, result, scale);

            // Добавляем эрозию к результатам и маске, удаляющую края, и пересчитываем контура
            // поскольку края гранул подсвечиваются близлежащими гранулами и сепаратором,
            // их нельзя считать представительными
            Cv2.Erode(result, result, Cv2.GetStructuringElement(MorphShapes.Ellipse, new Size(3, 3)));
            Cv2.Erode(mask, mask, Cv2.GetStructuringElement(MorphShapes.Ellipse, new Size(3, 3)));
            Cv2.FindContours(mask, out var newContours, out _, RetrievalModes.List, ContourApproximationModes.ApproxSimple);

            return (avgImage, roi, mask, contours, result, separatorType, sampleBrand);
        }
        catch (Exception ex)
        {
            _logger.LogError(ex, "Исключение при проведении анализа");
            return null;
        }
    }
    private async Task<(Mat avgImage, Rect roi, Mat mask, IEnumerable<Point[]> contours, string separatorType, string sampleBrand)?> Segment(IEnumerable<ImageData> imagesData)
    {
        var separatorPrefix = "SEPARATOR";
        var samplePrefix = "SAMPLE";

        if (imagesData.Count() != IMAGES_COUNT)
            return null;

        var avgImage = CalculateAverageImage(imagesData.Select(i => i.Path));
        if (avgImage is null)
            return null;

        var roi = CalculateSeparatorRoi(avgImage);
        if (roi is null)
        {
            avgImage.Release();
            return null;
        }

        avgImage = new Mat(avgImage, roi.Value);
        Cv2.Resize(avgImage, avgImage, new Size(ANALYSIS_IMAGE_SIZE, ANALYSIS_IMAGE_SIZE));

        var images = new List<(Mat mat, double visible, double uv365, double uv254)>();
        foreach (var imageData in imagesData
            .OrderBy(i => i.VisibleIntensity)
            .ThenBy(i => i.Uv365Intensity)
            .ThenBy(i => i.Uv254Intensity))
            images.Add((
                LoadAndPrepare(imageData.Path, roi.Value, ANALYSIS_IMAGE_SIZE),
                imageData.VisibleIntensity,
                imageData.Uv365Intensity,
                imageData.Uv254Intensity
            ));

        var separatorMask = new Mat(new Size(ANALYSIS_IMAGE_SIZE, ANALYSIS_IMAGE_SIZE), MatType.CV_8UC1, new Scalar(0));
        Cv2.Circle(separatorMask, ANALYSIS_IMAGE_SIZE / 2, ANALYSIS_IMAGE_SIZE / 2, ANALYSIS_IMAGE_SIZE / 2, new Scalar(255), -1);

        var result = new Mat(new Size(ANALYSIS_IMAGE_SIZE, ANALYSIS_IMAGE_SIZE), MatType.CV_16UC1, new Scalar(0));
        var keyMap = new Dictionary<string, ushort>();
        ushort lastKey = 0;

        for (var row = 0; row < ANALYSIS_IMAGE_SIZE; row++)
            for (var col = 0; col < ANALYSIS_IMAGE_SIZE; col++)
            {
                if (separatorMask.At<byte>(row, col) != 255)
                    continue;
                var vector = new float[FEATURES_LENGTH];
                foreach (var (i, image) in images.Enumerate())
                {
                    vector[i * (CHANNELS_COUNT + 3) + 0] = (float)image.visible;
                    vector[i * (CHANNELS_COUNT + 3) + 1] = (float)image.uv365;
                    vector[i * (CHANNELS_COUNT + 3) + 2] = (float)image.uv254;
                    var color = image.mat.At<Vec3b>(row, col);
                    vector[i * (CHANNELS_COUNT + 3) + 3] = color.Item0 / 255.0f;
                    vector[i * (CHANNELS_COUNT + 3) + 4] = color.Item1 / 255.0f;
                    vector[i * (CHANNELS_COUNT + 3) + 5] = color.Item2 / 255.0f;
                }
                var prediction = _clusterizationEngine?.Predict(new ClusterizationData
                {
                    Label = "",
                    Features = vector
                });
                if (prediction is null)
                    continue;
                if (!keyMap.ContainsKey(prediction.Label))
                    keyMap.Add(prediction.Label, ++lastKey);
                result.Set(row, col, keyMap[prediction.Label]);
            }

        foreach (var image in images)
            image.mat.Release();

        var valueMap = keyMap.ToDictionary(kv => kv.Value, kv => kv.Key);
        var frequency = new Dictionary<string, int>();

        for (var row = 0; row < ANALYSIS_IMAGE_SIZE; row++)
            for (var col = 0; col < ANALYSIS_IMAGE_SIZE; col++)
            {
                var value = result.At<ushort>(row, col);
                if (value == 0)
                    continue;

                if (!frequency.ContainsKey(valueMap[value]))
                    frequency.Add(valueMap[value], 0);
                frequency[valueMap[value]] += 1;
                if (valueMap[value].StartsWith(separatorPrefix))
                    result.Set<ushort>(row, col, 0);
                if (valueMap[value].StartsWith(samplePrefix))
                    result.Set(row, col, ushort.MaxValue);
            }

        result.ConvertTo(result, MatType.CV_8UC1, 1 / 255.0);

        Cv2.MorphologyEx(result, result, MorphTypes.Close, Cv2.GetStructuringElement(MorphShapes.Ellipse, new Size(3, 3)));
        Cv2.MorphologyEx(result, result, MorphTypes.Open, Cv2.GetStructuringElement(MorphShapes.Ellipse, new Size(3, 3)));
        Cv2.Dilate(result, result, Cv2.GetStructuringElement(MorphShapes.Ellipse, new Size(3, 3)));
        Cv2.GaussianBlur(result, result, new Size(5, 5), 1);
        Cv2.Threshold(result, result, 127, 255, ThresholdTypes.Binary);

        Cv2.FindContours(result, out var contours, out _, RetrievalModes.List, ContourApproximationModes.ApproxSimple);

        result = new Mat(new Size(ANALYSIS_IMAGE_SIZE, ANALYSIS_IMAGE_SIZE), MatType.CV_8UC1, new Scalar(0));
        Cv2.DrawContours(result, contours.Where(ContourIsValid), -1, new Scalar(255), -1);

        _logger.LogInformation("Найдено гранул: {}", contours.Where(ContourIsValid).Count());
        _logger.LogInformation("Детекторованные обекты:\n\t{}", string.Join("\n\t", frequency.Select(f => $"{f.Key}: {f.Value}")));

        return (
            avgImage,
            roi.Value,
            result,
            contours.Where(ContourIsValid),
            frequency
                .Where(kv => kv.Key.StartsWith(separatorPrefix))
                .MaxBy(kv => kv.Value)
                .Key
                .Substring(separatorPrefix.Count() + 1),
            frequency
                .Where(kv => kv.Key.StartsWith(samplePrefix))
                .MaxBy(kv => kv.Value)
                .Key
                .Substring(samplePrefix.Count() + 1)
        );
    }

    public static (Mat avgImage, Rect roi, List<float[]> vectors)? GetFeatureVectors(IEnumerable<ImageData> imagesData)
    {
        if (imagesData.Count() != IMAGES_COUNT)
            return null;

        var avgImage = CalculateAverageImage(imagesData.Select(i => i.Path));
        if (avgImage is null)
            return null;

        var roi = CalculateSeparatorRoi(avgImage);
        if (roi is null)
            return null;

        var images = new List<(Mat mat, double visible, double uv365, double uv254)>();
        foreach (var imageData in imagesData
            .OrderBy(i => i.VisibleIntensity)
            .ThenBy(i => i.Uv365Intensity)
            .ThenBy(i => i.Uv254Intensity))
            images.Add((
                LoadAndPrepare(imageData.Path, roi.Value),
                imageData.VisibleIntensity,
                imageData.Uv365Intensity,
                imageData.Uv254Intensity
            ));

        var separatorMask = new Mat(new Size(EFFECTIVE_IMAGE_SIZE, EFFECTIVE_IMAGE_SIZE), MatType.CV_8UC1, new Scalar(0));
        Cv2.Circle(separatorMask, EFFECTIVE_IMAGE_SIZE / 2, EFFECTIVE_IMAGE_SIZE / 2, EFFECTIVE_IMAGE_SIZE / 2, new Scalar(255), -1);

        var vectors = new List<float[]>();

        for (var row = 0; row < EFFECTIVE_IMAGE_SIZE; row++)
            for (var col = 0; col < EFFECTIVE_IMAGE_SIZE; col++)
            {
                if (separatorMask.At<byte>(row, col) != 255)
                    continue;
                var vector = new float[FEATURES_LENGTH];

                foreach (var (i, image) in images.Enumerate())
                {
                    vector[i * (CHANNELS_COUNT + 3) + 0] = (float)image.visible;
                    vector[i * (CHANNELS_COUNT + 3) + 1] = (float)image.uv365;
                    vector[i * (CHANNELS_COUNT + 3) + 2] = (float)image.uv254;

                    var color = image.mat.At<Vec3b>(row, col).ToVec3f() / 255.0f;

                    vector[i * (CHANNELS_COUNT + 3) + 3] = color.Item0;
                    vector[i * (CHANNELS_COUNT + 3) + 4] = color.Item1;
                    vector[i * (CHANNELS_COUNT + 3) + 5] = color.Item2;
                }

                vectors.Add(vector);
            }

        foreach (var (mat, _, _, _) in images)
            mat.Release();

        return (avgImage, roi.Value, vectors);
    }
    private static Mat LoadAndPrepare(string imagePath, Rect roi, int size = EFFECTIVE_IMAGE_SIZE)
    {
        var image = Cv2.ImRead(imagePath);
        image = new Mat(image, roi);
        Cv2.GaussianBlur(image, image, new Size(5, 5), 1);
        Cv2.Resize(image, image, new Size(size, size));
        return image;
    }
    public static Mat? CalculateAverageImage(IEnumerable<string> imagePaths)
    {
        if (imagePaths.Count() == 0)
            return null;

        Mat avgImage = new Mat();

        foreach (var imagePath in imagePaths)
        {
            var image = Cv2.ImRead(imagePath);

            if (image.Channels() != 3)
                return null;

            if (avgImage.Empty())
                avgImage = new Mat(image.Size(), MatType.CV_16UC3, new Scalar(0, 0, 0));

            if (image.Size().Width != avgImage.Size().Width ||
                image.Size().Height != avgImage.Size().Height)
                return null;

            Cv2.Add(image, avgImage, avgImage, dtype: MatType.CV_16U);

            image.Release();
        }

        avgImage.ConvertTo(avgImage, MatType.CV_8UC(avgImage.Channels()), 1f / imagePaths.Count());

        return avgImage;
    }
    public static Rect? CalculateSeparatorRoi(Mat avgImage)
    {
        if (avgImage.Empty())
            return null;
        return new Rect(SEPARATOR_X, SEPARATOR_Y, SEPARATOR_WIDTH, SEPARATOR_HEIGHT);
    }
    // Этот метод конечно хорош, но проблема с тем, что при калибровке на чёрном продукте
    // он выдаёт очень обрезает часть полезного изображения из-за неравномерности освещения
    // public static Rect? CalculateSeparatorRoi(Mat avgImage)
    // {
    //     if (avgImage.Empty())
    //         return null;

    //     Mat avgImageClone = new Mat();
    //     Mat mask = new Mat();
    //     Rect? roi = null;

    //     avgImage.CopyTo(avgImageClone);

    //     avgImageClone = avgImageClone.CvtColor(ColorConversionCodes.BGR2GRAY);
    //     Cv2.Threshold(avgImageClone, mask, SEPARATOR_ROI_THRESHOLD, 255, ThresholdTypes.Binary);
    //     Cv2.FindContours(mask, out var contours, out _, RetrievalModes.List, ContourApproximationModes.ApproxSimple);

    //     var maxContour = contours.MaxBy(c => Cv2.ContourArea(c));
    //     if (maxContour is null)
    //         return null;

    //     roi = Cv2.BoundingRect(maxContour);

    //     avgImageClone.Release();
    //     mask.Release();

    //     return roi;
    // }
    private static bool ContourIsValid(Point[] contour)
    {
        var bbox = Cv2.MinAreaRect(contour);
        var area = Cv2.ContourArea(contour);
        return CONTOUR_MIN_SIZE < bbox.Size.Width && CONTOUR_MIN_SIZE < bbox.Size.Height && CONTOUR_MIN_AREA < area &&
               CONTOUR_MAX_SIZE > bbox.Size.Width && CONTOUR_MAX_SIZE > bbox.Size.Height && CONTOUR_MAX_AREA > area;
    }

    // Формула расстояния Махаланобиса: d = sqrt((x - mu)^T * Sigma^-1 * (x - mu))
    private static double CalculateMahalanobis(Vector<double> x, Vector<double> mean, Matrix<double> invCovariance)
    {
        var diff = x - mean;
        var distanceSq = diff * invCovariance * diff;
        return Math.Sqrt(distanceSq);
    }

    // Ручной расчет матрицы ковариации (для экономии памяти и контроля)
    private static Matrix<double> MatrixCovariance(Matrix<double> matrix, Vector<double> meanRows)
    {
        int rows = matrix.RowCount;
        int cols = matrix.ColumnCount;
        var cov = new DenseMatrix(cols, cols);

        for (int i = 0; i < cols; i++)
            for (int j = i; j < cols; j++)
            {
                double sum = 0;
                for (int k = 0; k < rows; k++)
                {
                    sum += (matrix[k, i] - meanRows[i]) * (matrix[k, j] - meanRows[j]);
                }
                double val = sum / (rows - 1);
                cov[i, j] = val;
                cov[j, i] = val;
            }
        return cov;
    }
}