diff --git a/geometry/graham_scan_algorithm.cpp b/geometry/graham_scan_algorithm.cpp new file mode 100644 index 000000000..8ed144c11 --- /dev/null +++ b/geometry/graham_scan_algorithm.cpp @@ -0,0 +1,79 @@ +/****************************************************************************** +* @file +* @brief Implementation of the [Convex Hull](https://en.wikipedia.org/wiki/Convex_hull) +* implementation using [Graham Scan](https://en.wikipedia.org/wiki/Graham_scan) +* @details +* In geometry, the convex hull or convex envelope or convex closure of a shape +* is the smallest convex set that contains it. The convex hull may be defined +* either as the intersection of all convex sets containing a given subset of a +* Euclidean space, or equivalently as the set of all convex combinations of +* points in the subset. For a bounded subset of the plane, the convex hull may +* be visualized as the shape enclosed by a rubber band stretched around the subset. +* +* The worst case time complexity of Jarvis’s Algorithm is O(n^2). Using Graham’s +* scan algorithm, we can find Convex Hull in O(nLogn) time. +* +* ### Implementation +* +* Sort points +* We first find the bottom-most point. The idea is to pre-process +* points be sorting them with respect to the bottom-most point. Once the points +* are sorted, they form a simple closed path. +* The sorting criteria is to use the orientation to compare angles without actually +* computing them (See the compare() function below) because computation of actual +* angles would be inefficient since trigonometric functions are not simple to evaluate. +* +* Accept or Reject Points +* Once we have the closed path, the next step is to traverse the path and +* remove concave points on this path using orientation. The first two points in +* sorted array are always part of Convex Hull. For remaining points, we keep track +* of recent three points, and find the angle formed by them. Let the three points +* be prev(p), curr(c) and next(n). If orientation of these points (considering them +* in same order) is not counterclockwise, we discard c, otherwise we keep it. +* +* @author [Lajat Manekar](https://github.com/Lazeeez) +* +*******************************************************************************/ +#include /// for IO Operations +#include /// for std::assert +#include /// for std::vector +#include /// for all the functions used + +/******************************************************************************* + * @brief Self-test implementations + * @returns void + *******************************************************************************/ +void test() { + std::vector points = {{0, 3}, + {1, 1}, + {2, 2}, + {4, 4}, + {0, 0}, + {1, 2}, + {3, 1}, + {3, 3}}; + std::vector expected_result = {{0, 3}, + {4, 4}, + {3, 1}, + {0, 0}}; + std::vector derived_result; + std::vector res; + + derived_result = geometry::grahamscan::convexHull(points, points.size()); + + std::cout << "Test#1..."; + for (int i = 0; i < expected_result.size(); i++) { + assert(derived_result[i].x == expected_result[i].x); + assert(derived_result[i].y == expected_result[i].y); + } + std::cout << "Passed" << std::endl; +} + +/******************************************************************************* + * @brief Main function + * @returns 0 on exit + *******************************************************************************/ +int main() { + test(); + return 0; +} \ No newline at end of file diff --git a/geometry/graham_scan_functions.h b/geometry/graham_scan_functions.h new file mode 100644 index 000000000..291d84661 --- /dev/null +++ b/geometry/graham_scan_functions.h @@ -0,0 +1,193 @@ +/****************************************************************************** +* @file +* @brief Implementation of the [Convex Hull](https://en.wikipedia.org/wiki/Convex_hull) +* implementation using [Graham Scan](https://en.wikipedia.org/wiki/Graham_scan) +* @details +* In geometry, the convex hull or convex envelope or convex closure of a shape +* is the smallest convex set that contains it. The convex hull may be defined +* either as the intersection of all convex sets containing a given subset of a +* Euclidean space, or equivalently as the set of all convex combinations of +* points in the subset. For a bounded subset of the plane, the convex hull may +* be visualized as the shape enclosed by a rubber band stretched around the subset. +* +* The worst case time complexity of Jarvis’s Algorithm is O(n^2). Using Graham’s +* scan algorithm, we can find Convex Hull in O(nLogn) time. +* +* ### Implementation +* +* Sort points +* We first find the bottom-most point. The idea is to pre-process +* points be sorting them with respect to the bottom-most point. Once the points +* are sorted, they form a simple closed path. +* The sorting criteria is to use the orientation to compare angles without actually +* computing them (See the compare() function below) because computation of actual +* angles would be inefficient since trigonometric functions are not simple to evaluate. +* +* Accept or Reject Points +* Once we have the closed path, the next step is to traverse the path and +* remove concave points on this path using orientation. The first two points in +* sorted array are always part of Convex Hull. For remaining points, we keep track +* of recent three points, and find the angle formed by them. Let the three points +* be prev(p), curr(c) and next(n). If orientation of these points (considering them +* in same order) is not counterclockwise, we discard c, otherwise we keep it. +* +* @author [Lajat Manekar](https://github.com/Lazeeez) +* +*******************************************************************************/ +#include /// for IO Operations +#include /// for std::stack +#include /// for std::vector +#include /// for std::swap +#include /// for mathematics and datatype conversion + +/****************************************************************************** + * @namespace geometry::grahamscan + * @brief geometric algorithms + *******************************************************************************/ +namespace geometry::grahamscan { + + /****************************************************************************** + * @struct Point + * @brief for X and Y co-ordinates of the co-ordinate. + *******************************************************************************/ + struct Point { + int x, y; + }; + + // A global point needed for sorting points with reference + // to the first point Used in compare function of qsort() + Point p0; + + /****************************************************************************** + * @brief A utility function to find next to top in a stack. + * @param S Stack to be used for the process. + * @returns @param Point Co-ordinates of the Point + *******************************************************************************/ + Point nextToTop(std::stack &S) { + Point p = S.top(); + S.pop(); + Point res = S.top(); + S.push(p); + return res; + } + + /****************************************************************************** + * @brief A utility function to return square of distance between p1 and p2. + * @param p1 Co-ordinates of Point 1 . + * @param p2 Co-ordinates of Point 2 . + * @returns @param int distance between p1 and p2. + *******************************************************************************/ + int distSq(Point p1, Point p2) { + return (p1.x - p2.x) * (p1.x - p2.x) + (p1.y - p2.y) * (p1.y - p2.y); + } + + /****************************************************************************** + * @brief To find orientation of ordered triplet (p, q, r). + * @param p Co-ordinates of Point p . + * @param q Co-ordinates of Point q . + * @param r Co-ordinates of Point r . + * @returns @param int 0 --> p, q and r are collinear, 1 --> Clockwise, + * 2 --> Counterclockwise + *******************************************************************************/ + int orientation(Point p, Point q, Point r) { + int val = (q.y - p.y) * (r.x - q.x) - (q.x - p.x) * (r.y - q.y); + + if (val == 0) return 0; // collinear + return (val > 0) ? 1 : 2; // clock or counter-clock wise + } + + /****************************************************************************** + * @brief A function used by library function qsort() to sort an array of + * points with respect to the first point + * @param vp1 Co-ordinates of Point 1 . + * @param vp2 Co-ordinates of Point 2 . + * @returns @param int distance between p1 and p2. + *******************************************************************************/ + int compare(const void *vp1, const void *vp2) { + auto *p1 = static_cast(vp1); + auto *p2 = static_cast(vp2); + + // Find orientation + int o = orientation(p0, *p1, *p2); + if (o == 0) { + return (distSq(p0, *p2) >= distSq(p0, *p1)) ? -1 : 1; + } + + return (o == 2) ? -1 : 1; + } + + /****************************************************************************** + * @brief Prints convex hull of a set of n points. + * @param points vector of Point with co-ordinates. + * @param size Size of the vector. + * @returns @param vector vector of Conver Hull. + *******************************************************************************/ + std::vector convexHull(std::vector points, uint64_t size) { + + // Find the bottom-most point + int ymin = points[0].y, min = 0; + for (int i = 1; i < size; i++) { + int y = points[i].y; + + // Pick the bottom-most or chose the left-most point in case of tie + if ((y < ymin) || (ymin == y && points[i].x < points[min].x)) { + ymin = points[i].y, min = i; + } + } + + // Place the bottom-most point at first position + std::swap(points[0], points[min]); + + // Sort n-1 points with respect to the first point. A point p1 comes + // before p2 in sorted output if p2 has larger polar angle + // (in counterclockwise direction) than p1. + p0 = points[0]; + qsort(&points[1], size - 1, sizeof(Point), compare); + + // If two or more points make same angle with p0, Remove all but the one + // that is farthest from p0 Remember that, in above sorting, our criteria + // was to keep the farthest point at the end when more than one points have + // same angle. + int m = 1; // Initialize size of modified array + for (int i = 1; i < size; i++) { + // Keep removing i while angle of i and i+1 is same with respect to p0 + while (i < size - 1 && + orientation(p0, points[i], points[i + 1]) == 0) { + i++; + } + + points[m] = points[i]; + m++; // Update size of modified array + } + + // If modified array of points has less than 3 points, convex hull is not possible + if (m < 3) return {}; + + // Create an empty stack and push first three points to it. + std::stack S; + S.push(points[0]); + S.push(points[1]); + S.push(points[2]); + + // Process remaining n-3 points + for (int i = 3; i < m; i++) { + // Keep removing top while the angle formed by + // points next-to-top, top, and points[i] makes + // a non-left turn + while (S.size() > 1 && orientation(nextToTop(S), S.top(), points[i]) != 2) { + S.pop(); + } + S.push(points[i]); + } + + std::vector result; + // Now stack has the output points, push them into the resultant vector + while (!S.empty()) { + Point p = S.top(); + result.push_back(p); + S.pop(); + } + + return result; // return resultant vector with Convex Hull co-ordinates. + } +} // namespace geometry::grahamscan \ No newline at end of file