計測自動制御学会東北支部 第 342 回研究集会 (2023.5.26)
資料番号 342-4

力覚信号に基づく粒子フィルタを用いた
多数のパラメータで表現される道具形状の推定
Estimation of a tool shape with a number of parameters
using a particle ﬁlter based on force signals○沓澤京 ∗ ，林部充宏 ∗○ Kyo Kutsuzawa∗ , Mitsuhiro Hayashibe∗*東北大学*Tohoku University
キーワード : 力覚信号処理 (force-signal processing)，内在的接触計測 (intrinsic contact sensing)，粒子
フィルタ (particle ﬁlter)，道具使用 (robotic tool use)

連絡先 : 〒 980-8579 宮城県仙台市青葉区荒巻字青葉 6-6-01
東北大学 大学院工学研究科 ロボティクス専攻 沓澤京

Tel.: 022-795-5682，E-mail: kutsuzawa@tohoku.ac.jp

1.

序論の論文 4) が挙げられ、それ以降も「力とモーメ

人間が道具越しに環境に触れたとき、たとえ目をつむっていても、道具のどこにどの程度の力が加わったかを把握することができる。当然ながら道具には神経も感覚器もないので、人間は手元に伝わってきた力の情報から、道具上での接触状態を推定していることになる。身体感ントから剛体上の接触位置を推定する課題」を主として、様々に研究されている。この技術は環境との接触を伴う制御に有用である。例えば
接触状態を制御する位置と力のハイブリッド制
御 5, 6) では、接触点を原点とする座標系を扱うために接触位置を知る必要がある。
力覚信号から道具上の接触位置などを求める

覚を外部の物体まで拡張するこの能力は、物体を巧みに操って作業する上で非常に有用と考えられる。
力触覚を基に測定点から離れた箇所の接触状態を推定する技術は、intrinsic tactile sensing1)
あるいは intrinsic contact sensing2, 3) と呼ばれる。決まった和訳はないが、intrinsic contact

sensing を直訳するなら「内在的接触計測」などとなるだろう。最も古い研究事例には Salisbury

ために、従来手法のほとんどでは道具の形状が既知である必要があった 2, 4) 。しかし、力覚セン
サに対する道具の形状が常に正確にわかるとは限らない。開かれた環境で臨機応変な道具使用
や物体操作をおこなうときには道具もしくは物
体の形状はふつう未知である。道具形状をまったく必要としない手法 3, 7) も提案されているが、推定には接触力の擾乱が必要であり、加えて推定の収束が原理的に遅いという弱点がある。

–1–

道具の形状が与えられない条件で素早い接触

慣性力の影響は事前に取り除いてあるとする。

位置推定を実現するためには、道具形状を同時すると、接触点の位置 c は以下のように求まる。

に推定することが有用である。粒子フィルタを利用することで、接触位置と道具形状とを同時に推定する手法が提案された 8) 。粒子フィルタ
は確率分布を点群で近似して扱う手法であり、不確実性を強く伴う推定において有用である。
一方で、推定したい変数の数（状態変数の次元数）が増えるに伴い、近似に必要な粒子数が指数関数的に増加してしまう。そのため、先行研
究 8) の手法では道具形状が少数の変数で表現さc=

位置は一意に推定できないことがわかる。代わりに、接触点を通り接触力と平行な直線が、接触位置の候補直線として求まる。既存手法 4) では基本的に、道具の形状が既知として、候補直線と道具表面との交差点が接触位置であるとして解を絞り込む。

そこで本研究では、確率モデルの修正によっc=て、道具の形状が極めて多量の変数で表現される

multaneous localization and mapping（SLAM）
9, 10, 11) で使われる粒子フィルタと同様の工夫によって、粒子群を形状パラメータの空間に散布する必要がないように確率モデルを分解する
ことで、1000 次元以上のパラメータをわずか 30個の粒子で推定できるようになる。本稿では、平面状でのシミュレーションによって検証をおこなう。平面上ではあるが、道具形状をグリッド形式で表現すると、そのパラメータ数は容易に100 を超える。

(3)

ここで、α は未知定数である。式 (3) から、接触れる場合にしか適用されていなかった。

場合にも利用可能な粒子フィルタを構築する。si-

F ×M+ αF .
∥F ∥2

F ×M+ αF ,∥F ∥2s.t. g(c) = 0.

(4)

ここで g は、道具表面において値がゼロになる
ような関数である。
上記技術は、環境に接触する対象物の表面に分布型センサを敷き詰める必要がなく、単一の6 軸力覚センサでセンサから離れた箇所の接触を計測できることから、外部物体の接触を扱う
ときに有用である。この技術はさまざまに改良がされており、道具の運動速度情報 12) や力の微分値 13) などを使って精度を向上することができる。また、多点接触への拡張 14, 15, 16) やロボット制御への応用 14, 17, 18) もされている。
力の瞬時値のみを利用するのでは、原理的に式 (3) の情報までしかわからないため、性能に

2.

基本原理

力覚センサが取り付けられた剛体の道具上に、環境から力が印加されたときを考える。このとき、接触位置を c、接触力を f とする。また、力覚センサで計測される力とモーメントをそれぞれ F 、M とする。道具について、以下のよう限界がある。過去の測定値も活用することで原理的な性能向上が期待できる 4) 。接触位置の変化よりも接触力ベクトルの変動のほうが十分速いと仮定することで、非凸包の物体への拡張 19)
や、道具形状を必要としない接触位置推定 3, 7)が可能になる。

な力のつり合いが成り立つと仮定する。

F = f,(1)

M = c × f.

(2)

3.

手法

提案手法は、無香粒子フィルタ 20, 21) を使って、接触位置 c と道具形状を表すパラメータ s

環境からの接触は 1 点のみで生じ、接触点におとを同時推定する。本研究では s をグリッドマッ

いてトルクが加わらないとする。また、重力やプとして表現する。グリッドマップの各セルで–2–

形状パラメータ
1.0

存在確率

粒子

力覚センサ

cy

形状の表現方法
接触力 f0.0

形状パラメータ道具
赤色: 真の形状
黄色: 形状あり
紫色: 形状なし

cx

存在確率

接触位置 c1.0

0.0

Fig. 1 提案手法の概要。道具の形状をグリッドマップで表現する。粒子フィルタにて粒子群は接触位置の空間にのみ分布され、各粒子は各々で形状パラメータの推定値をもつ。
は、道具表面が存在するかしないかがスカラ値の確率分布は未知である。そこで、以下のようとして表現される。これによってさまざまな形に式変形して、提案分布と呼ばれる、設計者が状が表現できるようになるが、一方で、次元数任意に決める確率分布 q(x0:k |y0:k ) の期待値にが極めて大きくなる（2 次元平面の各軸を 100 分

置き換える。

割すると、パラメータの次元数は 1002 = 10000にもなる）。そのため、単純に x ≜ (c, s) として x を推定しようとすると、いわゆる次元の呪いによって指数関数的な数の粒子が必要になってしまう。本研究では、Rao–Blackwellized 粒子
フィルタ 9) と同様の工夫によって、粒子群を接触位置 c の空間にのみ散布すれば良いように確率モデルを修正して、s の次元数によらず少数の粒子で計算できるようにする。手法の概念図を図 1 に示す。

3.1

粒子フィルタ

Ex0:k ∼p(x0:k |y0:k ) [ϕ(x0:k )]Z= ϕ(x0:k )p(x0:k |y0:k )dx0:kZp(x0:k |y0:k )= ϕ(x0:k )q(x0:k |y0:k )dx0:kq(x0:k |y0:k )Zp(y0:k |x0:k )p(x0:k )q(x0:k |y0:k )dx0:k= ϕ(x0:k )p(y0:k )q(x0:k |y0:k )R|x0:k )p(x0:k )ϕ(x0:k ) p(y0:kq(x0:k |y0:k )dx0:kq(x0:k |y0:k )=R p(y0:k |x0:k )p(x0:k )q(x0:k |y0:k )dx0:kq(x0:k |y0:k )R
ϕ(x0:k )wk q(x0:k |y0:k )dx0:kR=
wk q(x0:k |y0:k )dx0:k
Ex0:k ∼q(x0:k |y0:k ) [ϕ(x0:k )wk ]=(6)
Ex0:k ∼q(x0:k |y0:k ) [wk ]まずは一般的な粒子フィルタの確率モデルに

ここで、wk は確率分布 p(x0:k |y0:k ) と q(x0:k |y0:k )ついて考える。時刻 k における観測変数を yk 、との違いを補正する重み係数であり、以下のよ

状態変数を xk とする。また、時刻 0 から k まで

うに表される。

の状態変数の時系列を以下のように表記する。

x0:k ≜ (x0 , x1 , . . . , xk )

wk ≜
(5)

本研究では、状態変数を xk ≜ (ck , sk ) とし、観測変数を yk ≜ (Fk , Mk ) とする。推定の目標p(y0:k |x0:k )p(x0:k )q(x0:k |y0:k )

(7)

さらに、

q(x0:k |y0:k ) = q(xk |x0:k−1 , y0:k )q(x0:k−1 |y0:k )は、期待値 Ex0:k ∼p(x0:k |y0:k ) [ϕ(x0:k )] を求めることである。ただし ϕ は任意の写像とする。この期待値は p(x0:k |y0:k ) に関するものであるが、こ–3–= q(xk |x0:k−1 , y0:k )q(x0:k−1 |y0:k−1 )(8)のように定義することで、wk を再帰的な形に変

この仮定は、形状パラメータ sk が過去の状態変形できる。

数および観測変数から一意に定まることを意味

p(yk |xk )p(xk |xk−1 )
wk =
wk−1
q(xk |x0:k−1 , y0:k )

(i)

(9)

する。これによって、粒子 xk は以下のように得られるようになる。

一般的な粒子フィルタでは、上記の期待値を(i)

以下のように N 個の粒子群の平均で近似する。

Ex0:k ∼q(x0:k |y0:k ) [ϕ(x0:k )wk ]≈Ex0:k ∼q(x0:k |y0:k ) [wk ]

(10)
(i)
ここで x0:k は提案分布からサンプリングされる(i)(i)
粒子であり、wk は粒子 x0:k から計算される。

(i)

(i)

(16)

(i)

ck は高々3 次元であり、q(ck |c0:k−1 , M0:k , F0:k )は無香カルマンフィルタで実用的に実装できる。
(i)

このとき wk は以下のようになる。
(i)

wk =

(11)

(i)

(i)

(i)

p(Mk , Fk |ck )p(ck |ck−1 , sk )(i) (i)q(ck |c0:k−1 , M0:k , F0:k )

(i)

wk−1
(17)

(12)

ただし、個々の確率分布は以下のように定義する。

このとき、提案分布 q(x0:k |y0:k ) が元の確率分

(i)

(i)

(i)

(i)

(i)

(i)

(i)

(i)

p(ck |ck−1 , sk ) = p(ck |ck−1 )p(ck |sk ),(19)近似できる。無香粒子フィルタでは、提案分布(i)q(xk |x0:k−1 , y0:k ) を無香カルマンフィルタで実
装する。

(i)

p(Mk , Fk |ck ) = p(Mk |Fk , ck )p(Fk ), (18)

布 p(x0:k |y0:k ) に近いほど期待値を高い精度で

さらに
(i)

3.2

(15)

(i)

(i)

(i)
(i)
x0:k ∼ q(x0:k |y0:k ),(i)(i) (i)p(yk |xk )p(xk |xk−1 ) (i)
(i)
wk =
wk−1
(i) (i)q(xk |x0:k−1 , y0:k )(i)sk = ψ(sk−1 , ck−1 , Fk−1 )

PN −1

(i)
(i)
i=0 ϕ(x0:k )wk
PN −1 (i)i=0 wk(i)ck ∼ q(ck |c0:k−1 , M0:k , F0:k ),(i)p(Mk |Fk , ck ) = N (Mk ; ck × Fk , ΣM ),確率モデルの変形(20)(i)
上記の方法では、sk が高次元なため、x0:k を

p(Fk ) = const,

直接サンプリングすると非常に多数の粒子が必

(i)

(i)

(i)

(21)
(i)

(22)

(i)

(23)

p(ck |ck−1 ) = N (ck ; ck−1 , Σc ),要になる。そこで本研究では、以下のように提

(i)

(i)

p(ck |sk ) = Grids(i) (ck ).

案分布を変形する。

kであるとする。ここで N (x; µ, Σ) は平均 µ、分q(xk |x0:k−1 , y0:k )2 I と Σc =散 Σ の正規分布を表す。ΣM = σM

≜ q(ck , sk |c0:k−1 , s0:k−1 , M0:k , F0:k )= q(ck |c0:k−1 , M0:k , F0:k )q(sk |sk−1 , ck−1 , Fk−1 )(13)

σc2 I はモーメントの観測値と接触位置との分散を表す。また、Grids (c) は、s で表現される格子のうち位置 c を含むセルを表す。

加えて、形状パラメータ sk に関する確率分布に
以下のような仮定をおく。

q(sk |sk−1 , ck−1 , Fk−1 )1 if s = ψ(s , c , F ),k
k−1 k−1k−1=
0 otherwise

3.3

形状パラメータの推定

形状パラメータは、式 (14) の写像 ψ によって計算される。本研究では、二つの指標を使って(14)

ψ を表現する。概念図を図 2 に示す。

–4–

力の測定値 F
セル s(ik,j)

推定接触位置 c(i)Table 1

提案手法のパラメータ項目
粒子数
セルの大きさセルの数
接触位置の標準偏差
モーメントの標準偏差
距離の閾値
角度の閾値
形状パラメータの増加量
形状パラメータの減少量

Pos(s(ik,)j )dth
指標1: ς(i)を+Δs1
指標2: ς(i)を−Δs2

Fig. 2 形状推定の概念図。各粒子は接触位置の推定値と測定値とを基に、2 つの指標によって各セルでの形状表面の存在確率を増減させる。

記号
N
dim s
σc
σM
dth
θth
∆s1
∆s2

値
30
0.5 cm80 × 80
10−4 cm10−3 Nm
2 cm
π
12 rad10−2
5 × 10−3

タの減少量であり、θth は角度の閾値である。この指標は、接触位置の候補直線からの角度 θth(i)

一つ目の指標は、接触位置 ck−1 が道具表面未満であるセル sk,j の値を減少させる。

上にのみ存在することを表すものである。これは、次式によって表現される。

1(i)sk,j =
(i)

(i)

,

(24)

4.

シミュレーション検証

4.1

設定

1 + exp(−ςk,j )提案手法を、平面上の接触を模擬したシミュ(i)

ςk,j =ςk−1,j + ∆s1 ,
(i)

(i)

if ∥Pos(sk,j ) − ck−1 ∥ < dthレーションで検証した。本検証で利用した提案(25)

手法のパラメータを表 1 に記す。

ここで ςk,j は形状パラメータの中間表現であり、道具の形状は平面上において直線状を取ると

j はセルを特定するインデックスである。∆s1 はする。この設定は、例えば刃物においては通常

形状パラメータの増加量であり、dth は距離の閾

は刃先の曲線状にのみ接触が生じる状況を模擬値である。Pos(x) はセル x の中心位置である。

したものである。形状は力覚センサに固定され

(i)

(i)

この指標は、ck−1 からの距離が dth 未満である(i)

セル sk,j の値を増加させる。
二つ目の指標は、接触位置の候補直線（式 (3)）
と道具表面との交差点は 1 点のみとするものである。この仮定は必ずしも正しくないが、広い範囲の形状パラメータを変更するため、推定の高速化が期待できる。この指標は、次式によった座標系において次式で表される：

(x, y) 0.1 m <=x<= 0.3 m, y = f (x) m(28)
ここで f は形状を指定する関数である。
接触位置 c = (cx , cy ) は、2 秒おきに次式に従って変化するとする。

て表現される。
(i)sk,j =

1

,

(26)

(i)
1 + exp(−ςk,j )(i)(i)
ςk,j =ςk−1,j − ∆s2 ,(i)(i)
if Pos(sk,j ) ∈ cone(ck−1 , Fk−1 , θth )(i)(i)
and ∥Pos(sk,j ) − ck−1 ∥ >
(27)
= dth

cx ∼ U (0.1 m, 0.3 m),(29)

cy = f (x) m(30)

ここで U (a, b) は a <= x < b の値域をとる一様乱数である。接触力 f = (fx , fy ) は、75 s を境にして値の取り方が変化するとした。75 s まで
は、次式のように接触力の方向が時間的に変動

ここで、cone(a, v, θ) は頂点を a、軸を v 、頂角を θ とする双円錐である。∆s2 は形状パラメー–5–ずれの形状においても、提案手法は形状既知のする。

場合にこそ劣るものの、比較手法よりも顕著に

fx = A sin θf (t),(31)

fy = A cos θf (t),πθf (t) = sin(4πt) + θ⊥
6(32)

形状推定をおこなったためであり、形状推定を(33)

おこなうことの利点を示している。ただし、図t は時刻であり、θ⊥ は ck において道具表面と垂直な方向の角度を表す。また、振幅 A は 2 秒おきに（すなわち、接触位置が変化するごとに）次式によってランダムに変化する。

A ∼ U (1.0 N, 3.0 N)

正確な推定ができた。これは、75 s までの間に3 に示されたように、推定された形状はわずかであるが真値からのずれがあるため、形状既知の手法と比べると誤差が大きい。
提案手法と比較手法とでの接触位置推定の時系列を図 5 および 6 に示す。提案手法では接触(34)

力の変動がなくなった以降も、75 s までと同程度に推定ができている。一方で、比較手法にお75 s 以降は、接触力の変動がなくなり、2 秒おいては、75 s 以降に接触位置の誤差が急増してきに次式のように定まるとした。

いる。実際の作業中においては接触力が常に変4.2

fx = A sin θf ,(35)

fy = A cos θf ,
 π π
θf ∼ U − ,+ θ⊥ .
6 6(36)

動するとは限らないため、そのような場合は形状推定をおこなう提案手法のほうが有効であると期待される。

(37)

5.

結果

形状推定の様子を図 3 に示す。ここで、赤色の曲線は形状の真値を表す。いずれの形状においても、形状パラメータが徐々に真値に収束している様子が観察できる。また、形状に応じて異なる形状パラメータへと収束していることから、提案手法によって形状推定ができたといえる。わずか 30 個の粒子で 6400 個ものパラメー
タで表現される形状が推定できたことは、提案手法の有効性を端的に示すといえる。
接触位置の推定結果を、既存手法と比較した。

結論

本研究では、確率モデルの工夫によって、高次元パラメータで表現される道具形状と道具上の接触位置とを同時に推定する技術を提案した。
提案手法によって、わずか 30 個の粒子で 6400
個もの形状パラメータが推定できることが示された。また、形状推定を経ることによって、形状情報をまったく使わない場合よりも正確な接触位置推定ができることも示された。今後は、3 次元空間への拡張やより複雑な形状への適用、ロボット制御への応用などの展開を想定している。

比較手法として、逐次最小二乗法を用いて形状情報を一切利用することなく接触位置を推定できる手法 3) を利用した。ただし、逐次最小二乗

参考文献

法の忘却係数は α = 0.9 とした。また、参考のため、形状既知における推定手法として式 (4)を用いた手法 4) を用いた。結果を図 4 に示す。

1) A. Bicchi and D. Pavllo, “Intrinsic tactile
sensing for artiﬁcial hands,” in Robotics
Research - IV, 1988, pp. 83–90.

ここでは、接触力の変動がなくなる 75 s 以降の2) A. Bicchi, “Intrinsic contact sensing for

接触位置の推定誤差について、各形状、各手法

soft ﬁngers,” in Proceedings., IEEE In-

で 5 回の試行を実行した結果を示している。いternational Conference on Robotics and–6–t=5 st = 10 s

t = 20 s

t = 30 s

t = 40 s

t = 50 s

t = 75 s1

t = 100 s

形状の存在確率

t=0 s

y [m]

0.2

−0.2
0.0

x [m]

0.4

0

接触位置の誤差 [m]

7.5
5.0

25

50

75

100

0

25

50

75

100

0.2
0.1
0.0

Fx

力 [N]

比較手法

提案手法

形状既知

比較手法

提案手法

形状既知

比較手法

提案手法

形状既知

比較手法

提案手法

形状既知

比較手法

5
提案手法

0.0

0
0.3

2.5
形状既知

接触位置の推定誤差 [cm]

10.0

Fy

M

1.0
0.5

0
0

Fig. 4 75 s 以降の接触位置推定誤差の比較。
1
提案手法はいずれの形状においても、形状情報を一切使わない比較手法よりも顕著に正確な推定ができた。

cx の推定値
cy の推定値

0

15.0
12.5

cx の真値cy の真値

1

25

50

75

0.0
100

モーメント [Nm]形状

接触位置 [m]Fig. 3 さまざまな形状における形状推定の様子。形状パラメータが徐々に真値（赤色で表示）に
1
収束していった。

時刻 [s]

Fig. 5 提案手法での接触位置の推定結果。
75 s
1で接触力の変動がなくなって以降も、変わらず
ほぼ正確に推定できている。

–7–

接触位置 [m]cx の真値cy の真値1tion, vol. 3, no. 1, pp. 43–53, feb 1987.

07) R. Koike, S. Sakaino, and T. Tsuji, “Con0

25

50

75

100

tact point estimation in tactile interface

0.3

using particle ﬁlter,” in IECON 2017 -

0.2

43rd Annual Conference of the IEEE In-

0.1
0.0

dustrial Electronics Society, vol. 24, no. 1.
0
Fx

5
力 [N]

25

50

75

Fy

100
M

1.0
0.5

0
0

25

50

75

0.0
100

IEEE, oct 2017, pp. 8267–8272.

モーメント [Nm]

接触位置の誤差 [m]IEEE Journal on Robotics and Automa-

cx の推定値
cy の推定値

8) K. Kutsuzawa, S. Sakaino, and T. Tsuji,“Simultaneous estimation of contact position and tool shape using an unscented

時刻 [s]

particle ﬁlter,” IEEJ Journal of Industry
Fig. 6 比較手法での接触位置の推定結果。
75 s
1で接触力の変動がなくなって以降は、推定誤差が急増している。

Applications, vol. 9, no. 5, pp. 505–514,
2020.
9) K. P. Murphy, “Bayesian map learning in

Automation.

IEEE Comput. Soc. Press,

1990, pp. 968–973.

dynamic environments,” in Advances in
Neural Information Processing Systems,3) T. Tsuji, T. Seki, and S. Sakaino, “In-

1999, pp. 1015–1021.

trinsic contact sensing for touch interface10) G. Grisetti, C. Stachniss, and W. Bur-

with movable structure,” vol. 64, no. 9,

gard, “Improving grid-based slam with

pp. 7342–7349, 2017.

rao-blackwellized particle ﬁlters by adap-

4) J. Salisbury, “Interpretation of contact
geometries from force measurements,”
in Proceedings. 1984 IEEE International
Conference on Robotics and Automation,
vol. 1.

Institute of Electrical and Elec-tronics Engineers, 1984, pp. 240–247.

tive proposals and selective resampling,”
in Proceedings of the 2005 IEEE International Conference on Robotics and Automation, vol. 2005, no. April.

IEEE,

2005, pp. 2432–2437.
11) ——, “Improved techniques for grid map-

5) M. H. Raibert and J. J. Craig, “Hybrid
position/force control of manipulators,”
Journal of Dynamic Systems, Measurement, and Control, vol. 103, no. 2, pp.
126–133, jun 1981.

ping with rao-blackwellized particle ﬁlters,” IEEE Transactions on Robotics,
vol. 23, no. 1, pp. 34–46, 2007.
12) S. Muto and K. Shimokura, “Contact
point detection using force and velocity6) O. Khatib, “A uniﬁed approach for motion and force control of robot manipulators: The operational space formulation,”–8–

information complementarily.” Journal of
the Robotics Society of Japan, vol. 11,
no. 3, pp. 429–435, 1993.

13) H. Kitamura, S. Sakaino, and T. Tsuji,19) N. Kurita, S. Sakaino, and T. Tsuji,“Contact point calculation on a haptic“Whole-body force sensation by force sen-interface utilizing diﬀerentiated force,”

sor with end-eﬀector of arbitrary shape,”

IEEJ Journal of Industry Applications,in 2012 IEEE/RSJ International Confer-vol. 6, no. 2, pp. 151–159, 2017.

ence on Intelligent Robots and Systems.

14) R.

Featherstone,

S.

Thiebaut,andIEEE, 10 2012, pp. 5428–5433.

O. Khatib, “A general contact model for20) R. van der Merwe, A. Doucet, N. de Fre-

dynamically-decoupled force/motion con-itas, and E. Wan, “The unscented particle

trol,” in Proceedings 1999 IEEE Inter-

ﬁlter,” 2000.

national Conference on Robotics and Automation (Cat. No.99CH36288C). IEEE,
1999, pp. 3281–3286.

21) ——, “The unscented particle ﬁlter,” inAdvances in Neural Information Processing Systems, vol. 13.

15) K. Kutsuzawa, S. Sakaino, and T. Tsuji,“Estimation of individual contact force
when two contact points exist during
robotic tool use,” in The Abstracts of
the international conference on advanced
mechatronics : toward evolutionary fusion of IT and mechatronics : ICAM, vol.
2015.6, no. 0. Japan Society of Mechanical Engineers, 2015, pp. 46–47.
16) ——, “Estimation of individual force atthree contact points on an end-eﬀector
by a six-axis force/torque sensor,” in
IECON 2016 - 42nd Annual Conferenceof the IEEE Industrial Electronics Society. IEEE, 10 2016, pp. 6409–6414.
17) J. Park and O. Khatib, “Robot multiple
contact control,” Robotica, vol. 26, no. 5,
pp. 667–677, 9 2008.
18) K. Kutsuzawa, S. Sakaino, and T. Tsuji,“A control system for a tool use robot:
Drawing a circle by educing functionsof a compass,” Journal of Robotics andMechatronics, vol. 29, no. 2, pp. 395–405,
4 2017.

–9–

MIT Press, 2000.