77
88# functions for cylinder rotation
99def cyl2cart (r , theta , z ):
10- return (r * np .cos (theta ), r * np .sin (theta ), z )
10+ return (r * np .cos (theta ), r * np .sin (theta ), z )
1111
1212
1313def roll (R , zi , zf , RT ):
@@ -20,23 +20,26 @@ def roll(R, zi, zf, RT):
2020 x , y , z = cyl2cart (r , t , z )
2121
2222 # Euler rotation
23- rot = np .dot (
24- RT .rot ().to_array (),
25- np .array ([x .ravel (), y .ravel (), z .ravel ()])
26- )
23+ rot = np .dot (RT .rot ().to_array (), np .array ([x .ravel (), y .ravel (), z .ravel ()]))
2724
2825 x_rt = rot [0 , :].reshape (x .shape ) + RTw .trans ().to_array ()[0 ]
2926 y_rt = rot [1 , :].reshape (y .shape ) + RTw .trans ().to_array ()[1 ]
3027 z_rt = rot [2 , :].reshape (z .shape ) + RTw .trans ().to_array ()[2 ]
3128
32- return x_rt , y_rt , z_rt ,
29+ return (
30+ x_rt ,
31+ y_rt ,
32+ z_rt ,
33+ )
34+
3335
3436# model
3537model = biorbd .Model ("WrappingObjectExample.bioMod" )
3638
37- # wrapping RT matrix
38- RTw = biorbd .WrappingCylinder (
39- model .muscle (0 ).pathModifier ().object (0 )).RT (model , np .array ([0.1 , 0.1 , 0.1 , 0.1 , 0.1 , 0.1 ]))
39+ # wrapping RT matrix
40+ RTw = biorbd .WrappingCylinder (model .muscle (0 ).pathModifier ().object (0 )).RT (
41+ model , np .array ([0.1 , 0.1 , 0.1 , 0.1 , 0.1 , 0.1 ])
42+ )
4043RTw_trans = RTw .transpose ()
4144
4245# Point in global base
@@ -47,27 +50,29 @@ def roll(R, zi, zf, RT):
4750
4851# Plot the surface
4952fig = plt .figure ()
50- ax = fig .add_subplot (111 , projection = '3d' )
53+ ax = fig .add_subplot (111 , projection = "3d" )
5154x , y , z = roll (0.05 , - 0.2 , 0.2 , RTw )
5255ax .plot_surface (x , y , z )
53- plt .plot ((Po [0 ], Po_wrap [0 ]), (Po [1 ], Po_wrap [1 ]), (Po [2 ], Po_wrap [2 ]), marker = 'o' )
54- plt .plot ((Pi [0 ], Pi_wrap [0 ]), (Pi [1 ], Pi_wrap [1 ]), (Pi [2 ], Pi_wrap [2 ]), marker = 'o' )
56+ plt .plot ((Po [0 ], Po_wrap [0 ]), (Po [1 ], Po_wrap [1 ]), (Po [2 ], Po_wrap [2 ]), marker = "o" )
57+ plt .plot ((Pi [0 ], Pi_wrap [0 ]), (Pi [1 ], Pi_wrap [1 ]), (Pi [2 ], Pi_wrap [2 ]), marker = "o" )
5558
5659# To set Axe3D to (almost) equal
57- max_range = np .array ([x .max ()- x .min (), y .max ()- y .min (), z .max ()- z .min ()]).max ()
58- Xb = 0.5 * max_range * np .mgrid [- 1 :2 :2 , - 1 :2 :2 , - 1 :2 :2 ][0 ].flatten () + 0.5 * (x .max ()+ x .min ())
59- Yb = 0.5 * max_range * np .mgrid [- 1 :2 :2 , - 1 :2 :2 , - 1 :2 :2 ][1 ].flatten () + 0.5 * (y .max ()+ y .min ())
60- Zb = 0.5 * max_range * np .mgrid [- 1 :2 :2 , - 1 :2 :2 , - 1 :2 :2 ][2 ].flatten () + 0.5 * (z .max ()+ z .min ())
60+ max_range = np .array ([x .max () - x .min (), y .max () - y .min (), z .max () - z .min ()]).max ()
61+ Xb = 0.5 * max_range * np .mgrid [- 1 :2 :2 , - 1 :2 :2 , - 1 :2 :2 ][0 ].flatten () + 0.5 * (x .max () + x .min ())
62+ Yb = 0.5 * max_range * np .mgrid [- 1 :2 :2 , - 1 :2 :2 , - 1 :2 :2 ][1 ].flatten () + 0.5 * (y .max () + y .min ())
63+ Zb = 0.5 * max_range * np .mgrid [- 1 :2 :2 , - 1 :2 :2 , - 1 :2 :2 ][2 ].flatten () + 0.5 * (z .max () + z .min ())
6164for xb , yb , zb in zip (Xb , Yb , Zb ):
62- ax .plot ([xb ], [yb ], [zb ], 'w' )
65+ ax .plot ([xb ], [yb ], [zb ], "w" )
6366plt .show ()
6467
6568# Compute muscle length
6669r = biorbd .WrappingCylinder (model .muscle (0 ).pathModifier ().object (0 )).radius ()
6770l_muscle_bd = model .muscle (0 ).musculoTendonLength (model , np .array ([0.1 , 0.1 , 0.1 , 0.1 , 0.1 , 0.1 ]))
6871
6972
70- l_wt_arc = np .sqrt ((Po_wrap [0 ]- Po [0 ])** 2 + (Po_wrap [1 ]- Po [1 ])** 2 + (Po_wrap [2 ]- Po [2 ])** 2 ) + np .sqrt ((Pi_wrap [0 ]- Pi [0 ])** 2 + (Pi_wrap [1 ]- Pi [1 ])** 2 + (Pi_wrap [2 ]- Pi [2 ])** 2 )
73+ l_wt_arc = np .sqrt ((Po_wrap [0 ] - Po [0 ]) ** 2 + (Po_wrap [1 ] - Po [1 ]) ** 2 + (Po_wrap [2 ] - Po [2 ]) ** 2 ) + np .sqrt (
74+ (Pi_wrap [0 ] - Pi [0 ]) ** 2 + (Pi_wrap [1 ] - Pi [1 ]) ** 2 + (Pi_wrap [2 ] - Pi [2 ]) ** 2
75+ )
7176
7277# Pi_wrap and Po_wrap provide in wrapping base to compute arc length
7378Pi_wrap = biorbd .Vector3d (0.28997869688863276 , 0.34739632020407846 , 0.6636920365713757 )
@@ -78,10 +83,14 @@ def roll(R, zi, zf, RT):
7883Po_wrap .applyRT (RTw_trans )
7984Po_wrap = Po_wrap .to_array ()
8085
81- arc = math .acos (
82- ((Po_wrap [0 ] * Pi_wrap [0 ]) + (Pi_wrap [1 ] * Po_wrap [1 ])) / (
83- math .sqrt (((Pi_wrap [1 ]** 2 )+ (Pi_wrap [0 ]** 2 )) * ((Po_wrap [1 ]** 2 ) + (Po_wrap [0 ]** 2 ))))) * r
86+ arc = (
87+ math .acos (
88+ ((Po_wrap [0 ] * Pi_wrap [0 ]) + (Pi_wrap [1 ] * Po_wrap [1 ]))
89+ / (math .sqrt (((Pi_wrap [1 ] ** 2 ) + (Pi_wrap [0 ] ** 2 )) * ((Po_wrap [1 ] ** 2 ) + (Po_wrap [0 ] ** 2 ))))
90+ )
91+ * r
92+ )
8493
85- l_arc = math .sqrt (arc ** 2 + (Po_wrap [2 ]- Pi_wrap [2 ])** 2 )
94+ l_arc = math .sqrt (arc ** 2 + (Po_wrap [2 ] - Pi_wrap [2 ]) ** 2 )
8695
87- l_m = l_wt_arc + l_arc # l_m = 0.582246096069153
96+ l_m = l_wt_arc + l_arc # l_m = 0.582246096069153
0 commit comments