sudo apt-get install libboost-dev

#include <iostream>
#include <fstream>
#include <vector>
#include <string>
#include <boost/geometry.hpp>
#include <boost/geometry/geometries/point.hpp>
#include <boost/geometry/geometries/box.hpp>
#include <boost/geometry/geometries/register/point.hpp>
#include <boost/geometry/index/rtree.hpp>

namespace bg = boost::geometry;
namespace bgi = boost::geometry::index;

// Define the point type
typedef bg::model::point<double, 2, bg::cs::cartesian> point;
typedef std::pair<point, unsigned> value;

// Define the GeoJSON feature type
struct GeoJSONFeature {
    std::string type;
    std::map<std::string, std::string> properties;
    point geometry;
};

int main() {
    // Read CSV file (File A)
    std::ifstream csvFile("file_a.csv");
    std::vector<point> pointsA;

    while (csvFile) {
        double lon, lat;
        char comma;
        csvFile >> lon >> comma >> lat;

        if (csvFile)
            pointsA.emplace_back(lon, lat);
    }

    // Read GeoJSON file (File B) and build spatial index
    std::ifstream geojsonFile("file_b.geojson");
    std::vector<GeoJSONFeature> featuresB;
    bgi::rtree<value, bgi::rstar<16>> rtreeB;

    while (geojsonFile) {
        GeoJSONFeature feature;
        geojsonFile >> feature.type;

        if (feature.type == "Feature") {
            geojsonFile >> std::ws; // Skip whitespace

            char c;
            while (geojsonFile.get(c) && c != '}') {
                if (c == '"') {
                    std::string key, value;
                    std::getline(geojsonFile, key, '"');
                    geojsonFile >> c >> value;
                    feature.properties[key] = value;
                }
            }

            geojsonFile >> std::ws; // Skip whitespace
            geojsonFile.ignore(2); // Ignore "type" and ","
            geojsonFile >> feature.geometry;

            featuresB.push_back(feature);
            rtreeB.insert({ feature.geometry, featuresB.size() - 1 });
        }
    }

    // Define query distance (adjust as needed)
    double queryDistance = 0.1;

    // Find nearby points and output to result GeoJSON file
    std::ofstream resultGeoJSON("result.geojson");

    for (const auto& pointA : pointsA) {
        // Create a bounding box for the query
        bg::model::box<point> queryBox;
        bg::envelope(pointA, queryBox);
        bg::buffer(queryBox, queryBox, queryDistance);

        // Query the spatial index for nearby features
        std::vector<value> result;
        rtreeB.query(bgi::intersects(queryBox), std::back_inserter(result));

        // Output nearby features to result GeoJSON file
        for (const auto& match : result) {
            const GeoJSONFeature& nearbyFeature = featuresB[match.second];

            resultGeoJSON << "{\n";
            resultGeoJSON << "  \"type\": \"Feature\",\n";
            resultGeoJSON << "  \"properties\": {\n";

            for (const auto& property : nearbyFeature.properties) {
                resultGeoJSON << "    \"" << property.first << "\": \"" << property.second << "\",\n";
            }

            resultGeoJSON << "  },\n";
            resultGeoJSON << "  \"geometry\": {\n";
            resultGeoJSON << "    \"type\": \"Point\",\n";
            resultGeoJSON << "    \"coordinates\": [" << bg::get<0>(nearbyFeature.geometry) << ", " << bg::get<1>(nearbyFeature.geometry) << "]\n";
            resultGeoJSON << "  }\n";
            resultGeoJSON << "}\n";
        }
    }

    std::cout << "Processing completed. Results saved in result.geojson." << std::endl;

    return 0;
}