class VoronoiXY
	def initialize() ### TIG
	  @model = Sketchup.active_model
	  @state=0
	  self.voronoi()
	end
	def points_to_array()
	  ###
	  @model = Sketchup.active_model
	  ss = @model.selection
	  if ss.empty? 
		UI.messagebox("No Selection !")
		return nil
	  end
	  ### only process cpoints [TIG]
	  ss.to_a.each{|e|ss.remove(e) unless e.kind_of?(Sketchup::ConstructionPoint)}
	  if ss.empty? 
		UI.messagebox("No Guide-Points in the Selection !");
		return nil;
	  end
	  if not ss[2]
		UI.messagebox("Insufficient Guide-Points in the Selection [needs >= 3] !")
		return nil
	  end
	  ###
	  points_array = []   
	  ss.to_a.each{|cpt| ### TIG'd
		x = cpt.position[0].to_f
		y = cpt.position[1].to_f
		z = cpt.position[2].to_f
		points_array << [x, y ,z]
	  }
	  points_array.uniq!
	  return points_array
	end #of def
	def voronoi()
	###
	entities = @model.active_entities
	points_table = []
	points_table = self.points_to_array()
	###
	triangles = []
	if points_table[0]
	  triangles = Delauney3::triangulate(points_table)
	else
	  puts "No triangles can be made."
	  return nil
	end
	@model.start_operation("Triangulate Points")
	### group for tri's output
	point_groups=[]
	for i in 0..points_table.length-1
		point_groups[i]=[]
		triangles.each do |tri|
			if tri.include? i
				point_groups[i]<<tri
			end
		end
	end
	point_groups.each do |pg|
		pg.flatten!
		pg.uniq!
	end

	@model.commit_operation
	conhull=[]
	qhull(Sketchup.active_model.selection).each do |p|
		for i in 0..points_table.length-1
			if p.position==points_table[i]
				conhull<<i
			end
		end
	end
	for i in 0..points_table.length-1
		connpoints=point_groups[i]
		if conhull.include?(i)
			j=conhull.index(i)
			if j == 0
				connpoints<<conhull[1]
				connpoints<<conhull[conhull.length-1]				
			else
				if j == conhull.length-1
					connpoints<<conhull[j-1]
					connpoints<<conhull[0]
				else
					connpoints<<conhull[j-1]
					connpoints<<conhull[j+1]
				end
			end
		end
		connpoints.uniq!
		# UI.messagebox connpoints
		connpoints.delete(i)
		# UI.messagebox connpoints
		sorted=sortradial(connpoints,i,points_table)
		# UI.messagebox sorted
		cellpoints=voropoints(i,sorted,points_table,conhull)
		cellpoints=removedupilcates(cellpoints)
		# UI.messagebox cellpoints
		grent=Sketchup.active_model.entities.add_group().entities
		grent.add_face(cellpoints)
		#UI.messagebox "i = "+ i.to_s
	end

	end#def
	def removedupilcates(points)
	uniqpoints=[]
	temp=[]
	n=points.length
	for i in 0..n-1
		for j in i+1..n-1
			if points[i]==points[j]
				temp<<points[j]
			end
		end
	end
	uniqpoints=points-temp
	uniqpoints
	end
	def qhull(points)
		vector=Geom::Vector3d.new 1,0,0 ######################################################################### to create voronoi in other planes this vector has to be changed (if all the points are coplanar then this vector can be constructed as a vector on this plane , parallel to XY plane)
		pminy=[points[0]]
		for i in 1..points.length-1
			if points[i].position.y < pminy.last.position.y
				pminy<<points[i]
			end
		end
		vertices=[pminy.last]
		temp=[]
		angle=180.degrees
		for i in 0..points.length-1
			if (vertices.last.position!=points[i].position)&&(vector.angle_between(vertices.last.position.vector_to(points[i].position)) < angle)
				temp<<points[i]
				angle = vector.angle_between(vertices.last.position.vector_to(points[i].position))
			end
		end
		vertices<<temp.last
		angle=180.degrees
		for k in 1..points.length-1
			vector=vertices[vertices.length-2].position.vector_to(vertices.last.position)
			for i in 0..points.length-1
				if (vertices.last.position!=points[i].position)&&(vector.angle_between(vertices.last.position.vector_to(points[i].position))<angle)
					temp<<points[i]
					angle = vector.angle_between(vertices.last.position.vector_to(points[i].position))
				end
			end
			vertices<<temp.last
			angle=vector.angle_between(vertices[vertices.length-2].position.vector_to(vertices.last.position))+180.degrees
			if vertices.first==vertices.last
				break
			end	
		end
		vertices.delete_at(vertices.length-1)
		vertices
	end
	def sortradial(pointsind,bpointindex,pointstable)
		angles=[]
		sortedp=[]
		vec=Geom::Vector3d.new(1,0,0) ######################################################################### to create voronoi in other planes this vector has to be changed (if all the points are coplanar then this vector can be constructed as a vector on this plane , parallel to XY plane)
		pointsind.each do |p|
			if vec.cross(pointstable[bpointindex].vector_to(pointstable[p])).z>=0
				angles<<vec.angle_between(pointstable[bpointindex].vector_to(pointstable[p]))
			else
				angles<<360.degrees-vec.angle_between(pointstable[bpointindex].vector_to(pointstable[p]))
			end
		end
		angles.sort!
# UI.messagebox angles
		for i in 0..angles.length-1
			for j in 0..pointsind.length-1
				if vec.cross(pointstable[bpointindex].vector_to(pointstable[pointsind[j]])).z>=0
					if vec.angle_between(pointstable[bpointindex].vector_to(pointstable[pointsind[j]]))==angles[i]
						sortedp<<pointsind[j]
# UI.messagebox "sorted point type 0 "+pointsind[j].to_s
					break
					end
				else
					if 360.degrees-vec.angle_between(pointstable[bpointindex].vector_to(pointstable[pointsind[j]]))==angles[i]
						sortedp<<pointsind[j]
# UI.messagebox "sorted point type 1 "+pointsind[j].to_s
					break
					end
				end
				# if vec.angle_between(pointstable[bpointindex].vector_to(pointstable[pointsind[j]]))==angles[i]
					# sortedp<<pointsind[j]
# UI.messagebox "sorted point type 2 "+pointsind[j].to_s
					# break
				# end
			end
		end
		sortedp
	end
	def voropoints(bpntindx,sorpntindx,pointstable,conhulindx)
		polygon=[]
		conhulindx.each do |i|
			polygon<<pointstable[i]
		end
		points=[]
		vec=Geom::Vector3d.new(0,0,1) ######################################################################### to create voronoi in other planes this vector has to be changed (if all the points are coplanar , then this vector is a normal vector of the plane)
		for i in 0..sorpntindx.length-1
			if i<sorpntindx.length-1
				j=i+1
			else
				j=0
			end
			test=[]
			if (conhulindx.include? bpntindx)
				test<<0
			else
				test<<1
			end
			if (conhulindx.include? sorpntindx[i])
				test<<0
			else
				test<<1
			end
			if (conhulindx.include? sorpntindx[j])
				test<<0
			else	
				test<<1
			end
			linevector=[pointstable[bpntindx].vector_to(pointstable[sorpntindx[i]]),pointstable[bpntindx].vector_to(pointstable[sorpntindx[j]])]
			midpoints=[Geom.linear_combination(0.5,pointstable[bpntindx],0.5,pointstable[sorpntindx[i]]),Geom.linear_combination(0.5,pointstable[bpntindx],0.5,pointstable[sorpntindx[j]])]
			midlines=[[midpoints[0],vec.cross(linevector[0])],[midpoints[1],vec.cross(linevector[1])]]
			int=Geom.intersect_line_line([midlines[0][0],midlines[0][1]], [midlines[1][0],midlines[1][1]])
			if (test[0]==0)&&(test[1]==0)&&(test[2]==0)
				intprim=int+int.vector_to(pointstable[bpntindx])+int.vector_to(pointstable[bpntindx])
				points<<intprim.project_to_line([midlines[0][0],midlines[0][1]])
				points<<intprim
				points<<intprim.project_to_line([midlines[1][0],midlines[1][1]])
			else
				points<<int
			end
		end
		adjpoints=[]
		for i in 0..points.length-1
			if i==points.length-1
				j=0
			else
				j=i+1
			end
			if Geom.point_in_polygon_2D(points[i], polygon, true)
				adjpoints<<points[i]
Sketchup.active_model.entities.add_cpoint(adjpoints.last)
# UI.messagebox "pointtype 0"
				if Geom.point_in_polygon_2D(points[j], polygon, true)==false
					adjpoints<<Geom::Point3d.new(intersect_p_p_ps(points[i],points[j],polygon))
Sketchup.active_model.entities.add_cpoint(adjpoints.last)
# UI.messagebox "pointtype 1"
					if Geom.point_in_polygon_2D(pointstable[bpntindx], polygon, false)==false
						adjpoints<<Geom::Point3d.new(pointstable[bpntindx])
Sketchup.active_model.entities.add_cpoint(adjpoints.last)
# UI.messagebox "pointtype 2"
					end
				end
			else
				if Geom.point_in_polygon_2D(points[j], polygon, true)
					adjpoints<<Geom::Point3d.new(intersect_p_p_ps(points[i],points[j],polygon))
Sketchup.active_model.entities.add_cpoint(adjpoints.last)
# UI.messagebox "pointtype 3"
				else
					# if Geom.point_in_polygon_2D(pointstable[bpntindx], polygon, false)
						# temp=intersections_p_p_ps(points[i],points[j],polygon)
						# temp.each do |p|
							# adjpoints<<Geom::Point3d.new(p)
						# end
					# end
				end
			end
		end
		if adjpoints.length==0
			temp=intersections_ps_ps(points,polygon)
			if temp.length==2
				adjpoints<<Geom::Point3d.new(temp[0])
				adjpoints<<Geom::Point3d.new(temp[1])
				adjpoints<<Geom::Point3d.new(pointstable[bpntindx])
			end
		end
		adjpoints.reject! { |c| c.nil? }
		# points
		adjpoints
	end
	def intersect_p_p_ps(p1,p2,polygon)
		point=[]
		p1p2line=[p1,p1.vector_to(p2)]
		for i in 0..polygon.length-1
			if i==polygon.length-1
				j=0
			else
				j=i+1
			end
			ijline=[polygon[i],polygon[i].vector_to(polygon[j])]
			int=Geom.intersect_line_line(p1p2line, ijline)
			if (int)&&(point_is_between(p1,p2,int)&&point_is_between(polygon[i],polygon[j],int))
				point<<int
				break
			end
		end
		point[0]
	end
	def intersections_p_p_ps(p1,p2,polygon)
		point=[]
		p1p2line=[p1,p1.vector_to(p2)]
		for i in 0..polygon.length-1
			if i==polygon.length-1
				j=0
			else
				j=i+1
			end
			ijline=[polygon[i],polygon[i].vector_to(polygon[j])]
			int=Geom.intersect_line_line(p1p2line, ijline)
			if point_is_between(p1,p2,int)&&point_is_between(polygon[i],polygon[j],int)
				point<<int
			end
		end
		point
	end
	def intersections_ps_ps(points,polygon)
		point=[]
		for i in 0..polygon.length-1
			if i==polygon.length-1
				j=0
			else
				j=i+1
			end
			ijline=[polygon[i],polygon[i].vector_to(polygon[j])]
			for k in 0..points.length-1
				if k==points.length-1
					l=0
				else
					l=k+1
				end
				klline=[points[k],points[k].vector_to(points[l])]
				int=Geom.intersect_line_line(klline, ijline)
				if (int!=nil)&&(point_is_between(points[k],points[l],int))&&(point_is_between(polygon[i],polygon[j],int))
					point<<int
				end
			end
		end
		point.reject! { |c| c.nil? }
		point
	end
	def point_is_between(p1,p2,pt)
		if (((pt[0]<=p1[0])&&(pt[0]>=p2[0]))||((pt[0]>=p1[0])&&(pt[0]<=p2[0])))&&(((pt[1]<=p1[1])&&(pt[1]>=p2[1]))||((pt[1]>=p1[1])&&(pt[1]<=p2[1])))&&(((pt[2]<=p1[2])&&(pt[2]>=p2[2]))||((pt[2]>=p1[2])&&(pt[2]<=p2[2])))
			test=true
		else
			test=false
		end
		test
	end	
end#class