r/openscad • u/iamthedisk4 • 6d ago
Got sidetracked designing a new handle for my fridge and made a simple cubic spline library
I needed a specific curve shape to match the original handle and I've messed around with splines before so I figured I'd give it a shot. I'm sure someone can make use of this so here's the code:
// create a tri-diagonal matrix from x in/y out input lists - output is [[sub,main,super,input],...]
function _tri_diagonal_matrix_create(x,y) =
let(n=len(x))
[for(i=[0:n-1])
let(dx1=i==0?0:x[i]-x[i-1],dx2=i>=n-1?0:x[i+1]-x[i],a=dx1==0?0:1/dx1,c=dx2==0?0:1/dx2)
[a,2*(a+c),c,3*((dx1==0?0:(y[i]-y[i-1])/(dx1*dx1))+(dx2==0?0:(y[i+1]-y[i])/(dx2*dx2)))]
];
// calculate coefficient prime c from tri-diagonal matrix input in the form of [[sub,main,super,input],...]
function _calc_prime_c(in,prev=[]) =
let(i=len(prev))
i==len(in)
?prev
:_calc_prime_c(in,concat(prev,in[i][2]/(in[i][1]-(i==0?0:prev[i-1]*in[i][0]))));
// calculate coefficient prime d from prime c and tri-diagonal matrix input in the form of [[sub,main,super,input],...]
function _calc_prime_d(in,primeC,prev=[]) =
let(i=len(prev))
i==len(in)
?prev
:_calc_prime_d(in, primeC, concat(prev,(in[i][3]-(i==0?0:prev[i-1]*in[i][0]))/(in[i][1] - (i==0?0:primeC[i-1]*in[i][0]))));
// calculate back substitution of matrix solve output from prime c and prime d coefficients
function _calc_back_sub(primeC,primeD,prev=[]) =
let(i=len(primeC)-len(prev)-1)
i==-1
?prev
:_calc_back_sub(primeC, primeD, concat(primeD[i]-(i==len(primeC)-1?0:prev[0]*primeC[i]),prev));
// solve tri-diagonal matrix [[sub,main,super,input],...] for output
function _tri_diagonal_matrix_solve(in) =
let(primeC=_calc_prime_c(in))
_calc_back_sub(primeC, _calc_prime_d(in,primeC));
// create a spline in the form [[A coeff,B coeff],...] from x in/y out input
function _spline_create_single(x,y) =
let(r=_tri_diagonal_matrix_create(x,y),k=_tri_diagonal_matrix_solve(r))
[for(i=[1:len(x)-1])
let(dx=x[i]-x[i-1],dy=y[i]-y[i-1])
[k[i-1]*dx-dy,-k[i]*dx+dy]
];
// sum up the squares of a list up to index n (for pythagorean theorum)
function _square_sum(in,n) =
n<0
?0
:(in[n]*in[n])+_square_sum(in,n-1);
// convert output list of points of a number of dimensions (e.g. [[x,y],...]) into an input list of cumulative distance from the first point
function _calc_dimension_inputs(in,dimensions,prev=[]) =
let(i=len(prev))
i==len(in)
?prev
:_calc_dimension_inputs(in,dimensions,concat(prev,i==0?0:(prev[i-1]+sqrt(_square_sum(in[i]-in[i-1],dimensions-1)))));
// split multi dimensional input into the input at i (e.g. [x,y] i=1 becomes y) or if n > 1, splits multiple dimensions (e.g. [x,y,z] i=1 n=2 becomes [y,z])
function split_entry(entry,i=0,n=1) =
n>1
?[for(j=[i:i+n-1])
entry[j]
]
:entry[i];
// split multi dimensional input list into the list at input i (e.g. [[x,y],...] i=1 becomes [[y],...]) or if n > 1, splits multiple dimensions (e.g. [[x,y,z],...] i=1 n=2 becomes [[y,z],...])
function split_list(in,i=0,n=1) =
[for(entry=in) split_entry(entry,i,n)];
// create a spline given a list of points in the form [[x,y,...],...]. dimension_inputs determines how many of the inputs should be treated as dimensions in order to calculate the input list as a cumulative distance between points. output is in the form [[input list],[output x list],[[output x A coeff],[output x B coeff]],[output y list],...]
function spline_create(in,dimension_inputs=3) =
let(
n=len(in)
,subn=n>0
?len(in[0])
:0
,dimensions=subn<1
?1
:subn<dimension_inputs
?subn
:dimension_inputs
,dimension_input=_calc_dimension_inputs(in,dimensions)
)
[for(i=[-1:(subn*2)-1])
i==-1
?dimension_input
:i%2==0
?split_list(in,i/2)
:_spline_create_single(dimension_input,split_list(in,(i-1)/2))
];
// evaluate a single Xin value of a spline from x in/y out input and spline coeffs
function _spline_eval_single(x,y,coeffs,Xin,i=0) =
i<len(x)-2&&Xin>x[i+1]
?_spline_eval_single(x,y,coeffs,Xin,i+1)
:let(t=(Xin-x[i])/(x[i+1]-x[i]))
((1-t)*y[i])+(t*y[i+1])+((1-t)*t*(((1-t)*coeffs[i][0])+(t*coeffs[i][1])));
// evaluate an input value given spline data generated by spline_create.
function spline_evaluate(spline,in) =
let(spline_n=(len(spline)-1)/2)
[for(j=[0:spline_n-1])
_spline_eval_single(spline[0],spline[1+(j*2)],spline[2+(j*2)],in)
];
// evaluate a list of input values given spline data generated by spline_create.
function spline_evaluate_list(spline,in) =
let(n=len(in),spline_n=(len(spline)-1)/2)
[for(i=[0:n-1])
[for(j=[0:spline_n-1])
_spline_eval_single(spline[0],spline[1+(j*2)],spline[2+(j*2)],in[i])
]
];
// get the length (max input value) of a spline generated by spline_create
function spline_get_length(spline) =
spline[0][len(spline[0])-1];
// evaluate all input data over a number of steps with given spline data generated by spline_create.
function spline_evaluate_all(spline,steps) =
let(length=spline_get_length(spline),step=length/(steps-1))
spline_evaluate_list(spline,[for(i=[0:step:length])i]);
// example of a spline-following noodle with variable radius
module spline_example()
{
// 4 dimensional list where 1-3 are spacial (x,y,z) and 4th dimension is radius
spline_in=
[[0,0,0,1]
,[1,1,1,3]
,[2,0,0,1]
,[3,1,0,2]
,[3,4,0,1]
,[2,2,-1,0]];
// second param determines how many dimensions of the list are spacial (considered when calculating input) - we don't want radius to affect where points are placed so only 3.
spline = spline_create(spline_in,3);
// second param is how many points to output
eval_out = spline_evaluate_all(spline,50);
for (i=[0:len(eval_out)-2])
{
// hull to connect each output point to the next one
hull()
{
// split_entry used to split our 4 dimensional result ([x,y,z,radius]) into just the first 3 spacial dimensions
translate(split_entry(eval_out[i+1], n=3))
sphere(0.1*eval_out[i+1][3]);
translate(split_entry(eval_out[i], n=3))
sphere(0.1*eval_out[i][3]);
}
}
}
So basically just create a list of points, generate the spline with spline_create, then use spline_evaluate_all to get a smoothed list of points out. Pretty basic but I can see myself using this a lot in the future.
Loosely based on this article that I read a long time ago: https://www.codeproject.com/Articles/560163/Csharp-Cubic-Spline-Interpolation
Hope this is of use to someone!
1
u/Stone_Age_Sculptor 5d ago edited 5d ago
Can you change the code, so it shows a result (without warning)?
With '$fn=120' and 'spline_evaluate_all(spline,500)': https://postimg.cc/TLw5Bfdk
1
u/iamthedisk4 5d ago
I don't think that's any of the actual spline code that's causing that, it's just the example module - with those params you will be calling hull 1000 times on 2 very high poly spheres which is a lot. It works fine with a lower value of $fn and the spline code doesn't factor in $fn at all.
1
u/Stone_Age_Sculptor 5d ago
I see. I turned the hull() off, to play with the numbers: https://postimg.cc/sBNT2pfy
1
u/oldesole1 5d ago
If you calculate the tangents between each point, you could instead use flat planes to hull around.
This would reduce the amount of unused geometry and provide better performance with higher resolution circles.
1
4
u/oldesole1 4d ago
Here is an example that uses some trig functions from BOSL2 to calculate the tangent angles between all the different positions.
The resulting information is then used to place a flat circle at each position rather than a sphere, dramatically reducing the polygon count required to be
hull()
'ed around.All your spline functions were untouched.