// regularisation.cpp

#include <opencv2/opencv.hpp>
#include <iostream>
#include <fstream>
#include <sstream>
#include <vector>
#include <cmath>
#include <algorithm>
#include <map>

// Type definitions
using Point2D = cv::Point2d;
using Path = std::vector<Point2D>;
using Paths = std::vector<Path>;

// Function to read CSV and organize points into paths
Paths read_csv(const std::string& csv_path) {
    Paths path_XYs;
    std::ifstream file(csv_path);
    if (!file.is_open()) {
        std::cerr << "Failed to open " << csv_path << std::endl;
        return path_XYs;
    }

    std::string line;
    std::map<int, std::vector<Point2D>> paths_map;

    while (std::getline(file, line)) {
        std::stringstream ss(line);
        std::string item;
        std::vector<double> values;
        while (std::getline(ss, item, ',')) {
            try {
                values.push_back(std::stod(item));
            } catch (const std::invalid_argument& e) {
                // Handle non-numeric data
                std::cerr << "Invalid data encountered: " << item << std::endl;
                values.push_back(0.0);
            }
        }
        if (values.size() >= 3) {
            int path_id = static_cast<int>(values[0]);
            Point2D pt = { values[1], values[2] };
            paths_map[path_id].push_back(pt);
        }
    }

    // Organize paths based on unique path IDs
    for (const auto& pair : paths_map) {
        path_XYs.push_back(pair.second);
    }

    return path_XYs;
}

// Function to calculate angle between three points
double calculate_angle(const Point2D& p1, const Point2D& p2, const Point2D& p3) {
    cv::Point2d a = p1 - p2;
    cv::Point2d b = p3 - p2;

    double dot_product = a.dot(b);
    double cross_product = a.x * b.y - a.y * b.x;
    double angle = std::atan2(cross_product, dot_product) * 180.0 / CV_PI;
    if (angle < 0) angle += 360.0;
    return angle;
}

// Function to check if points form a straight line
bool is_straight_line(const Path& XY, double tolerance = 0.05) {
    if (XY.size() < 3) return true;

    double total_length = 0.0;
    for (size_t i = 1; i < XY.size(); ++i) {
        double dx = XY[i].x - XY[i - 1].x;
        double dy = XY[i].y - XY[i - 1].y;
        total_length += std::sqrt(dx * dx + dy * dy);
    }
    double start_to_end = cv::norm(XY.back() - XY.front());
    return (std::abs(total_length - start_to_end) / total_length) < tolerance;
}

// Function to fit a circle using Kasa's method
bool kasa_circle_fitting(const Path& XY, double& xc, double& yc, double& r) {
    int n = XY.size();
    if (n < 3) return false;

    cv::Mat A(n, 3, CV_64F);
    cv::Mat B(n, 1, CV_64F);
    for (int i = 0; i < n; ++i) {
        A.at<double>(i, 0) = XY[i].x;
        A.at<double>(i, 1) = XY[i].y;
        A.at<double>(i, 2) = 1.0;
        B.at<double>(i, 0) = XY[i].x * XY[i].x + XY[i].y * XY[i].y;
    }

    cv::Mat C;
    if (!cv::solve(A, B, C, cv::DECOMP_SVD)) {
        std::cerr << "Circle fitting failed." << std::endl;
        return false;
    }
    xc = C.at<double>(0, 0) / 2.0;
    yc = C.at<double>(1, 0) / 2.0;
    r = std::sqrt(C.at<double>(2, 0) + xc * xc + yc * yc);

    return true;
}

// Function to regularize a circle
Path regularize_circle(const Path& XY) {
    double xc, yc, r;
    if (!kasa_circle_fitting(XY, xc, yc, r)) {
        return XY;
    }

    Path regularized;
    int num_points = 100;
    for (int i = 0; i < num_points; ++i) {
        double theta = 2 * CV_PI * i / num_points;
        regularized.emplace_back(cv::Point2d(xc + r * std::cos(theta), yc + r * std::sin(theta)));
    }
    regularized.emplace_back(regularized.front());
    return regularized;
}

// Function to regularize an ellipse
Path regularize_ellipse(const Path& XY) {
    if (XY.size() < 5) {
        std::cerr << "Not enough points to fit an ellipse." << std::endl;
        return XY;
    }

    // Convert Path to cv::Mat
    cv::Mat points_mat(XY.size(), 1, CV_32FC2);
    for (size_t i = 0; i < XY.size(); ++i) {
        points_mat.at<cv::Point2f>(i, 0) = cv::Point2f(static_cast<float>(XY[i].x), static_cast<float>(XY[i].y));
    }

    cv::RotatedRect ellipse = cv::fitEllipse(points_mat);
    Path regularized;
    int num_points = 200;
    for (int i = 0; i < num_points; ++i) {
        double theta = 2 * CV_PI * i / num_points;
        double x = ellipse.center.x + (ellipse.size.width / 2.0) * std::cos(theta);
        double y = ellipse.center.y + (ellipse.size.height / 2.0) * std::sin(theta);
        regularized.emplace_back(cv::Point2d(x, y));
    }
    regularized.emplace_back(regularized.front());
    return regularized;
}

// Function to regularize a rectangle
Path regularize_rectangle(const Path& XY) {
    double x_min = XY[0].x, y_min = XY[0].y;
    double x_max = XY[0].x, y_max = XY[0].y;

    for (const auto& pt : XY) {
        if (pt.x < x_min) x_min = pt.x;
        if (pt.x > x_max) x_max = pt.x;
        if (pt.y < y_min) y_min = pt.y;
        if (pt.y > y_max) y_max = pt.y;
    }

    Path regularized = {
        cv::Point2d(x_min, y_min),
        cv::Point2d(x_max, y_min),
        cv::Point2d(x_max, y_max),
        cv::Point2d(x_min, y_max),
        cv::Point2d(x_min, y_min)
    };
    return regularized;
}

// Function to regularize a square
Path regularize_square(const Path& XY) {
    double x_min = XY[0].x, y_min = XY[0].y;
    double x_max = XY[0].x, y_max = XY[0].y;

    for (const auto& pt : XY) {
        if (pt.x < x_min) x_min = pt.x;
        if (pt.x > x_max) x_max = pt.x;
        if (pt.y < y_min) y_min = pt.y;
        if (pt.y > y_max) y_max = pt.y;
    }

    double side_length = std::min(x_max - x_min, y_max - y_min);

    Path regularized = {
        cv::Point2d(x_min, y_min),
        cv::Point2d(x_min + side_length, y_min),
        cv::Point2d(x_min + side_length, y_min + side_length),
        cv::Point2d(x_min, y_min + side_length),
        cv::Point2d(x_min, y_min)
    };
    return regularized;
}

// Function to regularize a star
Path regularize_star(const Path& XY, int num_points = 5) {
    // Calculate centroid
    cv::Point2d centroid(0.0, 0.0);
    for (const auto& pt : XY) {
        centroid += pt;
    }
    centroid *= (1.0 / XY.size());

    // Calculate outer and inner radii
    double outer_radius = 0.0, inner_radius = std::numeric_limits<double>::max();
    for (const auto& pt : XY) {
        double dist = cv::norm(pt - centroid);
        if (dist > outer_radius) outer_radius = dist;
        if (dist < inner_radius) inner_radius = dist;
    }

    Path regularized;
    double angle_step = CV_PI / num_points;
    for (int i = 0; i < 2 * num_points; ++i) {
        double angle = i * angle_step;
        double r = (i % 2 == 0) ? outer_radius : inner_radius;
        double x = centroid.x + r * std::cos(angle);
        double y = centroid.y + r * std::sin(angle);
        regularized.emplace_back(cv::Point2d(x, y));
    }
    regularized.emplace_back(regularized.front());
    return regularized;
}

// Function to regularize a straight line
Path regularize_straight_line(const Path& XY) {
    return { XY.front(), XY.back() };
}

// Function to check if contour is a square
bool is_square(const std::vector<cv::Point>& contour, double tolerance = 0.05) {
    if (contour.size() != 4) return false;

    cv::RotatedRect rect = cv::minAreaRect(contour);
    double width = rect.size.width;
    double height = rect.size.height;
    double aspect_ratio = std::min(width, height) / std::max(width, height);

    double contour_area = cv::contourArea(contour);
    double rect_area = width * height;
    double area_ratio = contour_area / rect_area;

    bool is_square_shape = std::abs(1.0 - aspect_ratio) < tolerance;
    bool is_filled = area_ratio > (1.0 - tolerance);

    return is_square_shape && is_filled;
}

// Function to check if contour is a star
bool is_star(const std::vector<cv::Point>& approx) {
    int num_vertices = approx.size();
    if (num_vertices > 8) {
        std::vector<double> angles;
        for (int i = 0; i < num_vertices; ++i) {
            cv::Point2d p1 = cv::Point2d(approx[i].x, approx[i].y);
            cv::Point2d p2 = cv::Point2d(approx[(i + 1) % num_vertices].x, approx[(i + 1) % num_vertices].y);
            cv::Point2d p3 = cv::Point2d(approx[(i + 2) % num_vertices].x, approx[(i + 2) % num_vertices].y);
            double angle = calculate_angle(p1, p2, p3);
            angles.push_back(angle);
        }

        int sharp_angles = 0;
        for (double angle : angles) {
            if (angle < 60.0 || angle > 300.0) {
                sharp_angles++;
            }
        }
        if (sharp_angles >= 5) {
            return true;
        }
    }
    return false;
}

// Function to detect shape type
std::string detect_shape(const Path& XY) {
    // Convert Path to contour (vector<cv::Point>)
    std::vector<cv::Point> contour;
    for (const auto& pt : XY) {
        contour.emplace_back(cv::Point(static_cast<int>(std::round(pt.x)),
                                       static_cast<int>(std::round(pt.y))));
    }

    double peri = cv::arcLength(contour, true);
    std::vector<cv::Point> approx;
    cv::approxPolyDP(contour, approx, 0.04 * peri, true);
    double area = cv::contourArea(contour);

    if (area < 100.0) {
        return " ";
    }

    // Check for straight line
    // For straight line, we handle separately
    // So here, assume it's not a straight line

    if (approx.size() == 4) {
        if (is_square(contour)) {
            return "Square";
        }
        else {
            return "Rectangle";
        }
    }
    else if (approx.size() > 5) {
        if (is_star(approx)) {
            return "Star";
        }

        double circularity = 4 * CV_PI * area / (peri * peri);

        if (circularity > 0.95) {
            return "Circle";
        }
        else if (circularity >= 0.8 && circularity < 0.95) {
            // Fit ellipse to check aspect ratio
            if (contour.size() < 5) return " ";
            cv::RotatedRect ellipse = cv::fitEllipse(contour);
            double major_axis = std::max(ellipse.size.width, ellipse.size.height);
            double minor_axis = std::min(ellipse.size.width, ellipse.size.height);
            if ((major_axis / minor_axis) > 1.05) {
                return "Ellipse";
            }
            else {
                return "Circle";
            }
        }
    }

    return " ";
}

// Function to detect and regularize shape
std::pair<std::string, Path> detect_and_regularize_shape(const Path& XY) {
    if (is_straight_line(XY)) {
        return { "Straight Line", regularize_straight_line(XY) };
    }

    std::string shape = detect_shape(XY);
    Path regularized_XY = XY;

    if (shape == "Circle") {
        regularized_XY = regularize_circle(XY);
    }
    else if (shape == "Ellipse") {
        regularized_XY = regularize_ellipse(XY);
    }
    else if (shape == "Rectangle") {
        regularized_XY = regularize_rectangle(XY);
    }
    else if (shape == "Square") {
        regularized_XY = regularize_square(XY);
    }
    else if (shape == "Star") {
        regularized_XY = regularize_star(XY);
    }

    return { shape, regularized_XY };
}

// Function to calculate symmetry lines
std::vector<std::pair<Point2D, Point2D>> calculate_symmetry_lines(const Path& XY, const std::string& shape_type) {
    std::vector<std::pair<Point2D, Point2D>> lines;
    if (shape_type == "Circle" || shape_type == "Square") {
        // Vertical and Horizontal lines through centroid
        cv::Point2d centroid(0.0, 0.0);
        for (const auto& pt : XY) {
            centroid += pt;
        }
        centroid *= (1.0 / XY.size());

        lines.emplace_back(std::make_pair(cv::Point2d(centroid.x, -10000), cv::Point2d(centroid.x, 10000)));
        lines.emplace_back(std::make_pair(cv::Point2d(-10000, centroid.y), cv::Point2d(10000, centroid.y)));
    }
    else if (shape_type == "Ellipse" || shape_type == "Rectangle") {
        // Horizontal line through centroid
        cv::Point2d centroid(0.0, 0.0);
        for (const auto& pt : XY) {
            centroid += pt;
        }
        centroid *= (1.0 / XY.size());

        lines.emplace_back(std::make_pair(cv::Point2d(-10000, centroid.y), cv::Point2d(10000, centroid.y)));
    }
    else if (shape_type == "Star") {
        // Multiple symmetry lines based on number of points
        int num_lines = 5;
        cv::Point2d centroid(0.0, 0.0);
        for (const auto& pt : XY) {
            centroid += pt;
        }
        centroid *= (1.0 / XY.size());

        for (int i = 0; i < num_lines; ++i) {
            double angle = 2 * CV_PI * i / num_lines;
            double x = centroid.x + std::cos(angle) * 10000;
            double y = centroid.y + std::sin(angle) * 10000;
            lines.emplace_back(std::make_pair(cv::Point2d(centroid.x, centroid.y), cv::Point2d(x, y)));
        }
    }
    return lines;
}

// Function to draw symmetry lines on an image
void draw_symmetry_lines(cv::Mat& img, const std::vector<std::pair<Point2D, Point2D>>& lines, const cv::Scalar& color = cv::Scalar(0, 0, 0)) {
    for (const auto& line : lines) {
        cv::line(img,
                 cv::Point(static_cast<int>(std::round(line.first.x)),
                           static_cast<int>(std::round(line.first.y))),
                 cv::Point(static_cast<int>(std::round(line.second.x)),
                           static_cast<int>(std::round(line.second.y))),
                 color, 1, cv::LINE_AA);
    }
}

// Function to write regularized shapes to CSV
void polylines2csv(const Paths& paths_XYs, const std::string& csv_path) {
    std::ofstream file(csv_path);
    if (!file.is_open()) {
        std::cerr << "Failed to open " << csv_path << " for writing." << std::endl;
        return;
    }

    for (size_t i = 0; i < paths_XYs.size(); ++i) {
        const Path& path = paths_XYs[i];
        for (size_t j = 0; j < path.size(); ++j) {
            file << i << "," << path[j].x << "," << path[j].y << "\n";
        }
    }
}

// Function to plot shapes using OpenCV
void plot_shapes(const Paths& path_XYs, const Paths& regularized_XYs) {
    // Determine the size of the image
    double x_min = std::numeric_limits<double>::max(), y_min = std::numeric_limits<double>::max();
    double x_max = std::numeric_limits<double>::lowest(), y_max = std::numeric_limits<double>::lowest();

    for (const auto& path : path_XYs) {
        for (const auto& pt : path) {
            if (pt.x < x_min) x_min = pt.x;
            if (pt.y < y_min) y_min = pt.y;
            if (pt.x > x_max) x_max = pt.x;
            if (pt.y > y_max) y_max = pt.y;
        }
    }

    int padding = 50;
    int width = static_cast<int>(std::ceil(x_max - x_min)) + 2 * padding;
    int height = static_cast<int>(std::ceil(y_max - y_min)) + 2 * padding;

    // Create two blank images for original and regularized shapes
    cv::Mat original_img(height, width, CV_8UC3, cv::Scalar(255, 255, 255));
    cv::Mat regularized_img(height, width, CV_8UC3, cv::Scalar(255, 255, 255));

    std::vector<cv::Scalar> colors = { cv::Scalar(0,0,255), cv::Scalar(0,255,0),
                                      cv::Scalar(255,0,0), cv::Scalar(255,255,0),
                                      cv::Scalar(255,0,255), cv::Scalar(0,255,255),
                                      cv::Scalar(0,0,0) };

    std::cout << "\033[1mShapes Detected and Regularized:\033[0m\n";

    for (size_t i = 0; i < path_XYs.size(); ++i) {
        const Path& orig_path = path_XYs[i];
        const Path& reg_path = regularized_XYs[i];
        cv::Scalar color = colors[i % colors.size()];

        // Shift points to fit in image with padding
        std::vector<cv::Point> orig_contour;
        for (const auto& pt : orig_path) {
            orig_contour.emplace_back(cv::Point(static_cast<int>(std::round(pt.x - x_min)) + padding,
                                              static_cast<int>(std::round(pt.y - y_min)) + padding));
        }
        // Draw original shape
        if (orig_contour.size() > 1) {
            cv::polylines(original_img, orig_contour, false, color, 2);
        }

        // Regularized shape
        std::vector<cv::Point> reg_contour;
        for (const auto& pt : reg_path) {
            reg_contour.emplace_back(cv::Point(static_cast<int>(std::round(pt.x - x_min)) + padding,
                                             static_cast<int>(std::round(pt.y - y_min)) + padding));
        }
        if (reg_contour.size() > 1) {
            cv::polylines(regularized_img, reg_contour, false, color, 2);
        }

        // Detect shape and regularize
        auto [shape_type, regularized_shape] = detect_and_regularize_shape(orig_path);

        // Print detected shape
        std::cout << "- " << shape_type << "\n";

        // Calculate and draw symmetry lines
        std::vector<std::pair<Point2D, Point2D>> symmetry_lines = calculate_symmetry_lines(regularized_shape, shape_type);
        draw_symmetry_lines(regularized_img, symmetry_lines, cv::Scalar(0, 0, 0));

        // Print number of symmetry lines
        if (shape_type == "Ellipse" || shape_type == "Rectangle") {
            std::cout << "  - " << shape_type << " has 1 line of symmetry.\n";
        }
        else if (shape_type == "Circle" || shape_type == "Square") {
            std::cout << "  - " << shape_type << " has 2 lines of symmetry.\n";
        }
        else if (shape_type == "Star") {
            std::cout << "  - " << shape_type << " has 5 lines of symmetry.\n";
        }
        else if (shape_type == "Straight Line") {
            std::cout << "  - " << shape_type << " has infinite lines of symmetry.\n";
        }
        else {
            std::cout << "  - No specific symmetry lines defined for " << shape_type << ".\n";
        }
    }

    // Display images
    cv::imshow("Original Shapes", original_img);
    cv::imshow("Regularized Shapes with Symmetry Lines", regularized_img);
    cv::waitKey(0);
}

// Main function
int main() {
    std::string input_csv = "/Users/kartik/Desktop/CLG/MINOR/problems/isolated.csv"; // Update the path as needed
    std::string output_csv = "regularized_shapes.csv";

    // Read CSV
    Paths path_XYs = read_csv(input_csv);
    if (path_XYs.empty()) {
        std::cerr << "No paths found in the CSV file." << std::endl;
        return -1;
    }

    // Regularize shapes
    Paths regularized_XYs;
    std::vector<std::string> detected_shapes; // To store detected shapes

    for (const auto& path : path_XYs) {
        auto [shape, reg_XY] = detect_and_regularize_shape(path);
        regularized_XYs.push_back(reg_XY);
        detected_shapes.push_back(shape);
    }

    // Plot shapes and detect symmetry
    plot_shapes(path_XYs, regularized_XYs);

    // Write regularized shapes to CSV
    polylines2csv(regularized_XYs, output_csv);
    std::cout << "Regularized shapes have been written to " << output_csv << std::endl;

    return 0;
}



it's govomg blank output
base) kartik@Kartiks-MacBook-Air MINOR % g++ -std=c++17 -o shape_regularizer regularisation.cpp `pkg-config --cflags --libs opencv4`
(base) kartik@Kartiks-MacBook-Air MINOR %