@@ -816,43 +816,47 @@ function BSplineInterpolation(
816
816
u, t = munge_data (u, t)
817
817
n = length (t)
818
818
n < d + 1 && error (" BSplineInterpolation needs at least d + 1, i.e. $(d+ 1 ) points." )
819
-
819
+
820
820
# Initialize parameter vector
821
821
param_vec = zero (t)
822
822
param_vec[1 ] = zero (eltype (t))
823
823
param_vec[end ] = one (eltype (t))
824
-
824
+
825
825
# Initialize knot vector
826
826
knot_vec = zeros (eltype (t), n + d + 1 )
827
827
828
828
# Compute parameter vector based on type
829
829
if pVecType == :Uniform
830
830
for i in 2 : (n - 1 )
831
- param_vec[i] = param_vec[1 ] + (i - 1 ) * (param_vec[end ] - param_vec[1 ]) / (n - 1 )
831
+ param_vec[i] = param_vec[1 ] +
832
+ (i - 1 ) * (param_vec[end ] - param_vec[1 ]) / (n - 1 )
832
833
end
833
834
elseif pVecType == :ArcLen
834
835
# Only compute arc lengths when needed, using diff and broadcasting
835
- distances = sqrt .((diff (t)). ^ 2 .+ (diff (u)). ^ 2 )
836
+ distances = sqrt .((diff (t)) .^ 2 .+ (diff (u)) .^ 2 )
836
837
cumulative_distances = cumsum (distances)
837
838
total_distance = cumulative_distances[end ]
838
839
for i in 2 : (n - 1 )
839
- param_vec[i] = param_vec[1 ] + cumulative_distances[i - 1 ] / total_distance * (param_vec[end ] - param_vec[1 ])
840
+ param_vec[i] = param_vec[1 ] +
841
+ cumulative_distances[i - 1 ] / total_distance *
842
+ (param_vec[end ] - param_vec[1 ])
840
843
end
841
844
end
842
845
843
846
# Set boundary knots using vectorized assignment
844
- knot_vec[1 : d + 1 ] .= param_vec[1 ]
845
- knot_vec[end - d : end ] .= param_vec[end ]
847
+ knot_vec[1 : (d + 1 ) ] .= param_vec[1 ]
848
+ knot_vec[( end - d) : end ] .= param_vec[end ]
846
849
847
850
# Compute cumulative parameter sum for internal knots using cumsum
848
- param_cumsum = cumsum (param_vec[2 : end - 1 ])
851
+ param_cumsum = cumsum (param_vec[2 : ( end - 1 ) ])
849
852
850
853
if knotVecType == :Uniform
851
854
# uniformly spaced knot vector
852
855
# this method is not recommended because, if it is used with the chord length method for global interpolation,
853
856
# the system of linear equations would be singular.
854
857
for i in (d + 2 ): n
855
- knot_vec[i] = knot_vec[1 ] + (i - d - 1 ) // (n - d) * (knot_vec[end ] - knot_vec[1 ])
858
+ knot_vec[i] = knot_vec[1 ] +
859
+ (i - d - 1 ) // (n - d) * (knot_vec[end ] - knot_vec[1 ])
856
860
end
857
861
elseif knotVecType == :Average
858
862
# average spaced knot vector
@@ -865,7 +869,7 @@ function BSplineInterpolation(
865
869
index += 1
866
870
end
867
871
end
868
-
872
+
869
873
# control points
870
874
spline_coeffs = zeros (eltype (t), n, n)
871
875
spline_coefficients! (spline_coeffs, d, knot_vec, param_vec)
@@ -888,12 +892,12 @@ function BSplineInterpolation(
888
892
u, t = munge_data (u, t)
889
893
n = length (t)
890
894
n < d + 1 && error (" BSplineInterpolation needs at least d + 1, i.e. $(d+ 1 ) points." )
891
-
895
+
892
896
# Initialize parameter vector
893
897
param_vec = zero (t)
894
898
param_vec[1 ] = zero (eltype (t))
895
899
param_vec[end ] = one (eltype (t))
896
-
900
+
897
901
# Initialize knot vector
898
902
knot_vec = zeros (eltype (t), n + d + 1 )
899
903
@@ -902,33 +906,38 @@ function BSplineInterpolation(
902
906
903
907
if pVecType == :Uniform
904
908
for i in 2 : (n - 1 )
905
- param_vec[i] = param_vec[1 ] + (i - 1 ) * (param_vec[end ] - param_vec[1 ]) / (n - 1 )
909
+ param_vec[i] = param_vec[1 ] +
910
+ (i - 1 ) * (param_vec[end ] - param_vec[1 ]) / (n - 1 )
906
911
end
907
912
elseif pVecType == :ArcLen
908
913
# Only compute arc lengths when needed, using diff and broadcasting
909
914
time_diffs = diff (t)
910
- spatial_diffs = [sqrt (sum ((u[array_axes... , i+ 1 ] - u[array_axes... , i]). ^ 2 )) for i in 1 : (n- 1 )]
911
- distances = sqrt .(time_diffs.^ 2 .+ spatial_diffs.^ 2 )
915
+ spatial_diffs = [sqrt (sum ((u[array_axes... , i + 1 ] - u[array_axes... , i]) .^ 2 ))
916
+ for i in 1 : (n - 1 )]
917
+ distances = sqrt .(time_diffs .^ 2 .+ spatial_diffs .^ 2 )
912
918
cumulative_distances = cumsum (distances)
913
919
total_distance = cumulative_distances[end ]
914
920
for i in 2 : (n - 1 )
915
- param_vec[i] = param_vec[1 ] + cumulative_distances[i - 1 ] / total_distance * (param_vec[end ] - param_vec[1 ])
921
+ param_vec[i] = param_vec[1 ] +
922
+ cumulative_distances[i - 1 ] / total_distance *
923
+ (param_vec[end ] - param_vec[1 ])
916
924
end
917
925
end
918
926
919
927
# Set boundary knots using vectorized assignment
920
- knot_vec[1 : d + 1 ] .= param_vec[1 ]
921
- knot_vec[end - d : end ] .= param_vec[end ]
928
+ knot_vec[1 : (d + 1 ) ] .= param_vec[1 ]
929
+ knot_vec[( end - d) : end ] .= param_vec[end ]
922
930
923
931
# Compute cumulative parameter sum for internal knots using cumsum
924
- param_cumsum = cumsum (param_vec[2 : end - 1 ])
932
+ param_cumsum = cumsum (param_vec[2 : ( end - 1 ) ])
925
933
926
934
if knotVecType == :Uniform
927
935
# uniformly spaced knot vector
928
936
# this method is not recommended because, if it is used with the chord length method for global interpolation,
929
937
# the system of linear equations would be singular.
930
938
for i in (d + 2 ): n
931
- knot_vec[i] = knot_vec[1 ] + (i - d - 1 ) // (n - d) * (knot_vec[end ] - knot_vec[1 ])
939
+ knot_vec[i] = knot_vec[1 ] +
940
+ (i - d - 1 ) // (n - d) * (knot_vec[end ] - knot_vec[1 ])
932
941
end
933
942
elseif knotVecType == :Average
934
943
# average spaced knot vector
@@ -941,7 +950,7 @@ function BSplineInterpolation(
941
950
index += 1
942
951
end
943
952
end
944
-
953
+
945
954
# control points
946
955
spline_coeffs = zeros (eltype (t), n, n)
947
956
spline_coefficients! (spline_coeffs, d, knot_vec, param_vec)
@@ -1045,43 +1054,47 @@ function BSplineApprox(
1045
1054
u, t = munge_data (u, t)
1046
1055
n = length (t)
1047
1056
h < d + 1 && error (" BSplineApprox needs at least d + 1, i.e. $(d+ 1 ) control points." )
1048
-
1057
+
1049
1058
# Initialize parameter vector
1050
1059
param_vec = zero (t)
1051
1060
param_vec[1 ] = zero (eltype (t))
1052
1061
param_vec[end ] = one (eltype (t))
1053
-
1062
+
1054
1063
# Initialize knot vector
1055
1064
knot_vec = zeros (eltype (t), h + d + 1 )
1056
1065
1057
1066
# Compute parameter vector based on type
1058
1067
if pVecType == :Uniform
1059
1068
for i in 2 : (n - 1 )
1060
- param_vec[i] = param_vec[1 ] + (i - 1 ) * (param_vec[end ] - param_vec[1 ]) / (n - 1 )
1069
+ param_vec[i] = param_vec[1 ] +
1070
+ (i - 1 ) * (param_vec[end ] - param_vec[1 ]) / (n - 1 )
1061
1071
end
1062
1072
elseif pVecType == :ArcLen
1063
1073
# Only compute arc lengths when needed, using diff and broadcasting
1064
- distances = sqrt .((diff (t)). ^ 2 .+ (diff (u)). ^ 2 )
1074
+ distances = sqrt .((diff (t)) .^ 2 .+ (diff (u)) .^ 2 )
1065
1075
cumulative_distances = cumsum (distances)
1066
1076
total_distance = cumulative_distances[end ]
1067
1077
for i in 2 : (n - 1 )
1068
- param_vec[i] = param_vec[1 ] + cumulative_distances[i - 1 ] / total_distance * (param_vec[end ] - param_vec[1 ])
1078
+ param_vec[i] = param_vec[1 ] +
1079
+ cumulative_distances[i - 1 ] / total_distance *
1080
+ (param_vec[end ] - param_vec[1 ])
1069
1081
end
1070
1082
end
1071
1083
1072
1084
# Set boundary knots using vectorized assignment
1073
- knot_vec[1 : d + 1 ] .= param_vec[1 ]
1074
- knot_vec[end - d : end ] .= param_vec[end ]
1085
+ knot_vec[1 : (d + 1 ) ] .= param_vec[1 ]
1086
+ knot_vec[( end - d) : end ] .= param_vec[end ]
1075
1087
1076
1088
# Compute cumulative parameter sum for internal knots using cumsum
1077
- param_cumsum = cumsum (param_vec[2 : end - 1 ])
1089
+ param_cumsum = cumsum (param_vec[2 : ( end - 1 ) ])
1078
1090
1079
1091
if knotVecType == :Uniform
1080
1092
# uniformly spaced knot vector
1081
1093
# this method is not recommended because, if it is used with the chord length method for global interpolation,
1082
1094
# the system of linear equations would be singular.
1083
1095
for i in (d + 2 ): h
1084
- knot_vec[i] = knot_vec[1 ] + (i - d - 1 ) // (h - d) * (knot_vec[end ] - knot_vec[1 ])
1096
+ knot_vec[i] = knot_vec[1 ] +
1097
+ (i - d - 1 ) // (h - d) * (knot_vec[end ] - knot_vec[1 ])
1085
1098
end
1086
1099
elseif knotVecType == :Average
1087
1100
# average spaced knot vector using improved distribution
@@ -1096,11 +1109,11 @@ function BSplineApprox(
1096
1109
knot_vec[d + 2 ] = param_interp (0.5 )
1097
1110
else
1098
1111
internal_positions = range (0 , 1 , length = num_internal_knots)
1099
- knot_vec[(d+ 2 ): h] .= param_interp .(internal_positions)
1112
+ knot_vec[(d + 2 ): h] .= param_interp .(internal_positions)
1100
1113
end
1101
1114
end
1102
1115
end
1103
-
1116
+
1104
1117
# control points
1105
1118
control_points = zeros (eltype (u), h)
1106
1119
control_points[1 ] = u[1 ]
@@ -1143,12 +1156,12 @@ function BSplineApprox(
1143
1156
u, t = munge_data (u, t)
1144
1157
n = length (t)
1145
1158
h < d + 1 && error (" BSplineApprox needs at least d + 1, i.e. $(d+ 1 ) control points." )
1146
-
1159
+
1147
1160
# Initialize parameter vector
1148
1161
param_vec = zero (t)
1149
1162
param_vec[1 ] = zero (eltype (t))
1150
1163
param_vec[end ] = one (eltype (t))
1151
-
1164
+
1152
1165
# Initialize knot vector
1153
1166
knot_vec = zeros (eltype (t), h + d + 1 )
1154
1167
@@ -1157,33 +1170,38 @@ function BSplineApprox(
1157
1170
1158
1171
if pVecType == :Uniform
1159
1172
for i in 2 : (n - 1 )
1160
- param_vec[i] = param_vec[1 ] + (i - 1 ) * (param_vec[end ] - param_vec[1 ]) / (n - 1 )
1173
+ param_vec[i] = param_vec[1 ] +
1174
+ (i - 1 ) * (param_vec[end ] - param_vec[1 ]) / (n - 1 )
1161
1175
end
1162
1176
elseif pVecType == :ArcLen
1163
1177
# Only compute arc lengths when needed, using diff and broadcasting
1164
1178
time_diffs = diff (t)
1165
- spatial_diffs = [sqrt (sum ((u[array_axes... , i+ 1 ] - u[array_axes... , i]). ^ 2 )) for i in 1 : (n- 1 )]
1166
- distances = sqrt .(time_diffs.^ 2 .+ spatial_diffs.^ 2 )
1179
+ spatial_diffs = [sqrt (sum ((u[array_axes... , i + 1 ] - u[array_axes... , i]) .^ 2 ))
1180
+ for i in 1 : (n - 1 )]
1181
+ distances = sqrt .(time_diffs .^ 2 .+ spatial_diffs .^ 2 )
1167
1182
cumulative_distances = cumsum (distances)
1168
1183
total_distance = cumulative_distances[end ]
1169
1184
for i in 2 : (n - 1 )
1170
- param_vec[i] = param_vec[1 ] + cumulative_distances[i - 1 ] / total_distance * (param_vec[end ] - param_vec[1 ])
1185
+ param_vec[i] = param_vec[1 ] +
1186
+ cumulative_distances[i - 1 ] / total_distance *
1187
+ (param_vec[end ] - param_vec[1 ])
1171
1188
end
1172
1189
end
1173
1190
1174
1191
# Set boundary knots using vectorized assignment
1175
- knot_vec[1 : d + 1 ] .= param_vec[1 ]
1176
- knot_vec[end - d : end ] .= param_vec[end ]
1192
+ knot_vec[1 : (d + 1 ) ] .= param_vec[1 ]
1193
+ knot_vec[( end - d) : end ] .= param_vec[end ]
1177
1194
1178
1195
# Compute cumulative parameter sum for internal knots using cumsum
1179
- param_cumsum = cumsum (param_vec[2 : end - 1 ])
1196
+ param_cumsum = cumsum (param_vec[2 : ( end - 1 ) ])
1180
1197
1181
1198
if knotVecType == :Uniform
1182
1199
# uniformly spaced knot vector
1183
1200
# this method is not recommended because, if it is used with the chord length method for global interpolation,
1184
1201
# the system of linear equations would be singular.
1185
1202
for i in (d + 2 ): h
1186
- knot_vec[i] = knot_vec[1 ] + (i - d - 1 ) // (h - d) * (knot_vec[end ] - knot_vec[1 ])
1203
+ knot_vec[i] = knot_vec[1 ] +
1204
+ (i - d - 1 ) // (h - d) * (knot_vec[end ] - knot_vec[1 ])
1187
1205
end
1188
1206
elseif knotVecType == :Average
1189
1207
# average spaced knot vector using improved distribution
@@ -1198,11 +1216,11 @@ function BSplineApprox(
1198
1216
knot_vec[d + 2 ] = param_interp (0.5 )
1199
1217
else
1200
1218
internal_positions = range (0 , 1 , length = num_internal_knots)
1201
- knot_vec[(d+ 2 ): h] .= param_interp .(internal_positions)
1219
+ knot_vec[(d + 2 ): h] .= param_interp .(internal_positions)
1202
1220
end
1203
1221
end
1204
1222
end
1205
-
1223
+
1206
1224
# control points
1207
1225
control_points = zeros (eltype (u), size (u)[1 : (end - 1 )]. .. , h)
1208
1226
control_points[array_axes... , 1 ] = u[array_axes... , 1 ]
@@ -1221,7 +1239,8 @@ function BSplineApprox(
1221
1239
for i in 2 : (h - 1 )
1222
1240
residual_sum = zeros (eltype (spline_coeffs), size (u)[1 : (end - 1 )]. .. )
1223
1241
for k in 2 : (n - 1 )
1224
- residual_sum = residual_sum + spline_coeffs[k, i] * data_residual[array_axes... , k]
1242
+ residual_sum = residual_sum +
1243
+ spline_coeffs[k, i] * data_residual[array_axes... , k]
1225
1244
end
1226
1245
coeff_matrix[array_axes... , i - 1 ] = residual_sum
1227
1246
end
0 commit comments