Skip to main content

Advanced Geometry

Explore advanced geometric algorithms for complex spatial problems and computational geometry.

Sweep Line Algorithm​

Event-Driven Sweep Line​

struct Event {
    double x;
    int type; // 0: start, 1: end
    int id;
    
    Event(double x, int type, int id) : x(x), type(type), id(id) {}
    
    bool operator<(const Event& other) const {
        if (abs(x - other.x) < 1e-9) {
            return type < other.type; // Start events before end events
        }
        return x < other.x;
    }
};

// Sweep line for line segment intersection
vector<Point> sweepLineIntersection(vector<LineSegment>& segments) {
    vector<Event> events;
    vector<Point> intersections;
    
    // Create events for each segment
    for (int i = 0; i < segments.size(); i++) {
        double x1 = segments[i].p1.x;
        double x2 = segments[i].p2.x;
        
        if (x1 > x2) {
            swap(x1, x2);
        }
        
        events.push_back(Event(x1, 0, i)); // Start event
        events.push_back(Event(x2, 1, i)); // End event
    }
    
    sort(events.begin(), events.end());
    
    set<int> activeSegments;
    
    for (const Event& event : events) {
        if (event.type == 0) {
            // Start event - add segment to active set
            activeSegments.insert(event.id);
            
            // Check intersections with other active segments
            for (int id : activeSegments) {
                if (id != event.id) {
                    Point intersection = lineSegmentIntersection(segments[event.id], segments[id]);
                    if (intersection.x != INFINITY) {
                        intersections.push_back(intersection);
                    }
                }
            }
        } else {
            // End event - remove segment from active set
            activeSegments.erase(event.id);
        }
    }
    
    return intersections;
}

Sweep Line for Rectangle Union​

struct Rectangle {
    double x1, y1, x2, y2;
    
    Rectangle(double x1, double y1, double x2, double y2) 
        : x1(x1), y1(y1), x2(x2), y2(y2) {}
};

// Calculate union area of rectangles using sweep line
double rectangleUnionArea(vector<Rectangle>& rectangles) {
    vector<Event> events;
    
    // Create events for each rectangle
    for (int i = 0; i < rectangles.size(); i++) {
        events.push_back(Event(rectangles[i].x1, 0, i)); // Start
        events.push_back(Event(rectangles[i].x2, 1, i)); // End
    }
    
    sort(events.begin(), events.end());
    
    double totalArea = 0;
    double prevX = 0;
    set<int> activeRectangles;
    
    for (const Event& event : events) {
        if (!activeRectangles.empty()) {
            double width = event.x - prevX;
            double height = calculateActiveHeight(rectangles, activeRectangles);
            totalArea += width * height;
        }
        
        if (event.type == 0) {
            activeRectangles.insert(event.id);
        } else {
            activeRectangles.erase(event.id);
        }
        
        prevX = event.x;
    }
    
    return totalArea;
}

double calculateActiveHeight(vector<Rectangle>& rectangles, set<int>& active) {
    vector<pair<double, int>> yEvents;
    
    for (int id : active) {
        yEvents.push_back({rectangles[id].y1, 1}); // Start
        yEvents.push_back({rectangles[id].y2, -1}); // End
    }
    
    sort(yEvents.begin(), yEvents.end());
    
    double height = 0;
    int count = 0;
    double prevY = 0;
    
    for (const auto& event : yEvents) {
        if (count > 0) {
            height += event.first - prevY;
        }
        count += event.second;
        prevY = event.first;
    }
    
    return height;
}

Closest Pair of Points​

Divide and Conquer Approach​

// Find closest pair of points using divide and conquer
pair<Point, Point> closestPair(vector<Point>& points) {
    int n = points.size();
    if (n < 2) return {Point(), Point()};
    
    // Sort points by x-coordinate
    sort(points.begin(), points.end(), [](const Point& a, const Point& b) {
        return a.x < b.x;
    });
    
    return closestPairRecursive(points, 0, n - 1);
}

pair<Point, Point> closestPairRecursive(vector<Point>& points, int left, int right) {
    if (right - left <= 3) {
        // Base case: brute force for small number of points
        return bruteForceClosestPair(points, left, right);
    }
    
    int mid = (left + right) / 2;
    Point midPoint = points[mid];
    
    // Recursively find closest pair in left and right halves
    pair<Point, Point> leftPair = closestPairRecursive(points, left, mid);
    pair<Point, Point> rightPair = closestPairRecursive(points, mid + 1, right);
    
    double leftDist = leftPair.first.distanceTo(leftPair.second);
    double rightDist = rightPair.first.distanceTo(rightPair.second);
    
    double minDist = min(leftDist, rightDist);
    pair<Point, Point> minPair = (leftDist < rightDist) ? leftPair : rightPair;
    
    // Check for closer pairs across the dividing line
    vector<Point> strip;
    for (int i = left; i <= right; i++) {
        if (abs(points[i].x - midPoint.x) < minDist) {
            strip.push_back(points[i]);
        }
    }
    
    // Sort strip by y-coordinate
    sort(strip.begin(), strip.end(), [](const Point& a, const Point& b) {
        return a.y < b.y;
    });
    
    // Check for closer pairs in strip
    for (int i = 0; i < strip.size(); i++) {
        for (int j = i + 1; j < strip.size() && (strip[j].y - strip[i].y) < minDist; j++) {
            double dist = strip[i].distanceTo(strip[j]);
            if (dist < minDist) {
                minDist = dist;
                minPair = {strip[i], strip[j]};
            }
        }
    }
    
    return minPair;
}

pair<Point, Point> bruteForceClosestPair(vector<Point>& points, int left, int right) {
    double minDist = INFINITY;
    pair<Point, Point> minPair;
    
    for (int i = left; i <= right; i++) {
        for (int j = i + 1; j <= right; j++) {
            double dist = points[i].distanceTo(points[j]);
            if (dist < minDist) {
                minDist = dist;
                minPair = {points[i], points[j]};
            }
        }
    }
    
    return minPair;
}

Line Segment Intersection​

Bentley-Ottmann Algorithm​

struct Segment {
    Point p1, p2;
    int id;
    
    Segment(const Point& p1, const Point& p2, int id) : p1(p1), p2(p2), id(id) {}
    
    double getYAtX(double x) const {
        if (abs(p2.x - p1.x) < 1e-9) return p1.y;
        return p1.y + (p2.y - p1.y) * (x - p1.x) / (p2.x - p1.x);
    }
};

// Bentley-Ottmann algorithm for line segment intersection
vector<Point> bentleyOttmannIntersection(vector<Segment>& segments) {
    vector<Event> events;
    vector<Point> intersections;
    
    // Create events for each segment
    for (int i = 0; i < segments.size(); i++) {
        double x1 = segments[i].p1.x;
        double x2 = segments[i].p2.x;
        
        if (x1 > x2) {
            swap(x1, x2);
        }
        
        events.push_back(Event(x1, 0, i)); // Start
        events.push_back(Event(x2, 1, i)); // End
    }
    
    sort(events.begin(), events.end());
    
    set<Segment> activeSegments;
    
    for (const Event& event : events) {
        if (event.type == 0) {
            // Start event
            Segment current = segments[event.id];
            activeSegments.insert(current);
            
            // Check intersections with adjacent segments
            auto it = activeSegments.find(current);
            if (it != activeSegments.begin()) {
                auto prev = it;
                prev--;
                Point intersection = lineSegmentIntersection(current, *prev);
                if (intersection.x != INFINITY) {
                    intersections.push_back(intersection);
                }
            }
            
            auto next = it;
            next++;
            if (next != activeSegments.end()) {
                Point intersection = lineSegmentIntersection(current, *next);
                if (intersection.x != INFINITY) {
                    intersections.push_back(intersection);
                }
            }
        } else {
            // End event
            Segment current = segments[event.id];
            auto it = activeSegments.find(current);
            
            if (it != activeSegments.begin() && it != activeSegments.end()) {
                auto prev = it;
                prev--;
                auto next = it;
                next++;
                
                if (next != activeSegments.end()) {
                    Point intersection = lineSegmentIntersection(*prev, *next);
                    if (intersection.x != INFINITY) {
                        intersections.push_back(intersection);
                    }
                }
            }
            
            activeSegments.erase(current);
        }
    }
    
    return intersections;
}

Polygon Operations​

Point in Polygon​

// Check if point is inside polygon using ray casting
bool pointInPolygon(const Point& point, const vector<Point>& polygon) {
    int n = polygon.size();
    bool inside = false;
    
    for (int i = 0, j = n - 1; i < n; j = i++) {
        if (((polygon[i].y > point.y) != (polygon[j].y > point.y)) &&
            (point.x < (polygon[j].x - polygon[i].x) * (point.y - polygon[i].y) / 
             (polygon[j].y - polygon[i].y) + polygon[i].x)) {
            inside = !inside;
        }
    }
    
    return inside;
}

// Check if point is on polygon boundary
bool pointOnPolygonBoundary(const Point& point, const vector<Point>& polygon) {
    int n = polygon.size();
    
    for (int i = 0; i < n; i++) {
        int j = (i + 1) % n;
        LineSegment edge(polygon[i], polygon[j]);
        
        if (edge.containsPoint(point)) {
            return true;
        }
    }
    
    return false;
}

Polygon Intersection​

// Check if two polygons intersect
bool polygonsIntersect(const vector<Point>& poly1, const vector<Point>& poly2) {
    // Check if any vertex of poly1 is inside poly2
    for (const Point& p : poly1) {
        if (pointInPolygon(p, poly2)) {
            return true;
        }
    }
    
    // Check if any vertex of poly2 is inside poly1
    for (const Point& p : poly2) {
        if (pointInPolygon(p, poly1)) {
            return true;
        }
    }
    
    // Check if any edges intersect
    for (int i = 0; i < poly1.size(); i++) {
        int j = (i + 1) % poly1.size();
        LineSegment edge1(poly1[i], poly1[j]);
        
        for (int k = 0; k < poly2.size(); k++) {
            int l = (k + 1) % poly2.size();
            LineSegment edge2(poly2[k], poly2[l]);
            
            if (edge1.intersects(edge2)) {
                return true;
            }
        }
    }
    
    return false;
}

Polygon Triangulation​

// Triangulate polygon using ear clipping
vector<vector<Point>> triangulatePolygon(vector<Point>& polygon) {
    vector<vector<Point>> triangles;
    
    while (polygon.size() > 3) {
        bool earFound = false;
        
        for (int i = 0; i < polygon.size(); i++) {
            int prev = (i - 1 + polygon.size()) % polygon.size();
            int next = (i + 1) % polygon.size();
            
            if (isEar(polygon, prev, i, next)) {
                // Add triangle
                triangles.push_back({polygon[prev], polygon[i], polygon[next]});
                
                // Remove ear vertex
                polygon.erase(polygon.begin() + i);
                earFound = true;
                break;
            }
        }
        
        if (!earFound) {
            break; // No ear found, polygon might be degenerate
        }
    }
    
    // Add final triangle
    if (polygon.size() == 3) {
        triangles.push_back(polygon);
    }
    
    return triangles;
}

bool isEar(const vector<Point>& polygon, int prev, int curr, int next) {
    Point p1 = polygon[prev];
    Point p2 = polygon[curr];
    Point p3 = polygon[next];
    
    // Check if triangle is convex
    if ((p2 - p1).cross(p3 - p1) <= 0) {
        return false;
    }
    
    // Check if any other vertex is inside the triangle
    for (int i = 0; i < polygon.size(); i++) {
        if (i == prev || i == curr || i == next) continue;
        
        Point p = polygon[i];
        if (pointInTriangle(p, p1, p2, p3)) {
            return false;
        }
    }
    
    return true;
}

bool pointInTriangle(const Point& p, const Point& a, const Point& b, const Point& c) {
    double denom = (b.y - c.y) * (a.x - c.x) + (c.x - b.x) * (a.y - c.y);
    double alpha = ((b.y - c.y) * (p.x - c.x) + (c.x - b.x) * (p.y - c.y)) / denom;
    double beta = ((c.y - a.y) * (p.x - c.x) + (a.x - c.x) * (p.y - c.y)) / denom;
    double gamma = 1 - alpha - beta;
    
    return alpha >= 0 && beta >= 0 && gamma >= 0;
}

3D Geometry Basics​

3D Point Operations​

struct Point3D {
    double x, y, z;
    
    Point3D(double x = 0, double y = 0, double z = 0) : x(x), y(y), z(z) {}
    
    Point3D operator+(const Point3D& other) const {
        return Point3D(x + other.x, y + other.y, z + other.z);
    }
    
    Point3D operator-(const Point3D& other) const {
        return Point3D(x - other.x, y - other.y, z - other.z);
    }
    
    Point3D operator*(double scalar) const {
        return Point3D(x * scalar, y * scalar, z * scalar);
    }
    
    // Dot product
    double dot(const Point3D& other) const {
        return x * other.x + y * other.y + z * other.z;
    }
    
    // Cross product
    Point3D cross(const Point3D& other) const {
        return Point3D(
            y * other.z - z * other.y,
            z * other.x - x * other.z,
            x * other.y - y * other.x
        );
    }
    
    // Distance to another point
    double distanceTo(const Point3D& other) const {
        return sqrt((x - other.x) * (x - other.x) + 
                   (y - other.y) * (y - other.y) + 
                   (z - other.z) * (z - other.z));
    }
    
    // Magnitude
    double magnitude() const {
        return sqrt(x * x + y * y + z * z);
    }
};

3D Plane Operations​

struct Plane {
    double a, b, c, d; // ax + by + cz + d = 0
    
    Plane(double a = 0, double b = 0, double c = 0, double d = 0) 
        : a(a), b(b), c(c), d(d) {}
    
    // Plane from three points
    Plane(const Point3D& p1, const Point3D& p2, const Point3D& p3) {
        Point3D v1 = p2 - p1;
        Point3D v2 = p3 - p1;
        Point3D normal = v1.cross(v2);
        
        a = normal.x;
        b = normal.y;
        c = normal.z;
        d = -(a * p1.x + b * p1.y + c * p1.z);
    }
    
    // Distance from point to plane
    double distanceToPoint(const Point3D& p) const {
        return abs(a * p.x + b * p.y + c * p.z + d) / sqrt(a * a + b * b + c * c);
    }
    
    // Check if point is on plane
    bool containsPoint(const Point3D& p) const {
        return abs(a * p.x + b * p.y + c * p.z + d) < 1e-9;
    }
};

Performance Analysis​

Time Complexity​

  • Sweep line intersection: O(n log n + k) where k is number of intersections
  • Closest pair: O(n log n)
  • Point in polygon: O(n)
  • Polygon intersection: O(n × m)
  • 3D operations: O(1) for basic operations

Space Complexity​

  • Sweep line: O(n)
  • Closest pair: O(n)
  • Polygon operations: O(n)
  • 3D operations: O(1)

Common Patterns​

  1. Sweep line for intersection problems
  2. Divide and conquer for optimization problems
  3. Event-driven processing for geometric algorithms
  4. Cross product for orientation and area
  5. Ray casting for point-in-polygon tests

Applications​

  • Computer graphics: Rendering and collision detection
  • GIS: Geographic information systems
  • Robotics: Path planning and obstacle avoidance
  • Game development: Physics and collision detection
  • Computational geometry: Algorithm design and analysis