#!/usr/bin/env python # -*- coding: utf-8 -*- """Example b-spline basis functions and curves. :Date: 2020-01-02 .. module:: function_b_spline :platform: *nix, Windows :synopsis: Example b-spline basis functions and curves. .. moduleauthor:: Daniel Weschke """ from mpl_toolkits.mplot3d import Axes3D assert Axes3D import pylab from pylib.data import seq, unique_list from pylib.function import ( sample_half_open, sample_half_open_seq, b_spline_basis, b_spline_knots, b_spline_curve_with_knots) from pylib.helper import timeit def b_spline_basis_functions(knots, pmax=None): m = len(knots) - 1 # m = 3 if pmax is None or pmax > m: pmax = m u = seq(knots[0], knots[-1], (knots[-1]-knots[0])/50) Nij = [[[b_spline_basis(knots, i, j)(uj) for uj in u] for i in range(m-j)] for j in range(pmax+1)] Nj = [[sum(i) for i in zip(*Nij[:][j])] for j in range(pmax+1)] return u, Nij, Nj # @timeit('plot_matrix_seperate_single_sum') def plot_matrix_seperate_single_sum(u, Nij, Nj, pmax=None, figure_id=''): """ :: i: knot span p: degree sum of all i=0 i=1 i=2 ... basis function +-----+-----+-----+-----+-----+ p=0 | | | | | | +-----+-----+-----+-----+-----+ p=1 | | | | | | +-----+-----+-----+-----+-----+ ... | | | | | | """ p = len(Nij[0]) - 1 # degree if pmax is None or pmax > p: pmax = p m = len(Nij[0]) # knot spans pylab.figure(figure_id) [([(pylab.subplot(pmax+1, m+1, i*(m+1)+j+1), pylab.plot(u, Nij[i][j])) for j in range(m-i)], pylab.subplot(pmax+1, m+1, (i+1)*(m+1)), pylab.plot(u, Nj[i])) for i in range(pmax+1)] # @timeit('plot_matrix_joint_single_sum') def plot_matrix_joint_single_sum(u, Nij, Nj, pmax=None, figure_id=''): """ :: p: degree p=0 p=1 p=2 ... +-----+-----+-----+----- singe basis functions | | | | +-----+-----+-----+----- sum of all basis functions | | | | +-----+-----+-----+----- """ m = len(Nij[0]) # knot spans p = m - 1 # degree n = len(Nij) - 1 # degrees calculated if pmax is None or pmax > p: pmax = p if n < pmax: pmax = n pylab.figure(figure_id) [(pylab.subplot(2, pmax+1, 0*(pmax+1)+j+1),[pylab.plot(u, Nij[j][i]) for i in range(p+1-j)], pylab.subplot(2, pmax+1, 1*(pmax+1)+j+1).plot(u, Nj[j])) for j in range(pmax+1)] # @timeit('plot_b_spline_curve_with_knots') def plot_b_spline_curve_with_knots( Cu=None, C=None, U=None, CU=None, P=None, knots=True, figure_id='', ax=None): if ax is None: fig = pylab.figure(figure_id) ax = fig.add_subplot(111, projection='3d') # control points if P is not None: x, y, z = zip(*P) ax.plot(x, y, z, 'C1', linewidth=1) ax.plot(x, y, z, 'oC1', fillstyle='none') # B-spline curve if Cu is None and C is not None and U is not None: Cu = sample_half_open(C, U[0], U[-1]) if Cu is not None: ax.plot(*zip(*Cu)) if knots: # internal knots of clamped B-spline curve if CU is None and C is not None and U is not None: CU = sample_half_open_seq(C, unique_list(U)) if CU is not None: ax.plot(*zip(*CU), '.C2') return ax def example1(): # degree of basis functions # p = 0 # knots u0 = 0 u1 = 1 u2 = 2 u3 = 3 # knot vector U = (u0, u1, u2, u3) # knot spans 0, 1 and 2 are [0,1), [1,2), [2,3) # m = len(U) - 1 = 3 # basis functions of degree 0 are # N_0,0(u) = 1 on [0,1) and 0 elsewhere, # N_1,0(u) = 1 on [1,2) and 0 elsewhere, and # N_2,0(u) = 1 on [2,3) and 0 elsewhere. u, Nij, Nj = b_spline_basis_functions(U) plot_matrix_seperate_single_sum(u, Nij, Nj, figure_id='N1') def example2(): U = (0, 0.25, 0.5, 0.75, 1) # m = 4 u, Nij, Nj = b_spline_basis_functions(U) plot_matrix_seperate_single_sum(u, Nij, Nj, pmax=3, figure_id='N2') def example3(): U = (0, 0, 0, 0.3, 0.5, 0.5, 0.6, 1, 1, 1) u, Nij, Nj = b_spline_basis_functions(U) plot_matrix_joint_single_sum(u, Nij, Nj, pmax=3, figure_id='N3') def example4(): # clamped B-spline curves U = (0, 0, 0, 0, 0.14, 0.28, 0.42, 0.57, 0.71, 0.85, 1, 1, 1, 1) m = len(U) - 1 # 13 P = [[0,1], [.1,2], [.2,2], [.3,1], [.4,1], [.5,0], [.6,0], [.7,1], [.6,1], [.5,2]] P = [[0,1,0], [.1,2,0], [.2,2,0], [.3,1,1], [.4,1,0], [.5,0,0], [.6,0,0], [.7,1,0], [.6,1,0], [.5,2,0]] n = len(P) - 1 # 9 p = m - n - 1 C = b_spline_curve_with_knots(p, P, U) # Cu = sample_half_open(C, U[0], U[-1]) plot_b_spline_curve_with_knots(C=C, U=U, P=P, figure_id='N4') def example5(): P = [(48.0888108620491, 41.3794408583987, -0.372884555563078), (48.2156978454883, 41.4376819952444, -0.372981234249012), (48.3421249992392, 41.4949853259561, -0.373021507495587), (48.6647993346146, 41.6362394842807, -0.372960532155267), (48.8631931562637, 41.7184356658063, -0.372941685915525), (49.2627469867943, 41.8770863299231, -0.372965702682335), (49.4634735626356, 41.9533782153992, -0.37296977238557), (49.8669592828075, 42.1000733004711, -0.372963214968813), (50.0698102094774, 42.1705172272789, -0.372961547924453), (50.4773734605376, 42.3060584254141, -0.372966255857575), (50.6820075976044, 42.371132473348, -0.372964768740875), (51.0927404446479, 42.4960106864461, -0.37295550426811), (51.2989587087642, 42.5558597303533, -0.37296082550931), (51.7131511901618, 42.6711711778783, -0.372997134411391), (51.9206405900576, 42.7264781422536, -0.37297347298921), (52.3354542080309, 42.8306222566775, -0.372830008862485), (52.5444045251438, 42.8799268641912, -0.372905635604393), (52.7547311828766, 42.9295775362814, -0.373148784408123)] Uu = (0.0228378692889125, 0.0649893839528174, 0.130414237399236, 0.195679174703882, 0.260898113556269, 0.326065686455556, 0.391167888702703, 0.456323551413636, 0.521155975846261) mult = (4, 2, 2, 2, 2, 2, 2, 2, 4) U = [Uu[i] for i in range(len(Uu)) for j in range(mult[i])] def basis(): with timeit('b_spline_basis_functions'): u, Nij, Nj = b_spline_basis_functions(U, pmax=11) plot_matrix_joint_single_sum(u, Nij, Nj, figure_id='N5 N') def curve(knots_cp=True): n = len(P) - 1 m = len(U) - 1 p = m - n - 1 with timeit('b_spline_curve_with_knots_sampled'): C = b_spline_curve_with_knots(p, P, U) Cu = sample_half_open(C, U[0], U[-1]) plot_b_spline_curve_with_knots( Cu=Cu, CU=sample_half_open_seq(C, Uu) if knots_cp else None, P=P if knots_cp else None, figure_id='N5 C') pylab.gca().set_zlim3d([-1, 0]) #basis() curve() #curve(False) def example6(): # # given data, data to compare with own generated data # P = [(48.0888108620491, 41.3794408583987, -0.372884555563078), (48.2156978454883, 41.4376819952444, -0.372981234249012), (48.3421249992392, 41.4949853259561, -0.373021507495587), (48.6647993346146, 41.6362394842807, -0.372960532155267), (48.8631931562637, 41.7184356658063, -0.372941685915525), (49.2627469867943, 41.8770863299231, -0.372965702682335), (49.4634735626356, 41.9533782153992, -0.37296977238557), (49.8669592828075, 42.1000733004711, -0.372963214968813), (50.0698102094774, 42.1705172272789, -0.372961547924453), (50.4773734605376, 42.3060584254141, -0.372966255857575), (50.6820075976044, 42.371132473348, -0.372964768740875), (51.0927404446479, 42.4960106864461, -0.37295550426811), (51.2989587087642, 42.5558597303533, -0.37296082550931), (51.7131511901618, 42.6711711778783, -0.372997134411391), (51.9206405900576, 42.7264781422536, -0.37297347298921), (52.3354542080309, 42.8306222566775, -0.372830008862485), (52.5444045251438, 42.8799268641912, -0.372905635604393), (52.7547311828766, 42.9295775362814, -0.373148784408123)] Uu = (0.0228378692889125, 0.0649893839528174, 0.130414237399236, 0.195679174703882, 0.260898113556269, 0.326065686455556, 0.391167888702703, 0.456323551413636, 0.521155975846261) mult = (4, 2, 2, 2, 2, 2, 2, 2, 4) U = [Uu[i] for i in range(len(Uu)) for j in range(mult[i])] n = len(P) - 1 # number of control points - 1 m = len(U) -1 # number of knot spans p = m - n - 1 # degree # print(p) C = b_spline_curve_with_knots(p, P, U) Cu = sample_half_open(C, U[0], U[-1]) ax = plot_b_spline_curve_with_knots(Cu=Cu, figure_id='N6') # # now generate own knots # p = 3 # degree m = n + p + 1 # number of knot spans U = b_spline_knots(n, p) # U = (0, 0, 0, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 14, 14, 14) # print(len(U), U) C = b_spline_curve_with_knots(p, P, U) Cu = sample_half_open(C, U[0], U[-1]) plot_b_spline_curve_with_knots(Cu=Cu, C=C, U=U, ax=ax) #pylab.gca().set_zlim3d([-1, 0]) def example7(): P = [(-69.99999999999, 6.23746531921, 119.99999999999), (-69.99999999999, 5.72567326271, 120.0), (-70.0, 5.2192123782, 120.0), (-70.0, 4.718082665681, 120.0), (-70.0, 4.222284125153, 120.0), (-70.0, 3.731816756617, 120.0), (-70.0, 3.246680560071, 120.0), (-70.0, 2.766875535516, 120.0), (-70.0, 2.292401682953, 120.0), (-70.0, 1.82325900238, 120.0), (-70.0, 1.359447493799, 120.0), (-70.0, 0.900967157208, 120.0), (-70.0, 0.447817992609, 120.0), (-69.99999999999, -1.642767914458, 119.99999999999), (-69.95823302652, -3.284128624962, 119.96005936915), (-69.87482030195, -4.920567121601, 119.88028704462), (-69.7502001661, -6.548567099244, 119.76113869291), (-69.58513326705, -8.164660004507, 119.60341679526), (-69.38069590753, -9.765475979336, 119.40825949383), (-69.13826152984, -11.34779669959, 119.17711816806), (-68.85947036564, -12.90861011093, 118.91172378914), (-68.54618721656, -14.44516705927, 118.61404198996), (-68.20044740283, -15.95503981856, 118.28621691953), (-67.82439084537, -17.43618251348, 117.93050381741), (-67.42018430873, -18.88699343893, 117.54919035943), (-66.38758251206, -22.29350540161, 116.57795367665), (-65.73413448053, -24.21919416036, 115.96554908806), (-65.03556391343, -26.07956263873, 115.31370192302), (-64.29648522575, -27.87744521951, 114.62724902124), (-63.52042791056, -29.61397650413, 113.91024118834), (-62.71277001846, -31.2871331673, 113.16847249966), (-61.87901674025, -32.89613508726, 112.4077395496), (-61.02226890253, -34.44424317029, 111.63151855295), (-60.14386048877, -35.93622309031, 110.841629897), (-59.24623340527, -37.37387103089, 110.04091803203), (-58.33313678246, -38.756819007, 109.23340236149), (-57.40516397009, -40.08965381341, 108.42014822445), (-55.01642244455, -43.34968716212, 106.34685334397), (-53.53813210335, -45.21002358427, 105.08230356996), (-52.03458081272, -46.96355473328, 103.81650090025), (-50.50977333246, -48.61777575741, 102.55504546095), (-48.96692322622, -50.17949288551, 101.30260097008), (-47.40801072014, -51.65517742458, 100.06277532323), (-45.83395649691, -53.05076391735, 98.838417665337), (-44.24537288755, -54.37123987567, 97.632189028576), (-42.64319994423, -55.62048722024, 96.446976549361), (-41.02853017217, -56.80167054344, 95.285780596074), (-39.40157208613, -57.91802439479, 94.15110945761), (-37.76099941739, -58.97317114245, 93.044667943398), (-33.51657488475, -61.52870256173, 90.285109332808), (-30.89322834847, -62.94313830978, 88.677041723595), (-28.23238814008, -64.22251443102, 87.15321189152), (-25.55002271787, -65.36600940426, 85.732620462342), (-22.80465037824, -66.39367032451, 84.407541386373), (-20.07652968241, -67.27260606721, 83.231827543651), (-17.27816204332, -68.0519427367, 82.167075623263), (-14.44812618974, -68.69718956814, 81.25538726298), (-11.58521735637, -69.2186257426, 80.51004437835), (-8.697808521544, -69.60868640648, 79.943540692506), (-5.807821595725, -69.86923851767, 79.563401742733), (-2.906442029202, -70.00000000001, 79.372539331945), (0.372271370757, -69.99999999999, 79.372539331937), (0.762465017147, -69.99999999999, 79.372539331936), (1.17058093917, -69.99999999999, 79.372539331935), (1.596619136825, -69.99999999999, 79.372539331934), (2.040579610114, -69.99999999999, 79.372539331932), (2.502462359036, -69.99999999999, 79.372539331931), (2.982267383591, -69.99999999998, 79.37253933193), (3.479994683778, -69.99999999998, 79.372539331929), (3.995644259599, -69.99999999998, 79.372539331927), (4.529216111052, -69.99999999998, 79.372539331926), (5.080710238139, -69.99999999998, 79.372539331925), (5.650126640858, -69.99999999997, 79.372539331923), (6.23746531921, -69.99999999997, 79.372539331922)] Uu = (-3.932977778518, 1.445001333941e-14, 14.427668850884, 34.626262491172, 65.717581403178, 114.40570375563, 120.64191760942) mult = (14, 12, 12, 12, 12, 12, 14) U = [Uu[i] for i in range(len(Uu)) for j in range(mult[i])] n = len(P) - 1 # number of control points - 1 m = len(U) -1 # number of knot spans p = m - n - 1 # degree # print(p) with timeit('b_spline_curve_with_knots'): C = b_spline_curve_with_knots(p, P, U) with timeit('sampling'): Cu = sample_half_open(C, U[0], U[-1]) with timeit('plotting'): plot_b_spline_curve_with_knots(Cu=Cu, figure_id='N7') # b_spline_curve_with_knots took 0.005 ms, 0.004 ms, 0.002 ms # sampling took 20542.242 ms, 21172.614 ms, 21136.981 ms # plotting took 29.950 ms, 28.201 ms, 26.364 ms if __name__ == "__main__": example7() pylab.show()