// Description: implementation of D-CEL class

// Revision history: 28/6/99/LY - first few lines written

//							25/9/88/LY - development ended

 

#include "dcel.h"

#include <algorithm>

#include <iostream>

#include <limits>

#include <string>



using namespace std;



class D_CEL;  



D_CEL::D_CEL()

{	

}



// The destructor need only delete the half-edge records allocated - vertex and face records are deleted by the

// destructor of the storage container that stores them.

D_CEL::~D_CEL()

{

	for (HalfEdgeList::iterator he_iter = half_edge_records.begin(); he_iter != half_edge_records.end(); he_iter++)

	{	delete he_iter->second;	}

}



// Precondition: vertices1, edges1, faces1 and vertices2, edges2, faces2 are valid input containers holding 

//	correcly formatted input D_CEL data. 

// Postcondition: *this contains valid records for the vertices, half-edges and faces of a complete d-cel

void D_CEL::build_d_cel(InputVertices& vertices1, InputEdges& edges1, InputFaces& faces1, char subdiv1,

				 InputVertices& vertices2, InputEdges& edges2, InputFaces& faces2, 

				 char subdiv2)

{

	InputVertices::iterator v_iter;	

	InputEdges::iterator e_iter;

	InputFaces::iterator f_iter;

	// allocate vertex records

	Vertex*	vertex_ptr;

	for (v_iter = vertices1.begin(); v_iter != vertices1.end(); v_iter++)

	{

		vertex_ptr = new Vertex(v_iter - vertices1.begin()); 

		vertex_ptr->coordinates = *v_iter;

		vertex_ptr->auxiliary_data.subdivision = endpoint_from_subdiv1;

		vertex_records.push_back(vertex_ptr);

	}

	for (v_iter = vertices2.begin(); v_iter != vertices2.end(); v_iter++)

	{

		vertex_ptr = new Vertex(v_iter - vertices2.begin() + vertices1.size());

		vertex_ptr->coordinates = *v_iter;

		vertex_ptr->auxiliary_data.subdivision = endpoint_from_subdiv2;

		vertex_records.push_back(vertex_ptr);

	}

	// allocate half-edge records

	HalfEdge* he;

	for (e_iter = edges1.begin(); e_iter != edges1.end(); e_iter++)

	{

		for (InputEdge::iterator cycle_iter = e_iter->begin(); cycle_iter != e_iter->end()-1; cycle_iter++)

		{

			he = new HalfEdge;

			he->origin = vertex_records[*cycle_iter]; 

			half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))] = he;

			vertex_records[*cycle_iter]->incident_edge = he;		// update record of origin vertex

		}

	}

	for (e_iter = edges2.begin(); e_iter != edges2.end(); e_iter++)

	{

		for (InputEdge::iterator cycle_iter = e_iter->begin(); cycle_iter != e_iter->end()-1; cycle_iter++)

		{

			he = new HalfEdge;

			he->origin = vertex_records[*(cycle_iter) +  vertices1.size()]; 

			half_edge_records[calc_key(*(cycle_iter) + vertices1.size(), *(cycle_iter + 1) + vertices1.size())] = he;

			//cout << calc_key(*cycle_iter + vertices1.size(), *(cycle_iter+1) + vertices1.size()) << " ";



			vertex_records[*(cycle_iter) + vertices1.size()]->incident_edge = he;		// update record of origin vertex

		}

	}

	//cout << endl << endl;



	// allocate face records

	Face *face_ptr, *unbounded_face_ptr;

	vector<Face*> face_record_ptrs1(faces1.size() + 1);

	vector<unsigned int>	face_record_indices1(faces1.size() + 1);



	char convert_buf[numeric_limits<long>::digits10 + 1] = {'\0'};



	convert_buf[0] = subdiv1;

	convert_buf[1] = '0';

	unbounded_face_ptr = new Face(convert_buf);

	face_record_ptrs1[faces1.size()] = unbounded_face_ptr;

	face_record_indices1[faces1.size()] = faces1.size();

	int face_index, face_no;

	for (f_iter = faces1.begin(), face_index = 0, face_no = 1; f_iter != faces1.end(); f_iter++, face_index++)

	{

		InputEdge::iterator ic_iter = f_iter->begin();

		if (*(ic_iter+1) != -2)

		{	convert_buf[0] = subdiv1;

			itoa(face_no++, convert_buf+1, 10); 

			face_ptr = new Face(convert_buf);

			face_ptr->outer_component = half_edge_records[calc_key(*ic_iter, *(ic_iter+1))];

			//cout << *ic_iter << *(ic_iter+1) << " ";

			// build inner components list

			for (ic_iter+=2; ic_iter != f_iter->end(); ic_iter+=2)

			{	face_ptr->inner_components.push_back(half_edge_records[calc_key(*ic_iter, *(ic_iter+1))]);	}

			face_records.push_back(face_ptr);

			face_record_ptrs1[face_index] = face_ptr;			

			face_record_indices1[face_index] = face_index ;

		}

		else

		{	if (*ic_iter == -2)	// connect to unbounded face

			{	

				face_record_indices1[face_index] = face_record_indices1[faces1.size()];	

				HalfEdge* hePtr = half_edge_records[calc_key(*( (edges1.begin()+face_index)->begin() ) , 

																		   *( (edges1.begin()+face_index)->begin()+1)  )];

				face_record_ptrs1[ face_record_indices1[faces1.size()] ]->inner_components.push_back(hePtr);

			}

			else

			{	face_record_indices1[face_index] = *ic_iter;	

			}

		}

	}

	//cout << endl << endl << "Faces allocated: " << face_no << endl;

	vector<Face*> face_record_ptrs2(faces2.size() + 1);

	vector<unsigned int> face_record_indices2(faces2.size() + 1);



	convert_buf[0] = subdiv2;

	convert_buf[1] = '0';

	unbounded_face_ptr = new Face(convert_buf);

	face_record_ptrs2[faces2.size()] = unbounded_face_ptr;

	face_record_indices2[faces2.size()] = faces2.size();

	for (f_iter = faces2.begin(), face_index = 0, face_no = 1; f_iter != faces2.end(); f_iter++, face_index++)

	{

		InputEdge::iterator ic_iter = f_iter->begin();

		if (*(ic_iter+1) != -2)

		{	convert_buf[0] = subdiv2;

			itoa(face_no++, convert_buf+1, 10);

			face_ptr = new Face(convert_buf);

			face_ptr->outer_component = half_edge_records[calc_key(*(ic_iter) + vertices1.size(), *(ic_iter + 1) + vertices1.size())];

			// build inner components list

			for (ic_iter+=2; ic_iter != f_iter->end(); ic_iter+=2)

			{	face_ptr->inner_components.push_back(half_edge_records[

				calc_key(*(ic_iter) + vertices1.size(), *(ic_iter + 1) + vertices1.size())]);	

			}

			face_records.push_back(face_ptr);

			face_record_ptrs2[face_index] = face_ptr;

			face_record_indices2[face_index] = face_index;

		}

		else

		{	if (*(ic_iter) == -2)

			{	face_record_indices2[face_index] = face_record_indices2[faces2.size()];	

				HalfEdge* hePtr = half_edge_records[calc_key( *( (edges2.begin()+face_index)->begin() )+vertices1.size(), 

																			 *( (edges2.begin()+face_index)->begin()+1)+vertices1.size() )];

				//cout << calc_key(*( (edges2.begin()+face_index)->begin() )+vertices1.size() , 

				//					  *( (edges2.begin()+face_index)->begin()+1)+vertices1.size() ) << "IH ";



				face_record_ptrs2[ face_record_indices2[faces2.size()] ]->inner_components.push_back(hePtr);

			}

			else

			{	face_record_indices2[f_iter - faces2.begin()] = *ic_iter;

			}

		}

	}

	// finish half-edge records by setting references to now allocated half-edge and face records

	for (e_iter = edges1.begin(); e_iter != edges1.end(); e_iter++)

	{

		InputEdge::iterator cycle_iter = e_iter->begin();

		// set remaining data members in record for first half-edge in cycle

		half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->prev = 

			half_edge_records[calc_key(*(e_iter->end()-2), *(e_iter->end()-1))];

		half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->next =

			half_edge_records[calc_key(*(cycle_iter+1), *(cycle_iter+2))];

		half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->twin =

			half_edge_records[calc_key(*(cycle_iter+1), *cycle_iter)];

		half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->incident_face =

			face_record_ptrs1[ face_record_indices1[e_iter - edges1.begin()] ];

		

		// set remaining data members in records for in-between half-edges

		for (cycle_iter++; cycle_iter != e_iter->end()-2; cycle_iter++)

		{

			half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->prev =

				half_edge_records[calc_key(*(cycle_iter-1), *cycle_iter)];	

			half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->next = 

				half_edge_records[calc_key(*(cycle_iter+1), *(cycle_iter+2))];

			half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->twin = 

				half_edge_records[calc_key(*(cycle_iter+1), *cycle_iter)];

			half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->incident_face =

				face_record_ptrs1[ face_record_indices1[e_iter - edges1.begin()] ];

		}

		// set remaining data members in records for last half-edge in cycle

		half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->prev = 

			half_edge_records[calc_key(*(cycle_iter-1), *cycle_iter)];

		half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->next =

			half_edge_records[calc_key(*(e_iter->begin()), *(e_iter->begin()+1))];

		half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->twin =

			half_edge_records[calc_key(*(cycle_iter+1),*cycle_iter)];

		half_edge_records[calc_key(*cycle_iter, *(cycle_iter+1))]->incident_face =

			face_record_ptrs1[ face_record_indices1[e_iter - edges1.begin()] ];

	}

	// repeat for edges of second subdivision

	for (e_iter = edges2.begin(); e_iter != edges2.end(); e_iter++)

	{

		InputEdge::iterator cycle_iter = e_iter->begin();

		// set remaining data members in record for first half-edge in cycle

		half_edge_records[calc_key(*cycle_iter + vertices1.size(), *(cycle_iter + 1) + vertices1.size())]->prev = 

			half_edge_records[calc_key(*(e_iter->end() - 2) + vertices1.size(), *(e_iter->end() - 1) + vertices1.size())];

		half_edge_records[calc_key(*(cycle_iter) + vertices1.size(), *(cycle_iter + 1) + vertices1.size())]->next =

			half_edge_records[calc_key(*(cycle_iter + 1) + vertices1.size(), *(cycle_iter + 2) + vertices1.size())];

		half_edge_records[calc_key(*(cycle_iter) + vertices1.size(), *(cycle_iter + 1) + vertices1.size())]->twin =

			half_edge_records[calc_key(*(cycle_iter + 1) + vertices1.size(), *(cycle_iter) + vertices1.size())];

		half_edge_records[calc_key(*(cycle_iter) + vertices1.size(), *(cycle_iter + 1) + vertices1.size())]->incident_face =

			face_record_ptrs2[ face_record_indices2[ e_iter - edges2.begin() ] ];

		// set remaining data members in records for in-between half-edges

		for (cycle_iter++; cycle_iter != e_iter->end()-2; cycle_iter++)

		{

			half_edge_records[calc_key(*cycle_iter + vertices1.size(), *(cycle_iter+1) + vertices1.size())]->prev =

				half_edge_records[calc_key(*(cycle_iter-1) + vertices1.size(), *(cycle_iter) + vertices1.size())];	

			half_edge_records[calc_key(*cycle_iter + vertices1.size(), *(cycle_iter+1) + vertices1.size())]->next = 

				half_edge_records[calc_key(*(cycle_iter+1) + vertices1.size(), *(cycle_iter+2) + vertices1.size())];

			half_edge_records[calc_key(*cycle_iter + vertices1.size(), *(cycle_iter+1) + vertices1.size())]->twin = 

				half_edge_records[calc_key(*(cycle_iter+1) + vertices1.size(), *cycle_iter + vertices1.size())];

			half_edge_records[calc_key(*cycle_iter + vertices1.size(), *(cycle_iter+1) + vertices1.size())]->incident_face =

				face_record_ptrs2[ face_record_indices2[ e_iter - edges2.begin()] ];

		}

		// set remaining data members in records for last half-edge in cycle

		half_edge_records[calc_key(*cycle_iter + vertices1.size(), *(cycle_iter+1) + vertices1.size())]->prev = 

			half_edge_records[calc_key(*(cycle_iter-1) + vertices1.size(), *cycle_iter + vertices1.size())];

		half_edge_records[calc_key(*cycle_iter + vertices1.size(), *(cycle_iter+1) + vertices1.size())]->next =

			half_edge_records[calc_key(*(e_iter->begin()) + vertices1.size(), *(e_iter->begin()+1) + vertices1.size())];

		half_edge_records[calc_key(*cycle_iter + vertices1.size(), *(cycle_iter+1) + vertices1.size())]->twin =

			half_edge_records[calc_key(*(cycle_iter+1) + vertices1.size(),*cycle_iter + vertices1.size())];

		half_edge_records[calc_key(*cycle_iter + vertices1.size(), *(cycle_iter+1) + vertices1.size())]->incident_face =

			face_record_ptrs2[ face_record_indices2[ e_iter - edges2.begin()] ];

	}

}



// Precondition: *this must have been constructed as combination of two valid D_CELs A and B 

// Postcondition: *this is a valid D_CEL and the overlay of A and B	

void D_CEL::overlay()

{

	// initialise event queue

	for (HalfEdgeList::const_iterator c_he_iter = half_edge_records.begin(); 

		  c_he_iter != half_edge_records.end(); c_he_iter++)

	{	event_queue.insert(make_pair(c_he_iter->second->origin->coordinates, 

											  QueueElementValue(c_he_iter->second, c_he_iter->second->twin)));

//	cout << c_he_iter->second->origin->vertex << "-"

//		  << c_he_iter->second->next->origin->vertex << "-" 

//		  << c_he_iter->second->origin->coordinates.first << "-"

//		  << c_he_iter->second->origin->coordinates.second << " " ;

	}

//	cout << "EQ SIZE: " << event_queue.size() << endl;

//	for (EventQueue::const_iterator iter = event_queue.begin(); iter != event_queue.end(); iter++)

//	{

//			cout << iter->first.first << "-"

//		  << iter->first.second << " "; 

//	}

	plane_sweep();		// carry out plane sweep, validating D_CEL and gathering information as prescribed

	build_bc_graph_and_label();	// build boundary cycle graph by traversing overlay D_CEL

}



// Precondition: vertex_records, half_edge_records and face_records contain valid d-cel component records. Vertices

//					  and faces may overlap, and edges may intersect but are not allowed to overlap.

// Postcondition: *this is a valid d-cel in that vertex_records, half_edge_records and face_records each contain

//						only valid entries, none of which overlap or intersect. 

void D_CEL::plane_sweep()

{

	EventQueue::iterator	ep_iter = event_queue.begin();

	// perform plane sweep in accordance with the description of the map overlay algorithm

	Vertex	*vPtr, *v1Ptr, *v2Ptr;

	long org_size, vertex_no;

	org_size = vertex_no = vertex_records.size();

	

	while (!event_queue.empty())

	{	Coordinates event_point = (event_queue.begin())->first;

		vPtr = new Vertex;		

		v1Ptr = v2Ptr = 0; 



		vPtr->coordinates = event_point;	// insert coordinates of event point into new vertex

		vPtr->incident_edge = ep_iter->second.first; 

		vPtr->vertex = vertex_no++;

		

		while ((ep_iter = event_queue.begin())->first == event_point)	// while still at the current event point

		{ 	//cout << ep_iter->second.second->origin->vertex << ","

		  	//	  <<  ep_iter->second.second->next->origin->vertex << " ";

/*		cout << "Event Queue" << endl;

		for (EventQueue::const_iterator iter = event_queue.begin(); iter != event_queue.end(); iter++)

		{

			cout << iter->first.first << "-"

		  << iter->first.second << " "; 

		}

*/

			EventType	event = event_type(ep_iter->second, vPtr->coordinates);

			if (event==error_event)	// bør evt. håndteres med exception

			{	event_queue.erase(ep_iter);

				continue;	}

			if (event==upper_endpoint || event==lower_endpoint)

			{	

				if (ep_iter->second.first->origin->coordinates != event_point) 

				{	event_queue.erase(ep_iter);

					continue;	}	// continue if event has become spurious (as a result of intersection validation)

				// store references to subdivision vertices that are to be replaced by *vPtr

				if (v1Ptr && v1Ptr != ep_iter->second.first->origin)	// has v1Ptr been taken in an earlier iteration ?

				{	v2Ptr = ep_iter->second.first->origin;					// yes, occupy v2Ptr instead

					vPtr->auxiliary_data.subdivision = intersection_between_vertices;	// signal to the labelling stage that this is a special type of vertex

					// vertex-vertex intersection occurred so proceed to validate d_cel at vPtr->coordinates

					validate_vertex_vertex(ep_iter->second.first, vPtr);

				}

 				else  

				{	v1Ptr = ep_iter->second.first->origin;

					vPtr->incident_edge = v1Ptr->incident_edge; // store reference to half-edge from one of the subdivisions

					vPtr->auxiliary_data.subdivision = ep_iter->second.first->origin->auxiliary_data.subdivision;

				}

			}

			if (event==edge_edge_intersection)

			{		// retrieve and store information about the faces to which *vPtr belongs

					//cout << ep_iter->second.first->origin->vertex << "-"

					//	  << ep_iter->second.first->twin->origin->vertex << ","

					//	  << ep_iter->second.second->origin->vertex << "-"

					//	  << ep_iter->second.second->twin->origin->vertex << endl;

					if (ep_iter->second.first->origin->auxiliary_data.subdivision == endpoint_from_subdiv1)

					{ 

						label_new_vertex(ep_iter, vPtr->auxiliary_data.s1_face, vPtr->auxiliary_data.s2_face);

					}

					else

					{  

						label_new_vertex(ep_iter, vPtr->auxiliary_data.s2_face, vPtr->auxiliary_data.s1_face);

					}

			}

			validate_d_cel(event, ep_iter, vPtr);	

			event_queue.erase(ep_iter);

		}

/*		cout << "Sweep status structure" << endl;

		for (SweepStatus::iterator s_iter = sweep_status.begin(); s_iter != sweep_status.end(); s_iter++)

		{	cout << "(" << s_iter->first.upper.first << "," << s_iter->first.upper.second << " "

				  << s_iter->first.lower.first << "," << s_iter->first.lower.second << ") ";	}

		cout << "Count for " << vPtr->vertex << " :" << sweep_status.size();

*/		// if *vPtr came from an input subdivision, then determine the face from the other subdivision containing it

		if (v1Ptr && !v2Ptr)	

		{	SweepStatus::iterator sw_iter = sweep_status.upper_bound(

																	SweepElementKey(Coordinates(event_point.first+1, event_point.second+1),

																						 Coordinates(event_point.first+1, event_point.second+2)));	

			while (sw_iter != sweep_status.end()) 

			{	{	if (vPtr->auxiliary_data.subdivision == endpoint_from_subdiv1)

					{	if (sw_iter->second->twin->origin->auxiliary_data.subdivision == endpoint_from_subdiv2)

						{	vPtr->auxiliary_data.s2_face = sw_iter->second->twin->incident_face->face;

							break;}

					}

					else

					{	if (sw_iter->second->twin->origin->auxiliary_data.subdivision == endpoint_from_subdiv1)

						{	vPtr->auxiliary_data.s1_face = sw_iter->second->twin->incident_face->face;	

							break;}

					}

					sw_iter++;

				} 

			}

		}

		pair<SweepStatus::iterator, SweepStatus::iterator>	sw_iters;

		// use equal_range to determine set of line segments whose upper endpoint coincide with the event point

		sw_iters = equal_range(sweep_status.begin(), sweep_status.end(), 

							   // wrap the event point in a dummy sweep status element to allow use with SweepCmp predicate	

							   make_pair(SweepElementKey(event_point, Coordinates(event_point.first,event_point.second+1)),static_cast<HalfEdge*>(0)),

							   // pass predicate function object for custom comparison by calling nested class constructor

 							   SweepCmp());

/*		s_iter = sw_iters.first;

		cout << "(" << s_iter->first.upper.first << "," << s_iter->first.upper.second << " "

			  << s_iter->first.lower.first << "," << s_iter->first.lower.second << ") ";	

*/				

		// vPtr now points to a newly allocated vertex with the coordinates of the event point

		// *v1Ptr and *v2Ptr have been disconnected from the D-CEL. They will be deleted from vertex_records

		// upon completion of the plane sweep.

		vertex_records.push_back(vPtr);

		Coordinates new_event;

		// the naming of sw_left, sw_right reflects that event points denoting intersection are always stored

		// leftmost edge first. This invariant makes later processing (validation, mainly) simpler.

		SweepStatus::iterator sw_left, sw_right;		

		if (sw_iters.first == sw_iters.second)	// if no upper endpoints or intersection are coincident with event point	

		{	Coordinates	new_event;

			if (sw_iters.first != sweep_status.begin() && 

				sw_iters.second != sweep_status.end())

			{	sw_left = --sw_iters.first;

				sw_right = sw_iters.second;

				// update vPtr with reference to (successor of) half-edge to its immediate left.

				vPtr->auxiliary_data.left_half_edge = sw_left->second;

					// test for intersection and record new event, if any

				if (find_new_event(sw_left, sw_right, event_point, new_event)) 

				{	event_queue.insert(make_pair(new_event, make_pair(sw_left->second, sw_right->second)));	

				}

			}

		}

		else

		{	

			if (sw_iters.first != sweep_status.begin())	// if there is an edge E that intersects the sweep line, and which is to the left 

			{	sw_right = sw_iters.first;						// of the leftmost edge E' whose upper endpoint == event_point then check E,E' for intersection

				sw_left = --sw_iters.first;

				// same as above

				vPtr->auxiliary_data.left_half_edge = sw_left->second; 

				if (find_new_event(sw_left,sw_right, event_point, new_event))

				{	event_queue.insert(make_pair(new_event, make_pair(sw_left->second, sw_right->second)));	

				}	

			}

			if (sw_iters.second != sweep_status.end()) // if there is an edge E that intersects the sweep line, and which is to the right

			{	sw_right = sw_iters.second;			   // of the rightmost edge E' whose upper endpoint == event_point then check E,E' for intersection	

				sw_left = --sw_iters.second;

				if (find_new_event(sw_left, sw_right, event_point, new_event))

				{	event_queue.insert(make_pair(new_event, make_pair(sw_left->second,sw_right->second)));	

				}

			}

		}

	}

	// delete vertices from input D_CELs (now replaced by overlay vertices)

	vertex_records.erase(vertex_records.begin(), vertex_records.begin() + org_size);

}



// Private member function for calculating key value for entries in HalfEdgeList containers

// Appropriate pairing function could have been used but would have been less efficient

inline HalfEdgeListKey D_CEL::calc_key (long origin, long destination) const

{	// shift origin (number of bits in a long / 2) places to the left and add destination

	return ((origin + 1) << (std::numeric_limits<HalfEdgeListKey>::digits >> 2)) + destination + 1;

}



// Private member function for determining the type of an event point

inline D_CEL::EventType D_CEL::event_type(const QueueElementValue& qev, const Coordinates& event_point) const

{

	if (qev.second != qev.first->twin)

	{	if (qev.first->origin->coordinates == event_point ||

			 qev.first->twin->origin->coordinates == event_point ||

			 qev.second->origin->coordinates == event_point ||

			 qev.second->twin->origin->coordinates == event_point)

		{	return edge_vertex_intersection;	}

		else

		{	return edge_edge_intersection;	}

	}

	if (EqKeyCmp()(qev.first->origin->coordinates, qev.second->origin->coordinates))

	{	return upper_endpoint;	}

	if (EqKeyCmp()(qev.second->origin->coordinates, qev.first->origin->coordinates))

	{	return lower_endpoint;	}

	return error_event;

}



// private member function which finds intersection point, if any, between two segments and 

// stores the coordinates in event_point. Return value is true if intersection is found, false otherwise.

inline bool D_CEL::find_new_event(SweepStatus::iterator sw_iter1, SweepStatus::iterator sw_iter2, Coordinates& event_point, Coordinates& new_event)

{	

	Coordinates		up1 = sw_iter1->first.upper, low1 = sw_iter1->first.lower, 

						up2 = sw_iter2->first.upper, low2 = sw_iter2->first.lower;

	// check for endpoint overlap

	if (up1 == up2 || up1 == low2 || low1 == up2 || low1 == low2)

	{	return false;	}

	double slope1, slope2, disp1, disp2;

	double intersect_x, intersect_y;

	if (up1.first != low1.first)

	{	slope1 = 1.0 / sw_iter1->first.slope;	// slope stored with key is dx/dy, so invert to get dy/dx

		if (up2.first != low2.first)

		{	slope2 = 1.0 / sw_iter2->first.slope;	

			if (slope1 == slope2) 

			{	return false;

			}

			disp1 = up1.second - slope1*up1.first;

			disp2 = up2.second - slope2*up2.first;

			intersect_x = -(disp1-disp2) / (slope1-slope2);

			intersect_y = disp1 + slope1*intersect_x;

		}

		else

		{

			intersect_x = up2.first;

			

			disp1 = up1.second - slope1*up1.first;

			intersect_y = disp1 + slope1*intersect_x;

		}

	}

	else

	{

		intersect_x = up1.first;

		if (up2.first != low2.first)

		{

			slope2 = 1.0 / (sw_iter2->first.slope > 0 ? sw_iter2->first.slope : -sw_iter2->first.slope);

			intersect_y = up2.second + slope2*(intersect_x > up2.first ? intersect_x - up2.first : up2.first - intersect_x);

		}

		else

		{	return false;

		}

	}

	if ( (up1.first-intersect_x)*(intersect_x-low1.first)>=0 

	  && (up2.first-intersect_x)*(intersect_x-low2.first)>=0

	  && (up1.second-intersect_y)*(intersect_y-low1.second)>=0

	  && (up2.second-intersect_y)*(intersect_y-low2.second)>=0)

	{	new_event.first = intersect_x;

		new_event.second = intersect_y;

		

		return (new_event.second > event_point.second || 

				 (new_event.second == event_point.second &&  new_event.first > event_point.first)) ;

	}

	return false;

}



// function returns left if line segment p0-p1 is to the left of p0-p2 (when viewed from p0), right if it's not

// and error if neither applies. Taken from Sedgewick's "Algorithms in C++" and also documented in "Introduction to Algorithms".

inline D_CEL::Orientation D_CEL::turn(const Coordinates& p0, const Coordinates& p1, const Coordinates& p2) const

{

	double dx1, dx2, dy1, dy2;

	dx1 = p1.first - p0.first; 

	dx2 = p2.first - p0.first;

	dy1 = p1.second - p0.second;

	dy2 = p2.second - p0.second;

	if (dx1*dy2 > dy1*dx2) return left;	//turn left

	if (dx1*dy2 < dy1*dx2) return right;	//turn right

	if ((dx1*dx2 < 0) || (dy1*dy2 < 0)) return right;

	if ((dx1*dx1+dy1*dy1) < (dx2*dx2+dy2*dy2)) return left;

	return error;

}



inline void D_CEL::label_new_vertex(EventQueue::iterator ep_iter, string& face1, string& face2)

{

		// determine which of ep_iter->second.first and ep_iter->second.second has the smaller slope 

		// note that the origin vertices of the intersecting half-edges referred to by ep_iter are guaranteed to also be their upper endpoint

		if ((calc_slope(ep_iter->second.first->origin->coordinates, ep_iter->second.first->twin->origin->coordinates)) <

			 (calc_slope(ep_iter->second.second->origin->coordinates, ep_iter->second.second->twin->origin->coordinates)))

		{	face1 = ep_iter->second.first->incident_face->face;

			// use ep_iter->second.first to find out whether vPtr is a possible leftmost vertex of face above or beneath it

			if (ep_iter->second.first->origin->coordinates.first < ep_iter->second.first->twin->origin->coordinates.first)

			{	// vPtr may be leftmost vertex of face beneath it so we'll be needing right-left half-edge of ep_iter->second.second pair

				if (ep_iter->second.second->origin->coordinates.first < ep_iter->second.second->twin->origin->coordinates.first)

				{	face2 = ep_iter->second.second->twin->incident_face->face;	}

				else 

				{	face2 = ep_iter->second.second->incident_face->face;	}

			}

			else

			{	// vPtr may be leftmost vertex of face above it so we'll be needing left-right half-edge of ep_iter->second.second pair

				if (ep_iter->second.second->origin->coordinates.first > ep_iter->second.second->twin->origin->coordinates.first)

					{	face2 = ep_iter->second.second->twin->incident_face->face;	}

				else

					{	face2 = ep_iter->second.second->incident_face->face;	}

			}

		}

		else

		{	// use ep_iter->second.second to find out whether vPtr is a possible leftmost vertex of face above or beneath it

			if (ep_iter->second.second->origin->coordinates.first < ep_iter->second.second->twin->origin->coordinates.first)

			{	// vPtr may be leftmost vertex of face beneath it so we'll be needing right-left half-edge of ep_iter->second.first pair

				if (ep_iter->second.first->origin->coordinates.first < ep_iter->second.first->twin->origin->coordinates.first)

				{	face1 = ep_iter->second.first->twin->incident_face->face;	}

				else 

				{	face1 = ep_iter->second.first->incident_face->face;	}

			}

			else

			{	// vPtr may be leftmost vertex of face above it so we'll be needing left-right half-edge of ep_iter->second.first pair

				if (ep_iter->second.first->origin->coordinates.first > ep_iter->second.first->twin->origin->coordinates.first)

				{	face1 = ep_iter->second.first->twin->incident_face->face;	}

				else

				{	face1 = ep_iter->second.first->incident_face->face;	}

			}

			face2 = ep_iter->second.second->incident_face->face;

		}

}



// Precondition: vPtr marks the point of intersection between two vertices from distinct subdivisions.

// vPtr->incident_edge is a reference to a half-edge from one subdivision, and old_hePtr_first->incident_edge

// is a reference to half-edge from the other.

// Postcondition: *this is validated at the point of intersection.

void D_CEL::validate_vertex_vertex(HalfEdge* old_hePtr_first, Vertex *vPtr)

{

	static SweepStatus	incident_lower_endpoints, incident_upper_endpoints;

	if (vPtr->auxiliary_data.upper1.empty() && vPtr->auxiliary_data.lower1.empty())	// if this is the first time execution reaches this point for *vPtr

	{

		// first, build 4 ordered sequences for use in the labelling stage and also in the next/prev pointer update

		// code a little further below. Note that the sequences are stored with *vPtr so that they are accessible later.

		HalfEdge	*tmp_hePtr, *start_hePtr;

		//tmp_hePtr = vPtr->incident_edge->twin;

		tmp_hePtr = start_hePtr =	old_hePtr_first->origin->auxiliary_data.subdivision == endpoint_from_subdiv1 ? 

											old_hePtr_first->twin : vPtr->incident_edge->twin;

		do	// index by pair(other endpoint, vPtr->coordinates)

		{	if (tmp_hePtr->origin->coordinates.second < vPtr->coordinates.second)

			{	vPtr->auxiliary_data.lower1.insert(make_pair(SweepElementKey(tmp_hePtr->origin->coordinates, tmp_hePtr->twin->origin->coordinates),

																			tmp_hePtr));	}

			else

			{	vPtr->auxiliary_data.upper1.insert(make_pair(SweepElementKey(tmp_hePtr->twin->origin->coordinates, tmp_hePtr->origin->coordinates),

																			tmp_hePtr));	}

		tmp_hePtr = tmp_hePtr->next->twin;

		} while (tmp_hePtr != start_hePtr);

		// tmp_hePtr = old_hePtr_first->twin;

		tmp_hePtr = start_hePtr =	old_hePtr_first->origin->auxiliary_data.subdivision == endpoint_from_subdiv2 ? 

											old_hePtr_first->twin : vPtr->incident_edge->twin;

		do	// index by the endpoint != vPtr->coordinates

		{	if (tmp_hePtr->origin->coordinates.second < vPtr->coordinates.second)

			{	vPtr->auxiliary_data.lower2.insert(make_pair(SweepElementKey(tmp_hePtr->origin->coordinates, tmp_hePtr->twin->origin->coordinates),

																			tmp_hePtr));	}

			else

			{	vPtr->auxiliary_data.upper2.insert(make_pair(SweepElementKey(tmp_hePtr->twin->origin->coordinates, tmp_hePtr->origin->coordinates),

																			tmp_hePtr));	}

		tmp_hePtr = tmp_hePtr->next->twin;

		} while (tmp_hePtr != start_hePtr);



		incident_lower_endpoints.clear();

		incident_upper_endpoints.clear();

		// use the information about half-edges from both subdivisions incident to *vPtr by an upper or lower endpoint

		// to build two ordered sequences so that the neighbours of *old_hePtr_first can be determined and next/prev pointers updated

		merge(vPtr->auxiliary_data.lower1.begin(), vPtr->auxiliary_data.lower1.end(), 

				vPtr->auxiliary_data.lower2.begin(), vPtr->auxiliary_data.lower2.end(),

				inserter(incident_lower_endpoints, incident_lower_endpoints.begin()));

		merge(vPtr->auxiliary_data.upper1.begin(), vPtr->auxiliary_data.upper1.end(),

				vPtr->auxiliary_data.upper2.begin(), vPtr->auxiliary_data.upper2.end(),

				inserter(incident_upper_endpoints, incident_upper_endpoints.begin()));

	}

	// find *old_hePtr_first's neighbours and update next/prev pointers where appropriate

	pair<SweepStatus::iterator, SweepStatus::iterator>	neighbours;

	// if *old_hePtr_first has vPtr as its lower endpoint

	if (old_hePtr_first->twin->origin->coordinates.second < vPtr->coordinates.second)

	{	neighbours = incident_lower_endpoints.equal_range(SweepElementKey(old_hePtr_first->twin->origin->coordinates, old_hePtr_first->origin->coordinates));

		if (neighbours.first == incident_lower_endpoints.begin())	// leftmost lower endpoint reached ?

		{	if (!incident_upper_endpoints.empty())	// any upper endpoints incident on *vPtr ?

			{	old_hePtr_first->prev = incident_upper_endpoints.begin()->second;	

				incident_upper_endpoints.begin()->second->next = old_hePtr_first;	}	// connect to leftmost

			else

			{	old_hePtr_first->prev = (--incident_lower_endpoints.end())->second;	

				(--incident_lower_endpoints.end())->second->next = old_hePtr_first;	} // connect to rightmost lower endpoint

		}

		else

		{	--neighbours.first;

			old_hePtr_first->prev = neighbours.first->second;

			neighbours.first->second->next = old_hePtr_first;

		}

		if (neighbours.second == incident_lower_endpoints.end())	// rightmost lower endpoint reached ?

		{	if (!incident_upper_endpoints.empty())	// any upper endpoints incident on *vPtr ?

			{	old_hePtr_first->twin->next = (--incident_upper_endpoints.end())->second->twin;

				(--incident_upper_endpoints.end())->second->twin->prev = old_hePtr_first->twin;

			}

			else

			{	old_hePtr_first->twin->next = incident_lower_endpoints.begin()->second->twin;

				incident_lower_endpoints.begin()->second->twin->prev = old_hePtr_first->twin;

			}

		}

		else

		{	old_hePtr_first->twin->next = neighbours.second->second->twin;

			neighbours.second->second->twin->prev = old_hePtr_first->twin;

		}

	}

	else // *old_hePtr_first has *vPtr as its upper endpoint

	{	neighbours = incident_upper_endpoints.equal_range(SweepElementKey(old_hePtr_first->origin->coordinates, old_hePtr_first->twin->origin->coordinates));	

		if (neighbours.first == incident_upper_endpoints.begin())	// leftmost upper endpoint reached ?

		{	if (!incident_lower_endpoints.empty())	// any lower endpoints incident on *vPtr ?

			{	old_hePtr_first->twin->next = incident_lower_endpoints.begin()->second->twin;	

				incident_lower_endpoints.begin()->second->twin->prev = old_hePtr_first->twin;	}	// connect to leftmost

			else

			{	old_hePtr_first->twin->next = (--incident_upper_endpoints.end())->second->twin;	

				(--incident_lower_endpoints.end())->second->twin->prev = old_hePtr_first->twin;	} // connect to rightmost lower endpoint

		}

		else

		{	--neighbours.first; 

			old_hePtr_first->twin->next = neighbours.first->second->twin;

			neighbours.first->second->twin->prev = old_hePtr_first->twin;

		}

		if (neighbours.second == incident_upper_endpoints.end())	// rightmost upper endpoint reached ?

		{	if (!incident_lower_endpoints.empty())	// any lower endpoints incident on *vPtr ?

			{	old_hePtr_first->prev = (--incident_lower_endpoints.end())->second;

				(--incident_lower_endpoints.end())->second->next = old_hePtr_first;

			}

			else

			{	old_hePtr_first->prev = incident_upper_endpoints.begin()->second;

				incident_upper_endpoints.begin()->second->prev = old_hePtr_first;

			}

		}

		else

		{	old_hePtr_first->prev = neighbours.second->second;

			neighbours.second->second->next = old_hePtr_first;

		}

	}

	return;

}



// Precondition: vPtr marks the point of intersection between a half-edge and a vertex.

//					  *new_hePtr_first and *new_hePtr_first_twin are new half-edges, the origin vertices of which are

//					  *vPtr and *(old_hePtr_first->twin->origin), respectively.

// Postcondition: *this is validated at the point of intersection

void D_CEL::validate_edge_vertex(HalfEdge* old_hePtr_first, HalfEdge* old_hePtr_second, 

											HalfEdge* new_hePtr_first, HalfEdge* new_hePtr_first_twin, Vertex* vPtr)

{

	HalfEdge*	tmp_hePtr;

	tmp_hePtr = old_hePtr_second;

	// find first half-edge in counter-clockwise direction from old_hePtr_first

	do

	{

		if ((turn(tmp_hePtr->origin->coordinates, vPtr->coordinates, tmp_hePtr->twin->prev->origin->coordinates) == right) ||

			 (turn(tmp_hePtr->origin->coordinates, vPtr->coordinates, old_hePtr_first->origin->coordinates) == left &&

			  turn(old_hePtr_first->origin->coordinates, vPtr->coordinates, tmp_hePtr->twin->prev->origin->coordinates) == left))

		{	break;	}

		tmp_hePtr = tmp_hePtr->twin->prev;

	} while (tmp_hePtr != old_hePtr_second);

	// update next and prev pointers of existing half-edges

	if (turn(old_hePtr_first->origin->coordinates, vPtr->coordinates, tmp_hePtr->twin->prev->origin->coordinates) == right)

	{	old_hePtr_first->twin->prev = new_hePtr_first_twin;

		new_hePtr_first_twin->next = old_hePtr_first->twin;

	}

	else

	{	old_hePtr_first->twin->prev = tmp_hePtr->twin->prev;

		tmp_hePtr->twin->prev->next = old_hePtr_first->twin;

	}

	old_hePtr_first->next = tmp_hePtr->twin;

	tmp_hePtr->twin->prev = old_hePtr_first;

	// find first half-edge in counter-clockwise direction from new_hePtr_first

	tmp_hePtr = old_hePtr_second;

	do

	{

		if ((turn(tmp_hePtr->origin->coordinates, vPtr->coordinates, tmp_hePtr->next->twin->origin->coordinates) == left) ||

			 (turn(tmp_hePtr->origin->coordinates, vPtr->coordinates, old_hePtr_first->twin->origin->coordinates) == right &&

			  turn(old_hePtr_first->twin->origin->coordinates, vPtr->coordinates, tmp_hePtr->next->twin->origin->coordinates) == right))

		{	break;	}

	tmp_hePtr = tmp_hePtr->next->twin;

	} while (tmp_hePtr != old_hePtr_second);

	// update pointers

	if (turn(new_hePtr_first_twin->origin->coordinates, vPtr->coordinates, tmp_hePtr->next->twin->origin->coordinates) == left)

		{	new_hePtr_first_twin->next = old_hePtr_first->twin;

			old_hePtr_first->twin->prev = new_hePtr_first_twin;

		}

	else

		{	new_hePtr_first_twin->next = tmp_hePtr->next;

			tmp_hePtr->next->prev = new_hePtr_first_twin;

		}

	new_hePtr_first->prev = tmp_hePtr;

	tmp_hePtr->next = new_hePtr_first;

}



// Precondition: *vPtr marks an intersection between half-edges *old_hePtr_first and *old_hePtr_second.

//					  *new_hePtr_first and *new_hePtr_first_twin are new half-edges.

// Postcondition: *this is validated at the point of intersection.	

void D_CEL::validate_edge_edge(HalfEdge* old_hePtr_first, HalfEdge* old_hePtr_second, 

										 HalfEdge* new_hePtr_first, HalfEdge* new_hePtr_first_twin, Vertex* vPtr)

{

	HalfEdge	*new_hePtr_second = new HalfEdge, *new_hePtr_second_twin = new HalfEdge;

	// repeat the above manoeuvres for the other half-edge involved in the intersection

	// 1st new half-edge is connected to origin of existing half-edge and given the properties of that half-edge

	new_hePtr_second->origin = vPtr;

	new_hePtr_second->incident_face = old_hePtr_second->incident_face;

	new_hePtr_second->twin = new_hePtr_second_twin;

	// 2nd new half-edge is connected to new vertex and given properties of twin of existing half-edge

	new_hePtr_second_twin->origin = old_hePtr_second->twin->origin;

	new_hePtr_second_twin->incident_face = old_hePtr_second->twin->incident_face;

	new_hePtr_second_twin->twin = new_hePtr_second;

	// iterators to existing half-edge and its twin are retrieved before update

	HalfEdgeList::iterator he_iter_second = half_edge_records.find(calc_key(old_hePtr_second->origin->vertex, old_hePtr_second->twin->origin->vertex));



	HalfEdgeList::iterator he_iter_second_twin = half_edge_records.find(calc_key(old_hePtr_second->twin->origin->vertex, 

																		  old_hePtr_second->origin->vertex));

	// update next and prev pointers of old and new half-edges involved in the intersection

	// note that due to the way in which intersection event points are recorded, it can safely be assumed that old_hePtr_first

	// is the leftmost edge involved in the intersection

	new_hePtr_second->next = old_hePtr_second->next;

	new_hePtr_second->next->prev = new_hePtr_second;			



	new_hePtr_second_twin->prev = old_hePtr_second->twin->prev;

	new_hePtr_second_twin->prev->next = new_hePtr_second_twin;

	new_hePtr_first->prev = old_hePtr_second;

	old_hePtr_second->next = new_hePtr_first;	

	

	old_hePtr_second->twin->prev = old_hePtr_first;

	old_hePtr_first->next = old_hePtr_second->twin;

		

	old_hePtr_first->twin->prev = new_hePtr_second_twin;

	new_hePtr_second_twin->next = old_hePtr_first->twin;

	

	new_hePtr_second->prev = new_hePtr_first_twin;

	new_hePtr_first_twin->next = new_hePtr_second;

	// update existing half-edge with intersection vertex

	old_hePtr_second->twin->origin = vPtr;

	// old half-edge records must be deleted before insertion with updated key values

	half_edge_records.erase(he_iter_second);

	half_edge_records.erase(he_iter_second_twin);

	// insert records for old (but now updated) half-edges

	half_edge_records.insert(make_pair(calc_key(old_hePtr_second->origin->vertex, old_hePtr_second->twin->origin->vertex), old_hePtr_second));

	half_edge_records.insert(make_pair(calc_key(old_hePtr_second->twin->origin->vertex, old_hePtr_second->origin->vertex), old_hePtr_second->twin));

	// insert records for new half-edges

	half_edge_records.insert(make_pair(calc_key(new_hePtr_second->origin->vertex, new_hePtr_second_twin->origin->vertex), new_hePtr_second));

	half_edge_records.insert(make_pair(calc_key(new_hePtr_second_twin->origin->vertex, new_hePtr_second->origin->vertex), new_hePtr_second_twin));

	// update event queue with lower endpoints of newly created half-edges 

	event_queue.insert(make_pair(new_hePtr_second_twin->origin->coordinates,

							 QueueElementValue(new_hePtr_second_twin, new_hePtr_second)));



	// add new half-edges to sweep status structure

	sweep_status.insert(make_pair(SweepElementKey(new_hePtr_second->origin->coordinates,

												  new_hePtr_second_twin->origin->coordinates), new_hePtr_second));

}	



// Precondition: event designates the type of intersection to be validated. *vPtr marks the point of intersection.

//					  ep_iter points to an event queue entry which references the half-edges involved in the intersection.

// Postcondition: *this is validated at the point of intersection.	

void D_CEL::validate_d_cel(EventType event, EventQueue::iterator ep_iter, Vertex* vPtr)

{

	HalfEdgeList::iterator he_iter_first, he_iter_second;

	HalfEdge	*old_hePtr_first = ep_iter->second.first, *old_hePtr_second = ep_iter->second.second;



	he_iter_first = half_edge_records.find(calc_key(old_hePtr_first->origin->vertex, old_hePtr_first->twin->origin->vertex));

	if (event == edge_edge_intersection || event == edge_vertex_intersection)

	{	// shared code for edge-edge and edge-vertex type intersections...

		HalfEdge	*new_hePtr_first = new HalfEdge, *new_hePtr_first_twin = new HalfEdge;

		if (old_hePtr_first->twin->origin == vPtr)	// if intersection has already been validated

		{	return;	}

		// remove intersecting half-edges from sweep status structure

		sweep_status.erase(SweepElementKey(old_hePtr_first->origin->coordinates, old_hePtr_first->twin->origin->coordinates));	

		sweep_status.erase(SweepElementKey(old_hePtr_second->origin->coordinates, old_hePtr_second->twin->origin->coordinates));

		// 1st new half-edge is connected to origin of existing half-edge and assigned its member data

		new_hePtr_first->origin = vPtr;

		new_hePtr_first->incident_face = old_hePtr_first->incident_face;

		new_hePtr_first->twin = new_hePtr_first_twin;

		// 2nd new half-edge is connected to new vertex and assigned some of its member data

		new_hePtr_first_twin->origin = old_hePtr_first->twin->origin;

		new_hePtr_first_twin->incident_face = old_hePtr_first->twin->incident_face;

		new_hePtr_first_twin->twin = new_hePtr_first;

		// update the prev and next pointers for which information is available at this time

		new_hePtr_first->next = old_hePtr_first->next;

		new_hePtr_first->next->prev = new_hePtr_first;



		new_hePtr_first_twin->prev = old_hePtr_first->twin->prev;

		new_hePtr_first_twin->prev->next = new_hePtr_first_twin;

		// iterator to twin of existing half-edge is retrieved before update

		HalfEdgeList::iterator he_iter_first_twin = half_edge_records.find(calc_key(old_hePtr_first->twin->origin->vertex, 

																											 old_hePtr_first->origin->vertex));



		switch(event)	// call code specific to each of the two similar types of intersection

		{	case edge_vertex_intersection	:	validate_edge_vertex(old_hePtr_first, old_hePtr_second, new_hePtr_first, new_hePtr_first_twin, vPtr);	

														break;

			case edge_edge_intersection	:	validate_edge_edge(old_hePtr_first, old_hePtr_second, new_hePtr_first, new_hePtr_first_twin, vPtr);

														break;

		}

		// more shared code for intersections of type edge_vertex and edge_edge...

		// update existing half-edge with new endpoint

		old_hePtr_first->twin->origin = vPtr;	

		// old half-edge records must be deleted before insertion with updated key values

		half_edge_records.erase(he_iter_first);

		half_edge_records.erase(he_iter_first_twin);

		// insert records for old (but now updated) half-edges

		half_edge_records.insert(make_pair(calc_key(old_hePtr_first->origin->vertex, old_hePtr_first->twin->origin->vertex), old_hePtr_first));

		half_edge_records.insert(make_pair(calc_key(old_hePtr_first->twin->origin->vertex, old_hePtr_first->origin->vertex), old_hePtr_first->twin));

		// insert records for first pair of new half-edges (second pair, if any, inserted by validate_edge_edge function)

		half_edge_records.insert(make_pair(calc_key(new_hePtr_first->origin->vertex, new_hePtr_first_twin->origin->vertex), new_hePtr_first));

		half_edge_records.insert(make_pair(calc_key(new_hePtr_first_twin->origin->vertex, new_hePtr_first->origin->vertex), new_hePtr_first_twin));

		// update event queue with lower endpoint of newly created half-edge duo

		event_queue.insert(make_pair(new_hePtr_first_twin->origin->coordinates, 

								 QueueElementValue(new_hePtr_first_twin, new_hePtr_first)));

		// add new half-edge duo to sweep status structure

		sweep_status.insert(make_pair(SweepElementKey(new_hePtr_first->origin->coordinates, 

													  new_hePtr_first_twin->origin->coordinates), new_hePtr_first));

		return;

	}

	he_iter_second = half_edge_records.find(calc_key(old_hePtr_second->origin->vertex, old_hePtr_second->twin->origin->vertex));

	// shared code for lower_endpoint and upper_endpoint cases

	old_hePtr_first->origin = vPtr;

	// note that iterators into a map structure such as half_edge_records are NOT invalidated by insertion or deletion

	half_edge_records.erase(he_iter_first);

	half_edge_records.erase(he_iter_second);

	// insert updated half-edges with key values calculated from new endpoints

	half_edge_records[calc_key(old_hePtr_first->origin->vertex, old_hePtr_second->origin->vertex)] = old_hePtr_first;

	half_edge_records[calc_key(old_hePtr_second->origin->vertex, old_hePtr_first->origin->vertex)] = old_hePtr_second;



	if (event == upper_endpoint)

	{		//cout << old_hePtr_second->origin->vertex << "/" << old_hePtr_second->twin->origin->vertex << " ";



		sweep_status.insert(make_pair(SweepElementKey(old_hePtr_first->origin->coordinates,

												      old_hePtr_second->origin->coordinates), old_hePtr_first));

	}

	if (event == lower_endpoint)

	{	sweep_status.erase(SweepElementKey(old_hePtr_second->origin->coordinates, old_hePtr_first->origin->coordinates));

	}

}



// Precondition: vertex_records, half_edge_records, and face_records form a valid D_CEL D

// Postcondition: a boundary cycle graph representing the faces and components of D has been computed and created

void D_CEL::build_bc_graph_and_label()

{

	vector<BCGraphNode*>	ob_cycles, hb_cycles;



	BCGraphNode		*new_node;

	HalfEdge			*hePtr, *leftmost_edgePtr;



	for (HalfEdgeList::iterator he_iter = half_edge_records.begin(); he_iter != half_edge_records.end(); he_iter++)

	{

		if (he_iter->second->BCG_nodePtr != 0) // if half-edge has been processed as part of earlier traversal 

		{	continue;	}

		new_node = new BCGraphNode;

		leftmost_edgePtr = hePtr = he_iter->second;

		int count = 0;

		do		// traverse cycle

		{

			hePtr->BCG_nodePtr = new_node;

			if (hePtr->origin->coordinates < leftmost_edgePtr->origin->coordinates) // if origin vertex is further to the left than current candidate

			{	leftmost_edgePtr = hePtr;	}

			hePtr = hePtr->next;

			count++;

		} while (hePtr != he_iter->second);

//		cout << "Edges in cycle:" << count << endl;

		new_node->leftmost_edge = leftmost_edgePtr;

		// check angle at leftmost vertex, and insert into ob_cycles if < 180 deg. or into hb_cycles if >

		Orientation turn_res  = turn(leftmost_edgePtr->next->origin->coordinates, leftmost_edgePtr->origin->coordinates, 

										leftmost_edgePtr->prev->origin->coordinates);

		if (turn_res == right)	

		{	

			hb_cycles.push_back(new_node);

		}

		if (turn_res == left)

		{	

			ob_cycles.push_back(new_node);

		}			

	}

	// connect leftmost half-edge of each cycle to the node representing that cycle 

	BCGraphNode*	unbounded_face = new BCGraphNode;

	for (vector<BCGraphNode*>::iterator hb_iter = hb_cycles.begin(); hb_iter != hb_cycles.end(); hb_iter++)

	{	HalfEdge*	lhePtr = (*hb_iter)->leftmost_edge->origin->auxiliary_data.left_half_edge;

		if (lhePtr == 0) // if cycle abuts on unbounded face

		{	unbounded_face->neighbours.push_back(*hb_iter);	}

		else

		{	lhePtr->BCG_nodePtr->neighbours.push_back(*hb_iter);	}

	} 

	ob_cycles.push_back(unbounded_face);

	// release memory allocated for existing face records

	Face*	fPtr;

	FaceList::size_type	num_faces = face_records.size();

	for (vector<BCGraphNode*>::iterator ob_iter = ob_cycles.begin(); ob_iter != ob_cycles.end(); ob_iter++)

	{	fPtr = new Face; 

		fPtr->outer_component = (*ob_iter)->leftmost_edge;

		if (fPtr->outer_component)	// if *fPtr is not the unbounded face, then assign new label to its face field

		{	VertexAux*	aux_data = &(*ob_iter)->leftmost_edge->origin->auxiliary_data;

			// cout << aux_data->s1_face << "," << aux_data->s2_face << " ";

			// the subdivision data member indicates where to get the labelling information

			switch (aux_data->subdivision)

			{

				case intersection_between_halfedges	:

								fPtr->face = aux_data->s1_face;

								fPtr->face += aux_data->s2_face;

								break;

				case endpoint_from_subdiv1	:	

								// the leftmost half-edge comes from input subdivision 1, so take first part of label from that

									// if at the lowermost half-edge in a "fan" connected to a single, shared leftmost vertex

									// SKAL FORBEDRES - check virker kun hvis vertex ikke er resultat af intersektion

								if ((*ob_iter)->leftmost_edge->twin->origin->auxiliary_data.subdivision == endpoint_from_subdiv2)

								{	fPtr->face = (*ob_iter)->leftmost_edge->prev->incident_face->face;	}

								else

								{	fPtr->face = (*ob_iter)->leftmost_edge->incident_face->face;	}

								fPtr->face += aux_data->s2_face;

								break;

				case endpoint_from_subdiv2	:	

								// the leftmost half-edge comes from input subdivision 2, so take second part of label from that

								fPtr->face = aux_data->s1_face;	

									// if at the lowermost half-edge in a "fan" connected to a single, shared leftmost vertex

									// SKAL FORBEDRES - check virker kun hvis vertex ikke er resultat af intersektion

								if ((*ob_iter)->leftmost_edge->twin->origin->auxiliary_data.subdivision == endpoint_from_subdiv1)

								{	fPtr->face += (*ob_iter)->leftmost_edge->prev->incident_face->face;	}	

								else												

								{	fPtr->face += (*ob_iter)->leftmost_edge->incident_face->face;	}

								break;

				case intersection_between_vertices	:	

								// handle case of vertex-vertex intersection (which disallows direct look-up of 

								// face labelling information since potentially O(n) half-edges from both subdivisions may

								// be incident on the overlay vertex marking the point of intersection)

								if ((*ob_iter)->leftmost_edge->next->origin->coordinates.second < 

									 (*ob_iter)->leftmost_edge->origin->coordinates.second) // if incident by lower endpoint

								{	SweepStatus::iterator edge_iter = aux_data->lower1.lower_bound(

																					SweepElementKey((*ob_iter)->leftmost_edge->next->origin->coordinates,

																										 (*ob_iter)->leftmost_edge->origin->coordinates));

									if (edge_iter == aux_data->lower1.end() ||

											(edge_iter == aux_data->lower1.begin() &&

											 edge_iter->second->twin != (*ob_iter)->leftmost_edge))	 

									{	if (!aux_data->upper1.empty())

										{	edge_iter = aux_data->upper1.begin();	}

										else

										{	edge_iter = --(aux_data->lower1.end());	}

									}

									else if (edge_iter->second->twin != (*ob_iter)->leftmost_edge) // if leftmost edge from subdivision 2

									{	--edge_iter;	}

									

									fPtr->face = edge_iter->second->twin->incident_face->face;	

									edge_iter = aux_data->lower2.lower_bound(

																					SweepElementKey((*ob_iter)->leftmost_edge->next->origin->coordinates,

																										 (*ob_iter)->leftmost_edge->origin->coordinates));

									if (edge_iter == aux_data->lower2.end() || 

											(edge_iter == aux_data->lower2.begin() && 

											 edge_iter->second->twin != (*ob_iter)->leftmost_edge))

									{	if (!aux_data->upper2.empty())

										{	edge_iter = aux_data->upper2.begin();	}

										else

										{	edge_iter = --(aux_data->lower2.end());	}

									}

									else if (edge_iter->second->twin != (*ob_iter)->leftmost_edge )	// if leftmost edge from subdivision 1

									{	--edge_iter;	}



									fPtr->face += edge_iter->second->twin->incident_face->face;

								}

								else

								{	// label face with the s1 face that contains it

									SweepStatus::iterator edge_iter = aux_data->upper1.lower_bound(

																					SweepElementKey((*ob_iter)->leftmost_edge->origin->coordinates,

																										 (*ob_iter)->leftmost_edge->next->origin->coordinates));

									if (edge_iter == aux_data->upper1.end() ||

											(edge_iter == aux_data->upper1.begin() &&

											 edge_iter->second->twin != (*ob_iter)->leftmost_edge))	 

									{	if (!aux_data->lower1.empty())

										{	edge_iter = aux_data->lower1.begin();	}

										else

										{	edge_iter = --(aux_data->upper1.end());	}

									}

									else if (edge_iter->second->twin != (*ob_iter)->leftmost_edge) // if leftmost edge from subdivision 2

									{	--edge_iter;	}

									fPtr->face = edge_iter->second->twin->incident_face->face;	

									// label face with the s2 face that contains it

									edge_iter = aux_data->upper2.lower_bound(

																					SweepElementKey((*ob_iter)->leftmost_edge->origin->coordinates,

																										 (*ob_iter)->leftmost_edge->next->origin->coordinates));

									// if there are no subdivision 2 edges to the left of edge_iter->second

									if (edge_iter == aux_data->upper2.end() || 

											(edge_iter == aux_data->upper2.begin() && 

											 edge_iter->second->twin != (*ob_iter)->leftmost_edge))

									{	if (!aux_data->lower2.empty())

										{	edge_iter = aux_data->lower2.begin();	}

										else

										{	edge_iter = --(aux_data->upper2.end());	}

									}

									else if (edge_iter->second->twin != (*ob_iter)->leftmost_edge )	// if leftmost edge from subdivision 1

									{	--edge_iter;	}

									

									fPtr->face += edge_iter->second->twin->incident_face->face; 

								}	

								break;

			}

			

		}

		build_inner_components((*ob_iter)->neighbours.begin(), (*ob_iter)->neighbours.end(), fPtr);

		set_incident_face(fPtr->outer_component, fPtr);

		face_records.push_back(fPtr);

	}	

	face_records.erase(face_records.begin(), face_records.begin() + num_faces);

}



// Precondition: hePtr is a pointer to a valid half-edge which is part of a face cycle. fPtr is a pointer to a valid face record.

// Postcondition: For every half-edge he in the cycle to which *hePtr belongs, it holds that he.incident_face == fPtr. 

void D_CEL::set_incident_face(HalfEdge* hePtr, Face* fPtr)

{

	if (!hePtr) return;

	HalfEdge*	tmp_hePtr = hePtr;

	do

	{

		tmp_hePtr->incident_face = fPtr;

		tmp_hePtr = tmp_hePtr->next;

	} while (tmp_hePtr != hePtr);

}

// recursive utility function for building inner components list of a face record from information

// represented by a BCGraph connected component

void D_CEL::build_inner_components(vector<BCGraphNode*>::iterator first, vector<BCGraphNode*>::iterator last, Face* fPtr)

{

	while (first != last)

	{	fPtr->inner_components.push_back((*first)->leftmost_edge);

		set_incident_face((*first)->leftmost_edge, fPtr);	// half-edges of hole boundary are incident to face being constructed

		build_inner_components((*first)->neighbours.begin(), (*first)->neighbours.end(), fPtr);

		first++;

	}

	return;

}



inline double D_CEL::calc_slope(const Coordinates& p0, const Coordinates& p1) const

{	// using operator ?: performs better than fabs(), I suspect 

	double res;

	return (res = (p1.first - p0.first) / (p1.second - p0.second)) >= 0 ? res : -res;

}