plotter.hpp

Go to the documentation of this file.

//
// Created by domenico on 11/27/25.
//
 
#ifndef MONTECARLO_1_PLOTTER_HPP
#define MONTECARLO_1_PLOTTER_HPP
 
#include <fstream>
#include <iostream>
#include <string>
#include <sstream>
#include <vector>
#include <cstdlib> // std::system
#include <filesystem>
#include "../domains/integration_domain.hpp"
 
namespace mc::utils {
 
inline void closeGnuplotWindows() {
    std::system("pkill -f gnuplot > /dev/null 2>&1");
}
 
inline std::string formatTitle(std::string name) {
    std::string target = "_samples";
    size_t pos = name.find(target);
    if (pos != std::string::npos) name.erase(pos, target.length());
    for (auto &c : name) if (c == '_') c = ' ';
    return name;
}
 
template <size_t dim>
inline void createGnuplotScript(const std::string& tempRawDataFile,
                                const mc::domains::IntegrationDomain<dim>& domain,
                                size_t currentSamples) {
 
    if (dim > 3) return; // Cannot plot 4D+ geometry easily
 
    // Prepare unique filenames
    std::string baseName = tempRawDataFile;
    size_t lastDot = baseName.find_last_of(".");
    if (lastDot != std::string::npos) baseName = baseName.substr(0, lastDot);
    size_t lastSlash = baseName.find_last_of("/\\");
    if (lastSlash != std::string::npos) baseName = baseName.substr(lastSlash + 1);
 
    std::string uniqueID = baseName + "_geom_" + std::to_string(currentSamples);
    std::string uniqueDataFile = uniqueID + ".dat";
    std::string scriptName = "plot_" + uniqueID + ".gp";
 
    // Read temp raw data and write unique file with Inside/Outside status
    std::ifstream inFile(tempRawDataFile);
    if (!inFile.is_open()) return;
 
    std::vector<mc::geom::Point<dim>> points;
    std::string line;
    while (std::getline(inFile, line)) {
        if (line.empty()) continue;
        std::stringstream ss(line);
        mc::geom::Point<dim> p;
        for (size_t i = 0; i < dim; ++i) ss >> p[i];
        points.push_back(p);
    }
    inFile.close();
 
    std::ofstream outFile(uniqueDataFile);
    for (const auto& p : points) {
        for (size_t i = 0; i < dim; ++i) outFile << p[i] << " ";
        bool isInside = domain.isInside(p);
        outFile << (isInside ? 1 : 0) << "\n";
    }
    outFile.close();
 
    // Create Gnuplot script
    std::ofstream gp(scriptName);
    if (!gp.is_open()) return;
 
    gp << "set grid\n";
    gp << "set title 'Domain Geometry (" << dim << "D) - N=" << currentSamples << "'\n";
    gp << "set key outside\n";
 
    // Set Axis Ranges
    auto bounds = domain.getBounds();
    double margin = 0.1;
    auto setRange = [&](int axisIdx, std::string axisName) {
        if (axisIdx < dim) {
            double min = bounds[axisIdx].first;
            double max = bounds[axisIdx].second;
            double span = max - min;
            gp << "set " << axisName << "range [" << (min - span * margin) << ":" << (max + span * margin) << "]\n";
        }
    };
    setRange(0, "x");
    if (dim >= 2) setRange(1, "y");
    if (dim >= 3) setRange(2, "z");
 
    int statusCol = dim + 1; // The added status column index (1-based)
 
    if (dim == 1) {
        gp << "set xlabel 'X'\n unset ytics\n set yrange [-1:1]\n";
        gp << "plot '" << uniqueDataFile << "' u 1:($" << statusCol << "==0?0:1/0) w p pt 7 ps 0.4 lc rgb 'blue' t 'Out', \\\n";
        gp << "     '" << uniqueDataFile << "' u 1:($" << statusCol << "==1?0:1/0) w p pt 7 ps 0.4 lc rgb 'red' t 'In'\n";
    }
    else if (dim == 2) {
        gp << "set size square\n set xlabel 'X'\n set ylabel 'Y'\n";
        gp << "plot '" << uniqueDataFile << "' u 1:($" << statusCol << "==0?$2:1/0) w p pt 7 ps 0.4 lc rgb 'blue' t 'Out', \\\n";
        gp << "     '" << uniqueDataFile << "' u 1:($" << statusCol << "==1?$2:1/0) w p pt 7 ps 0.4 lc rgb 'red' t 'In'\n";
    }
    else if (dim == 3) {
        gp << "set view equal xyz\n set xlabel 'X'\n set ylabel 'Y'\n set zlabel 'Z'\n set view 60, 30\n";
        gp << "splot '" << uniqueDataFile << "' u 1:2:($" << statusCol << "==0?$3:1/0) w p pt 7 ps 0.4 lc rgb 'blue' t 'Out', \\\n";
        gp << "      '" << uniqueDataFile << "' u 1:2:($" << statusCol << "==1?$3:1/0) w p pt 7 ps 0.4 lc rgb 'red' t 'In'\n";
    }
 
    gp << "pause mouse close\n"; // Keep window open
    gp.close();
 
    // Execute Gnuplot in background
    std::string command = "gnuplot " + scriptName + " > /dev/null 2>&1 &";
    std::system(command.c_str());
}
 
template <size_t dim, typename Func>
inline void createFunctionGnuplotScript(const std::string& tempRawDataFile,
                                        const mc::domains::IntegrationDomain<dim>& domain,
                                        const Func& func,
                                        size_t currentSamples) {
    if (dim > 3) return;
 
    std::string baseName = tempRawDataFile;
    size_t lastDot = baseName.find_last_of(".");
    if (lastDot != std::string::npos) baseName = baseName.substr(0, lastDot);
    size_t lastSlash = baseName.find_last_of("/\\");
    if (lastSlash != std::string::npos) baseName = baseName.substr(lastSlash + 1);
 
    // Unique filename for function values
    std::string uniqueID = baseName + "_func_" + std::to_string(currentSamples);
    std::string uniqueDataFile = uniqueID + ".dat";
    std::string scriptName = "plot_" + uniqueID + ".gp";
 
    // Read raw data, evaluate f(x) ONLY for points inside the domain, and write to file
    std::ifstream inFile(tempRawDataFile);
    if (!inFile.is_open()) return;
 
    std::ofstream outFile(uniqueDataFile);
 
    std::string line;
    while (std::getline(inFile, line)) {
        if (line.empty()) continue;
        std::stringstream ss(line);
        mc::geom::Point<dim> p;
        for (size_t i = 0; i < dim; ++i) ss >> p[i];
 
        if (domain.isInside(p)) {
            double val = func(p);
            for(size_t i=0; i<dim; ++i) outFile << p[i] << " ";
            outFile << val << "\n";
        }
    }
    inFile.close();
    outFile.close();
 
    // Generate Gnuplot script
    std::ofstream gp(scriptName);
    if (!gp.is_open()) return;
 
    gp << "set grid\n";
    gp << "set title 'Function Value (" << dim << "D) - N=" << currentSamples << "'\n";
 
    // Configure Color Palette logic
    if (dim == 3) {
        // Only 3D Domain (4D total visualization) uses heatmap
        gp << "set palette rgbformulae 33,13,10\n";
        gp << "set colorbox\n";
    } else {
        // 1D and 2D Domains use solid green
        gp << "unset colorbox\n";
    }
 
    int valCol = dim + 1; // Column containing function value
 
    if (dim == 1) {
        gp << "set xlabel 'X'\n set ylabel 'f(X)'\n";
        // 1D Domain -> 2D Plot (x, f(x)) -> Solid Green
        gp << "plot '" << uniqueDataFile << "' u 1:" << valCol << " w p pt 7 ps 0.4 lc rgb 'green' t 'f(x)'\n";
    }
    else if (dim == 2) {
        gp << "set size square\n set view 60, 30\n";
        gp << "set xlabel 'X'\n set ylabel 'Y'\n set zlabel 'f(X,Y)'\n";
        // 2D Domain -> 3D Plot (x, y, z=f(x,y)) -> Solid Green
        gp << "splot '" << uniqueDataFile << "' u 1:2:" << valCol << " w p pt 7 ps 0.4 lc rgb 'green' t 'f(x,y)'\n";
    }
    else if (dim == 3) {
        gp << "set view equal xyz\n set view 60, 30\n";
        gp << "set xlabel 'X'\n set ylabel 'Y'\n set zlabel 'Z'\n";
        // 3D Domain -> 4D Visualization (x, y, z, color=f(x,y,z)) -> Palette
        gp << "splot '" << uniqueDataFile << "' u 1:2:3:" << valCol << " w p pt 7 ps 0.4 lc palette t 'f(x,y,z)'\n";
    }
 
    gp << "pause mouse close\n";
    gp.close();
 
    std::string command = "gnuplot " + scriptName + " > /dev/null 2>&1 &";
    std::system(command.c_str());
}
 
template <typename Func>
inline void saveFunctionGrid(const std::string& filename, const Func& func,
                             double x_min, double x_max, double y_min, double y_max,
                             int resolution = 100) {
    // Extract directory from filename
    std::string dir = "./pso_frames";
    if (filename.find("ga_") != std::string::npos) {
        dir = "./ga_frames";
    }
    
    // Create directory if it doesn't exist
    try {
        std::filesystem::create_directories(dir);
    } catch (const std::exception&) {}
    
    std::string fullpath = dir + "/" + filename;
    std::ofstream out(fullpath);
    if (!out.is_open()) return;
 
    double dx = (x_max - x_min) / resolution;
    double dy = (y_max - y_min) / resolution;
 
    for (int i = 0; i <= resolution; ++i) {
        double x = x_min + i * dx;
        for (int j = 0; j <= resolution; ++j) {
            double y = y_min + j * dy;
            std::vector<double> p = {x, y};
            double val = func(p);
            // Write format: X Y Value
            out << x << " " << y << " " << val << "\n";
        }
        out << "\n"; // Blank line required for Gnuplot pm3d mode
    }
    out.close();
}
 
template <typename ParticleT>
inline void saveSwarmFrame(const std::string& basename, size_t iteration, const std::vector<ParticleT>& swarm) {
    // Determine directory based on basename
    std::string dir = "./pso_frames";
    if (basename.find("ga_") != std::string::npos) {
        dir = "./ga_frames";
    }
    
    // Create directory if it doesn't exist
    try {
        std::filesystem::create_directories(dir);
    } catch (const std::exception&) {}
    
    std::string filename = dir + "/" + basename + "_iter_" + std::to_string(iteration) + ".dat";
    std::ofstream out(filename);
    if (!out.is_open()) return;
 
    for (const auto& p : swarm) {
        if (p.position.size() >= 3) {
            // Write X Y Z
            out << p.position[0] << " " << p.position[1] << " " << p.position[2] << "\n";
        } else if (p.position.size() == 2) {
            // Fallback for 2D (Z=0)
            out << p.position[0] << " " << p.position[1] << " 0.0\n";
        }
    }
    out.close();
}
 
inline void createPSOAnimationScript(const std::string& scriptName,
                                     const std::string& gridFile,
                                     const std::string& swarmBasename,
                                     size_t max_iter,
                                     const std::string& title) {
    // Determine directory based on swarmBasename
    std::string dir = "./pso_frames";
    if (swarmBasename.find("ga_") != std::string::npos) {
        dir = "./ga_frames";
    }
    
    std::ofstream gp(scriptName);
    if (!gp.is_open()) return;
 
    gp << "set title '" << title << "'\n";
    gp << "set view map\n"; // Top-down view (Heatmap)
    gp << "set dgrid3d\n";  // Enable 3D grid data processing
    gp << "set pm3d interpolate 0,0\n"; // Smooth color interpolation
    gp << "set palette rgbformulae 33,13,10\n"; // Rainbow palette
    gp << "unset key\n";
    gp << "set size square\n";
 
    // Gnuplot Loop for Animation
    gp << "do for [i=0:" << (max_iter-1) << "] {\n";
    gp << "    set title sprintf('" << title << " - Iter: %d', i)\n";
    gp << "    plot '" << dir << "/" << gridFile << "' with image, \\\n";
    gp << "         sprintf('" << dir << "/" << swarmBasename << "_iter_%d.dat', i) u 1:2 with points pt 7 ps 1.5 lc rgb 'white'\n";
    gp << "    pause 0.1\n"; // Delay between frames (0.1 seconds)
    gp << "}\n";
    gp << "pause mouse close\n"; // Keep window open until clicked
    gp.close();
 
    // Execute Gnuplot in background
    std::string command = "gnuplot " + scriptName + " > /dev/null 2>&1 &";
    std::system(command.c_str());
}
 
template <typename Func>
inline void saveFunctionSlices3D(const std::string& filename, const Func& func,
                                 double min, double max, int resolution = 50) {
    // Extract directory from filename
    std::string dir = "./pso_frames";
    if (filename.find("ga_") != std::string::npos) {
        dir = "./ga_frames";
    }
    
    // Create directory if it doesn't exist
    try {
        std::filesystem::create_directories(dir);
    } catch (const std::exception&) {}
    
    std::string fullpath = dir + "/" + filename;
    std::ofstream out(fullpath);
    if (!out.is_open()) return;
 
    double step = (max - min) / resolution;
 
    // 1. XY Plane (at Z = min) - "Floor"
    for (int i = 0; i <= resolution; ++i) {
        double x = min + i * step;
        for (int j = 0; j <= resolution; ++j) {
            double y = min + j * step;
            double z = min; // Bottom
            double val = func({x, y, z});
            out << x << " " << y << " " << z << " " << val << "\n";
        }
        out << "\n";
    }
    out << "\n\n";
 
    // 2. XZ Plane (at Y = min) - "Back Right Wall" (FIXED: was max)
    // Changing to Y=min moves the wall to the back, opening the view.
    for (int i = 0; i <= resolution; ++i) {
        double x = min + i * step;
        for (int k = 0; k <= resolution; ++k) {
            double z = min + k * step;
            double y = min; // Fixed at back (min)
            double val = func({x, y, z});
            out << x << " " << y << " " << z << " " << val << "\n";
        }
        out << "\n";
    }
    out << "\n\n";
 
    // 3. YZ Plane (at X = min) - "Back Left Wall"
    for (int j = 0; j <= resolution; ++j) {
        double y = min + j * step;
        for (int k = 0; k <= resolution; ++k) {
            double z = min + k * step;
            double x = min; // Fixed at left (min)
            double val = func({x, y, z});
            out << x << " " << y << " " << z << " " << val << "\n";
        }
        out << "\n";
    }
    out.close();
}
 
inline void createPSOAnimationScript3D(const std::string& scriptName,
                                       const std::string& slicesFile,
                                       const std::string& swarmBasename,
                                       size_t max_iter,
                                       const std::string& title,
                                       double min_bound, double max_bound) {
    // Determine subdirectory based on filename
    std::string dir = "pso_frames";
    if (swarmBasename.find("ga_") != std::string::npos) {
        dir = "ga_frames";
    }
    
    // Ensure subdirectory exists
    try {
        std::filesystem::create_directories(dir);
    } catch (const std::exception& e) {
        std::cerr << "Error creating directory " << dir << ": " << e.what() << std::endl;
    }
    
    // Prepend directory to file paths
    std::string fullSlicesFile = dir + "/" + slicesFile;
    std::string fullSwarmBasename = dir + "/" + swarmBasename;
    
    std::ofstream gp(scriptName);
    if (!gp.is_open()) return;
 
    gp << "set title '" << title << "'\n";
 
    // View from the "open" side of the cube towards the back walls
    gp << "set view 60, 120\n";
 
    gp << "set grid\n";
    gp << "set xtics\n";
    gp << "set ytics\n";
    gp << "set ztics\n";
 
    // Use 'explicit' to prevent pm3d from interfering with points
    gp << "set pm3d explicit\n";
    gp << "set palette rgbformulae 33,13,10\n";
    gp << "unset colorbox\n";
 
    // Fixed range
    gp << "set xrange [" << min_bound << ":" << max_bound << "]\n";
    gp << "set yrange [" << min_bound << ":" << max_bound << "]\n";
    gp << "set zrange [" << min_bound << ":" << max_bound << "]\n";
 
    // Animation Loop
    gp << "do for [i=0:" << (max_iter-1) << "] {\n";
    gp << "    set title sprintf('" << title << " - Iter: %d', i)\n";
 
    // Draw slices first (background), then particles (foreground)
    gp << "    splot '" << fullSlicesFile << "' u 1:2:3:4 with pm3d nocontour title '', \\\n";
    gp << "          sprintf('" << fullSwarmBasename << "_iter_%d.dat', i) u 1:2:3 with points pt 7 ps 1.5 lc rgb 'white' title 'Particles'\n";
 
    gp << "    pause 0.1\n";
    gp << "}\n";
    gp << "pause mouse close\n";
    gp.close();
 
    // Execute visibly
    std::cout << "Executing Gnuplot 3D script: " << scriptName << std::endl;
    std::string command = "gnuplot -p " + scriptName;
    int ret = std::system(command.c_str());
    if (ret != 0) {
        std::cerr << "GNUPLOT ERROR (Code " << ret << "): Ensure gnuplot is installed." << std::endl;
    }
}
 
inline void createDroneVisualizationScript(const std::string& scriptName,
                                          const std::string& geometryFile,
                                          const std::string& title = "Drone Arm Domain Geometry") {
    std::ofstream gp(scriptName);
    if (!gp.is_open()) {
        std::cerr << "Could not create gnuplot script: " << scriptName << std::endl;
        return;
    }
 
    gp << "#!/usr/bin/gnuplot\n";
    gp << "# Visualization script for drone arm domain geometry\n";
    gp << "# Auto-generated by drone_optimization\n\n";
    
    gp << "set title '" << title << "'\n";
    gp << "set xlabel 'X'\n";
    gp << "set ylabel 'Y'\n";
    gp << "set zlabel 'Z'\n";
    
    gp << "# Equal aspect ratio to prevent distortion\n";
    gp << "set view equal xyz\n";
    gp << "set view 60, 30\n";
    gp << "set grid\n\n";
    
    gp << "# Data file location\n";
    gp << "datafile = '" << geometryFile << "'\n\n";
    
    gp << "# Plot arm (blue) and motor (red) with smaller points for volume effect\n";
    // Use stringcolumn(1) since column 1 contains strings (arm/motor)
    gp << "splot datafile u (stringcolumn(1) eq \"arm\" ? $2 : 1/0):3:4 with points pt 7 ps 0.5 lc rgb '#3366cc' title 'Arm', \\\n";
    gp << "      datafile u (stringcolumn(1) eq \"motor\" ? $2 : 1/0):3:4 with points pt 7 ps 0.5 lc rgb '#dc3912' title 'Motor'\n\n";
    
    gp << "pause mouse close\n";
    gp.close();
    
    std::cout << "Gnuplot visualization script created: " << scriptName << std::endl;
    std::cout << "To visualize: gnuplot -persist " << scriptName << std::endl;
}
 
} // namespace mc::utils
 
#endif //MONTECARLO_1_PLOTTER_HPP