(57)【要約】【課題】体内の各部位の干渉回避演算を従来とは異なる方法で実現する。【解決手段】運動解析装置は、体内の各部位を表わす第１〜第３の三次元モデルを格納する記憶部を備える。三次元モデルの各々は、複数の質点で表わされている。複数の質点の各々は、三次元上の座標値と質量とに対応付けられており、かつ、三次元モデル内の関連付けられている他の質点との距離に依存する力学的関係よって当該他の質点と結合している。第３の部位モデルの複数の質点は、第１，第２の部位モデルを構成する複数の質点のいずれかと関連付けられている。運動解析装置の演算部は、各質点について、当該質点と隣接する質点との間の距離に応じて、当該質点が当該隣接する質点から受ける受動的な力を計算し、各質点の座標値を順次計算する。【選択図】図５
2019192053.tif 000002
【特許請求の範囲】
【請求項１】
動物の運動を解析するための運動解析装置であって、
前記動物の体内の各部位を表わす複数の三次元モデルを格納するための記憶部を備え、
前記複数の三次元モデルの各々は、複数の質点で表わされ、
前記複数の質点の各々は、三次元上の座標値と質量とに対応付けられており、かつ、三次元モデル内の関連付けられている他の質点との距離に依存する力学的関係よって当該他の質点と結合しており、
前記複数の三次元モデルは、前記動物の体内の第１〜第３の部位を表わす第１〜第３の部位モデルを含み、
前記第３の部位モデルを構成する複数の質点の一部または全部は、前記第１の部位モデルを構成する複数の質点と、前記第２の部位モデルを構成する複数の質点との内のいずれかの質点に関連付けられており、
前記第１〜第３の部位モデルの複数の質点の各々についての座標値を計算するための演算部を備え、
前記演算部は、
前記第１〜第３の部位モデルの複数の質点の各々について、当該質点と関連付けられている他の質点との間の距離に応じて、当該質点が当該他の質点から受ける力を計算し、
前記第１〜第３の部位モデルの複数の質点の各々について、当該質点に対応付けられている質量と、当該質点について計算された前記力とに依存する運動方程式によって当該質点の座標値を順次計算する、運動解析装置。
【請求項２】
前記第３の部位モデルは、前記第３の部位の第１の表面を構成する複数の質点と、前記第３の部位の第２の表面を構成する複数の質点とで構成されている、請求項１に記載の運動解析装置。
【請求項３】
前記第１の表面を構成する複数の質点の一部または全部は、前記第１の部位モデルを構成する少なくとも１つの質点と関連付けられており、
前記第２の表面を構成する複数の質点の一部または全部は、前記第２の部位モデルを構成する少なくとも１つの質点と関連付けられている、請求項２に記載の運動解析装置。
【請求項４】
前記第３の部位モデルは、前記動物の皮膚を表わし、
前記第１の表面を構成する複数の質点の一部または全部は、前記第１の部位モデルを構成する複数の質点と、前記第２の部位モデルを構成する複数の質点との内のいずれかの質点に関連付けられている、請求項２に記載の運動解析装置。
【請求項５】
前記第１の部位モデルは、前記動物の体内の第１の筋を表わし、
前記第２の部位モデルは、前記動物の体内の第２の筋を表を表わす、請求項１〜４のいずれか１項に記載の運動解析装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本開示は、動物の運動を解析するための運動解析装置に関する。
【背景技術】
【０００２】
近年、人や動物等の筋の運動を解析するための技術が開発されている。当該運動は、人の筋の力学的機能をコンピュータ上で表現した筋モデルを用いて解析される。筋の運動解析技術に関し、特開２０１７−０３７５５３号公報（特許文献１）は、バネマスダンパモデル方式を採用した人体の筋モデルを開示している。
【０００３】
より複雑な運動解析を実現するためには、体内の部位間の干渉を考慮する必要がある。三次元オブジェクト同士の干渉を回避するための技術に関し、特開２００６−３０２０９８号公報（特許文献２）がある。特許文献２は、「仮想空間に配置されたオブジェクト同士の干渉を簡易な演算により演算可能な三次元オブジェクト制御装置」を開示している。特許文献２に開示される三次元オブジェクト制御装置は、体内の部位間の干渉回避を前提としていない。
【先行技術文献】
【特許文献】
【０００４】
【特許文献１】特開２０１７−０３７５５３号公報
【特許文献２】特開２００６−３０２０９８号公報
【発明の概要】
【発明が解決しようとする課題】
【０００５】
体内の各部位は、他の部位と常に接触しており、より複雑な運動解析を実現するためには、各部位の干渉を常に考慮する必要がある。体内の各部位について干渉判定処理および干渉回避処理を行うと、処理に膨大な時間を要する。そのため、これらの処理を常に行うことは現実的ではない。そのため、体内の各部位の干渉回避を従来とは異なる手法で実現する必要がある。
【０００６】
本開示は上述のような問題点を解決するためになされたものであって、ある局面における目的は、体内の各部位の干渉回避演算を従来とは異なる方法で実現することが可能な運動解析装置を提供することである。
【課題を解決するための手段】
【０００７】
ある局面に従うと、動物の運動を解析するための運動解析装置は、上記動物の体内の各部位を表わす複数の三次元モデルを格納するための記憶部を備える。上記複数の三次元モデルの各々は、複数の質点で表わされる。上記複数の質点の各々は、三次元上の座標値と質量とに対応付けられており、かつ、三次元モデル内の関連付けられている他の質点との距離に依存する力学的関係よって当該他の質点と結合している。上記複数の三次元モデルは、上記動物の体内の第１〜第３の部位を表わす第１〜第３の部位モデルを含む。上記第３の部位モデルを構成する複数の質点の一部または全部は、上記第１の部位モデルを構成する複数の質点と、上記第２の部位モデルを構成する複数の質点との内のいずれかの質点に関連付けられている。上記第１〜第３の部位モデルの複数の質点の各々についての座標値を計算するための演算部を備える。上記演算部は、上記第１〜第３の部位モデルの複数の質点の各々について、当該質点と関連付けられている他の質点との間の距離に応じて、当該質点が当該他の質点から受ける力を計算し、上記第１〜第３の部位モデルの複数の質点の各々について、当該質点に対応付けられている質量と、当該質点について計算された上記力とに依存する運動方程式によって当該質点の座標値を順次計算する。
【０００８】
好ましくは、上記第３の部位モデルは、上記第３の部位の第１の表面を構成する複数の質点と、上記第３の部位の第２の表面を構成する複数の質点とで構成されている。
【０００９】
好ましくは、上記第１の表面を構成する複数の質点の一部または全部は、上記第１の部位モデルを構成する少なくとも１つの質点と関連付けられている。上記第２の表面を構成する複数の質点の一部または全部は、上記第２の部位モデルを構成する少なくとも１つの質点と関連付けられている。
【００１０】
好ましくは、上記第３の部位モデルは、上記動物の皮膚を表わす。上記第１の表面を構成する複数の質点の一部または全部は、上記第１の部位モデルを構成する複数の質点と、上記第２の部位モデルを構成する複数の質点との内のいずれかの質点に関連付けられている。
【００１１】
好ましくは、上記第１の部位モデルは、上記動物の体内の第１の筋を表わす。上記第２の部位モデルは、上記動物の体内の第２の筋を表を表わす。
【発明の効果】
【００１２】
ある局面において、体内の各部位の干渉回避を従来とは異なる方法で実現することができる。
【００１３】
本発明の上記および他の目的、特徴、局面および利点は、添付の図面と関連して理解される本発明に関する次の詳細な説明から明らかとなるであろう。
【図面の簡単な説明】
【００１４】
【図１】筋モデルの一例を概略的に示す図である。
【図２】筋モデルの各質点に関する質点情報１１２を示す図である。
【図３】運動解析装置の機能構成の一例を示すブロック図である。
【図４】各質点の座標計算時におけるＧＰＵ（Graphics Processing Unit）へのスレッドの配置を示す概念図である。
【図５】筋モデル同士の干渉回避処理を説明するための図である。
【図６】結合情報のデータ構造の一例を示す図である。
【図７】筋モデルと皮膚モデルとの干渉回避処理を説明するための図である。
【図８】筋骨格モデルを概略的に示す図である。
【図９】運動解析装置の主要なハードウェア構成を示すブロック図である。
【図１０】本提案手法のシミュレーション結果の一例を示す図である。
【発明を実施するための形態】
【００１５】
以下、図面を参照しつつ、本発明に従う各実施の形態について説明する。以下の説明では、同一の部品および構成要素には同一の符号を付してある。それらの名称および機能も同じである。したがって、これらについての詳細な説明は繰り返さない。また、以下で説明する変形例は、適宜選択的に組み合わされてもよい。
【００１６】
［１  筋モデルの座標系］
図１を参照して、実施の形態に従う筋モデル５０を表わすための座標系について説明する。図１は、筋モデル５０の一例を概略的に示す図である。
【００１７】
筋モデル５０は、複数の質点で表わされる。質点とは、筋を構成する各部分を表わす要素のことをいう。各質点には、三次元上の座標値と質量とが関連付けられている。
【００１８】
図１に示されるように、筋モデル５０を構成する各質点は、三次元の格子座標系（ｉ、ｊ、ｋ）で表わされる。当該座標系の「ｋ」は、筋の長手方向（すなわち、筋の始端から終端に向かう方向）に直交する断面の番号を表わす。筋モデル５０の各断面に含まれる各質点の番号は、２次元格子座標（ｉ、ｊ）で管理される。
【００１９】
ｉ方向の質点数は、Ｎｉ個とする。ｊ方向の質点数は、Ｎｊ個とする。ｋ方向の質点数は、Ｎｋ個とする。すなわち、1つの筋は、６面体のポリゴンを直列に（Ｎｋ−１）個並べてできる細長いポリゴンを、並列に（Ｎｉ−１）＊（Ｎｊ−１）個束ねることによって構成される。（ｉ、ｊ、ｋ）番目の質点の通し番号ｎは、「ｎ＝ｉ＋（Ｎｉ＊ｊ）＋（Ｎｉ＊Ｎｊ＊ｋ）」で表わされる。たとえば、Ｎｉ＝５、Ｎｊ＝７である場合、（ｉ、ｊ、ｋ）＝（０、２、１）における質点の通し番号は、４５（＝０＋（５＊２）＋（５＊７＊１））となる。Ｎｉ、Ｎｊ、Ｎｋの値が調整されることで、筋モデル５０は、細長い筋や平たい筋等の様々な形状で表現され得る。
【００２０】
なお、この形状表現は、筋だけでなく、骨と骨をつないでいる靭帯や、多頭筋を束ねている腱膜の表現にも用いられ得る。これにより、人体の各部位（たとえば、筋、腱、靭帯、腱膜等）の形状が一元的に管理される。また、このような形状表現により、後述するように同一ＧＰＵのカーネル関数で全ての人体部位の形状が計算され得る。
【００２１】
［２  筋モデルの運動解析］
（２．１  概要）
筋モデル５０（図１参照）には、隣り合う質点を力学的に結合させたバネマスダンパモデルが採用される。バネマスダンパモデルの特徴は、各質点の運動方程式が同じ計算式で表わされ、かつ、各質点のシミュレーションの１タイムステップ（Δｔ）の積分計算が他の質点とは独立に行うことができる点である。これにより、大量の質点の１タイムステップの積分計算が並列に処理され得る。時間積分の方法としては、たとえば、ベレの方法が用いられる。
【００２２】
（２．２  運動方程式）
図２および図３を参照して、運動解析装置１００による筋モデル５０（図１参照）の運動解析方法について説明する。図２は、筋モデル５０の各質点に関する質点情報１１２を示す図である。図３は、運動解析装置１００の機能構成の一例を示すブロック図である。
【００２３】
以下の式（１）に示されるように、ｎ番目の質点の運動方程式は、隣接する質点から受ける受動的な力「ｆPassive」と、能動的な筋収縮力「ｆActive」と、隣接する筋または骨から受ける外力「ｆCol」によって決定される。
【００２４】
【数１】
2019192053.tif 000003
【００２５】
式（１）の「ｍn」は、ｎ番目の質点の質量を表わす。各質点の質量「ｍn」は、たとえば、図２の質点情報１１２に規定されている。一例として、各質点の質量「ｍn」は、人体の構造を規定している人体モデルから予め得られている。式（１）の「ａn」は、ｎ番目の質点の加速度ベクトルを表わす。
【００２６】
式（１）は、たとえば、運動解析装置１００に含まれる演算部１５０によって計算される。演算部１５０は、後述するＣＰＵ（Central  Processing  Unit）１０２やＧＰＵ１０４に相当する。演算部１５０は、筋モデル５０を構成する複数の質点のうち、筋の表面を構成する質点が外力を受けたことに基づいて、各質点についての座標値を順次計算する。演算部１５０は、外力「ｆCol」を受けた質点については、当該質点に関連付けられている質量と、当該質点についての受動的な力「ｆPassive」と、外力「ｆCol」とに依存する運動方程式によって当該質点の座標値を順次計算する。演算部１５０は、外力「ｆCol」を受けていない質点については、当該質点に関連付けられている質量と、当該質点について計算された力「ｆPassive」とに依存する運動方程式によって当該質点の座標値を順次計算する。すなわち、演算部１５０は、外力を受けていない質点については、外力「ｆCol」をゼロに設定した上で式（１）の運動方程式を解く。なお、より詳細に各質点の座標値を計算する場合には、式（１）に示されるように能動的な筋収縮力「ｆActive」が運動方程式に追加される。
【００２７】
図３に示されるように、演算部１５０は、機能構成として、第１算出部１５１と、第２算出部１５２と、第３算出部１５３と、位置算出部１６０とを含む。第１算出部１５１は、式（１）の受動的な力「ｆPassive」を算出する。第２算出部１５２は、式（１）に示される筋収縮力「ｆActive」を算出する。第３算出部１５３は、式（１）に示される外力「ｆCol」を算出する。位置算出部１６０は、上記式（１）の運動方程式を解くことにより、筋モデル５０の各質点の座標値を算出する。
【００２８】
（ａ）第１算出部１５１
以下、隣接する質点から受ける力「ｆPassive」の算出方法について説明する。力「ｆPassive」は、第１算出部１５１によって算出される。第１算出部１５１は、対象の質点と隣接する質点との力学的な結合関係をバネに見立てて「ｆPassive」を計算する。すなわち、筋モデル５０の各質点は、隣接する質点との距離に依存する力学的関係よって隣接する質点と結合しており、第１算出部１５１は、対象の質点と隣接する質点との間の距離に応じて、対象の質点が隣接する質点から受ける力を算出する。
【００２９】
各質点には、１つ以上の質点が隣接する。安定的にシミュレーションを行うためには、質点同士の力学的結合関係（すなわち、バネのつなぎ方）に工夫が必要である。ｉ、ｊ、ｋ方向（図１参照）にバネが配置されただけでは、筋モデル５０は、形状を保てないことがあるためである。そのため、ｉ、ｊ、ｋ方向の全ての断面内において斜めのバネが配置されることが好ましい。これにより、運動解析装置１００は、筋モデル５０について安定して形状シミュレーションを行うことができる。たとえば、受動的な力「ｆPassive」は、以下の式（２）で表わされる。
【００３０】
【数２】
2019192053.tif 000004
【００３１】
式（２）の上付き文字「ｅ」は、東方向を表わす（図１参照）。上付き文字「ｗ」は、西方向を表わす。上付き文字「ｓ」は、南方向を表わす。上付き文字「ｎ」は、北方向を表わす。上付き文字「ｆ」は、前方向（図１の紙面下方向）を表わす。上付き文字「ｂ」は、後ろ方向（図１の紙面上方向）を表わす。たとえば、式（２）の力「ｆne」は、東北方向の質点（バネ）から受ける力を表わす。隣接する質点が存在しない場合には、その方向から受ける力はゼロとなる。たとえば、筋の始端に当たるｋ＝０の断面内の質点は、後方向（すなわち、骨側）からは力を受けないものとする。実際には、始端の質点は骨から力を受けることで骨に対して一定の位置を保っているが、このような拘束条件の実装方法については後述する。
【００３２】
以下の式（３）に示されるように、ｎ番目の質点が隣接するｍ番目の質点から受ける力「ｆnm」は、弾性力「ｆnmS」と粘性力「ｆnmD」との和で表わされる。力「ｆnm」は、式（２）の右辺の各項に相当する。
【００３３】
【数３】
2019192053.tif 000005
【００３４】
式（３）の弾性力「ｆnmS」は、以下の式（４）で表わされる。
【００３５】
【数４】
2019192053.tif 000006
【００３６】
式（４）の「ｋS」は、バネ定数である。「ｘnm」は、「ｘm−ｘn」を意味する。「ｘm」は、ｍ番目の質点の座標を表わす。「ｘn」は、ｎ番目の質点の座標を表わす。各質点の座標は、たとえば、図２の質点情報１１２に規定されている。「｜ｘnm｜」は、ｎ番目の質点とｍ番目の質点との間の距離を表わす。「Ｌnm」は、力を受けていないときにおけるｎ番目の質点とｍ番目の質点との間の距離（すなわち、自然長）を表わす。各質点間の自然長は、たとえば、図２の質点情報１１２に予め規定されている。
【００３７】
式（３）の粘性力「ｆnmD」は、以下の式（５）で表わされる。
【００３８】
【数５】
2019192053.tif 000007
【００３９】
式（５）の「ｋd」は、粘性抵抗係数を表わす。「ｖn」は、ｎ番目の質点の速度を表わす。「ｖm」は、ｍ番目の質点の速度を表わす。各質点の速度は、たとえば、当該質点の前回の位置から現在位置までの移動距離を移動時間で割ることにより算出される。
【００４０】
（ｂ）第２算出部１５２
上記式（１）に示される筋収縮力「ｆActive」の算出方法について説明する。筋収縮力「ｆActive」は、筋の能動的な活動によって生じる力である。筋収縮力「ｆActive」は、第２算出部１５２によって算出される。筋収縮力「ｆActive」は、ｋ方向にのみ働くものとし、その大きさは、筋の活動度、張力と長さとの関係「Ｆl（｜ｘnm｜）」、張力と速度との関係「Ｆv（｜ｖnm｜）」、および最大等尺性収縮筋力「Ｆmax」により決まる。より具体的には、筋収縮力「ｆActive」は、以下の式（６）で表わされる。
【００４１】
【数６】
2019192053.tif 000008
【００４２】
（ｃ）第３算出部１５３
上述したように、第３算出部１５３は、筋モデル５０に与える外力「ｆCol」を算出する。外力「ｆCol」は、たとえば、筋モデル５０が他の物体と衝突することによる反発力や干渉力として生じる。筋モデル５０に外力を与える物体は、たとえば、他の筋モデルおよび骨モデルの少なくとも一方を含む。外力「ｆCol」は、後述する計算方法で算出されてもよいし、運動解析装置１００の操作者によって任意の値で与えられてもよい。
【００４３】
外力「ｆCol」は、筋モデル５０の表面を構成する質点に与えられる。すなわち、外力「ｆCol」は、筋モデル５０の内部を構成する質点には与えられない。
【００４４】
（ｄ）位置算出部１６０
位置算出部１６０は、上記式（１）に示される運動方程式を各質点について解くことにより、各質点の加速度を算出する。位置算出部１６０は、各質点の初期位置や各質点の初期速度（典型的には、ゼロ）を反映した上で、各質点の加速度を時間で二回積分することにより各質点の座標を算出する。時間積分の方法としては、たとえば、ベレの方法が用いられる。
【００４５】
位置算出部１６０は、各質点の新たな座標値を算出したことに基づいて、図２の質点情報１１２における現在の位置を各質点の新たな座標値で書き換えるとともに、書き換え前の位置を質点情報１１２の履歴に格納する。
【００４６】
（２．３  ＧＰＵへの実装方法）
（ａ）概要
上記式（１）〜（６）に示されるように、各質点の新たな座標値は、現在の各質点の座標値および前回の各質点の座標値から算出される。すなわち、現在の各質点の座標値および前回の各質点の座標値が既知であれば、各質点の新たな座標値の計算は、他の質点の新たな座標値の計算結果に影響されない。この点に着目して、運動解析装置１００は、筋モデル５０の各質点についての運動方程式（上記式（１）参照）の計算を、一部または全部について並列に行う。これにより、運動解析装置１００は、筋モデル５０の運動を解析するのに要する処理時間を大幅に削減することができる。
【００４７】
（ｂ）スレッドの管理方法
運動解析装置１００は、たとえば、各質点の運動方程式の並列処理をＧＰＵで実現する。ＧＰＵは、複数のコア（たとえば、数百個のコア）を有する。ＧＰＵは、筋モデル５０の運動解析を実現するためのプログラムを多数のスレッドに分割し、各コアに対応するスレッドを割り当てることで、筋モデル５０の運動解析を並列処理する。
【００４８】
図４を参照して、ＧＰＵによるスレッドの管理方法について説明する。図４は、各質点の座標計算時におけるＧＰＵへのスレッドの配置を示す概念図である。
【００４９】
図４に示されるように、スレッドは、ＧＰＵ内において、「グリッド」および「ブロック」の２つの階層構造で管理される。各グリッド内には、複数のブロックが２次元で配置される。図４（Ａ）に示されるように、１つのブロックには、1つの筋が割り当てられる。好ましくは、1列目のブロックには左半身の各筋が割り当てられ、２列目のブロックには右半身の各筋が割り当てられる。
【００５０】
各ブロックには、複数のスレッドが２次元で配置される。図４（Ｂ）に示されるように、運動解析装置１００は、1つのスレッドに対しては、1つの筋繊維を割り当てる。人のほとんどの筋の形状は、Ｎｉ＝７、Ｎｊ＝１６、Ｎｋ＝３３の範囲内で表現され得るため、i方向の最大スレッド数が７に設定され、j方向の最大スレッド数が１６に設定される。
【００５１】
（ｃ）シェアードメモリの有効的な利用方法
ＧＰＵ並列処理を高速化するためには、アクセス速度の速いメモリを有効に利用する必要がある。ＤＲＡＭ（Dynamic Random Access Memory）内のグローバルメモリは、全てのブロックからがアクセス可能であるが、アクセス速度が遅いという欠点がある。一方で、ＳＭ（Streaming Multiprocessor）内にあるシェアードメモリは、同一ブロック内からしかアクセスすることができないが、アクセス速度が速いという利点がある。上記式（１）の運動方程式を計算する際に参照する機会の多い同一筋モデル内の情報（たとえば、位置や自然長等）がシェアードメモリに配置されることで処理が高速化されるため、運動解析装置１００は、１つのブロックに１つの筋が割り当てられるように設計される。
【００５２】
しかしながら、シェアードメモリの容量は、数十ＫＢと小さいため（たとえば、１６ＫＢ、３２ＫＢ、４８ＫＢの３種類）、1つの筋の情報を全てシェアードメモリに割り当てることは困難である。たとえば、筋のサイズがＮｉ＝７、Ｎｊ＝１６、Ｎｋ＝３３である場合であって、かつ、各質点がｆｌｏａｔ型（４バイト）の２４変数（内１８変数：自然長Ｌ、内３変数：現在の座標値、内３変数：過去の座標値）を有する場合には、総サイズが約３４６．５ＫＢ（＝３５４８１６Ｂ＝７＊１６＊３３＊２４＊４Ｂ）となりシェアードメモリの容量を超えてしまう。
【００５３】
このことを解決するため、シェアードメモリを節約する方法を以下に示す。ある質点の運動方程式を計算する際に必要な情報は、隣接する断面に含まれる情報があれば十分である。つまり、ｋ番目の断面に含まれる質点の運動方程式を計算する際には、ｋ−１番目、ｋ番目、ｋ＋１番目の３つの断面の情報があれば十分である。３つの断面に限れば、シェアードメモリに割り当てるデータ容量は、約３１．５ｋＢ（＝３２２５６Ｂ＝７＊１６＊３＊２４＊４Ｂ）に節約され得る。このような節約を実現するために、図４（Ｂ）に示されるように、運動解析装置１００は、１つのスレッドに１本の筋線維を割り当て、各断面の各質点の座標値をｋ＝０〜Ｎｋまで順に計算する。より具体的には、運動解析装置１００は、ｋ＝１の断面における質点の計算時には、ｋ＝０、ｋ＝１、ｋ＝２の３つの断面に含まれる質点の情報をシェアードメモリに配置し、その上で、ｋ＝１の断面における各質点の座標値を計算する。全てのスレッドが計算を終えた段階で、運動解析装置１００は、シェアードメモリをｋ＝１、ｋ＝２、ｋ＝３の断面に含まれる質点の情報に書き換え、次のｋ＝２の断面における質点の計算処理へと進む。
【００５４】
［３  筋モデルの干渉回避処理］
（３．１）概要
以下、筋モデル５０についての干渉計算について説明する。本実施の形態に従う運動解析装置１００は、バネマスダンパモデルに基づく上述の運動解析方法を、筋モデル５０と他の部位との干渉回避処理にも適用する。上述のように、バネマスダンパモデルの特徴は、各質点の運動方程式が同じ計算式で表わされ、かつ、各質点のシミュレーションの１タイムステップ（Δｔ）の積分計算が他の質点とは独立に行うことができる点である。これにより、大量の質点の１タイムステップの積分計算が並列に処理されるので、干渉回避処理についても同様に高速化することができる。
【００５５】
以下では、筋モデル５０についての干渉回避処理について説明を行うが、本明細書で説明する干渉回避処理は、筋以外の他の部位を表わす三次元モデルにも適用され得る。
【００５６】
（３．２）筋モデル同士の干渉回避処理
図５を参照して、干渉回避処理の一例について説明する。図５は、筋モデル同士の干渉回避処理を説明するための図である。
【００５７】
図５には、人体の第１の筋（第１の部位）を表わす筋モデル５０Ａ（第１の部位モデル）と、人体の第２の筋（第２の部位）を表わす筋モデル５０Ｂ（第２の部位モデル）とが示されている。図５では、説明の便宜のために、筋モデル５０Ａと筋モデル５０Ｂとが互いに離れている例が示されている、実際には、筋モデル５０Ａ，５０Ｂは、互いに接しているか、互いに近接している。
【００５８】
筋モデル５０Ａ，５０Ｂの干渉を回避するために、筋モデル５０Ａ，５０Ｂの間に緩衝オブジェクト６０が配置される。緩衝オブジェクト６０は、ユーザによって手動で配置されてもよいし、干渉する可能性のある部位間に自動で配置されてもよい。
【００５９】
緩衝オブジェクト６０は、人体の体内の部位に見立てることができる。本例では、緩衝オブジェクト６０は、人体の筋膜に見立てられる。
【００６０】
緩衝オブジェクト６０は、筋モデル５０Ａ，５０Ｂと同様に図２で説明した質点情報１１２によって規定される。すなわち、緩衝オブジェクト６０は、複数の質点で表現され、各質点には、質量、バネの自然長などが関連付けられる。バネの自然長は、初期姿勢における各質点間の距離とする。
【００６１】
緩衝オブジェクト６０のｊ方向の質点数を表わす「Ｎｊ」は、たとえば、２に設定される。つまり、緩衝オブジェクト６０は、第１の表面を構成する複数の質点（すなわち、ｊ＝０に相当する質点）と、第２の表面を構成する複数の質点（すなわち、ｊ＝１に相当する質点）とで構成される。緩衝オブジェクト６０を構成する各質点は、緩衝オブジェクト６０内の隣り合う質点とバネで結合している。
【００６２】
緩衝オブジェクト６０を構成する各質点は、緩衝オブジェクト６０内の質点だけでなく、緩衝オブジェクト６０に接する筋モデル５０Ａ，５０Ｂの各質点とも仮想的なバネで結合している。すなわち、緩衝オブジェクト６０を構成する一部または全部の質点は、筋モデル５０Ａ（第１の部位モデル）を構成する複数の質点と、筋モデル５０Ｂ（第２の部位モデル）を構成する複数の質点との内のいずれかの質点に関連付けられている。
【００６３】
図５の例では、緩衝オブジェクト６０の第１の表面を表わす各質点（すなわち、ｊ＝０に相当する質点）は、筋モデル５０Ａの質点に結合している。緩衝オブジェクト６０の第１の表面を構成する質点は、筋モデル５０Ａの筋膜として見なすことできる。一方で、緩衝オブジェクト６０の第２の表面を表わす各質点（すなわち、ｊ＝１に相当する質点）は、筋モデル５０Ｂの質点に結合している。緩衝オブジェクト６０の第２の表面を表わす質点は、筋モデル５０Ｂの筋膜として見なすことできる。
【００６４】
緩衝オブジェクト６０と結合する筋モデル５０Ａ，５０Ｂ内の質点は、たとえば、運動解析装置１００によって自動で登録される。より具体的には、運動解析装置１００は、緩衝オブジェクト６０内の第１の表面を構成する各質点（すなわち、ｊ＝０に相当する質点）について筋モデル５０Ａ内の全質点との距離を算出する。そして、運動解析装置１００は、緩衝オブジェクト６０内の第１の表面を構成する各質点を、筋モデル５０Ａを構成する各質点の内の最も距離が近い１つ以上の質点と結合する。同様に、運動解析装置１００は、緩衝オブジェクト６０内の第２の表面を構成する各質点（すなわち、ｊ＝１に相当する質点）について、筋モデル５０Ｂ内の全質点との距離を算出する。そして、運動解析装置１００は、緩衝オブジェクト６０内の第２の表面を構成する各質点を、筋モデル５０Ｂを構成する各質点の内の最も距離が近い１つ以上の質点と結合する。緩衝オブジェクト６０の各質点に結合される筋モデル５０Ａ，５０Ｂの質点の数は、たとえば、２点程度である。
【００６５】
筋モデル５０Ａ，５０Ｂの各質点と緩衝オブジェクト６０の各質点との結合結果は、図６に示される結合情報１１３として記憶される。図６は、結合情報１１３のデータ構造の一例を示す図である。
【００６６】
図６に示されるように、結合情報１１３には、緩衝オブジェクト６０の各質点について筋モデル５０Ａの質点が対応付けられている。一例として、緩衝オブジェクト６０の質点番号「ｎ１」の質点は、筋モデル５０Ａの質点番号「ｎ１０」の質点と対応付けられている。
【００６７】
なお、図６には、筋モデル５０Ａと緩衝オブジェクト６０との結合情報１１３が示されているが、実際には、他の筋モデル（たとえば、筋モデル５０Ｂ）についても、当該他の筋モデルと緩衝オブジェクト６０との結合情報１１３が存在する。
【００６８】
運動解析装置１００は、筋モデル５０Ａ，５０Ｂおよび緩衝オブジェクト６０を構成する各質点について、上述の「筋モデルの運動解析」で説明した計算方法で運動解析を行う。すなわち、運動解析装置１００は、筋モデル５０Ａ，５０Ｂおよび緩衝オブジェクト６０を構成する各質点について結合先の他の質点から受ける力を計算し、当該力と各質点に対応付けられている質量とに依存する運動方程式（すなわち、上記式（１））によって各質点の座標値を順次計算する。
【００６９】
より具体的には、運動解析装置１００は、筋モデル５０Ａ，５０Ｂを構成する各質点の内、緩衝オブジェクト６０と結合している各質点については、筋モデル５０Ａ，５０Ｂ内の結合先の各質点から受ける力と、緩衝オブジェクト６０内の結合先の各質点から受ける力との両方に応じて、筋モデル５０Ａ，５０Ｂ内の各質点の位置を計算する。一方で、運動解析装置１００は、筋モデル５０Ａ，５０Ｂを構成する各質点の内、緩衝オブジェクト６０と結合していない各質点については、筋モデル５０Ａ，５０Ｂ内の結合先の各質点から受ける力に応じて、筋モデル５０Ａ，５０Ｂ内の各質点の位置を計算する。
【００７０】
また、運動解析装置１００は、緩衝オブジェクト６０を構成する各質点については、緩衝オブジェクト６０内の結合先の各質点から受ける力と、筋モデル５０Ａ，５０Ｂ内の結合先の各質点から受ける力との両方に応じて、緩衝オブジェクト６０内の各質点の位置を計算する。
【００７１】
これにより、運動解析装置１００は、筋モデル５０Ａ，５０Ｂ内の各質点間の力のやり取りを計算できるだけでなく、緩衝オブジェクト６０を介した筋モデル５０Ａと筋モデル５０Ｂとの間の力のやり取りを計算することができ、筋モデル５０Ａと筋モデル５０Ｂとの干渉を回避することができる。
【００７２】
（３．３）筋モデルと皮膚モデルとの干渉回避処理
次に、図７を参照して、干渉回避処理の他の例について説明する。図７は、筋モデルと皮膚モデルとの干渉回避処理を説明するための図である。
【００７３】
図７には、人体の第１の筋（第１の部位）を表わす筋モデル５０Ａ（第１の部位モデル）と、人体の第２の筋（第２の部位）を表わす筋モデル５０Ｂ（第２の部位モデル）と、緩衝オブジェクト６０とが示されている。図７では、説明の便宜のために、筋モデル５０Ａ，５０Ｂと緩衝オブジェクト６０とが互いに離れている例が示されている、実際には、筋モデル５０Ａ，５０Ｂと緩衝オブジェクト６０とは、互いに接しているか、互いに近接している。
【００７４】
基本的には、図７に示される緩衝オブジェクト６０は、図５に示される緩衝オブジェクト６０と同じである。図５に示される緩衝オブジェクト６０においては、筋モデル５０Ａの質点が緩衝オブジェクト６０の第１の表面を構成する複数の質点（すなわち、ｊ＝０に相当する質点）と結合しており、筋モデル５０Ｂの質点が緩衝オブジェクト６０の第２の表面を構成する複数の質点（すなわち、ｊ＝１に相当する質点）と結合していた。これに対して、図７に示される緩衝オブジェクト６０においては、筋モデル５０Ａの質点と筋モデル５０Ｂの質点とが、共に、緩衝オブジェクト６０の第１の表面（すなわち、ｊ＝０に相当する質点）を構成する質点と結合している。このように配置された緩衝オブジェクト６０は、人体の皮膚（第３の部位）として見なすことできる。
【００７５】
緩衝オブジェクト６０と結合する筋モデル５０Ａ，５０Ｂ内の質点は、たとえば、運動解析装置１００によって自動で登録される。より具体的には、運動解析装置１００は、緩衝オブジェクト６０内の第１の表面を構成する各質点（すなわち、ｊ＝０に相当する質点）について筋モデル５０Ａ内および筋モデル５０Ｂ内の全質点との距離を算出する。そして、運動解析装置１００は、緩衝オブジェクト６０内の第１の表面を構成する各質点を、筋モデル５０Ａ，５０Ｂを構成する各質点の内の最も距離が近い１つ以上の質点と結合する。緩衝オブジェクト６０の各質点に結合される筋モデル５０Ａ，５０Ｂの質点の数は、たとえば、２点程度である。
【００７６】
なお、運動解析装置１００は、緩衝オブジェクト６０内の第１の表面を構成する各質点だけではなく、緩衝オブジェクト６０内の第２の表面を構成する各質点（すなわち、ｊ＝１に相当する質点）も、筋モデル５０Ａ，５０Ｂ内のいずれかの質点と結合させてもよい。この場合、緩衝オブジェクト６０内の全質点が筋モデル５０Ａ，５０Ｂ内の質点のいずれかと結合することになる。より具体的には、運動解析装置１００は、緩衝オブジェクト６０内の全質点について筋モデル５０Ａ，５０Ｂ内の各質点との距離を算出する。そして、運動解析装置１００は、緩衝オブジェクト６０内の全質点を、筋モデル５０Ａ，５０Ｂを構成する各質点の内の最も距離が近い１つ以上の質点と結合する。
【００７７】
運動解析装置１００は、筋モデル５０Ａ，５０Ｂおよび緩衝オブジェクト６０を構成する各質点について、上述の「筋モデルの運動解析」で説明した計算方法で運動解析を行う。すなわち、運動解析装置１００は、筋モデル５０Ａ，５０Ｂおよび緩衝オブジェクト６０を構成する各質点について結合先の他の質点から受ける力を計算し、当該力と各質点に対応付けられている質量とに依存する運動方程式（すなわち、上記式（１））によって各質点の座標値を順次計算する。
【００７８】
より具体的には、運動解析装置１００は、筋モデル５０Ａ，５０Ｂを構成する各質点の内、緩衝オブジェクト６０と結合している各質点については、緩衝オブジェクト６０内の結合先の各質点から受ける力と、筋モデル５０Ａ，５０Ｂ内の結合先の各質点から受ける力との両方に応じて、筋モデル５０Ａ，５０Ｂ内の各質点の位置を計算する。一方で、運動解析装置１００は、筋モデル５０Ａ，５０Ｂを構成する各質点の内、緩衝オブジェクト６０と結合していない各質点については、筋モデル５０Ａ，５０Ｂ内の結合先の各質点から受ける力に応じて、筋モデル５０Ａ，５０Ｂ内の各質点の位置を計算する。
【００７９】
また、運動解析装置１００は、緩衝オブジェクト６０を構成する各質点については、緩衝オブジェクト６０内の結合先の各質点から受ける力と、筋モデル５０Ａ，５０Ｂ内の結合先の各質点から受ける力との両方に応じて、緩衝オブジェクト６０内の各質点の位置を計算する。
【００８０】
これにより、運動解析装置１００は、筋モデル５０Ａ，５０Ｂ内の各質点間の力のやり取りを計算できるだけでなく、筋モデル５０Ａ，５０Ｂと緩衝オブジェクト６０との間の力のやり取りを計算することができ、筋モデル５０Ａ，５０Ｂと緩衝オブジェクト６０との干渉を回避することができる。
【００８１】
［４  筋モデルおよび骨モデルの統合］
筋モデル等の形状が変化するモデルは、形状が変化しない剛体モデルに組み合わせることができる。図８を参照して、剛体モデルと、筋モデルとを組み合わせた筋骨格モデルについて説明する。以下では、剛体モデルの一例である骨モデルを例に挙げて筋骨格モデルの説明を行うが、剛体モデルは、骨モデルに限定されない。図８は、筋骨格モデルを概略的に示す図である。
【００８２】
筋骨格モデルを扱うためには、運動解析装置１００は、骨モデルだけでなく、筋モデルの運動方程式（上記式（１）参照）も同時に扱う必要がある。運動解析装置１００は、２種の運動方程式を同時に扱うために、筋の端部に含まれる質点（以下、「端点」ともいう。）が骨に対して不動であるという拘束条件の元に、両運動方程式の計算を行う。端点とは、筋の始端（すなわち、ｋ＝０の断面）と終端（すなわち、ｋ＝Ｎｋの断面）とに含まれる質点のことをいう。
【００８３】
骨モデルと筋モデルとの間の力のやり取りは、端点で行われる。筋モデルが骨モデルに与える力は、図８（Ａ）に示される力ベクトル「ｆ1」および「ｆ2」である。
【００８４】
骨モデルと筋モデルとの間の力のやり取りは、骨モデルと筋モデルとの干渉点においても行われる。当該干渉点においては、大きさが等しい反対方向の反発力「ｆCol」のそれぞれが骨モデルおよび筋モデルのそれぞれに与えられる。以上をまとめると、筋モデルから受ける力を組み込んだ骨モデルの運動方程式は、図８に示される上腕骨を例に挙げると以下の式（７），（８）で表わされる。
【００８５】
【数７】
2019192053.tif 000009
【００８６】
【数８】
2019192053.tif 000010
【００８７】
ここで、「ｍ1」は、上腕の質量を表わす。「g」は、重力加速度を表わす。「Ｆs」は、肩関節の関節間力を表わす。「Ｆe」は、肘関節の関節間力を表わす。「ｆi,jEnd」は、ｊ番目の筋におけるi番目の筋線維が上腕骨に与える力を表わす。「ａg1」は、上腕重心位置の加速度を表わす。「Ｌs」は、上腕重心位置から肩関節までのベクトルを表わす。「Ｌe」は、上腕重心位置から肘関節までのベクトルを表わす。「ｌi,jEnd」は、上腕重心位置からｊ番目の筋のi番目の端点までのベクトルを表わす。「ｌCol」は、上腕重心位置から干渉点までのベクトルを表わす。Ｉ1は、上腕の慣性テンソルを表わす。ドット付の「ω1」は、上腕の角加速度ベクトルを表わす。ドット無しの「ω1」は、上腕の角速度ベクトルを表わす。
【００８８】
なお、端点は、筋モデルの上記式（１）の運動方程式で表わされるものの、骨に対して不動であるという拘束条件があるため、実質的には骨モデルのシミュレーション結果に従うべき質点である。つまり、端点については、運動解析装置１００は、骨モデルのシミュレーション結果を参照して１タイムステップ後の位置を更新する。一方、端点以外の質点については、運動解析装置１００は、筋モデルのシミュレーション結果をそのまま用いて１タイムステップ後の位置を更新する。このように、運動解析装置１００は、端点とそれ以外の質点とを別に扱うことで、骨モデルと筋モデルとを統合的に扱うことができる。
【００８９】
［５  運動解析装置１００のハードウェア構成］
図９を参照して、運動解析装置１００のハードウェア構成の一例について説明する。図９は、運動解析装置１００の主要なハードウェア構成を示すブロック図である。図９に示されるように、運動解析装置１００は、ＲＯＭ（Read  Only  Memory）１０１と、ＣＰＵ１０２と、ＲＡＭ（Random  Access  Memory）１０３と、ＧＰＵ１０４と、ネットワークインターフェース１０５と、ディスプレイインターフェース１０６と、記憶装置１１０とを含む。
【００９０】
ＲＯＭ１０１は、オペレーティングシステム、運動解析装置１００で実行される制御プログラム等を格納する。ＣＰＵ１０２は、オペレーティングシステムや運動解析プログラム１１４等の各種プログラムを実行することで、運動解析装置１００の動作を制御する。ＲＡＭ１０３は、ワーキングメモリとして機能し、プログラムの実行に必要な各種データを一時的に格納する。
【００９１】
ＧＰＵ１０４は、筋モデル５０Ａ，５０Ｂや緩衝オブジェクト６０の運動解析処理を高速化するための演算装置である。ＧＰＵ１０４は、筋モデル５０Ａを構成する各質点と、筋モデル５０Ｂを構成する各質点と、緩衝オブジェクト６０を構成する各質点との運動解析を上述の「２．３  ＧＰＵへの実装方法」で説明した方法で並列処理する。
【００９２】
ディスプレイインターフェース１０６は、運動解析装置１００とディスプレイ１０７とを接続するための端子である。ディスプレイ１０７は、たとえば、液晶ディスプレイ、有機ＥＬ（Electro  Luminescence）ディスプレイ、またはその他の表示機器等を含む。運動解析プログラム１１４が実行されることにより、ディスプレイ１０７は、筋モデルや骨モデルの運動解析結果を表示する。
【００９３】
記憶装置１１０は、たとえば、ハードディスクや外付けの記憶装置等の記憶媒体である。一例として、記憶装置１１０は、人体の各部位をコンピュータ上で表現するための質点情報１１２（図２参照）と、各部位間の結合関係を規定している結合情報１１３（図６参照）と、本実施の形態に従う各種の処理を実現するための運動解析プログラム１１４とを格納する。
【００９４】
運動解析プログラム１１４は、単体のプログラムとしてではなく、任意のプログラムの一部に組み込まれて提供されてもよい。この場合、任意のプログラムと協働して本実施の形態に従う処理が実現される。このような一部のモジュールを含まないプログラムであっても、本実施の形態に従う運動解析装置１００の趣旨を逸脱するものではない。さらに、本実施の形態に従う運動解析プログラム１１４によって提供される機能の一部または全部は、専用のハードウェアによって実現されてもよい。さらに、運動解析装置１００とサーバとが協働して、本実施の形態に従う処理を実現するようにしてもよい。さらに、少なくとも１つのサーバが本実施の形態に従う処理を実現する、所謂クラウドサービスの形態で運動解析装置１００が構成されてもよい。
【００９５】
［７  まとめ］
（７．１  利点）
（ａ）演算時間の短縮
上述の運動解析方法を用いることで、運動解析装置１００は、筋モデル５０Ａと筋モデル５０Ｂとの間の距離を緩衝オブジェクト６０を通じて一定に保つことができる。その結果、筋モデル５０Ａと筋モデル５０Ｂとが初期姿勢で互いに接していた場合には、その接していた状態を保つことになる。
【００９６】
従来の方法では、干渉判定を行う際には、三次元オブジェクトを三角形の集合（すなわち、ポリゴンオブジェクト）として扱うことが多かった。この場合には、各オブジェクトから一つずつ三角形を選択し、２つの三角形の干渉判定を行う必要がある。三角形の数の多い精巧なオブジェクトの場合には、この作業に膨大な時間を要する。また、干渉判定があった場合に、その干渉深度に応じて両オブジェクトに反対方向の力（すなわち、ペナルティフォース）を与えることで干渉を回避していた。
【００９７】
本提案手法は、従来の方法での多くのデメリットを解消することができる。本提案手法の１つ目のメリットは、三次元オブジェクト同士の干渉判定を行う必要がなく、計算時間を大幅に短縮できる点である。緩衝オブジェクト６０の質点と筋モデル５０Ａ，５０Ｂの質点との間で生じる力ベクトルは、両質点の位置関係だけで決まるため、ポリゴンオブジェクトという概念を利用する必要がない。
【００９８】
本提案手法の２つ目のメリットは、急激な力の変化が発生しやすいペナルティフォースを用いる必要がない点である。干渉が判定された場合にだけペナルティフォースを与える方法では、判定時に急激に大きな力が発生してしまう場合があり、シミュレーションが安定しない。一方で、本実施の形態に従う運動解析方法では、緩衝オブジェクト６０の質点と筋モデル５０Ａ，５０Ｂの質点との間に常に非ゼロの力が作用するため、急激な力の変化が起こりにくく、シミュレーションが安定する。
【００９９】
（ｂ）シミュレーション性能の改善
従来、大胸筋、三角筋、上腕二頭筋等の多数の筋が重なる上腕骨大結節稜付近においては、シミュレーションが不安定になることが多かった。特に、肩関節の外転角度を９０度程度まで大きくした場合には、シミュレーションが発散してしまうことがあった。一方で、本提案手法を用いると発散することはなく、肩関節の外転角度が９０度を超えても安定してシミュレーションを継続することが可能となった。図１０は、本提案手法のシミュレーション結果の一例を示す図である。図１０に示されるように、本提案手法では、肩関節の外転角度が９０度を超えている場合でも、安定してシミュレーションが継続している。
【０１００】
（７．２  利用例）
本提案手法は、人体ＣＧ（Computer Graphics）モデルのように、複数の変形可能オブジェクトが常に接している状況を再現する場合に有用である。特に、人体筋骨格モデルの場合には、各筋を個別にモデル化する必要があるため、オブジェクトの数、および、オブジェクト同士の干渉数が膨大になる。そのため、ペナルティフォース法による干渉回避方法では、いたる箇所でシミュレーションが不安定となるため、このことを解消できる本提案手法は、極めて画期的である。
【０１０１】
また、本提案手法は、筋だけではなく、内蔵・血管などの他の部位にも適用可能である。内蔵・血管も、筋と同様に、結合組織を介して他の組織と緩やかに結合している。そのため、干渉判定による離散的な制御よりも、緩衝オブジェクト６０を配置することによる連続的な制御の方が効果的である。
【０１０２】
今回開示された実施の形態は全ての点で例示であって制限的なものではないと考えられるべきである。本発明の範囲は上記した説明ではなくて特許請求の範囲によって示され、特許請求の範囲と均等の意味および範囲内での全ての変更が含まれることが意図される。
【符号の説明】
【０１０３】
５０，５０Ａ，５０Ｂ  筋モデル、６０  緩衝オブジェクト、１００  運動解析装置、１０１  ＲＯＭ、１０２  ＣＰＵ、１０３  ＲＡＭ、１０４  ＧＰＵ、１０５  ネットワークインターフェース、１０６  ディスプレイインターフェース、１０７  ディスプレイ、１１０  記憶装置、１１２  質点情報、１１３  結合情報、１１４  運動解析プログラム、１５０  演算部、１５１  第１算出部、１５２  第２算出部、１５３  第３算出部、１６０  位置算出部。
【図１】
2019192053.tif 000011
【図２】
2019192053.tif 000012
【図３】
2019192053.tif 000013
【図４】
2019192053.tif 000014
【図５】
2019192053.tif 000015
【図６】
2019192053.tif 000016
【図７】
2019192053.tif 000017
【図８】
2019192053.tif 000018
【図９】
2019192053.tif 000019
【図１０】
2019192053.tif 000020