#include "tsplib_graph.hpp"

#include <algorithm>
#include <cctype>
#include <cmath>
#include <fstream>
#include <limits>
#include <sstream>
#include <stdexcept>
#include <string>
#include <utility>
#include <vector>

// ─── internal helpers ─────────────────────────────────────────────────────────
namespace {

std::string trim(const std::string& s) {
    const char* ws = " \t\r\n";
    const size_t b = s.find_first_not_of(ws);
    if (b == std::string::npos) return {};
    return s.substr(b, s.find_last_not_of(ws) - b + 1);
}

std::string upper(std::string s) {
    for (char& c : s)
        c = static_cast<char>(std::toupper(static_cast<unsigned char>(c)));
    return s;
}

bool is_section_keyword(const std::string& u) {
    return u == "NODE_COORD_SECTION"   ||
           u == "EDGE_WEIGHT_SECTION"  ||
           u == "DISPLAY_DATA_SECTION" ||
           u == "FIXED_EDGES_SECTION"  ||
           u == "TOUR_SECTION"         ||
           u == "EOF";
}

// ─── GEO distance (TSPLIB spec, section 2.1) ─────────────────────────────────

// TSPLIB GEO coordinates are encoded as DDD.MM (degrees + minutes).
double geo_rad(double x) {
    const int    deg = static_cast<int>(x);
    const double min = x - deg;
    return M_PI * (deg + 5.0 * min / 3.0) / 180.0;
}

double geo_dist(double lat1, double lon1, double lat2, double lon2) {
    constexpr double RRR = 6378.388;
    const double q1  = std::cos(lon1 - lon2);
    const double q2  = std::cos(lat1 - lat2);
    const double q3  = std::cos(lat1 + lat2);
    // clamp avoids NaN from floating-point rounding outside [-1, 1]
    const double arg = std::clamp(0.5 * ((1.0 + q1) * q2 - (1.0 - q1) * q3),
                                  -1.0, 1.0);
    return RRR * std::acos(arg);
}

} // namespace

// ─── TSPLIBGraph::load ────────────────────────────────────────────────────────
TSPLIBGraph TSPLIBGraph::load(const std::string& filepath) {
    std::ifstream fin(filepath);
    if (!fin.is_open())
        throw std::runtime_error("TSPLIBGraph: cannot open '" + filepath + "'");

    // ── Pass 1: slurp all lines; tag section-header lines ─────────────────────
    struct SecInfo { std::string name; size_t idx; };
    std::vector<std::string> lines;
    std::vector<SecInfo>     secs;

    for (std::string raw; std::getline(fin, raw); ) {
        std::string t = trim(raw);
        const std::string u = upper(t);
        if (is_section_keyword(u))
            secs.push_back({u, lines.size()});
        lines.push_back(std::move(t));
    }

    // ── Pass 2: parse header key: value pairs ──────────────────────────────────
    const size_t header_end = secs.empty() ? lines.size() : secs[0].idx;

    std::string ewt;   // EDGE_WEIGHT_TYPE
    std::string ewf;   // EDGE_WEIGHT_FORMAT
    int         dim = 0;

    for (size_t i = 0; i < header_end; ++i) {
        const std::string& ln    = lines[i];
        const size_t       colon = ln.find(':');
        if (colon == std::string::npos) continue;

        const std::string key = upper(trim(ln.substr(0, colon)));
        const std::string val = trim(ln.substr(colon + 1));

        if      (key == "DIMENSION")          dim = std::stoi(val);
        else if (key == "EDGE_WEIGHT_TYPE")   ewt = upper(val);
        else if (key == "EDGE_WEIGHT_FORMAT") ewf = upper(val);
    }

    if (dim <= 0)
        throw std::runtime_error("TSPLIBGraph: DIMENSION not found or invalid");
    if (ewt.empty())
        throw std::runtime_error("TSPLIBGraph: EDGE_WEIGHT_TYPE not found");

    // ── Helpers: locate and extract section text ───────────────────────────────
    auto section_text = [&](const std::string& name)
        -> std::pair<bool, std::string>
    {
        for (const auto& s : secs) {
            if (s.name != name) continue;
            // collect lines until the next section keyword (or end of file)
            size_t end = lines.size();
            for (const auto& ns : secs)
                if (ns.idx > s.idx && ns.idx < end) end = ns.idx;
            std::string buf;
            for (size_t i = s.idx + 1; i < end; ++i) {
                buf += lines[i];
                buf += '\n';
            }
            return {true, std::move(buf)};
        }
        return {false, {}};
    };

    auto require_section = [&](const std::string& name) -> std::string {
        auto [found, text] = section_text(name);
        if (!found)
            throw std::runtime_error(
                "TSPLIBGraph: section '" + name + "' not found in '" + filepath + "'");
        return text;
    };

    // ── Allocate storage ───────────────────────────────────────────────────────
    TSPLIBGraph g;
    g.n_ = dim;
    g.data_.assign(static_cast<size_t>(dim) * (dim - 1) / 2, 0.0);

    const int n = dim;

    // ── Coordinate-based distance types ───────────────────────────────────────
    const bool is_3d = (ewt == "EUC_3D" || ewt == "MAX_3D" || ewt == "MAN_3D");
    const bool coord_type = (ewt == "EUC_2D"  || ewt == "EUC_3D"  ||
                             ewt == "MAX_2D"  || ewt == "MAX_3D"  ||
                             ewt == "MAN_2D"  || ewt == "MAN_3D"  ||
                             ewt == "CEIL_2D" || ewt == "GEO"     || ewt == "ATT");

    if (coord_type) {
        std::istringstream ss(require_section("NODE_COORD_SECTION"));

        std::vector<double> cx(n), cy(n), cz(n, 0.0);
        for (int i = 0; i < n; ++i) {
            int node_id;
            if (!(ss >> node_id >> cx[i] >> cy[i]))
                throw std::runtime_error(
                    "TSPLIBGraph: unexpected end of NODE_COORD_SECTION");
            if (is_3d && !(ss >> cz[i]))
                throw std::runtime_error(
                    "TSPLIBGraph: missing z-coordinate in NODE_COORD_SECTION");
        }

        for (int i = 1; i < n; ++i) {
            for (int j = 0; j < i; ++j) {
                double d = 0.0;
                if (ewt == "EUC_2D" || ewt == "GEO") {
                    // GEO deliberately uses Euclidean distance on the raw coordinates
                    // (NOT the TSPLIB great-circle formula) to match the distance
                    // convention of the benchmark literature we compare against; using
                    // true GEO distances would make those results incomparable. This is
                    // why the geo_dist/geo_rad helpers and the GEO branch below are unused.
                    const double dx = cx[i]-cx[j], dy = cy[i]-cy[j];
                    d = std::sqrt(dx*dx + dy*dy);
                } else if (ewt == "EUC_3D") {
                    const double dx = cx[i]-cx[j], dy = cy[i]-cy[j], dz = cz[i]-cz[j];
                    d = std::sqrt(dx*dx + dy*dy + dz*dz);
                } else if (ewt == "MAX_2D") {
                    d = std::max(std::abs(cx[i]-cx[j]), std::abs(cy[i]-cy[j]));
                } else if (ewt == "MAX_3D") {
                    d = std::max({std::abs(cx[i]-cx[j]),
                                  std::abs(cy[i]-cy[j]),
                                  std::abs(cz[i]-cz[j])});
                } else if (ewt == "MAN_2D") {
                    d = std::abs(cx[i]-cx[j]) + std::abs(cy[i]-cy[j]);
                } else if (ewt == "MAN_3D") {
                    d = std::abs(cx[i]-cx[j]) + std::abs(cy[i]-cy[j]) + std::abs(cz[i]-cz[j]);
                } else if (ewt == "CEIL_2D") {
                    const double dx = cx[i]-cx[j], dy = cy[i]-cy[j];
                    d = std::ceil(std::sqrt(dx*dx + dy*dy));
//                } else if (ewt == "GEO") {
//                    d = geo_dist(geo_rad(cx[i]), geo_rad(cy[i]),
//                                 geo_rad(cx[j]), geo_rad(cy[j]));
                } else { // ATT: pseudo-Euclidean
                    const double dx = cx[i]-cx[j], dy = cy[i]-cy[j];
                    d = std::sqrt((dx*dx + dy*dy) / 10.0);
                }
                g.data_[tri_index(i, j)] = d;
            }
        }
        return g;
    }

    // ── EXPLICIT ───────────────────────────────────────────────────────────────
    if (ewt == "EXPLICIT") {
        if (ewf.empty())
            throw std::runtime_error(
                "TSPLIBGraph: EDGE_WEIGHT_FORMAT required for EXPLICIT type");

        std::istringstream ss(require_section("EDGE_WEIGHT_SECTION"));

        auto read_val = [&]() -> double {
            double v;
            if (!(ss >> v))
                throw std::runtime_error(
                    "TSPLIBGraph: unexpected end of EDGE_WEIGHT_SECTION");
            return v;
        };

        if (ewf == "FULL_MATRIX") {
            for (int i = 0; i < n; ++i)
                for (int j = 0; j < n; ++j) {
                    const double v = read_val();
                    if (i != j) g.store(i, j, v);
                }

        } else if (ewf == "UPPER_ROW") {
            // row i: columns i+1 … n-1
            for (int i = 0; i < n-1; ++i)
                for (int j = i+1; j < n; ++j)
                    g.store(i, j, read_val());

        } else if (ewf == "LOWER_ROW") {
            // row i: columns 0 … i-1
            for (int i = 1; i < n; ++i)
                for (int j = 0; j < i; ++j)
                    g.store(i, j, read_val());

        } else if (ewf == "UPPER_DIAG_ROW") {
            // row i: columns i … n-1  (diagonal entries are ignored)
            for (int i = 0; i < n; ++i)
                for (int j = i; j < n; ++j) {
                    const double v = read_val();
                    if (i != j) g.store(i, j, v);
                }

        } else if (ewf == "LOWER_DIAG_ROW") {
            // row i: columns 0 … i  (diagonal entries are ignored)
            for (int i = 0; i < n; ++i)
                for (int j = 0; j <= i; ++j) {
                    const double v = read_val();
                    if (i != j) g.store(i, j, v);
                }

        } else if (ewf == "UPPER_COL") {
            // column j: rows 0 … j-1
            for (int j = 1; j < n; ++j)
                for (int i = 0; i < j; ++i)
                    g.store(i, j, read_val());

        } else if (ewf == "LOWER_COL") {
            // column j: rows j+1 … n-1
            for (int j = 0; j < n-1; ++j)
                for (int i = j+1; i < n; ++i)
                    g.store(i, j, read_val());

        } else if (ewf == "UPPER_DIAG_COL") {
            // column j: rows 0 … j  (diagonal entries are ignored)
            for (int j = 0; j < n; ++j)
                for (int i = 0; i <= j; ++i) {
                    const double v = read_val();
                    if (i != j) g.store(i, j, v);
                }

        } else if (ewf == "LOWER_DIAG_COL") {
            // column j: rows j … n-1  (diagonal entries are ignored)
            for (int j = 0; j < n; ++j)
                for (int i = j; i < n; ++i) {
                    const double v = read_val();
                    if (i != j) g.store(i, j, v);
                }

        } else {
            throw std::runtime_error(
                "TSPLIBGraph: unsupported EDGE_WEIGHT_FORMAT '" + ewf + "'");
        }
        return g;
    }

    throw std::runtime_error(
        "TSPLIBGraph: unsupported EDGE_WEIGHT_TYPE '" + ewt + "'");
}

// ─── TSPLIBGraph::load_txt ──────────────────────────────────────────────────
TSPLIBGraph TSPLIBGraph::load_txt(const std::string& filepath) {
    std::ifstream fin(filepath);
    if (!fin.is_open())
        throw std::runtime_error("TSPLIBGraph: cannot open '" + filepath + "'");

    int n, m;
    if (!(fin >> n >> m) || n <= 0)
        throw std::runtime_error("TSPLIBGraph: invalid header in '" + filepath + "'");

    TSPLIBGraph g;
    g.n_ = n;
    g.data_.assign(static_cast<size_t>(n) * (n - 1) / 2,
                   std::numeric_limits<double>::infinity());

    for (int i = 0; i < m; ++i) {
        int u, v;
        double c;
        if (!(fin >> u >> v >> c))
            throw std::runtime_error(
                "TSPLIBGraph: unexpected end of edge list at line " +
                std::to_string(i + 2) + " in '" + filepath + "'");
        if (u < 0 || u >= n || v < 0 || v >= n || u == v)
            throw std::runtime_error(
                "TSPLIBGraph: invalid edge (" + std::to_string(u) + ", " +
                std::to_string(v) + ") in '" + filepath + "'");
        g.data_[tri_index(u > v ? u : v, u > v ? v : u)] = c;
    }

    return g;
}