/*
 * Copyright (c) 2026 University of Insubria <https://www.uninsubria.it>
 *
 * This library is free software; you can redistribute it and/or modify it
 * under the terms of the GNU Lesser General Public License as published by
 * the Free Software Foundation; either version 2.1 of the License, or (at
 * your option) any later version. All we ask is that proper credit is given
 * for our work, which includes - but is not limited to - adding the above
 * copyright notice to the beginning of your source code files, and to any
 * copyright notice that you may distribute with programs based on this work.
 *
 * This library is distributed in the hope that it will be useful, but WITHOUT
 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
 * FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
 * License for more details.
 *
 * You should have received a copy of the GNU Lesser General Public License
 * along with this program; if not, write to the Free Software Foundation, Inc.
 * <https://fsf.org/>
 */

package qsardes.descriptor;

import Jama.Matrix;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.HashMap;
import java.util.HashSet;
import java.util.List;
import java.util.Map;
import java.util.Set;
import org.jgrapht.Graph;
import org.jgrapht.GraphPath;
import org.jgrapht.alg.connectivity.BlockCutpointGraph;
import org.jgrapht.alg.shortestpath.AllDirectedPaths;
import org.jgrapht.graph.DefaultDirectedGraph;
import org.openscience.cdk.graph.PathTools;
import org.openscience.cdk.interfaces.IAtomContainer;
import qsardes.molecule.AbstractMoleculeBuilder;
import qsardes.molecule.CustomEdge;
import qsardes.molecule.Molecule;
import qsardes.util.UtilsCalc;
import qsardes.util.UtilsConv;

/**
 * Detour matrix descriptors abstract class.
 *
 * @author Nicola Chirico
 */
public abstract class DescriptorDetourMatrix implements Descriptor
{
    protected DescriptorDetourMatrix()
    {
    }

    @Override
    public String toString()
    {
        return "Detour matrix";
    }

    @Override
    public DescriptorType get_descriptor_type()
    {
        return DescriptorType.valueOf("DESCRIPTOR");
    }

    @Override
    public List<String> get_descriptor_names()
    {
        return Arrays.asList(
                "SpMax_Dt", "SpDiam_Dt", "SpAD_Dt", "SpMAD_Dt",
                "EE_Dt",    "VE1_Dt",    "VE2_Dt",  "VE3_Dt",
                "VR1_Dt",   "VR2_Dt",    "VR3_Dt");

    }

    @Override
    public Map<String, String> get_descriptor_values(
            AbstractMoleculeBuilder mol_build)
    {
        prepare_descr_map();

        try {
            mol = new Molecule(mol_build);
            calculate_descr();
        } catch (Exception ex) {
            System.err.println(ex);
        }

        return descr_map;
    }

    protected abstract double calculate_sp_ve1_sum(UtilsCalc utils_calc,
            double[][] eig_vec_array);

    protected abstract double calculate_sp_vr1_sum(IAtomContainer iac,
            UtilsCalc utils_calc, double[][] eig_vec_array);

    protected abstract double get_sp_ve1_sum(double sp_ve1_sum);

    protected abstract double get_ve3(double sp_ve1_sum, int atom_count);

    protected abstract double get_vr3(double vr1_sum, int atom_count);

    private void prepare_descr_map()
    {
        final String NAN = String.valueOf(Double.NaN);

        get_descriptor_names().forEach(name -> descr_map.put(name, NAN));
    }

    /*
      Detour matrix is calculated by merging molecular fragments obtained
      by articulation points splitting, as suggested by Moriwaki et al. 2018
      (Mordred descriptor)
    */
    private void calculate_descr()
    {
        final int ATOM_COUNT = mol.get_atom_container().getAtomCount();

        detect_vert_paths_and_artpoints();

        var matr = new Matrix(new UtilsConv().convert_to_double(
                get_detour_matr()));

        if (utils_calc.are_all_matrix_values_valid(matr.getArray())) {
            var eig_decomp = matr.eig();

            double[] eig_values = eig_decomp.getRealEigenvalues();
            double[][] eig_vec_array = eig_decomp.getV().getArray();

            double eig_val_min = utils_calc.get_min(eig_values);
            double eig_val_max = utils_calc.get_max(eig_values);

            double sp_ad_sum = utils_calc.calc_ad_sum(eig_values);
            double sp_ve1_sum = calculate_sp_ve1_sum(utils_calc,
                        eig_vec_array);

            double ee_lse_exp_sum = utils_calc.calc_ee_lse_exp_sum(eig_values);

            double vr1_sum = calculate_sp_vr1_sum(mol.get_atom_container(),
                        utils_calc, eig_vec_array);

            descr_map.put("SpMax_Dt", String.valueOf(eig_val_max));

            descr_map.put("SpDiam_Dt",
                    String.valueOf(eig_val_max - eig_val_min));

            descr_map.put("SpAD_Dt", String.valueOf(sp_ad_sum));

            descr_map.put("SpMAD_Dt", String.valueOf(sp_ad_sum
                    / ATOM_COUNT));

            /*
               log(1 + b) -> log(a + b) = log(exp(log(a)) + exp(log(b)))
               eig_val_max reassembles the final part of the Log-Sum-Exp func.
             */
            try {
                descr_map.put("EE_Dt", String.valueOf(
                        StrictMath.log(StrictMath.exp(StrictMath.log(1.0))
                                + StrictMath.exp(eig_val_max
                                        + StrictMath.log(ee_lse_exp_sum)))));
            } catch (Exception ex) {
                descr_map.put("EE_Dt", String.valueOf(Double.NaN));
            }

            descr_map.put("VE1_Dt", String.valueOf(
                    get_sp_ve1_sum(sp_ve1_sum)));

            descr_map.put("VE2_Dt", String.valueOf(
                    get_sp_ve1_sum(sp_ve1_sum) / (double) ATOM_COUNT));

            descr_map.put("VE3_Dt", String.valueOf(
                    get_ve3(sp_ve1_sum, ATOM_COUNT)));

            descr_map.put("VR1_Dt", String.valueOf(vr1_sum));

            descr_map.put("VR2_Dt", String.valueOf(vr1_sum
                    / (double) ATOM_COUNT));

            descr_map.put("VR3_Dt", String.valueOf(
                        get_vr3(vr1_sum, ATOM_COUNT)));
        }
    }

    private void detect_vert_paths_and_artpoints()
    {
        vert_paths_list = new ArrayList<>();

        var blocks_graph = new BlockCutpointGraph<>(mol.get_graph());

        var block_iter = blocks_graph.getBlocks().iterator();

        while (block_iter.hasNext()) {
            var graph = convert_to_directed_graph(block_iter.next());

            List<Integer> vert_list = new ArrayList<>();

            vert_list.addAll(graph.vertexSet());

            vert_paths_list.add(new VertPaths(vert_list,
                    new AllDirectedPaths<>(graph)));
        }
    }

    private int[][] get_detour_matr()
    {
        final int BLOCKS_NUM = vert_paths_list.size();

        /*
          If the molecule cannot be splitted by articulation points, use a
          CDK PathTools.getAllPaths based algorithm because, in this case, is
          faster than running the whole merge algorithm.
        */
        if (BLOCKS_NUM <= 1) {
            return get_detour_matr_np();
        }

        int[][] detour_matr = get_blocks_detour_matr();

        Boolean[] is_block_visited = new Boolean[BLOCKS_NUM];
        Arrays.fill(is_block_visited, Boolean.FALSE);

        List<HashSet<Integer>> vert_set_list = get_blocks_vertex_list();

        Set<Integer> cumul_vert_set = new HashSet<>();

        cumul_vert_set.addAll(vert_set_list.get(0));
        is_block_visited[0] = Boolean.TRUE;

        while (Arrays.asList(is_block_visited).contains(Boolean.FALSE)) {
            for (int block = 1; block < BLOCKS_NUM; block++) {
                if (is_block_visited[block]) {
                    continue;
                }

                var cur_vert_set = vert_set_list.get(block);

                var joint_set = copy_set(cumul_vert_set);

                joint_set.retainAll(cur_vert_set);

                if (joint_set.size() < 1) {
                    continue;
                }

                if (joint_set.size() == 1) {
                    Integer joint = joint_set.iterator().next();

                    for (Integer cur_vert : cur_vert_set) {
                        if (cur_vert.equals(joint)) {
                            continue;
                        }

                        Integer cur_dist = detour_matr[cur_vert][joint];

                        for (Integer cum_vert : cumul_vert_set) {
                            if (cum_vert.equals(joint)) {
                                continue;
                            }

                            Integer dist = cur_dist
                                    + detour_matr[cum_vert][joint];

                            detour_matr[cur_vert][cum_vert] = dist;
                            detour_matr[cum_vert][cur_vert] = dist;
                        }

                        is_block_visited[block] = Boolean.TRUE;
                    }

                    cumul_vert_set.addAll(cur_vert_set);
                } else {
                    throw new UnsupportedOperationException(
                            "Cannot handle " + joint_set.size()
                            + " articulation points at the same time.");
                }
            }
        }

        return detour_matr;
    }

    private int[][] get_detour_matr_np()
    {
        var iac = mol.get_atom_container();

        final int ATOM_COUNT = iac.getAtomCount();

        int[][] detour_matr = new int[ATOM_COUNT][ATOM_COUNT];

        for (var row : detour_matr) {
            Arrays.fill(row, 0);
        }

        for (int i = 0; i < ATOM_COUNT - 1; i++) {
            var atom_i = iac.getAtom(i);

            for (int j = i + 1; j < ATOM_COUNT; j++) {
                var atom_j = iac.getAtom(j);

                var paths = PathTools.getAllPaths(iac, atom_i, atom_j);

                int max_len = 0;

                for (var path : paths) {
                    int len = path.size();

                    if (len > max_len) {
                        max_len = len;
                    }
                }

                max_len--;

                detour_matr[i][j] = max_len;
                detour_matr[j][i] = max_len;
            }
        }

        return detour_matr;
    }

    private List<HashSet<Integer>> get_blocks_vertex_list()
    {
        List<HashSet<Integer>> vert_list = new ArrayList<>();

        vert_paths_list.forEach(item -> {
            HashSet<Integer> list = new HashSet<>();

            item.vert_list().forEach(vert -> {
                list.add(vert);
            });

            vert_list.add(list);
        });

        return vert_list;
    }

    private Set<Integer> copy_set(Set<Integer> src_set)
    {
        Set<Integer> set = new HashSet<>();

        src_set.forEach(item -> set.add(item));

        return set;
    }

    private int[][] get_blocks_detour_matr()
    {
        final int ATOM_COUNT = mol.get_atom_container().getAtomCount();

        int[][] detour_matr = new int[ATOM_COUNT][ATOM_COUNT];

        for (var row : detour_matr) {
            Arrays.fill(row, 0);
        }

        for (var vert_paths : vert_paths_list) {
            final int VERT_NUM = vert_paths.vert_list.size();

            for (int i = 0; i < VERT_NUM; i++) {
                final Integer SRC_V = vert_paths.vert_list.get(i);

                for (int j = i + 1; j < VERT_NUM; j++) {
                    final Integer DST_V = vert_paths.vert_list.get(j);

                    var sub_paths = vert_paths.paths
                            .getAllPaths(SRC_V, DST_V, true, null);

                    int max_len = get_max_path_length(sub_paths);

                    detour_matr[SRC_V][DST_V] = max_len;
                    detour_matr[DST_V][SRC_V] = max_len;
                }
            }
        }

        return detour_matr;
    }

    private int get_max_path_length(
            List<GraphPath<Integer, CustomEdge>> sub_paths)
    {
        final int PATHS_NUM = sub_paths.size();

        int[] dist = new int[PATHS_NUM];

        for (int i = 0; i < PATHS_NUM; i++) {
            dist[i] = sub_paths.get(i).getLength();
        }

        return utils_calc.get_max(dist);
    }

    private DefaultDirectedGraph<Integer, CustomEdge>
            convert_to_directed_graph(Graph<Integer, CustomEdge> block)
    {
        DefaultDirectedGraph<Integer, CustomEdge> dir_graph
                = new DefaultDirectedGraph<>(CustomEdge.class);

        for (var item : block.vertexSet()) {
            dir_graph.addVertex(item);
        }

        for (var item : block.edgeSet()) {
            dir_graph.addEdge(item.get_source(), item.get_target(),
                    new CustomEdge());
        }

        return dir_graph;
    }

    private final record VertPaths(List<Integer> vert_list,
            AllDirectedPaths<Integer, CustomEdge> paths)
    {
    }

    private List<VertPaths> vert_paths_list;

    private final Map<String, String> descr_map = new HashMap<>();

    private Molecule mol;

    private final UtilsCalc utils_calc = new UtilsCalc();
}
