@@ -45,125 +45,6 @@ impl Hash for CurveManipulatorGroup {
4545 }
4646}
4747
48- #[ derive( Debug ) ]
49- pub struct CubicSplines {
50- pub x : [ f32 ; 4 ] ,
51- pub y : [ f32 ; 4 ] ,
52- }
53-
54- impl CubicSplines {
55- pub fn solve ( & self ) -> [ f32 ; 4 ] {
56- let ( x, y) = ( & self . x , & self . y ) ;
57-
58- // Build an augmented matrix to solve the system of equations using Gaussian elimination
59- let mut augmented_matrix = [
60- [
61- 2. / ( x[ 1 ] - x[ 0 ] ) ,
62- 1. / ( x[ 1 ] - x[ 0 ] ) ,
63- 0. ,
64- 0. ,
65- // |
66- 3. * ( y[ 1 ] - y[ 0 ] ) / ( ( x[ 1 ] - x[ 0 ] ) * ( x[ 1 ] - x[ 0 ] ) ) ,
67- ] ,
68- [
69- 1. / ( x[ 1 ] - x[ 0 ] ) ,
70- 2. * ( 1. / ( x[ 1 ] - x[ 0 ] ) + 1. / ( x[ 2 ] - x[ 1 ] ) ) ,
71- 1. / ( x[ 2 ] - x[ 1 ] ) ,
72- 0. ,
73- // |
74- 3. * ( ( y[ 1 ] - y[ 0 ] ) / ( ( x[ 1 ] - x[ 0 ] ) * ( x[ 1 ] - x[ 0 ] ) ) + ( y[ 2 ] - y[ 1 ] ) / ( ( x[ 2 ] - x[ 1 ] ) * ( x[ 2 ] - x[ 1 ] ) ) ) ,
75- ] ,
76- [
77- 0. ,
78- 1. / ( x[ 2 ] - x[ 1 ] ) ,
79- 2. * ( 1. / ( x[ 2 ] - x[ 1 ] ) + 1. / ( x[ 3 ] - x[ 2 ] ) ) ,
80- 1. / ( x[ 3 ] - x[ 2 ] ) ,
81- // |
82- 3. * ( ( y[ 2 ] - y[ 1 ] ) / ( ( x[ 2 ] - x[ 1 ] ) * ( x[ 2 ] - x[ 1 ] ) ) + ( y[ 3 ] - y[ 2 ] ) / ( ( x[ 3 ] - x[ 2 ] ) * ( x[ 3 ] - x[ 2 ] ) ) ) ,
83- ] ,
84- [
85- 0. ,
86- 0. ,
87- 1. / ( x[ 3 ] - x[ 2 ] ) ,
88- 2. / ( x[ 3 ] - x[ 2 ] ) ,
89- // |
90- 3. * ( y[ 3 ] - y[ 2 ] ) / ( ( x[ 3 ] - x[ 2 ] ) * ( x[ 3 ] - x[ 2 ] ) ) ,
91- ] ,
92- ] ;
93-
94- // Gaussian elimination: forward elimination
95- for row in 0 ..4 {
96- let pivot_row_index = ( row..4 )
97- . max_by ( |& a_row, & b_row| augmented_matrix[ a_row] [ row] . abs ( ) . partial_cmp ( & augmented_matrix[ b_row] [ row] . abs ( ) ) . unwrap_or ( std:: cmp:: Ordering :: Equal ) )
98- . unwrap ( ) ;
99-
100- // Swap the current row with the row that has the largest pivot element
101- augmented_matrix. swap ( row, pivot_row_index) ;
102-
103- // Eliminate the current column in all rows below the current one
104- for row_below_current in row + 1 ..4 {
105- assert ! ( augmented_matrix[ row] [ row] . abs( ) > f32 :: EPSILON ) ;
106-
107- let scale_factor = augmented_matrix[ row_below_current] [ row] / augmented_matrix[ row] [ row] ;
108- for col in row..5 {
109- augmented_matrix[ row_below_current] [ col] -= augmented_matrix[ row] [ col] * scale_factor
110- }
111- }
112- }
113-
114- // Gaussian elimination: back substitution
115- let mut solutions = [ 0. ; 4 ] ;
116- for col in ( 0 ..4 ) . rev ( ) {
117- assert ! ( augmented_matrix[ col] [ col] . abs( ) > f32 :: EPSILON ) ;
118-
119- solutions[ col] = augmented_matrix[ col] [ 4 ] / augmented_matrix[ col] [ col] ;
120-
121- for row in ( 0 ..col) . rev ( ) {
122- augmented_matrix[ row] [ 4 ] -= augmented_matrix[ row] [ col] * solutions[ col] ;
123- augmented_matrix[ row] [ col] = 0. ;
124- }
125- }
126-
127- solutions
128- }
129-
130- pub fn interpolate ( & self , input : f32 , solutions : & [ f32 ] ) -> f32 {
131- if input <= self . x [ 0 ] {
132- return self . y [ 0 ] ;
133- }
134- if input >= self . x [ self . x . len ( ) - 1 ] {
135- return self . y [ self . x . len ( ) - 1 ] ;
136- }
137-
138- // Find the segment that the input falls between
139- let mut segment = 1 ;
140- while self . x [ segment] < input {
141- segment += 1 ;
142- }
143- let segment_start = segment - 1 ;
144- let segment_end = segment;
145-
146- // Calculate the output value using quadratic interpolation
147- let input_value = self . x [ segment_start] ;
148- let input_value_prev = self . x [ segment_end] ;
149- let output_value = self . y [ segment_start] ;
150- let output_value_prev = self . y [ segment_end] ;
151- let solutions_value = solutions[ segment_start] ;
152- let solutions_value_prev = solutions[ segment_end] ;
153-
154- let output_delta = solutions_value_prev * ( input_value - input_value_prev) - ( output_value - output_value_prev) ;
155- let solution_delta = ( output_value - output_value_prev) - solutions_value * ( input_value - input_value_prev) ;
156-
157- let input_ratio = ( input - input_value_prev) / ( input_value - input_value_prev) ;
158- let prev_output_ratio = ( 1. - input_ratio) * output_value_prev;
159- let output_ratio = input_ratio * output_value;
160- let quadratic_ratio = input_ratio * ( 1. - input_ratio) * ( output_delta * ( 1. - input_ratio) + solution_delta * input_ratio) ;
161-
162- let result = prev_output_ratio + output_ratio + quadratic_ratio;
163- result. clamp ( 0. , 1. )
164- }
165- }
166-
16748pub struct ValueMapperNode < C > {
16849 lut : Vec < C > ,
16950}
0 commit comments