Adding bspline function in polygons.jl (evaluates npoints B-spline curve based on control points) #305

g-adel · 2024-03-28T21:48:23Z

B-Spline function with clamped/unclamped options using De Boor's Algorithm. It's implemented in a pretty efficient way. However, I'm sure there's plenty of ways to speed it up using vectorization/parallelization.
I think this was the objective of the polyfit() experimental function? I didn't want to replace it so I leave that decision to the maintainer.

But yeah I've been working with B-Splines lately and thought it would be useful to have for this library! I tried to make it similar to the coding style and inputs/outputs of the library functions however any changes in that regard might be a good idea.
I think the only commonly used feature of splines is making them closed (i.e. a loop) but I haven't gotten around to that yet.
The other kind of spline that would be crucial for the library is a Catmull-Rom spline. With several spline types, I was thinking maybe it can be a unified function curve() or spline() which takes as input the spline type, sorta like how Processing does it?

Let me know if you think any changes are necessary for the PR. This is my first significant PR actually so I'm guessing I'm missing a lot of stuff lol.

B-Spline function with clamped/unclamped options using [De Boor's Algorithm](https://en.wikipedia.org/wiki/De_Boor%27s_algorithm).

cormullion · 2024-03-29T09:43:57Z

src/polygons.jl

@@ -775,6 +775,76 @@ function polyfit(plist::Array{Point,1}, npoints = 30)
    return resultpoly
 end

+"""
+    bspline(controlPoints::Array{Point,1}, npoints; degree=3, clamped=true)
+


If this creates a polygon, the name could start with poly to help people find it: polyspline? polybspline and to be consistent with all the other ones. But bspline could be OK too...

controlpoints or control_points is more Julian than camel-lycontrollPoints :)

polybspline makes sense since we're discussing adding other spline types I guess. Also ugh I yeah you're right for some reason the Julia community hates capital letters lol

cormullion · 2024-03-29T09:45:29Z

src/polygons.jl

+    npoints <= 0 && error("Error: npoints must be greater than zero.")
+    degree <= 0 && error("Error: degree must be greater than zero.")
+    degree >= nCP && error("Error: degree cannot be greater than the number of control points.")
+    points = []


this should be points = Point[] - might be faster

cormullion · 2024-03-29T09:47:16Z

src/polygons.jl

+    degree <= 0 && error("Error: degree must be greater than zero.")
+    degree >= nCP && error("Error: degree cannot be greater than the number of control points.")
+    points = []
+    T = []


T = Float64[] ?

cormullion · 2024-03-29T09:48:01Z

src/polygons.jl

+        p: Degree of B-spline.
+    """
+    function deBoor(k,x,t,c,p)
+        d = []


d = Float64[] ?

Actually, not necessarily. It can be d = Point[] which will satisfy most use cases. If we want to keep it capable of dealing with Float and other arithmetic data types, we'd have to keep it as []

I suppose it could be a Union? But I'm guessing the call to this internal function will always be using Points? I know that Any[] arrays are not recommended for performance reasons...

I'm unsure of whether having a union would provide any performance benefits over an Any[]. I think the most obvious way to maintain generality and performance is to make sure that the dot notation over Point as a tuple. That way the function can take in a Point[] and do both x and y calculations separately. But then again we'd lose out a lot of the contiguous memory accesses and vector operations. Do you know if there a way to make Julia precompile for two different possible data types (i.e. tell it it's either going to be a Float64[] or Point[])?

cormullion · 2024-03-29T09:49:21Z

src/polygons.jl

+    for i=0:npoints
+        t=i/npoints
+        if !clamped
+            t = t*(T[nCP-degree]-T[degree+1])+T[degree+1]


t *= should be OK/better?

cormullion · 2024-03-29T09:49:58Z

src/polygons.jl

+        end
+        k=1 #index of knot interval
+        while k < nCP
+            if (t<T[k+1])


if t < T[k+1] is OK

cormullion · 2024-03-29T09:53:01Z

src/polygons.jl

+            end
+            k += 1
+        end
+        println("t: ",t," k: ",k)


You could have @debug "t: ",t," k: ",k if you needed to keep some debugging in.

just forgot to remove that line😅

cormullion · 2024-03-29T09:59:02Z

This is great. Thanks very much!

Other things to do:

decide on name
add to exports in Luxor.jl
write tests in a test file in tests/
add test file to tests/runtest.jl
add at least a mention in docs, more if you like
add to CHANGELOG

cormullion · 2024-03-29T11:05:17Z

As for

The other kind of spline that would be crucial for the library is a Catmull-Rom spline.

perhaps just use https://github.com/JeffreySarnoff/CatmullRom.jl ?

hyrodium · 2024-03-29T11:21:42Z

Thank you for the great PR!

It's implemented in a pretty efficient way. However, I'm sure there's plenty of ways to speed it up using vectorization/parallelization.

Yeah, the performance can be improved. Here are comparisons with my package, BasicBSpline.jl:

julia> using Luxor, BasicBSpline, BenchmarkTools

julia> ps = [Point(randn(), randn()) for _ in 1:8]  # control points
8-element Vector{Point}:
 Point(-0.7994182143543495, 0.9390793075737545)
 Point(-0.22753695196095303, -0.7882913054541987)
 Point(1.2845720190835004, 0.5623233646865434)
 Point(-1.8355529584037287, 0.8438730777812197)
 Point(0.5510575401065174, 1.402746732606732)
 Point(0.45539845640833243, -1.265031994706883)
 Point(-0.1908845986894419, -0.10172057455987063)
 Point(-1.386438336193102, -0.29019694915633065)

julia> P = BSplineSpace{3}(KnotVector([0,0,0,0,1,2,3,4,5,5,5,5]))
BSplineSpace{3, Int64, KnotVector{Int64}}(KnotVector([0, 0, 0, 0, 1, 2, 3, 4, 5, 5, 5, 5]))

julia> C = BSplineManifold(ps, P)
BSplineManifold{1, (3,), Point, Int64, Tuple{BSplineSpace{3, Int64, KnotVector{Int64}}}}(Point[Point(-0.7994182143543495, 0.9390793075737545), Point(-0.22753695196095303, -0.7882913054541987), Point(1.2845720190835004, 0.5623233646865434), Point(-1.8355529584037287, 0.8438730777812197), Point(0.5510575401065174, 1.402746732606732), Point(0.45539845640833243, -1.265031994706883), Point(-0.1908845986894419, -0.10172057455987063), Point(-1.386438336193102, -0.29019694915633065)], (BSplineSpace{3, Int64, KnotVector{Int64}}(KnotVector([0, 0, 0, 0, 1, 2, 3, 4, 5, 5, 5, 5])),))

julia> C.(range(domain(P), length=21))
21-element Vector{Point}:
 Point(-0.7994182143543495, 0.9390793075737545)
 Point(-0.35288449596438315, 0.051968983954802477)
 Point(0.0612559346484744, -0.18664398382581504)
 Point(0.34154483103380284, -0.029266470466635504)
 Point(0.38652394674118207, 0.2715946493338039)
 Point(0.1400172011610212, 0.5107818139356997)
 Point(-0.27302282295295666, 0.6725347139341815)
 Point(-0.6823613770061979, 0.7884423529831135)
 Point(-0.9177637124041496, 0.890093734736359)
 Point(-0.8581722393237627, 0.99468809056206)
 Point(-0.5792380030280005, 1.0618655626854696)
 Point(-0.205789207551331, 1.0368765210461206)
 Point(0.1373459430717789, 0.864971335583544)
 Point(0.35033275819609594, 0.5204195808298326)
 Point(0.4333106007333227, 0.09356764968732044)
 Point(0.41141234698439616, -0.29621886034909783)
 Point(0.309770873250253, -0.5295743517845274)
 Point(0.13797423332130898, -0.5305523328668529)
 Point(-0.1565688090541049, -0.3968827348150751)
 Point(-0.6419943126381785, -0.26971459459097424)
 Point(-1.386438336193102, -0.29019694915633065)

julia> Luxor.bspline(ps, 20)  # The same results as previous one (except for floating point errors)
21-element Vector{Any}:
 Point(-0.7994182143543495, 0.9390793075737545)
 Point(-0.35288449596438315, 0.0519689839548025)
 Point(0.061255934648474314, -0.18664398382581504)
 Point(0.3415448310338027, -0.029266470466635615)
 Point(0.3865239467411821, 0.27159464933380395)
 Point(0.14001720116102095, 0.5107818139356994)
 Point(-0.27302282295295643, 0.6725347139341815)
 Point(-0.6823613770061976, 0.7884423529831134)
 Point(-0.9177637124041498, 0.8900937347363591)
 Point(-0.8581722393237624, 0.9946880905620599)
 Point(-0.5792380030280007, 1.0618655626854698)
 Point(-0.2057892075513309, 1.0368765210461206)
 Point(0.13734594307177872, 0.8649713355835442)
 Point(0.3503327581960958, 0.5204195808298324)
 Point(0.4333106007333225, 0.09356764968732088)
 Point(0.41141234698439616, -0.29621886034909783)
 Point(0.30977087325025293, -0.5295743517845275)
 Point(0.1379742333213091, -0.5305523328668529)
 Point(-0.15656880905410495, -0.39688273481507497)
 Point(-0.6419943126381779, -0.26971459459097435)
 Point(-1.386438336193102, -0.29019694915633065)

julia> @benchmark C.(range(domain(P), length=21))
BenchmarkTools.Trial: 10000 samples with 23 evaluations.
 Range (min … max):  950.913 ns … 219.537 μs  ┊ GC (min … max): 0.00% … 98.23%
 Time  (median):       1.001 μs               ┊ GC (median):    0.00%
 Time  (mean ± σ):     1.156 μs ±   2.916 μs  ┊ GC (mean ± σ):  3.50% ±  1.39%

  ▂▆██▆▄▂    ▁▂▂▁▁▁▃▄▄▃▂▁▁    ▁ ▁▁▁▁▂▂▂▂▂▂▁ ▁▁                  ▂
  ████████████████████████████████████████████████▇▇▇▆▆▅▆▆▆▅▄▅▅ █
  951 ns        Histogram: log(frequency) by time       1.79 μs <

 Memory estimate: 576 bytes, allocs estimate: 4.

julia> @benchmark Luxor.bspline(ps, 20)  # I have commented out `println("t: ",t," k: ",k)` in this PR.
BenchmarkTools.Trial: 10000 samples with 1 evaluation.
 Range (min … max):  16.361 μs …  3.592 ms  ┊ GC (min … max): 0.00% … 97.29%
 Time  (median):     18.024 μs              ┊ GC (median):    0.00%
 Time  (mean ± σ):   22.123 μs ± 66.233 μs  ┊ GC (mean ± σ):  6.41% ±  2.16%

  ▂▇██▆▅▄▃▄▄▄▃▂▁▁     ▁▁▁▁▂▂▂▂▃▄▃▃▂▂▂▂▁▁                      ▂
  █████████████████████████████████████████▇▇▆▆▆▆▅▆▆▆▆▅▅▅▄▄▄▅ █
  16.4 μs      Histogram: log(frequency) by time      39.5 μs <

 Memory estimate: 24.31 KiB, allocs estimate: 908.

I'm not sure where the performance difference comes from, but one reason would be type-instability as commented in #305 (comment). (Note that BasicBSpline.jl does not use vectorization/parallelization explicitly.)

As a maintainer of BasicBSpline.jl, I'm happy to help updating the code with BasicBSpline.jl. But there are more B-spline packages in Julia ecosystem, so another package might be better choice. One possible concern with BasicBSpline.jl would be increased load-time.

julia> @time using Luxor
  0.567440 seconds (546.01 k allocations: 39.472 MiB, 11.75% gc time, 3.79% compilation time)

julia> @time using BasicBSpline
  0.347235 seconds (334.66 k allocations: 17.003 MiB)

In any cases, avoiding reinventing the wheel is advisable for better interoperability with other packages and enhanced performance.

hyrodium · 2024-03-29T13:00:02Z

A B-spline curve can be split into Bezier curves, so I think adding a method to generate a BezierPath instance would be useful.
Here is an example with BasicBSpline.jl.

julia> ps  # Defined randomly as in the previous comment
8-element Vector{Point}:
 Point(-0.7994182143543495, 0.9390793075737545)
 Point(-0.22753695196095303, -0.7882913054541987)
 Point(1.2845720190835004, 0.5623233646865434)
 Point(-1.8355529584037287, 0.8438730777812197)
 Point(0.5510575401065174, 1.402746732606732)
 Point(0.45539845640833243, -1.265031994706883)
 Point(-0.1908845986894419, -0.10172057455987063)
 Point(-1.386438336193102, -0.29019694915633065)

julia> P = BSplineSpace{3}(KnotVector([0,0,0,0,1,2,3,4,5,5,5,5]))
BSplineSpace{3, Int64, KnotVector{Int64}}(KnotVector([0, 0, 0, 0, 1, 2, 3, 4, 5, 5, 5, 5]))

julia> C = BSplineManifold(ps, P)  # Same curve definition as in the previous comment
BSplineManifold{1, (3,), Point, Int64, Tuple{BSplineSpace{3, Int64, KnotVector{Int64}}}}(Point[Point(-0.7994182143543495, 0.9390793075737545), Point(-0.22753695196095303, -0.7882913054541987), Point(1.2845720190835004, 0.5623233646865434), Point(-1.8355529584037287, 0.8438730777812197), Point(0.5510575401065174, 1.402746732606732), Point(0.45539845640833243, -1.265031994706883), Point(-0.1908845986894419, -0.10172057455987063), Point(-1.386438336193102, -0.29019694915633065)], (BSplineSpace{3, Int64, KnotVector{Int64}}(KnotVector([0, 0, 0, 0, 1, 2, 3, 4, 5, 5, 5, 5])),))

julia> P2 = BSplineSpace{3}(unique(knotvector(P)) * 4)
BSplineSpace{3, Int64, KnotVector{Int64}}(KnotVector([0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 5, 5, 5, 5]))

julia> C2 = refinement(C, P2)
BSplineManifold{1, (3,), Point, Int64, Tuple{BSplineSpace{3, Int64, KnotVector{Int64}}}}(Point[Point(-0.7994182143543495, 0.9390793075737545), Point(-0.22753695196095303, -0.7882913054541987), Point(0.5285175335612737, -0.11298397038382763), Point(0.3865239467411823, 0.27159464933380395), Point(0.3865239467411823, 0.27159464933380395), Point(0.24453035992109085, 0.6561732690514355), Point(-0.7955112992413191, 0.7500231734163277), Point(-0.9177637124041496, 0.8900937347363591), Point(-0.9177637124041496, 0.8900937347363591), Point(-1.0400161255669802, 1.0301642960563906), Point(-0.24447929273023133, 1.2164555143315614), Point(0.1373459430717789, 0.8649713355835442), Point(0.1373459430717789, 0.8649713355835442), Point(0.5191711788737892, 0.513487156835527), Point(0.4872848176410608, -0.37577241893567803), Point(0.309770873250253, -0.5295743517845274), Point(0.309770873250253, -0.5295743517845274), Point(0.13225692885944523, -0.6833762846333767), Point(-0.1908845986894419, -0.10172057455987063), Point(-1.386438336193102, -0.29019694915633065)], (BSplineSpace{3, Int64, KnotVector{Int64}}(KnotVector([0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 5, 5, 5, 5])),))

julia> BezierPath([BezierPathSegment(controlpoints(C2)[1+4i:4+4i]) for i in 0:4])
5-element BezierPath:
 [Point(-0.7994182143543495, 0.9390793075737545), Point(-0.22753695196095303, -0.7882913054541987), Point(0.5285175335612737, -0.11298397038382763), Point(0.3865239467411823, 0.27159464933380395)]
 [Point(0.3865239467411823, 0.27159464933380395), Point(0.24453035992109085, 0.6561732690514355), Point(-0.7955112992413191, 0.7500231734163277), Point(-0.9177637124041496, 0.8900937347363591)]
 [Point(-0.9177637124041496, 0.8900937347363591), Point(-1.0400161255669802, 1.0301642960563906), Point(-0.24447929273023133, 1.2164555143315614), Point(0.1373459430717789, 0.8649713355835442)]
 [Point(0.1373459430717789, 0.8649713355835442), Point(0.5191711788737892, 0.513487156835527), Point(0.4872848176410608, -0.37577241893567803), Point(0.309770873250253, -0.5295743517845274)]
 [Point(0.309770873250253, -0.5295743517845274), Point(0.13225692885944523, -0.6833762846333767), Point(-0.1908845986894419, -0.10172057455987063), Point(-1.386438336193102, -0.29019694915633065)]

The curves C and C2 can be visualized with Plots.jl. They are equivalent to each other as a curve, but the control points are different. C2 has more control points, and the positions are equal to the control points of Bezier curves.

using StaticArrays, Plots
# We need to convert the type of control points from `Luxor.Point` to `StaticArrays.SVector{2, Float64}` to plot the curves.
plot(BSplineManifold(SVector{2}.(Tuple.(controlpoints(C))), bsplinespaces(C)), xlims=(-2,2), ylims=(-2,2), label="original B-spline curve")
plot(BSplineManifold(SVector{2}.(Tuple.(controlpoints(C2))), bsplinespaces(C2)), xlims=(-2,2), ylims=(-2,2), label="B-spline curve with Bezier control points")

An example script to plot the BezierPath instance with Luxor.jl:

# Invert around the y-axis because Luxor and Plots have different coordinate systems.
# Magnify the points by `*100` to update the size of the curve.
ps2 = [Point(p.x, -p.y)*100 for p in controlpoints(C2)]
bezpath = BezierPath([BezierPathSegment(ps2[1+4i:4+4i]) for i in 0:4])
@svg begin
    sethue(1.0, 0.0, 0.0)
    drawbezierpath(bezpath, action = :stroke, close=true)
    sethue(0.5, 0.5, 0.5)
    circle.(ps2, 3, action = :fill)
    poly(ps2, action = :stroke, close=false)
end

cormullion · 2024-03-29T14:37:33Z

@hyrodium Thanks for your great analysis!

>julia> @time using Luxor
  0.567440 seconds (546.01 k allocations: 39.472 MiB, 11.75% gc time, 3.79% compilation time)

julia> @time using BasicBSpline
  0.347235 seconds (334.66 k allocations: 17.003 MiB)

– this is mostly StaticArrays.jl?

julia-1.10> @time_imports using BasicBSpline
               ┌ 5.5 ms SuiteSparse_jll.__init__() 
     26.5 ms  SuiteSparse_jll 76.05% compilation time
               ┌ 5.5 ms SparseArrays.CHOLMOD.__init__() 93.91% compilation time
    154.0 ms  SparseArrays 3.36% compilation time
     12.3 ms  IntervalSets
      0.5 ms  IntervalSets → IntervalSetsRandomExt
     10.3 ms  Preferences
      0.6 ms  PrecompileTools
      1.0 ms  StaticArraysCore
    225.0 ms  StaticArrays
    356.7 ms  BasicBSpline

In any cases, avoiding reinventing the wheel is advisable for better interoperability with other packages and enhanced performance.

Yes. Due to Luxor.jl's great age (10 years old this year!) it was never very good at interoperability with other packages - it predates most of the cool and useful ones. Even StaticArrays.jl was a risky and costly dependency back when it was first developed. Recently I've started trying use external packages more - for example, replacing some of the polygon operations with PolygonAlgorithms.jl. And with the new package dependencies/extensions, we can use more external packages without incurring the cost of loading. Perhaps BasicBSplines would work as an extension?

A B-spline curve can be split into Bezier curves, so I think adding a method to generate a BezierPath instance would be useful.

Yes. I've forgotten how they work, now, but I think more Bezier handling is always good.

g-adel

The changes can be found in this commit. I don't think I can edit this PR but you can merge the commit to this PR. I can create a new one which I think is not preferred?

g-adel · 2024-03-29T18:47:33Z

src/polygons.jl

@@ -775,6 +775,76 @@ function polyfit(plist::Array{Point,1}, npoints = 30)
    return resultpoly
 end

+"""
+    bspline(controlPoints::Array{Point,1}, npoints; degree=3, clamped=true)
+


polybspline makes sense since we're discussing adding other spline types I guess. Also ugh I yeah you're right for some reason the Julia community hates capital letters lol

g-adel · 2024-03-29T18:51:31Z

src/polygons.jl

+            end
+            k += 1
+        end
+        println("t: ",t," k: ",k)


just forgot to remove that line😅

g-adel · 2024-03-29T19:01:07Z

perhaps just use https://github.com/JeffreySarnoff/CatmullRom.jl ?

Maybe... I can try it out and see if it's easily compatible with Luxor. I don't think implementing the basic algorithm will take too long though. Depends on how extensive we want the feature set for the splines and whether you're ok with adding dependencies.

hyrodium · 2024-03-30T01:55:26Z

– this is mostly StaticArrays.jl?

Yes! There is a plan to replace the dependency with StaticArraysCore.jl, but I haven't finished yet.

Perhaps BasicBSplines would work as an extension?

I think adding a constructor method BezierPath(::BasicBSpline.BSplineManifold{1, (3,), Point}) in an extension is enough and we don't have to add a new function like polybspline here.

julia> C2  # The same definition as in my previous comment
BSplineManifold{1, (3,), Point, Int64, Tuple{BSplineSpace{3, Int64, KnotVector{Int64}}}}(Point[Point(-0.7994182143543495, 0.9390793075737545), Point(-0.22753695196095303, -0.7882913054541987), Point(0.5285175335612737, -0.11298397038382763), Point(0.3865239467411823, 0.27159464933380395), Point(0.3865239467411823, 0.27159464933380395), Point(0.24453035992109085, 0.6561732690514355), Point(-0.7955112992413191, 0.7500231734163277), Point(-0.9177637124041496, 0.8900937347363591), Point(-0.9177637124041496, 0.8900937347363591), Point(-1.0400161255669802, 1.0301642960563906), Point(-0.24447929273023133, 1.2164555143315614), Point(0.1373459430717789, 0.8649713355835442), Point(0.1373459430717789, 0.8649713355835442), Point(0.5191711788737892, 0.513487156835527), Point(0.4872848176410608, -0.37577241893567803), Point(0.309770873250253, -0.5295743517845274), Point(0.309770873250253, -0.5295743517845274), Point(0.13225692885944523, -0.6833762846333767), Point(-0.1908845986894419, -0.10172057455987063), Point(-1.386438336193102, -0.29019694915633065)], (BSplineSpace{3, Int64, KnotVector{Int64}}(KnotVector([0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 5, 5, 5, 5])),))

julia> path = BezierPath([BezierPathSegment(controlpoints(C2)[1+4i:4+4i]) for i in 0:4])  # This instance will be obtained via the method `BezierPath(::BasicBSpline.BSplineManifold{1, (3,), Point})`.
5-element BezierPath:
 [Point(-0.7994182143543495, 0.9390793075737545), Point(-0.22753695196095303, -0.7882913054541987), Point(0.5285175335612737, -0.11298397038382763), Point(0.3865239467411823, 0.27159464933380395)]
 [Point(0.3865239467411823, 0.27159464933380395), Point(0.24453035992109085, 0.6561732690514355), Point(-0.7955112992413191, 0.7500231734163277), Point(-0.9177637124041496, 0.8900937347363591)]
 [Point(-0.9177637124041496, 0.8900937347363591), Point(-1.0400161255669802, 1.0301642960563906), Point(-0.24447929273023133, 1.2164555143315614), Point(0.1373459430717789, 0.8649713355835442)]
 [Point(0.1373459430717789, 0.8649713355835442), Point(0.5191711788737892, 0.513487156835527), Point(0.4872848176410608, -0.37577241893567803), Point(0.309770873250253, -0.5295743517845274)]
 [Point(0.309770873250253, -0.5295743517845274), Point(0.13225692885944523, -0.6833762846333767), Point(-0.1908845986894419, -0.10172057455987063), Point(-1.386438336193102, -0.29019694915633065)]

julia> bezierpathtopoly(path; steps=3)  # A bezier path can be converted into a vector of points. I'm not sure how to include the last end point. (See `[1:20]` in the next evaluation)
20-element Vector{Point}:
 Point(-0.7994182143543495, 0.9390793075737545)
 Point(-0.3528844959643831, 0.05196898395480254)
 Point(0.06125593464847445, -0.18664398382581515)
 Point(0.341544831033803, -0.02926647046663583)
 Point(0.3865239467411823, 0.27159464933380395)
 Point(0.14001720116102112, 0.5107818139356995)
 Point(-0.27302282295295655, 0.6725347139341815)
 Point(-0.6823613770061978, 0.7884423529831135)
 Point(-0.9177637124041496, 0.8900937347363591)
 Point(-0.8581722393237626, 0.99468809056206)
 Point(-0.5792380030280007, 1.0618655626854698)
 Point(-0.205789207551331, 1.0368765210461206)
 Point(0.1373459430717789, 0.8649713355835442)
 Point(0.35033275819609594, 0.5204195808298328)
 Point(0.43331060073332284, 0.09356764968732056)
 Point(0.4114123469843964, -0.2962188603490976)
 Point(0.309770873250253, -0.5295743517845274)
 Point(0.13797423332130898, -0.5305523328668528)
 Point(-0.15656880905410486, -0.396882734815075)
 Point(-0.6419943126381785, -0.2697145945909743)

julia> C.(range(domain(bsplinespaces(C)[1]), length=21))[1:20]  # The same results as previous one (except for floating point errors)
20-element Vector{Point}:
 Point(-0.7994182143543495, 0.9390793075737545)
 Point(-0.35288449596438315, 0.051968983954802477)
 Point(0.0612559346484744, -0.18664398382581504)
 Point(0.34154483103380284, -0.029266470466635504)
 Point(0.38652394674118207, 0.2715946493338039)
 Point(0.1400172011610212, 0.5107818139356997)
 Point(-0.27302282295295666, 0.6725347139341815)
 Point(-0.6823613770061979, 0.7884423529831135)
 Point(-0.9177637124041496, 0.890093734736359)
 Point(-0.8581722393237627, 0.99468809056206)
 Point(-0.5792380030280005, 1.0618655626854696)
 Point(-0.205789207551331, 1.0368765210461206)
 Point(0.1373459430717789, 0.864971335583544)
 Point(0.35033275819609594, 0.5204195808298326)
 Point(0.4333106007333227, 0.09356764968732044)
 Point(0.41141234698439616, -0.29621886034909783)
 Point(0.309770873250253, -0.5295743517845274)
 Point(0.13797423332130898, -0.5305523328668529)
 Point(-0.1565688090541049, -0.3968827348150751)
 Point(-0.6419943126381785, -0.26971459459097424)

The changes can be found in this commit. I don't think I can edit this PR but you can merge the commit to this PR. I can create a new one which I think is not preferred?

I'm not sure I understand what you mean. You have already added the commit 2095bbb to this PR, right?

g-adel · 2024-03-30T09:48:58Z

@hyrodium First of all thank you for the feedback! I didn't see your BasicBSpline.jl package when I wrote this code, seems really good! I just wanted to write one myself and then thought I should submit it to Luxor. I am not the most experienced with package load times and pre-compilation so I'm unsure what is the best option for Luxor. I could spend some time to make polybspline run more efficiently to somewhat match that of BasicBSpline.

I'm not sure I understand what you mean. You have already added the commit 2095bbb to this PR, right?

Ah good, yes that's the one. As I said this my first ever PR I don't have that much experience with them😅

g-adel · 2024-04-02T11:12:33Z

Has a decision been made for going with this PR or instead using BasicBSpline.jl? I'm willing to spend some time to get the function up to speed with the package and doing the remainder documentation/package complementary tasks stated in this comment.
I have to say I was quite excited to make my first ever major code contribution and to one of my favorite packages, non other. However, the quality and speed of the package should take priority so if going with the BasicBSpline is significantly better, then of course that's what you should go with.

cormullion · 2024-04-02T15:05:37Z

Yes, Bsplines are cool. Happy to have this function added, once the PR is ready to merge. Also happy to have more package extensions linking to specialized packages for other curves - although someone needs to write them first... 😂

I think it's good to have different ways of doing things here; the end result is to put graphics on drawings, so various approaches with different strengths and weaknesses with regard to performance or flexibility (assuming they're mathematically correct) is fine by me.

[Contributing to open source occasionally leads to disappointment as well as excitement - I've spent quite a few hours of my life on PRs that went nowhere. You get used to it! 😀]

Significantly improved the performance of the algorithm (mainly through eliminating a lot of type instabilities). Also added package requirements listed in JuliaGraphics#305 (comment) except for ChangeLog

g-adel · 2024-04-13T15:09:41Z

Ok the PR is now updated with significantly improved performance of the algorithm (mainly through eliminating a lot of type instabilities). Also added package requirements listed in #305 (comment) except for ChangeLog.
I've added a suite of tests that I think if passed would mean the spline is very likely to be still working. The documentation I left quite light and similar to polyfit however I did not clarify the difference between them (perhaps polyfit should be deprecated as polybspline advertises the same purpose?) As for the Changelog, perhaps @cormullion could add it as I'm not exactly sure what to write there.

Following the test made by @hyrodium here

julia> @benchmark polybspline(ps,20)
BenchmarkTools.Trial: 10000 samples with 10 evaluations.
 Range (min … max):  1.300 μs … 61.330 μs  ┊ GC (min … max): 0.00% … 93.01%
 Time  (median):     1.400 μs              ┊ GC (median):    0.00%
 Time  (mean ± σ):   1.565 μs ±  1.412 μs  ┊ GC (mean ± σ):  2.03% ±  2.28%

   ▄█     
  ▂██▆▃▂▂▃▃▃▂▂▂▁▂▂▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁ ▁
  1.3 μs         Histogram: frequency by time        3.45 μs <

 Memory estimate: 3.05 KiB, allocs estimate: 22.
 


 julia> @benchmark C.(range(domain(P), length=21))
BenchmarkTools.Trial: 10000 samples with 25 evaluations.
 Range (min … max):  928.000 ns … 90.380 μs  ┊ GC (min … max): 0.00% … 97.73%
 Time  (median):     968.000 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):     1.042 μs ±  1.541 μs  ┊ GC (mean ± σ):  2.47% ±  1.69%

  ▃█▇▄▂▂▂  ▁▄▃▁                                                ▂
  ████████▇████▇▇▆▆▆▅▃▄▁▄▄▃▃▁▁▃▆▆▇▆▆▆▅▄▆▅▆▅▅▅▆▄▅▄▅▅▅▅▃▃▁▃▄▄▅▅▅ █
  928 ns        Histogram: log(frequency) by time      2.12 μs <

 Memory estimate: 576 bytes, allocs estimate: 4.

I couldn't get the performance to match that of BasicBSpline which is currently ~150% the speed of polybspline. I've tried all the known suspects for code speed in Julia but I don't I can improve upon it without StaticArrays.jl or the new 1.11 Memory type.

cormullion · 2024-04-15T18:12:47Z

@j-adel Thanks for your contribution!

Adding bspline function (evaluates npoints based on control points)

45df34b

B-Spline function with clamped/unclamped options using [De Boor's Algorithm](https://en.wikipedia.org/wiki/De_Boor%27s_algorithm).

cormullion reviewed Mar 29, 2024

View reviewed changes

Update polybspline (prev. bspline) after first review

2095bbb

g-adel commented Mar 29, 2024

View reviewed changes

Changed d to Point[]

ff5db4b

Update to polybspline() for performance and other requirements

6a814bb

Significantly improved the performance of the algorithm (mainly through eliminating a lot of type instabilities). Also added package requirements listed in JuliaGraphics#305 (comment) except for ChangeLog

cormullion merged commit f0c6897 into JuliaGraphics:master Apr 15, 2024
5 of 8 checks passed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Adding bspline function in polygons.jl (evaluates npoints B-spline curve based on control points) #305

Adding bspline function in polygons.jl (evaluates npoints B-spline curve based on control points) #305

g-adel commented Mar 28, 2024

cormullion Mar 29, 2024

g-adel Mar 29, 2024

cormullion Mar 29, 2024

cormullion Mar 29, 2024

cormullion Mar 29, 2024

g-adel Mar 29, 2024

cormullion Mar 30, 2024

g-adel Mar 30, 2024

cormullion Mar 29, 2024

cormullion Mar 29, 2024

cormullion Mar 29, 2024

g-adel Mar 29, 2024

cormullion commented Mar 29, 2024

cormullion commented Mar 29, 2024

hyrodium commented Mar 29, 2024 •

edited

Loading

hyrodium commented Mar 29, 2024 •

edited

Loading

cormullion commented Mar 29, 2024

g-adel left a comment

g-adel Mar 29, 2024

g-adel Mar 29, 2024

g-adel commented Mar 29, 2024 •

edited

Loading

hyrodium commented Mar 30, 2024 •

edited

Loading

g-adel commented Mar 30, 2024 •

edited

Loading

g-adel commented Apr 2, 2024

cormullion commented Apr 2, 2024 •

edited

Loading

g-adel commented Apr 13, 2024

cormullion commented Apr 15, 2024

Adding bspline function in polygons.jl (evaluates npoints B-spline curve based on control points) #305

Adding bspline function in polygons.jl (evaluates npoints B-spline curve based on control points) #305

Conversation

g-adel commented Mar 28, 2024

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

cormullion commented Mar 29, 2024

cormullion commented Mar 29, 2024

hyrodium commented Mar 29, 2024 • edited Loading

hyrodium commented Mar 29, 2024 • edited Loading

cormullion commented Mar 29, 2024

g-adel left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

g-adel commented Mar 29, 2024 • edited Loading

hyrodium commented Mar 30, 2024 • edited Loading

g-adel commented Mar 30, 2024 • edited Loading

g-adel commented Apr 2, 2024

cormullion commented Apr 2, 2024 • edited Loading

g-adel commented Apr 13, 2024

cormullion commented Apr 15, 2024

hyrodium commented Mar 29, 2024 •

edited

Loading

hyrodium commented Mar 29, 2024 •

edited

Loading

g-adel commented Mar 29, 2024 •

edited

Loading

hyrodium commented Mar 30, 2024 •

edited

Loading

g-adel commented Mar 30, 2024 •

edited

Loading

cormullion commented Apr 2, 2024 •

edited

Loading