ReSCue - Retrieval System for Curves

---

         

---

s1_interval_container.h

#ifndef S1_INTERVAL_CONTAINER_H
#define S1_INTERVAL_CONTAINER_H

#include "Point2.h"
#include "rigid_motion.h"
#include "angle.h"
#include "rescutil.h"
#include "rescuenumerics.h"
#include <utility>
#include <set>

namespace rescue
{

// ==========================================================================================
// ==========================================================================================
    

enum interval_bound_type {lower_bound,upper_bound};

template <class arithmetic>
class angle_interval_bound : public std::pair<unit_circle_angle<arithmetic>, interval_bound_type>
{

public:

    angle_interval_bound<arithmetic>(const unit_circle_angle<arithmetic>& uca, const interval_bound_type& ib) : 
            std::pair<unit_circle_angle<arithmetic>, interval_bound_type>(uca,ib)
            {}

    bool operator<(const angle_interval_bound<arithmetic>& a) const
    {
        if (first<a.first)
            return true;
        else if (first==a.first)
            return second<a.second;
        return false;
    }

};

    

// ==========================================================================================
// ==========================================================================================

template <class arithmetic>
class s1_interval_container
{

public:

    bool is_empty();

    typedef unit_circle_angle<arithmetic> angle_t;
    typedef Point2<arithmetic> point;
    typedef std::pair<angle_t, angle_t > angle_pair;
    typedef std::pair<point,point> point_pair;
    typedef std::list<point_pair> point_pair_list;
    typedef angle_interval_bound<arithmetic> interval_bound;
    typedef std::multiset<interval_bound> interval_multiset;
    typedef rigid_motion<arithmetic> motion;

    s1_interval_container();

    void push_back(const point_pair& pq, const arithmetic& eps);

    bool intersect(point_pair& pq);
    void del(point_pair& pq);

    angle_t which_angle();
    motion where();

private:

    angle_t intersection_pos;

    angle_pair get_angle_interval(const point_pair&, const point_pair&, const arithmetic&);

    static bool empty(const angle_pair&);

    static bool full(const angle_pair&);

    int intersects;
    point_pair_list point_pairs;
    interval_multiset interval_bounds;

};

// ==========================================================================================
// ==========================================================================================



// ==========================================================================================
template <class arithmetic>
s1_interval_container<arithmetic>::s1_interval_container()
// ==========================================================================================
{
    intersects = 0;
}

// ==========================================================================================
template <class arithmetic>
void s1_interval_container<arithmetic>::push_back(const point_pair& pq, const arithmetic& eps)
// ==========================================================================================
{

    intersects++;

    if (point_pairs.empty())
    {
        point_pairs.push_back(pq);
        return;
    }
        
    angle_t angle_zero(false);
    angle_t angle_2pi(true);

    point_pair pq_pred = point_pairs.back();
    point_pair pq_first = point_pairs.front();

    // compute angle intervals for new point pair.
    // we need two intervals on S^1, i.e. up to four intervals on the real line.
    angle_pair I1 = get_angle_interval(pq_first,pq,eps);
    angle_pair I2 = get_angle_interval(pq_pred,pq,eps);

    point_pairs.push_back(pq);
    
    if (empty(I1) || empty(I2))
        return;
    
    if (full(I1) && full(I2))
    {
        intersects--;
        return;
    }

    // convert I1 and I2 to (one to three) intervals on the real line
    if (I1.first<=I1.second)
    {
        if (I2.first<=I2.second)    
        {
            if (I2.second<I1.first || I1.second<I2.first)
                return;
            if (I1.first<=I2.first)
                interval_bounds.insert(interval_bound(I2.first,lower_bound));
            else
                interval_bounds.insert(interval_bound(I1.first,lower_bound));
            if (I1.second<=I2.second)
                interval_bounds.insert(interval_bound(I1.second,upper_bound));
            else
                interval_bounds.insert(interval_bound(I2.second,upper_bound));
        } else { //(I2.first > I2.second)    
            if (I2.first<=I1.second)
            {
                interval_bounds.insert(interval_bound(I2.first,lower_bound));
                interval_bounds.insert(interval_bound(I1.second,upper_bound));
            }
            if (I1.first<=I2.second)
            {
                interval_bounds.insert(interval_bound(I1.first,lower_bound));
                interval_bounds.insert(interval_bound(I2.second,upper_bound));
            }
        }
    } else { // (I1.first > I1.second)
        if (I2.first<=I2.second)    
        {
            if (I1.first<=I2.second)
            {
                interval_bounds.insert(interval_bound(I1.first,lower_bound));
                interval_bounds.insert(interval_bound(I2.second,upper_bound));
            }
            if (I2.first<=I1.second)
            {
                interval_bounds.insert(interval_bound(I2.first,lower_bound));
                interval_bounds.insert(interval_bound(I1.second,upper_bound));
            }
        } else { // (I1.first > I1.second) UND (I2.first > I2.second)
            if (I1.first<=I2.second)
            {
                interval_bounds.insert(interval_bound(I1.first,lower_bound));
                interval_bounds.insert(interval_bound(I2.second,upper_bound));
            } else if (I2.first<=I1.second)
            {
                interval_bounds.insert(interval_bound(I2.first,lower_bound));
                interval_bounds.insert(interval_bound(I1.second,upper_bound));
            }
            if (I1.first<=I2.first)
            {
                interval_bounds.insert(interval_bound(I2.first,lower_bound));
                interval_bounds.insert(interval_bound(angle_2pi,upper_bound));
            } else {
                interval_bounds.insert(interval_bound(I1.first,lower_bound));
                interval_bounds.insert(interval_bound(angle_2pi,upper_bound));
            }
            if (I1.second<=I2.second)
            {
                interval_bounds.insert(interval_bound(angle_zero,lower_bound));
                interval_bounds.insert(interval_bound(I1.second,upper_bound));
            } else {
                interval_bounds.insert(interval_bound(angle_zero,lower_bound));
                interval_bounds.insert(interval_bound(I2.second,upper_bound));
            }
        }
    }

}

// ==========================================================================================
template <class arithmetic>
bool s1_interval_container<arithmetic>::intersect(s1_interval_container<arithmetic>::point_pair& pq)
// ==========================================================================================
{
}

// ==========================================================================================
template <class arithmetic>
void s1_interval_container<arithmetic>::del(s1_interval_container<arithmetic>::point_pair& pq)
// ==========================================================================================
{
}

// ==========================================================================================
template <class arithmetic>
bool s1_interval_container<arithmetic>::is_empty()
// ==========================================================================================
{

    interval_multiset::iterator i_it;
    int open_intervals = 0;
    bool success = false;

    for (i_it = interval_bounds.begin(); i_it != interval_bounds.end(); i_it++)
    {
        std::cerr<<i_it->first.rad();
        if (i_it->second==lower_bound) {
            open_intervals++; 
            if (open_intervals>=intersects-1)
            {
                angle_t phi_1 = i_it->first;
                //phi_1.bisect();
                //i_it++;
                //angle_t phi_2 = i_it->first;
                //phi_2.bisect();
                intersection_pos = phi_1;
                //intersection_pos += phi_2;
                return true;
            }
            std::cerr<<"[";
        }
        else
            {open_intervals--; std::cerr<<"]";}
    }

    std::cerr<<"\n";
    return success;

}

// ==========================================================================================
template <class arithmetic>
s1_interval_container<arithmetic>::angle_pair  s1_interval_container<arithmetic>::get_angle_interval(const point_pair& pp1, const point_pair& pp2, const arithmetic& eps)
// ==========================================================================================
{

    point p, q, S0, S1;

    arithmetic tp_x = (pp1.first.x+pp2.first.x) / 2.0;
    arithmetic tp_y = (pp1.first.y+pp2.first.y) / 2.0;
    arithmetic tq_x = (pp1.second.x+pp2.second.x) / 2.0;
    arithmetic tq_y = (pp1.second.y+pp2.second.y) / 2.0;

    p.x = pp2.first.x-tp_x;
    p.y = pp2.first.y-tp_y;
    q.x = pp2.second.x-tq_x;
    q.y = pp2.second.y-tq_y;

    arithmetic l_p_sq = p.x*p.x+p.y*p.y;
    arithmetic l_q_sq = q.x*q.x+q.y*q.y;
    arithmetic l_p = sqrt(l_p_sq);
    arithmetic l_q = sqrt(l_q_sq);

    if (l_p==0.0)
        if (l_q<=eps)
            // interval [0,2pi]
            return angle_pair(angle_t(false),angle_t(true));
        else
            // empty interval
            return angle_pair(angle_t(true),angle_t(false));
    if (l_q==0.0)
        if (l_p<=eps)
            // interval [0,2pi]
            return angle_pair(angle_t(false),angle_t(true));
        else
            // empty interval
            return angle_pair(angle_t(true),angle_t(false));
    arithmetic diff = (l_p-l_q);
    if (diff>eps || (-diff)>eps)
        // empty interval
        return angle_pair(angle_t(true),angle_t(false));

    arithmetic cos_phi = (l_p_sq+l_q_sq-eps*eps)/(2.0*l_p*l_q);
    arithmetic sin_phi = sqrt(1-cos_phi*cos_phi);

    angle_t phi_0(sin_phi,cos_phi);
    angle_t phi_1((-1.0)*sin_phi,cos_phi);

    angle_t alpha(p.y/l_p,p.x/l_p);
    angle_t minus_beta((-1.0)*q.y/l_q,q.x/l_p);

    angle_t theta_0, theta_1;

    theta_0 = alpha;
    theta_0 += minus_beta;
    theta_0 += phi_0;

    theta_1 = alpha;
    theta_1 += minus_beta;
    theta_1 += phi_1;

    return angle_pair(theta_1,theta_0);

}

template <class arithmetic>
    s1_interval_container<arithmetic>::angle_t s1_interval_container<arithmetic>::which_angle()
// ==========================================================================================
{
    return intersection_pos;
}

template <class arithmetic>
    s1_interval_container<arithmetic>::motion s1_interval_container<arithmetic>::where()
// ==========================================================================================
{
    motion rho_phi = motion(intersection_pos);

    point pp = point_pairs.front().first;    
    point qq = point_pairs.front().second;
    qq *= rho_phi;

    arithmetic t_x = pp.x-qq.x;
    arithmetic t_y = pp.y-qq.y;

    motion there(t_x,t_y);
    there *= rho_phi;

    return there;
}

// ==========================================================================================
template <class arithmetic>
    bool s1_interval_container<arithmetic>::full(const s1_interval_container<arithmetic>::angle_pair& I)
// ==========================================================================================
{
    if (I.first==angle_t(false) && I.second==angle_t(true))
        return true;
    return false;
}

// ==========================================================================================
template <class arithmetic>
    bool s1_interval_container<arithmetic>::empty(const s1_interval_container<arithmetic>::angle_pair& I)
// ==========================================================================================
{
    if (I.first==angle_t(true) && I.second==angle_t(false))
        return true;
    return false;
}


// ==========================================================================================
// ==========================================================================================

}

#endif

ReSCue - Retrieval System for Curves
Universität Bonn, Institut für Informatik III, 2001