=begin
#-------------------------------------------------------------------------------------------------------------------------------------------------
#*************************************************************************************************
# Copyright 2009 - Designed September 2009 by Fredo6

# Permission to use this software for any purpose and without fee is hereby granted
# Distribution of this software for commercial purpose is subject to:
#  - the expressed, written consent of the author
#  - the inclusion of the present copyright notice in all copies.

# THIS SOFTWARE IS PROVIDED "AS IS" AND WITHOUT ANY EXPRESS OR
# IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED
# WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
#-----------------------------------------------------------------------------
# Name			:   CurviloftSpline.rb
# Original Date	:   4 Jan 2010 - version 1.0
# Description	:   Specific method for Spline junctions in the Curviloft script family
#-------------------------------------------------------------------------------------------------------------------------------------------------
#*************************************************************************************************
=end

module F6_Curviloft

class CVL_LoftAlgo 

#----------------------------------------------------------------------------
# SPLINE: analysis of Contours
#----------------------------------------------------------------------------

#Verify contours for Spline
def spline_check_contours(lst_contours, manual_selection=false)
	#checking if enough contours
	return false if lst_contours.length < 2
	
	#Ordering the plates
	if manual_selection
		@lorder = (0..lst_contours.length-1).to_a
	else
		@lorder = spline_sort_contours(lst_contours)
	end
	
	lst_plates = []
	@lorder.each { |i| lst_plates.push plate_create(lst_contours[i]) }	

	#Creating the cloture plate
	if lst_plates.length > 2
		@cloture_plate = lst_plates.first
		platelast = lst_plates.last
		d1 = platelast.pts[0].distance @cloture_plate.pts[0]
		d2 = platelast.pts[0].distance @cloture_plate.pts[-1]
	end	
	
	#Creating Plates and Links
	n = lst_plates.length
	@single_link = (n <= 2)
	for i in 0..n-2
		j = (i+1).modulo(n)
		link_create lst_plates[i], lst_plates[j]
	end	

	#Compute the directions for bezier chunks
	lst_plates.each do |plate|
		plate_vecdir plate
	end

	#Creating the Master link and cloture links
	@slave_links = @slave_links_simple = @lst_links.clone
	if @lst_links.length > 1
		@cloture_link = link_create lst_plates.last, @cloture_plate
		@master_link_simple = link_create lst_plates[0], lst_plates[-1]
		@master_link_loop = link_create lst_plates[0], @cloture_plate
		@slave_links_loop = @slave_links_simple + [@cloture_link]
	end

	#Computing the borders for highlight
	@lst_links.each do |link|
		link.hiborders = link.plates.collect { |plate| plate.pts + ((plate.loop) ? [plate.pts.first] : []) }
		link.contours = link.hiborders
	end
		
	#Preparing the cloture link
	if @cloture_link
		link_prepare @cloture_link
		@lst_links.pop
	end
	
	#Switching to Slave mode
	if @master_link_simple
		link_prepare @master_link_simple
		link_prepare @master_link_loop
		@lst_links.pop
		@lst_links.pop
		@state_spline_master = false
	end

	#Final step
	@lst_plates = lst_plates
	
	true
end

#Sort the contours along a topological order
def spline_sort_contours(lpts)
	lplates = []
	lpts.each { |pts| lplates.push plate_create(pts) }
	
	#Ordering and orienting plates
	hsh = {}
	lsp_R = spline_chain_plate lplates, hsh, false
	lsp_L = spline_chain_plate lplates, hsh, true
	lsp_R = lsp_R - lsp_L
	
	lsp = lsp_L.reverse
	lsp += lsp_R if lsp_R.length > 0
	
	ldif = lplates - lsp
	ldif.each { |plate| spline_insert_alone_plates(lsp, plate) if lsp.length > 1 }

	@lorder = []
	lsp.each do |plate|
		@lorder.push lplates.rindex(plate)
	end
	@lorder 
end

#Return the final ordering of plates
def spline_get_order
	@lorder
end

#Compute the chain of plates in a given direction
def spline_chain_plate(lplates, hsh_plate, iside)
	plateprev = lplates[0]
	hsh_plate[plateprev.object_id] = 0
	normprev = plateprev.vecnormal
	normprev = normprev.reverse if normprev && iside
	lsp_ordered = [plateprev]
	
	while true
		baryprev = plateprev.bary
		lsp = []
		lplates.each do |plate|
			next if plate == plateprev
			hhd = hsh_plate[plate.object_id]
			next if hhd && hhd < plateprev.bary.distance(plate.bary)
			vec = baryprev.vector_to plate.bary
			if plate.vecnormal && normprev && plate.vecnormal.valid?
				norm = plate.vecnormal
				norm = norm.reverse if norm % vec < 0
				if vec % normprev >= 0
					lsp.push [plate, norm, plate.bary.distance(baryprev)]
				end	
			else
				lsp.push [plate, vec, plate.bary.distance(baryprev)]
			end
		end
		
		#End of sequence
		break if lsp.length == 0
		
		#Sorting and selecting the next plate
		lsp.sort! { |a, b| a[2] <=> b[2] }
		newplate = lsp[0][0]
		hsh_plate[newplate.object_id] = newplate.bary.distance(baryprev)
		plateprev = newplate
		normprev = lsp[0][1]
		lsp_ordered.push plateprev
	end

	lsp_ordered
end

#Insert the missing plates if any
def spline_insert_alone_plates(lplates, plate)
	ld = []
	bary = plate.bary
	for i in 1..lplates.length - 1
		ld.push [i, bary.distance_to_line([lplates[i-1].bary, lplates[i].bary])]
	end
	ld.sort! { |a, b| a[1] <=> b[1] }
	lplates[ld[0][0], 0] = plate
end

#Compute the right link configuration depending on mode Loops and spline method
def spline_decide_on_link
	loop = hprop_get :option_global_loop
	method = hprop_get :spline_method
	global_spline = (method == :bspline || method == :fspline)

	if loop
		@slave_links = @slave_links_loop
		@master_link = @master_link_loop
		@state_global_loop = true
	else
		@slave_links = @slave_links_simple
		@master_link = @master_link_simple
		if @state_global_loop && @active_link == @slave_links_loop[-1]
			@active_link = nil
			@edit_link = false
		end
		@state_global_loop = false
	end
	
	if global_spline
		@lst_links = [@master_link]
	else
		@lst_links = @slave_links
	end	
end

#----------------------------------------------------------------------------------
# Top level methods for constructing the junctions
#----------------------------------------------------------------------------------

#Construct the Bezier curves for a link
def spline_junction_construct(link)
	#Master or slave mode
	
	#Computing the matrix
	unless @state_recalculate
		spline_matrix_analyze 
		@state_recalculate = true
	end	

	#Construction of the junction curves
	method = hprop_get(:spline_method)
	global_spline = (method == :bspline || method == :fspline)
	if global_spline && @slave_links.length == 1
		method = :catmull
		global_spline = false
	end	
	@type_numseg = nil
	
	if method == :droite
		lbz_pts = spline_droite_construct(link)
	elsif method == :ortho_elliptic
		lbz_pts = spline_profiling_construct(link, 'C')
		@type_numseg = :ortho
	elsif method == :bezier1
		lbz_pts = spline_simple_construct(link, 1)
	elsif method == :bezier2
		lbz_pts = spline_simple_construct(link, 2)
	elsif method == :ortho_spline
		lbz_pts = spline_profiling_construct(link, 'BZ')
		@type_numseg = :ortho
	elsif method == :orthospline2
		lbz_pts = spline_orthospline_construct(link, 2)
	elsif method == :catmull || method == :cubic
		lbz_pts = spline_multi_construct(link)
	elsif method == :bspline || method == :fspline
		@type_numseg = :master
		lbz_pts = spline_global_construct(link, method)
	end
	lbz_pts
end

def spline_which_numseg?
	@type_numseg
end

#Verify if the option requires a change between Master and Slaves
def spline_verify_master?
	return false if @method != :splineloft || @master_link_simple == nil
	spline_decide_on_link
	method = hprop_get :spline_method
	global_spline = (method == :bspline || method == :fspline)
	if global_spline && @slave_links.last.pairs == nil
		@lst_links = @slave_links
		@slave_links.each { |link| link_calculate link }
		spline_matrix_analyze
		@state_recalculate = true
		spline_decide_on_link
	end	
	
	return false unless (!global_spline && @state_spline_master) || (global_spline && !@state_spline_master)
	
	@state_spline_master = !@state_spline_master
	if @state_spline_master
		@active_link = nil
		@edit_link = false
	end
	
	true
end


#----------------------------------------------------------------------------------
# Construction of Spline junctions
#----------------------------------------------------------------------------------

#Construct individual junctions by Bezier curves
def spline_multi_construct(link)
	if @slave_links.length == 1
		return spline_double_construct(link)
	end	
	
	#Reconstructing the junction curves if needed
	unless @state_reconstruct
		spline_multi_junctions
		@state_reconstruct = true
	end	
	
	#Constructing the splines
	ilink = @lst_links.rindex(link)
	lbz = @spline_link_curves[ilink]
	link.pairs = @spline_link_pairs[ilink]
	@spline_link_curves[ilink]
end

#----------------------------------------------------------------------------------
# Construction of non-jointed Bezier junctions
#----------------------------------------------------------------------------------

#Construction of the orthospline curves (Bezier with 3 points)
def spline_droite_construct(link)
	num_bz = hprop_get :num_bz, link
	lbz_pts = []

	link.pairs.each_with_index do |pair, i|
		pt1 = pair[0].pt
		pt2 = pair[1].pt
		bz = G6.divide_segment pt1, pt2, num_bz
		lbz_pts.push bz
	end
		
	lbz_pts
end

#Construction of the elliptic curves (Based on 1/4 circle)
def spline_profiling_construct(link, type='C')
	plate1 = link.plates[0]
	plate2 = link.plates[1]

	vec1 = plate1.vecdir
	vec2 = plate2.vecdir

	#Planes are parallel
	normal1 = plate1.normal
	normal2 = plate2.normal
	para = normal1.parallel?(normal2)
		
	#Computing the Elliptic curves
	num_bz = hprop_get :num_bz_ortho, link
	lbz_pts = []
	prf_type = [type, num_bz]
	link.pairs.each_with_index do |pair, i|
		pt1 = pair[0].pt
		pt2 = pair[1].pt
		if para && pt1.on_line?([pt2, normal2])
			bz = G6.divide_segment pt1, pt2, num_bz
		else
			bz = Traductor::Profiling.compute_by_normals prf_type, pt1, normal1, pt2, normal2
		end	
		lbz_pts.push bz
	end

	lbz_pts
end

#Construction of the orthospline curves (Bezier with 3 points)
def spline_orthospline_construct(link, style)
	plate1 = link.plates[0]
	plate2 = link.plates[1]
	
	vec1 = plate1.vecdir
	vec2 = plate2.vecdir
	
	num_bz = hprop_get :num_bz, link
	lbz_pts = []
	
	#Correction of first and ante-last links when Loop
	if hprop_get(:option_global_loop)
		if link == @lst_links[0]
			vec = @lst_links[-2].plates[1].bary.vector_to link.plates[0].bary
			vec1 = Geom.linear_combination 0.5, vec, 0.5, vec1
		elsif link == @lst_links[-2]
			vec = link.plates[1].bary.vector_to  @lst_links[0].plates[0].bary
			vec2 = Geom.linear_combination 0.5, vec, 0.5, vec2
		elsif link == @lst_links[-1]
			vec2m = @lst_links[-2].plates[1].vecdir
			vec = link.plates[0].bary.vector_to link.plates[1].bary
			vec1 = Geom.linear_combination 0.5, vec, 0.5, vec2m
			vec0 = @lst_links[0].plates[0].vecdir
			vec2 = Geom.linear_combination 0.5, vec, 0.5, vec0	
		end	
	end
	
	#Computing the Bezier curves - Style 1
	if style == 1
		normal1 = plate1.normal
		normal2 = plate2.normal
		link.pairs.each_with_index do |pair, i|
			pt1 = pair[0].pt
			pt2 = pair[1].pt
			ptmid = spline_middle_bezier pt1, normal1, pt2, normal2
			pts = []
			pts = [pt1, ptmid, pt2]
			bz = G6::BezierCurve.compute pts, num_bz
			lbz_pts.push bz
		end
		
	else
		normal1 = plate1.normal
		normal2 = plate2.normal
		link.pairs.each_with_index do |pair, i|
			pt1 = pair[0].pt
			pt2 = pair[1].pt
			ptmid1, ptmid2 = spline_middle_bezier_points pt1, normal1, pt2, normal2
			unless ptmid2
				bz = G6::BezierCurve.compute [pt1, ptmid1, pt2], num_bz
			else
				bz1 = G6::BezierCurve.compute [pt1, ptmid1, pt2], num_bz
				bz2 = G6::BezierCurve.compute [pt1, ptmid2, pt2], num_bz
				bz = G6.curl_average_curve bz1, bz2
			end	
			lbz_pts.push bz
		end
	end
	
	lbz_pts
end

def spline_middle_bezier_points(pt1, normal1, pt2, normal2)
	return pt1 if pt1 == pt2
	
	#Normals are parallel
	if normal1.parallel?(normal2)
		ptproj1 = Geom.intersect_line_plane [pt1, normal1], [pt2, normal2]
		ptmid1 = Geom.linear_combination 0.5, pt1, 0.5, ptproj1
		ptproj2 = Geom.intersect_line_plane [pt2, normal2], [pt1, normal1]
		ptmid2 = Geom.linear_combination 0.5, pt2, 0.5, ptproj2
		ptmid = Geom.linear_combination(0.5, ptmid1, 0.5, ptmid2)
		return [ptmid1, ptmid2]
	end	
	
	#Closest intersection points between Normal lines 
	lpt = Geom.closest_points [pt1, normal1], [pt2, normal2]
end

def spline_middle_bezier(pt1, normal1, pt2, normal2)
	return pt1 if pt1 == pt2
	
	#Normals are parallel
	if normal1.parallel?(normal2)
		ptproj1 = Geom.intersect_line_plane [pt1, normal1], [pt2, normal2]
		ptmid1 = Geom.linear_combination 0.5, pt1, 0.5, ptproj1
		ptproj2 = Geom.intersect_line_plane [pt2, normal2], [pt1, normal1]
		ptmid2 = Geom.linear_combination 0.5, pt2, 0.5, ptproj2
		return Geom.linear_combination(0.5, ptmid1, 0.5, ptmid2)
	end	
	
	#Closest intersection points between Normal lines 
	lpt = Geom.closest_points [pt1, normal1], [pt2, normal2]
	Geom.linear_combination 0.5, lpt[0], 0.5, lpt[1]
end

#Construction by a set of independent Bezier curves
def spline_simple_construct(link, style)
	plate1 = link.plates[0]
	plate2 = link.plates[1]
	
	vec1 = plate1.vecdir
	vec2 = plate2.vecdir
	factor1 = hprop_get(:tension1, link)
	factor2 = hprop_get(:tension2, link)
	
	num_bz = hprop_get :num_bz, link
	lbz_pts = []
	
	#Correction of first and ante-last links when Loop
	if hprop_get(:option_global_loop)
		if link == @lst_links[0]
			vec = @lst_links[-2].plates[1].bary.vector_to link.plates[0].bary
			vec1 = Geom.linear_combination 0.5, vec, 0.5, vec1
		elsif link == @lst_links[-2]
			vec = link.plates[1].bary.vector_to  @lst_links[0].plates[0].bary
			vec2 = Geom.linear_combination 0.5, vec, 0.5, vec2
		elsif link == @lst_links[-1]
			vec2m = @lst_links[-2].plates[1].vecdir
			vec = link.plates[0].bary.vector_to link.plates[1].bary
			vec1 = Geom.linear_combination 0.5, vec, 0.5, vec2m
			vec0 = @lst_links[0].plates[0].vecdir
			vec2 = Geom.linear_combination 0.5, vec, 0.5, vec0	
		end	
	end
	
	#Computing the Bezier curves - Style 1
	if style == 1
		link.pairs.each_with_index do |pair, i|
			pt1 = pair[0].pt
			pt2 = pair[1].pt
			d = pt1.distance(pt2)
			pts = []
			pt11 = pt1.offset vec1, d * factor1
			pt21 = pt2.offset vec2.reverse, d * factor2
			pts = [pt1, pt11, pt21, pt2]
			bz = G6::BezierCurve.compute pts, num_bz
			lbz_pts.push bz
		end
		
	#Computing the Bezier curves - Style 2
	else
		vec1 = plate1.normal
		vec2 = plate2.normal
		link.pairs.each_with_index do |pair, i|
			pt1 = pair[0].pt
			pt2 = pair[1].pt
			vec = pt1.vector_to pt2
			ratio1 = 0.5
			ratio2 = 0.5
			v1 = vec1
			v2 = vec2
			d = pt1.distance(pt2)
			pts = []
			pt11 = pt1.offset v1, d * factor1
			pt21 = pt2.offset v2.reverse, d * factor2
			pts = [pt1, pt11, pt21, pt2]
			bz = G6::BezierCurve.compute pts, num_bz
			lbz_pts.push bz
		end
	
	end
	
	lbz_pts
end

#Construct a spline junction out of 2 points only
def spline_double_construct(link)
	plate1 = link.plates[0]
	plate2 = link.plates[1]
	
	vec1 = plate1.vecdir
	vec2 = plate2.vecdir.reverse
	
	#Computing the Bezier curves
	num_bz = hprop_get :num_bz, link
	lbz_pts = []
	link.pairs.each_with_index do |pair, i|
		pt1 = pair[0].pt
		pt2 = pair[1].pt
		
		#First spline
		line1 = [pt1, vec1]
		pt2b = pt2.project_to_line line1
		d12 = pt1.distance pt2b
		ptmirror = pt2b.offset vec1, d12
		if pt1 == ptmirror
			ptmirror = pt2b.offset vec1, -d12
		end	
		bz = G6::CatmullCurve.compute [pt1, pt2, ptmirror], num_bz+1
		bz1 = bz[0..num_bz]

		#Second spline
		line2 = [pt2, vec2]
		pt1b = pt1.project_to_line line2
		d12 = pt2.distance pt1b
		ptmirror = pt1b.offset vec2, d12
		if pt2 == ptmirror
			ptmirror = pt1b.offset vec2, -d12
		end	
		bz = G6::CatmullCurve.compute [pt2, pt1, ptmirror], num_bz+1
		bz2 = bz[0..num_bz].reverse

		#Weighted average
		bz = G6.curl_average_curve bz1, bz2
		
		lbz_pts.push bz
	end
	
	lbz_pts
end

#----------------------------------------------------------------------------------
# Construction of Global end-to-End junctions
#----------------------------------------------------------------------------------

#Construct global  junctions with a Unfirm BSpline or FSpline
def spline_global_construct(link, method)
	#Computing the control points for the BSpline
	lbz_ctrl_pts = spline_master_nurbs
	
	#Computing the BSpline
	numseg = hprop_get(:num_bz_global)

	if method == :fspline
		lbz_pts = lbz_ctrl_pts.collect { |lpt| G6::FSpline.compute lpt, numseg }
	else
		order = 3
		lbz_pts = lbz_ctrl_pts.collect { |lpt| G6::UniformBSpline.compute lpt, numseg, order }
	end	
end

#----------------------------------------------------------------------------------
# Construction of Master end-to-End junctions
#----------------------------------------------------------------------------------

#Construct the construct points for the Master
def spline_master_nurbs
	#Computing the extremity points of the matrix
	lbz_ctrl_pts = []
	ptsmat = @pts_matrix
	
	#Case where the contours are in loop
	if @master_link.biloop
		ptsmat = ptsmat + ptsmat + ptsmat
	end	
	
	#Processing the pairs
	ibeg = iend = nil
	@master_link.pairs.each do |pair|
		ptbeg = pair[0].pt
		ptend = pair[1].pt
		ibeg = spline_master_nurbs_locate_point(ptsmat, ibeg, ptbeg, 0)
		iend = spline_master_nurbs_locate_point(ptsmat, iend, ptend, -1)
		if ibeg < iend
			ifirst = ibeg
			ilast = iend + 1
		else	
			ifirst = iend
			ilast = ibeg + 1
		end	
		ptarg = spline_master_nurbs_process_pair ptsmat, pair, ptbeg, ptend, ifirst, ilast
		lbz_ctrl_pts.push ptarg
	end
	lbz_ctrl_pts
end

#Computing the sequence of control points
def spline_master_nurbs_process_pair(ptsmat, pair, ptbeg, ptend, ifirst, ilast)
	pts_beg = ptsmat[ifirst]
	pts_end = ptsmat[ilast]
	
	#Distances at extremities
	dbeg = spline_master_nurbs_curvi_distance ptsmat, ifirst, ilast, 0
	deltabeg = spline_master_nurbs_curvi_reach ptsmat, ifirst, ilast, ptbeg, 0
	ratiobeg = deltabeg / dbeg
	
	dend = spline_master_nurbs_curvi_distance ptsmat, ifirst, ilast, -1
	deltaend = spline_master_nurbs_curvi_reach ptsmat, ifirst, ilast, ptend, -1
	ratioend = deltaend / dend
	
	#Distances between Beg and End
	nlink = @slave_links.length	
	tdist = []
	for ipos in 0..nlink
		d = 0
		for i in ifirst+1..ilast
			d += ptsmat[i-1][ipos].distance ptsmat[i][ipos]
		end	
		tdist[ipos] = d
	end
	
	#Cumulative ecarts
	ecart = [0]
	for ipos in 1..nlink
		ecart[ipos] = ecart[ipos-1] + ptsmat[ifirst][ipos-1].distance(ptsmat[ifirst][ipos])
	end
	ecartot = ecart.last
	
	#Intermediate ratios
	newratio = []
	targdist = []
	for ipos in 0..nlink
		r = ecart[ipos] / ecartot
		newratio[ipos] = (1 - r) * ratiobeg + r * ratioend
		targdist[ipos] = newratio[ipos] * tdist[ipos]
	end
	
	#Computing the target points
	ptarg = []
	for ipos in 0..nlink
		ptarg[ipos] = spline_master_nurbs_interpol ptsmat, ifirst, ilast, ipos, targdist[ipos]
		if ptarg[ipos].x.nan?
			puts "NAN "
		end	
	end
	
	ptarg
end

def spline_master_nurbs_interpol(ptsmat, ifirst, ilast, ipos, targdist)
	ptfirst = ptsmat[ifirst][ipos]
	d = 0
	ptarg = nil
	for i in ifirst+1..ilast
		ptlast = ptsmat[i][ipos]
		dd = ptfirst.distance ptlast
		if targdist >= d
			u = targdist - d
			r = u / dd
			ptarg = Geom.linear_combination 1 - r, ptfirst, r, ptlast
			break
		end
		d += dd
		ptfirst = ptlast
	end	
	ptarg
end

#Finding interval in the point matrix
def spline_master_nurbs_curvi_distance(ptsmat, ifirst, ilast, ipos)
	ptfirst = ptsmat[ifirst][ipos]
	d = 0
	for i in ifirst+1..ilast
		ptlast = ptsmat[i][ipos]
		d += ptfirst.distance ptlast
		ptfirst = ptlast
	end	
	d
end
	
def spline_master_nurbs_curvi_reach(ptsmat, ifirst, ilast, pt, ipos)
	ptfirst = ptsmat[ifirst][ipos]
	d = 0
	for i in ifirst+1..ilast
		ptlast = ptsmat[i][ipos]
		if G6.point_within_segment?(pt, ptfirst, ptlast)
			d += pt.distance ptfirst
			break
		end	
		d += ptfirst.distance ptlast
		ptfirst = ptlast
	end	
	d
end

#Finding interval in the point matrix
def spline_master_nurbs_locate_point(ptsmat, ibeg0, pt, iborder)
	nmat = ptsmat.length - 1	
	ibeg0 = 0 unless ibeg0
	ibeg = nil
	for i in (ibeg0+1)..nmat
		ptsprev = ptsmat[i-1]
		ptsnext = ptsmat[i]
		if G6.point_within_segment?(pt, ptsprev[iborder], ptsnext[iborder])
			ibeg = i - 1
			break
		end	
	end	
	ibeg
end

#----------------------------------------------------------------------------------
# Spline Matrix management
#----------------------------------------------------------------------------------

#Build the matrix of all junctions
def spline_matrix_analyze
	nlink = @slave_links.length - 1
	return if nlink == 0
			
	#Buidling the initial matrix
	@spline_matrix = []
	link0 = @slave_links.first
	
	#First link
	pairs = link0.pairs
	@slave_links.first.pairs.each { |pair| @spline_matrix.push [pair] }
	
	#Adding the other links
	for ilink in 1..nlink
		@slave_links[ilink].pairs.each { |pair| spline_matrix_fit_pair ilink, pair }
	end
	
	#Complementing the NEXT missing knots
	for ipos in 0..nlink-1
		spline_matrix_column ipos, 0
	end
	
	#Complementing the PREVIOUS missing knots
	for ipos in 0..nlink-1
		spline_matrix_column nlink - ipos, 1
	end
	
	#Computing the matrix of points
	spline_matrix_of_points
	
	gloop = hprop_get :option_global_loop
	if gloop
		spline_complement_matrix_loop
	end	
end

def spline_loop_match(pts_first, pts_last)		
	n1 = pts_first.length - 1
	
	ipts_supp = []
	ifirst = 0
	ilast = 0
	for ilast in 0.. pts_last.length-1
		ptlast = pts_last[ilast]
		while ifirst < n1
			pt1 = pts_first[ifirst]
			pt2 = pts_first[ifirst+1]
			if G6.point_within_segment?(ptlast, pt1, pt2)
				ratio = pt1.distance(ptlast) / pt1.distance(pt2)
				ipts_supp[ilast] = [ifirst, ratio]
				break
			else
				ifirst += 1
			end	
		end	
	end
	
	ipts_supp
end

#Compute the supplementary points for closing loops
def spline_complement_matrix_loop
	#Forward - End of loop	
	pts_last = @spline_matrix.collect { |ll| ll[-1][1].pt }
	pts_first = @spline_matrix.collect { |ll| ll[0][0].pt }	
	n1 = pts_first.length
	pts_first += pts_first[0..-1]
	ipts_supp = spline_loop_match pts_first, pts_last
		
	pts1 = @spline_matrix.collect { |ll| ll[0][1].pt }	
	pts1 += pts1[0..-1]
	pts2 = @spline_matrix.collect { |ll| ll[1][1].pt }	
	pts2 += pts2[0..-1]
	@matrix_loop_forward = []
	ipts_supp.each do |ll|
		ipos, ratio = ll
		pt1 = Geom.linear_combination 1 - ratio, pts1[ipos], ratio, pts1[ipos+1]
		pt2 = Geom.linear_combination 1 - ratio, pts2[ipos], ratio, pts2[ipos+1]
		@matrix_loop_forward.push [pt1, pt2]
	end	

	#Backward - Beg of loop	
	pts_last = @spline_matrix.collect { |ll| ll[0][0].pt }
	pts_first = @spline_matrix.collect { |ll| ll[-1][1].pt }	
	n1 = pts_first.length
	pts_first += pts_first[0..-1]
	ipts_supp = spline_loop_match pts_first, pts_last
		
	pts1 = @spline_matrix.collect { |ll| ll[-1][0].pt }	
	pts1 += pts1[0..-1]
	pts2 = @spline_matrix.collect { |ll| ll[-2][0].pt }	
	pts2 += pts2[0..-1]
	@matrix_loop_backward = []
	ipts_supp.each do |ll|
		ipos, ratio = ll
		pt1 = Geom.linear_combination 1 - ratio, pts1[ipos], ratio, pts1[ipos+1]
		pt2 = Geom.linear_combination 1 - ratio, pts2[ipos], ratio, pts2[ipos+1]
		@matrix_loop_backward.push [pt2, pt1]
	end	
	
end

#Computing the Catmull curves
def spline_multi_junctions
	global_loop = hprop_get :option_global_loop
	numseg = hprop_get(:num_bz)
	lnumseg = @slave_links.collect { |link| hprop_get(:num_bz, link) }
	
	#Selecting the computation method: Catmull or Cubic Bezier
	method = hprop_get(:spline_method)
	modclass = (method == :catmull) ? G6::CatmullCurve : G6::CubicBezier
	callmeth = modclass.method 'compute'
	
	#Computing the interpolated spline curves
	if global_loop
		catcurves = spline_compute_curves_loop lnumseg, callmeth
	else
		catcurves = spline_compute_curves_simple lnumseg, callmeth
	end
	
	#Computing the visible junctions and the matrix of all junction points
	nlink = @slave_links.length - 1
	ns = @spline_matrix.length - 1
	@spline_link_curves = []
	@spline_link_pairs = []
	ibeg = 0
	for ilink in 0..nlink
		numseg = lnumseg[ilink]
		lbz = []
		lptmat = []
		iend = ibeg + numseg
		for i in 0..ns 
			ll = @spline_matrix[i]
			pair = ll[ilink]
			lptmat.push [pair[0].dist, pair[1].dist, pair[0].pt, pair[1].pt]
			next if pair[0].type == :virtual || pair[1].type == :virtual
			lbz.push [pair[0].dist, pair[1].dist, pair, catcurves[i][ibeg..iend]]
		end	
		lbz.sort! { |a, b| (a[0] == b[0]) ? a[1] <=> b[1] : a[0] <=> b[0] }
		lptmat.sort! { |a, b| (a[0] == b[0]) ? a[1] <=> b[1] : a[0] <=> b[0] }
		ibeg = iend
		@spline_link_pairs[ilink] = lbz.collect { |a| a[2] }
		@spline_link_curves[ilink] = lbz.collect { |a| a[3] }
	end	
end

#Compute the spline curves in Loop mode
def spline_compute_curves_loop(lnumseg, callmeth)
	catcurves = []
	
	#Specifying the segments for the fake additional sections
	nbeg0 = lnumseg[0]
	nbeg1 = lnumseg[1]
	nbeg1m = lnumseg[-1]
	nbeg2m = lnumseg[-2]
	lnumseg = [lnumseg[-2], lnumseg[-1]] + lnumseg + [lnumseg[0], lnumseg[1]]
	
	#Computing the spline curves, with the begin portions
	catcurves = []
	@spline_matrix.each_with_index do |ll, i|
		pts = @matrix_loop_backward[i]
		pts += ll.collect { |pair| pair[0].pt }
		pts += [ll[-1][1].pt]
		pts_end = @matrix_loop_forward[i]
		pts += pts_end
		bz = callmeth.call(pts, lnumseg)
		bz2 = bz[nbeg2m+nbeg1m..-(nbeg0+nbeg1+1)]
		catcurves.push bz2
	end	
	
	catcurves
end

#Compute the spline curves in Simple mode (no loop)
def spline_compute_curves_simple(lnumseg, callmeth)
	catcurves = []
	@spline_matrix.each_with_index do |ll, i|
		pts = ll.collect { |pair| pair[0].pt }
		pts.push ll.last[1].pt
		catcurves.push callmeth.call(pts, lnumseg)
	end	
	catcurves
end

#Sorting method for Spline matrix
def spline_matrix_sorting(a, b, n)
	for i in 0..n
		comp = a[i][0].dist <=> b[i][0].dist
		return comp if comp != 0
	end
	a[n][1].dist <=> b[n][1].dist
end

#Constructing the matrix of points for junctions
def spline_matrix_of_points
	nlink = @slave_links.length - 1
	ns = @spline_matrix.length - 1
	
	#Sorting the matrix
	splmat = @spline_matrix.sort { |a, b| spline_matrix_sorting(a, b, nlink) }
	@spline_matrix = splmat #####
	
	#Intermediate matrix by link
	pts_by_link = []
	for ilink in 0..nlink
		lptmat = []
		for i in 0..ns 
			ll = splmat[i]
			pair = ll[ilink]
			lptmat.push [pair[0].pt, pair[1].pt]
		end	
		pts_by_link[ilink] = lptmat
	end	
	
	#Constructing the matrix of points (by transposition)
	ptsmat = []
	for i in 0..ns
		ptsmat[i] = []
	end	
	pts_by_link.each_with_index do |lpt, ilink|
		lpt.each_with_index do |lspt, ipos|
			ptsmat[ipos][ilink] = lspt[0]
			ptsmat[ipos][ilink+1] = lspt[1] if ilink == nlink
		end
	end	
	@pts_matrix = ptsmat
end

def spline_matrix_print(spl=nil)
	spl = @spline_matrix unless spl
	for i in 0..spl.length-1
		ls = spl[i].collect { |pair| (pair) ? 1 : 0 }
		lm = spl[i]
		pt = (lm[0]) ? lm[-1][-1].pt : "nil"
		d = (lm[0]) ? Sketchup.format_length(lm[0][0].dist) : "nil"
		puts "#{i} - #{ls.inspect} - #{d} - #{pt}"
	end
end

#Create the initial matrix by fitting the pairs at poistion ipos
def spline_matrix_fit_pair(ipos, pair)
	pt1 = pair[0].pt
	lsll = @spline_matrix.find { |ll| ll[ipos] == nil && ll[ipos-1] && ll[ipos-1][1].pt == pt1 }
	if lsll
		lsll[ipos] = pair
	else
		lls = []
		lls[ipos] = pair
		@spline_matrix.push lls
	end	
end

#Treatment of the matrix, column by column
def spline_matrix_column(ipos, iborder)
	#Extracting the current column and column where to create the new pair
	newpos = (iborder == 0) ? ipos + 1 : ipos - 1
	iother = (iborder + 1).modulo(2)
	column = []
	for i in 0..@spline_matrix.length-1
		ll = @spline_matrix[i]
		column.push [ll[ipos], ll[newpos], i] if ll[ipos]
	end
	column.sort! do |a, b| 
		ua = a[0][iborder].dist
		ub = b[0][iborder].dist
		(ua == ub) ? a[0][iother].dist <=> b[0][iother].dist : ua <=> ub
	end	
	
	#Prolongating with new pairs if necessary
	for i in 0..column.length-1
		unless column[i][1]	
			pt = column[i][0][iother].pt
			newpair = spline_matrix_match_point @slave_links[newpos], pt, iborder, iother
			@spline_matrix[column[i][2]][newpos] = newpair
		end	
	end
end

#Insert a new pair to match the point on the subsequent link
def spline_matrix_match_point(link, pt, iborder, iother)
	#Locating the point in the current link
	pairs = link.pairs
	pairs = pairs + [pairs.first] if link.biloop
	inext = nil
	ptprev0 = ptnext0 = nil
	for i in 1..pairs.length-1
		ptprev0 = pairs[i-1][iborder].pt
		ptnext0 = pairs[i][iborder].pt
		if G6.point_within_segment?(pt, ptprev0, ptnext0)
			inext = i
			break
		end	
	end
	
	#Relative position of the point on the common border
	pairprev = pairs[inext-1]
	pairnext = pairs[inext]
	d = ptnext0.distance pt
	d0 = ptprev0.distance ptnext0
	ratio = d / d0
	zdist = pairnext[iborder].dist - d
	
	#Poistion of the new point on the other border
	ptprev1 = pairprev[iother].pt
	ptnext1 = pairnext[iother].pt
	newpt = Geom.linear_combination ratio, ptprev1, 1 - ratio, ptnext1
	newzdist = pairprev[iother].dist + ptprev1.distance(newpt)
	
	#Creating the new pair
	if iborder == 0
		newpair = pair_create pt, zdist, :virtual, newpt, newzdist, :virtual
	else	
		newpair = pair_create newpt, newzdist, :virtual, pt, zdist, :virtual
	end			
	newpair
end

end	#End Class CVL_LoftAlgo

end	#End Module F6_Curviloft
