Разделяй и властвуй, поиск пары ближайших точек

#include <vector>
#include <ostream>
#include <cmath>
#include "point.h"
#include <utility>
 
// Find the closest pair of points from given points.
std::pair<Point, Point> closest_pair(const std::vector<Point>& points);
std::pair<Point, Point> closest_pair1(std::vector<Point>& points);

#include "point.h"
#include <utility>
#include <cmath>
using namespace std;
 
 
namespace {
    static const double EPS = 1e-8;
 
    bool equal(double a, double b) {
        return abs(a - b) < EPS;
    }
}
 
double Point::distance(const Point& p) const {
    return sqrt((x - p.x) * (x - p.x) + (y - p.y) * (y - p.y));
}
 
bool Point::operator==(const Point& p) const {
    return equal(x, p.x) && equal(y, p.y);
}
 
bool Point::operator!=(const Point& p) const {
    return !equal(x, p.x) || !equal(y, p.y);
}
 
bool Point::operator<(const Point& p) const {
    return (x < p.x) || (equal(x, p.x) && y < p.y);
}
 
bool Point::operator<=(const Point& p) const {
    return (x < p.x) || (equal(x, p.x) && ((y < p.y) || equal(y, p.y)));
}
 
ostream& operator<<(ostream& os, const Point& p) {
    os << "(" << p.x << ", " << p.y << ")";
    return os;
}

#include <ostream>
#include <utility>
class Point {
public:
    Point() {}
    Point(float xx, float yy) : x(xx), y(yy) {}
 
    double distance(const Point& p) const;
 
    bool operator==(const Point& p) const;
    bool operator!=(const Point& p) const;
    bool operator<(const Point& p) const;
    bool operator<=(const Point& p) const;
 
public:
    double x{ 0.0 };
    double y{ 0.0 };
};
 
std::ostream& operator<<(std::ostream& os, const Point& p);

#include "Header.h"
#include <stdexcept>
#include <algorithm>
#include <ctime>
#include <iostream>
#include <chrono>
#include <utility>
using namespace std;
 
bool pointxcmp(const Point& p1, const Point& p2) // Используется для сортировки точек в порядке возрастания по координате x
{
    return p1.x < p2.x;
}
bool pointycmp(const Point& p1,const Point& p2) // Используется для сортировки точек в порядке возрастания по координате y
{
    return p1.y < p2.y;
}
 
std::pair<Point, Point> closest_pair( std::vector<Point>& points) {
    int i, j;
    double mindist = INFINITY;
    float d;
    int leftindex = 0;
    int rightindex = points.size() - 1;
    if (points.size() < 2) {
        throw invalid_argument("Not enough points");
    }
    auto minimum = make_pair(points[points.size() - 2], points[points.size() - 1]);
    for (i = leftindex; i <= rightindex; i++)
        for (j = i + 1; j <= rightindex; j++)
        {
            d = points[i].distance(points[j]);
            if (d < mindist) {
                mindist = d;
 
                minimum = make_pair(points[i], points[j]);
            }
        }
    return minimum;
}
pair<Point, Point> closest_pair_between(vector<Point>& PLeft, vector<Point>& PRight, float d)
{
 
    int i;
    double Xm;
    vector<Point> PStripe;
    Xm = (PLeft[PLeft.size() - 1].x + PRight[0].x) / 2; // граница между PLeft и PRight по x - координате
    for (i = 0; i < PLeft.size(); i++) {
        if (fabs(PLeft[i].x - Xm) < d)
            PStripe.push_back(PLeft[i]);
    }
    for (i = 0; i < PRight.size(); i++) {
        if (fabs(PRight[i].x - Xm) < d)
            PStripe.push_back(PRight[i]);
    }
    sort(PStripe.begin(), PStripe.end(), pointycmp);// Сортировать по координате y от маленького к большом
    pair<Point, Point> k;
    k.first = PStripe[0];
    k.second = PStripe[1];
    double s = k.second.distance(k.first);
    double t;
    for (i = 0; i < PStripe.size(); i++)
        for (int j = i + 1; j < PStripe.size() && (PStripe[j].y - PStripe[i].y) < d; j++) {
            if ((t = PStripe[i].distance(PStripe[j])) < s) {
                s = t;
                k.first = PStripe[i];
                k.second = PStripe[j];
            }
        }
    return k;
}
 
std::pair<Point, Point> closest_pair1( std::vector<Point>& points) {
    vector<Point> b;
    int leftindex = 0;
    int rightindex = points.size() - 1;
    vector<Point> PLeft, PRight, b1;
    int i, midindex;
    double pl, pr;
    double d3 = INFINITY, d;
 
    if (points.size() <= 3) {
        return closest_pair(points);
    }
    
    else {
        std::pair<Point, Point> tmp,tmp2,pb;
 
        sort(points.begin(), points.end(), pointxcmp); // Сортировать по координате x от малого к большему
 
        midindex = (leftindex + rightindex) / 2; // Находим среднюю позицию
        for (i = 0; i <=midindex; i++) // делим набор средних точек b на две части
                PLeft.push_back(b[i]);
        for ( ;i <points.size(); i++)
                PRight.push_back(b[i]);
        tmp = closest_pair1(PLeft);
        pl =  tmp.first.distance(tmp.second);
        tmp2 = closest_pair1(PRight);
        pr = tmp2.first.distance(tmp2.second);
        d = min(pl,pr); // текущее минимальное расстояние d 
        pb = closest_pair_between(PLeft, PRight, d); // пара ближайших точек между половинами PLeft и PRight.
    }
}
 
int main() {
    vector<Point> points;
    srand(time(NULL));
    int N = 5;
    for (int i = 0; i < N; i++) {
        points.push_back(Point(rand() % 20, rand() % 20));
        cout << points[i]<<" ";
    }
    cout << endl;
    auto t1 = std::chrono::high_resolution_clock::now();
    auto minimum = closest_pair1(points);
    auto t2 = std::chrono::high_resolution_clock::now();
    auto time = std::chrono::duration<double>(t2 - t1).count();
    
    cout << minimum.first<<endl;
    cout << minimum.second;
    std::cout << "Time: " << time << " sec." << std::endl;
 }

    vector<Point> b;
    ...
        for (i = 0; i <=midindex; i++) // делим набор средних точек b на две части
                PLeft.push_back(b[i]);
        for ( ;i <points.size(); i++)
                PRight.push_back(b[i]);