Space-efficient geometric algorithms and data structures By Ilya Katz and Hervé Brönnimann

intersection.hpp

Go to the documentation of this file.

#include <bits/stl_tree.h>
//#include "my_tree.hpp"
#include "explicit_tree_heap.hpp"
#include "cgl.hpp"
#include <queue>
#include <utility>

using namespace std;
using namespace inplaceds;
//using namespace inplaceds::tree;

typedef cgl<int> Geometry;
typedef Geometry::line_type line_type;
typedef Geometry::point_type point_type;
typedef Geometry::coord_type coord_type;
typedef DataNode<Geometry> node;

int SWEEP_LINE_COORD=-100; 

template <class Geometry>
struct PQComp
{
        typedef DataNode<Geometry>*     node_type;
        typedef typename Geometry::line_type    line_type;

        bool operator()(const node_type& a, const node_type& b){
                return !geom.less_x_object()(a->point,b->point);
        }

        Geometry geom;

};

template <class type>
struct kov
{

        const type& operator()(const type& a)const{return a;}
        type& operator()( type& a){return a;}
};

#define slope(x1,y1,x2,y2) ( (y2-y1)/(x2-x1) )

template <class Geometry>
struct BSTComp
{
        typedef typename Geometry::line_type    line_type;
        typedef typename Geometry::coord_type   coord_type;
        typedef DataNode<Geometry>                              node;

        bool operator()(const node& a, const node& b){
                return operator()(a.line, b.line);
        }

        bool operator()( const node* const& a, const node* const& b){
                return operator()(a->line, b->line);
        }

        bool operator()(const line_type& a, const line_type& b){

                //two lines are equal
                if(a[0][0]==b[0][0] && a[0][1]==b[0][1] &&
                        a[1][0]==b[1][0] && a[1][1]==b[1][1]){
                                return false;
                        }

                //else, need to check which line is has higher
                //y-coordinate at the current sweepline position

                coord_type currentY1 =
                        slope(a[0][0],a[0][1], a[1][0], a[1][1])*SWEEP_LINE_COORD + a[0][0],
                                         currentY2 =
                        slope(b[0][0],b[0][1], b[1][0], b[1][1])*SWEEP_LINE_COORD + b[0][0];

                return currentY1<currentY2;

        }

        Geometry geom;
};

typedef _Rb_tree<node*, node*, kov<node*>, BSTComp<Geometry> > tree_type;
typedef priority_queue<node*, vector<node*>, PQComp<Geometry> > priority_queue_type;

std::pair<bool, point_type> intersection(const line_type& l1, const line_type& l2);
void update_event(node*,node*,priority_queue_type&);
void intersection(vector<line_type>);

/*
 * temporary, for convience
 */

ostream& operator<<(ostream& os, node* x)
{
        os<<x->line<<" E:"<<int(x->myEvent);
        return os;
}

void intersection(vector<line_type> all_lines){

        Geometry geom;

        //initialize T and Q
        tree_type T;
        priority_queue_type Q;

        tree_type::iterator tree_it, tree_it_prev, tree_it_next;

        for(vector<line_type>::iterator it = all_lines.begin(); it!=all_lines.end(); it++){
                Q.push(new node(*it,LEFT_EVENT));
        }

        //move sweepline to e (leftmost event)


        node* current = Q.top();
        cerr<<"Q size "<<Q.size()<<endl;
        SWEEP_LINE_COORD = current->point[0];


        //the problem with this code is that
        //T.find is using the BSTComp functor that
        //compares using the current value of the sweep_line
                while(!Q.empty()){
                Q.pop();
                switch(current->myEvent){
                        case(LEFT_EVENT):
                                tree_it = T.insert_equal(current);              //insert line into T
                                tree_it_next=tree_it_prev=tree_it;
                                cerr<<"Inserted in T"<<*tree_it<<endl;

                                tree_it_next++;
                                if(!(tree_it_next==T.end())){
                                        update_event(*tree_it, *tree_it_next, Q);
                                }else{
                                        (*tree_it)->right_event();
                                        Q.push(*tree_it);
                                }

                                if(!(tree_it==T.begin())){
                                        tree_it_prev--;
                                        cerr<<"Updating previous:"<<*tree_it_prev<<endl;
                                        update_event(*tree_it_prev, *tree_it, Q);       //update previous line event
                                }
                                break;
                        case(RIGHT_EVENT):
                                tree_it = T.find(current);
                                if(tree_it==T.end()){
                                        cerr<<"Couldnt find "<<current;
                                        exit(-1);
                                }
                                tree_it_next =tree_it;
                                tree_it_next++;
                                cerr<<"Removing from T:"<<*tree_it<<endl;
                                T.erase(tree_it, tree_it_next);                 //remove line from T
                                cerr<<"Removed from T"<<endl;
                                if(!(tree_it==T.begin())){
                                        tree_it_prev = tree_it;
                                        tree_it_prev--;
                                        update_event(*tree_it_prev, *tree_it,Q);        //update previous line event
                                }
                                break;
                        case(INTERSECTION_EVENT):
                                cout<<"Got one!"<<endl;
                                cout<<current->line<<" and "<<current->intersecting->line<<endl;
                                tree_it = T.find(current);
                                cerr<<"T is now:"<<endl;
                                copy(T.begin(), T.end(), ostream_iterator<node*>(cerr, " "));
                                cerr<<endl;
                                if(tree_it==T.end()){
                                        cerr<<"Couldnt find "<<current;
                                        exit(-1);
                                }

                                tree_it_prev = tree_it_next = tree_it;



                                tree_it_next++;
                                if(tree_it_next!=T.end()){
                                        std::swap(*tree_it, *tree_it_next);                     //switch order of current and next - check if this is correct!
                                        cerr<<"Swapped"<<endl;
                                }



                                if(!(tree_it==T.begin())){
                                        cerr<<"Updating prev "<<endl;
                                        tree_it_prev--;
                                        update_event(*tree_it_prev, *tree_it,Q);                //update previous line event
                                        cerr<<"Updated prev "<<endl;
                                }
                                cerr<<(tree_it_next==T.end())<<(tree_it==T.end());
                                if(tree_it_next!=T.end()){
                                        cerr<<"Updating curr "<<endl;//<<*tree_it<<" "<<endl;
                                        update_event(*tree_it, *tree_it_next,Q);                //update current line event
                                        cerr<<"Updated curr "<<endl;
                                }

                                tree_it = tree_it_next;
                                tree_it_next++;
                                if(!(tree_it==T.end() || tree_it_next==T.end())){
                                        cerr<<"Updating next "<<endl;
                                        update_event(*tree_it, *(tree_it_next),Q);      //update next line event
                                        cerr<<"Updated next "<<endl;
                                }
                                cerr<<"ended intersection case"<<endl;
                                break;
                }
                current = Q.top();
                SWEEP_LINE_COORD = current->point[0];
        }

}

void update_event(node* current, node* next, priority_queue_type& Q){

        //if there is an intersection point and it is to the right of the sweepline,
        // enter it here, otherwise right-end point

        std::pair<bool, point_type> inter = intersection(current->line, next->line);

        if(!inter.first){ //right-end
                current->right_event();
        }
        else if(inter.second[0]>=SWEEP_LINE_COORD){ 
                current->intersection_event(inter.second, next);
        }
        Q.push(current);
}


std::pair<bool, point_type> intersection(const line_type& l1, const line_type& l2){

        std::pair<bool, point_type> ret(false, point_type());

        int m1 = slope(l1[0][0], l1[0][1], l1[1][0], l1[1][1]),
                        m2 = slope(l2[0][0], l2[0][1], l2[1][0], l2[1][1]);

        if(m1 == m2){  //lines are parallel
                return ret;
        }

        coord_type b1 = l1[0][1]-m1*l1[0][0],
                                 b2 = l2[0][1]-m2*l2[0][0];

        coord_type x = (b2-b1)/(m1-m2);


        //originally, the left end point is at position [0]
        //but, beacuse of encoding intersection, these positions
        //may have changed
        if(x<std::min(l1[0][0],l1[1][0]) || x>std::max(l1[0][0],l1[1][0]) ||
                 x<std::min(l2[0][0],l2[1][0]) || x>std::max(l2[0][0],l2[1][0])){
                return ret;
        }

        coord_type y1 = m1 * x + b1;
        coord_type y2 = m2 * x + b2;

        if(!(y1==y2)){
                return ret;
        }
        ret.second[0]=x;
        ret.second[1]=y1;
        ret.first=true;

        return ret;

}

---

Code Documentation generated Using Doxygen

Copyright © Ilya Katz and Hervé Brönnimann, 2005.