+void intersect_line_line_x(double a1, double b1, double d1, double a2, double b2, double d2, double *min_x, double *max_x)
+{
+ // find somewhere (x,y) where:
+ //
+ // A1x + B1y + D1 = 0 [I]
+ // A2x + B2y + D2 = 0 [II]
+ //
+ // [I] gives: y = -((A1/B1) x + D1/B1)
+ // Insert into [II]:
+ //
+ // A2x - B2 ((A1/B1) x + D1/B1) + D2 = 0
+ // A2x - A1 B2 / B1 x - D1 B2 / B1 + D2 = 0
+ // ((A2 B1 - A1 B2) / B1) x = D1 B2 / B1 - D2 [A2 B1 - A1 B2 = det]
+ // x = D1 B2 / det - D2 B1 / det = (D1 B2 - D2 B1) / det
+ double det = a2 * b1 - a1 * b2;
+ if (fabs(det) <= 1e-6) {
+ return;
+ }
+ double intersect_x = (d1 * b2 - d2 * b1) / det;
+
+ if (-det * (*max_x - *min_x) < 0) {
+ *min_x = min(*min_x, intersect_x);
+ *max_x = min(*max_x, intersect_x);
+ } else {
+ *min_x = max(*min_x, intersect_x);
+ *max_x = max(*max_x, intersect_x);
+ }
+}
+
+void intersect_line_line_y(double a1, double b1, double d1, double a2, double b2, double d2, double *min_y, double *max_y)
+{
+ // find somewhere (x,y) where:
+ //
+ // A1x + B1y + D1 = 0 [I]
+ // A2x + B2y + D2 = 0 [II]
+ //
+ // [I] gives: y = -((B1/A1) x + D1/A1)
+ // Insert into [II]:
+ //
+ // B2y - A2 ((B1/A1) x + D1/A1) + D2 = 0
+ // B2y - B1 A2 / A1 x - D1 A2 / A1 + D2 = 0
+ // ((A1 B2 - B1 A2) / A1) y = D1 A2 / A1 - D2 [A1 B2 - B1 A2 = det]
+ // y = D1 A2 / det - D2 A1 / det = (D1 B2 - D2 B1) / det
+ double det = a1 * b2 - a2 * b1;
+ if (fabs(det) <= 1e-6) {
+ return;
+ }
+ double intersect_y = (d1 * a2 - d2 * a1) / det;
+ if (det * (*max_y - *min_y) < 0) {
+ *min_y = min(*min_y, intersect_y);
+ *max_y = min(*max_y, intersect_y);
+ } else {
+ *min_y = max(*min_y, intersect_y);
+ *max_y = max(*max_y, intersect_y);
+ }
+}
+
+void plot_line_based_on_x(FILE *fp, double a, double b, double d, double x1, double x2)
+{
+ // see intersect_line_line()
+ double y1 = - (a * x1 + d) / b;
+ double y2 = - (a * x2 + d) / b;
+ fprintf(fp, "%f %f %f %f\n", x1, y1, x2, y2);
+}
+
+void plot_line_based_on_y(FILE *fp, double a, double b, double d, double y1, double y2)
+{
+ // see intersect_line_line()
+ double x1 = - (b * y1 + d) / a;
+ double x2 = - (b * y2 + d) / a;
+ fprintf(fp, "%f %f %f %f\n", x1, y1, x2, y2);
+}
+
+void plot_line(double a, double b, double d)
+{
+ if (fabs(a) > fabs(b)) {
+ double x1 = -2000.0, x2 = 2000.0;
+ if (b > 0) swap(x1, x2);
+ intersect_line_line_x(a, b, d, 0, 1.0, 2000, &x1, &x2);
+ intersect_line_line_x(a, b, d, 0, -1.0, 2000, &x1, &x2);
+ plot_line_based_on_x(stdout, a, b, d, x1, x2);
+ } else {
+ double y1 = -2000.0, y2 = 2000.0;
+ if (a < 0) swap(y1, y2);
+ intersect_line_line_y(a, b, d, 1.0, 0.0, 2000, &y1, &y2);
+ intersect_line_line_y(a, b, d, -1.0, 0.0, 2000, &y1, &y2);
+ plot_line_based_on_y(stdout, a, b, d, y1, y2);
+ }
+}
+
+void plot_bsp_tree(FILE *fp, BSPTreeNode *node, vector<pair<BSPTreeNode*, bool> > previous_visited)
+{
+ double a = node->a, b = node->b, d = node->d;
+
+ double x1, x2, y1, y2;
+
+ if (node->is_leaf) {
+ return;
+ }
+
+ if (fabs(a) > fabs(b)) {
+ x1 = -2000.0, x2 = 2000.0;
+ if (b > 0) swap(x1, x2);
+ intersect_line_line_x(a, b, d, 0, 1.0, 2000, &x1, &x2);
+ intersect_line_line_x(a, b, d, 0, -1.0, 2000, &x1, &x2);
+
+ // clip against all previous lines
+ for (unsigned i = 0; i < previous_visited.size(); ++i) {
+ BSPTreeNode* o = previous_visited[i].first;
+ if (previous_visited[i].second) {
+ intersect_line_line_x(a, b, d, o->a, o->b, o->d, &x1, &x2);
+ } else {
+ intersect_line_line_x(a, b, d, -o->a, -o->b, -o->d, &x1, &x2);
+ }