特開 2024-121236 光吸収係数算出装置および光吸収係数算出プログラム

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】公開特許公報(A)

(11)【公開番号】P2024121236

(43)【公開日】2024-09-06

(54)【発明の名称】光吸収係数算出装置および光吸収係数算出プログラム

(51)【国際特許分類】

   G01N 23/046 20180101AFI20240830BHJP

   G01N 23/083 20180101ALI20240830BHJP

【ＦＩ】

G01N23/046

G01N23/083

【審査請求】未請求

【請求項の数】4

【出願形態】ＯＬ

(21)【出願番号】P 2023028211

(22)【出願日】2023-02-27

(71)【出願人】

【識別番号】504255685

【氏名又は名称】国立大学法人京都工芸繊維大学

(74)【代理人】

【識別番号】100174780

【弁理士】

【氏名又は名称】小野 敦史

(72)【発明者】

【氏名】西川 幸宏

【テーマコード（参考）】

2G001

【Ｆターム（参考）】

2G001AA01

2G001BA11

2G001CA01

2G001DA09

2G001FA08

2G001JA08

(57)【要約】

【課題】高速に対象物の任意の部位の光吸収係数を算出することができる光吸収係数算出技術を提供する。
【解決手段】光吸収係数算出装置Ａは、特定領域中の特定座標をグローバル座標値に変換するグローバル座標変換作用素を決定するグローバル座標変換作用素決定部３と、グローバル座標値を相対回転させた後の座標に変換する回転変換作用素を決定する回転変換作用素決定部４と、特定座標にグローバル座標変換作用素と回転変換作用素とを作用させて得られる座標に基づいて、検出器座標値を決定する検出器座標値決定部５と、透過投影データに対してフィルタ関数との畳み込み演算を行う畳み込み演算部６と、回転角度毎に検出器座標値の透過投影データに基づいて算出される畳み込み演算の演算値を取得し、その演算値を積算することにより特定座標の光吸収係数を算出する積分部７と、を備えている。
【選択図】図１

【特許請求の範囲】

【請求項1】

光源から対象物に対して照射された光線の検出器の検出面上で得られた透過投影データに基づき、当該対象物の特定領域の光吸収係数を算出する光吸収係数算出装置であって、
前記透過投影データは、前記対象物を前記光線の光軸に対して相対回転させ、回転角度毎に得られたものであり、
前記特定領域中の座標である特定座標をグローバル座標系における座標に変換するための作用素であるグローバル座標変換作用素を決定するグローバル座標変換作用素決定部と、
前記グローバル座標系における座標を前記相対回転させた後の座標に変換するための作用素である回転変換作用素を決定する回転変換作用素決定部と、
前記特定座標に前記グローバル座標変換作用素を作用させた後に、前記回転変換作用素を作用させて得られる座標に基づいて、前記検出面上での座標である検出器座標値を決定する検出器座標値決定部と、
前記透過投影データに対してフィルタ関数との畳み込み演算を行う畳み込み演算部と、
前記回転角度毎に前記検出器座標値の前記透過投影データに基づいて算出される前記畳み込み演算の演算値を取得し、当該演算値を積算することにより前記特定座標の前記光吸収係数を算出する積分部と、を備える光吸収係数算出装置。

【請求項2】

前記回転変換作用素決定部は、想定する座標軸に対する前記相対回転の回転軸のずれを補正する作用素を含む前記回転変換作用素を決定する請求項１記載の光吸収係数算出装置。

【請求項3】

前記グローバル座標変換作用素と前記回転変換作用素とは同次座標系における正方行列であり、
前記検出器座標値決定部は、前記回転変換作用素の右側に前記グローバル座標変換作用素を掛けることにより求められる合成変換行列と、前記特定座標を同次座標系で表現したベクトルと、の積に基づいて前記検出器座標値を決定する請求項１または２記載の光吸収係数算出装置。

【請求項4】

光源から対象物に対して照射された光線の検出器の検出面上で得られた透過投影データに基づき、当該対象物の特定領域の光吸収係数を算出する光吸収係数算出プログラムであって、
前記透過投影データは、前記対象物を前記光線の光軸に対して相対回転させ、回転角度毎に得られたものであり、
前記特定領域中の座標である特定座標をグローバル座標系における座標に変換するための作用素であるグローバル座標変換作用素を決定するグローバル座標変換作用素決定機能と、
前記グローバル座標系における座標を前記相対回転させた後の座標に変換するための作用素である回転変換作用素を決定する回転変換作用素決定機能と、
前記特定座標に前記グローバル座標変換作用素を作用させた後に、前記回転変換作用素を作用させて得られる座標に基づいて、前記検出面上での座標である検出器座標値を決定する検出器座標値決定機能と、
前記透過投影データに対してフィルタ関数との畳み込み演算を行う畳み込み演算機能と、
前記回転角度毎に前記検出器座標値の前記透過投影データに基づいて算出される前記畳み込み演算の演算値を取得し、当該演算値を積算することにより前記特定座標の前記光吸収係数を算出する積分機能と、をコンピュータに実現させる光吸収係数算出プログラム。

【発明の詳細な説明】

【技術分野】

【0001】

本発明は、対象物に対して照射されたＸ線の透過投影データから対象物の各々の部位の光吸収係数を算出する技術に関する。

【背景技術】

【0002】

対象物の内部構造を調べるための技術としてＣＴ（Ｃｏｍｐｕｔｅｄ Ｔｏｍｏｇｒａｐｈｙ）が用いられている。ＣＴでは、Ｘ光源と検出器との間に対象物を配置し、対象物を相対回転させることにより複数方向から透過投影データを取得（撮影）している。この複数の透過投影データからは、対象物の部位の光吸収係数が求められ、その光吸収係数に基づいて断層画像等が生成されている。

【0003】

透過投影データから光吸収係数を算出する方法の一つとしてＦＢＰ（Ｆｉｌｔｅｒｅｄ Ｂａｃｋ Ｐｒｏｊｅｃｔｉｏｎ）法が知られている。ＦＢＰ法では、理論的には対象物の断面像ｆの座標（ｘ，ｙ）における光吸収係数は以下の式で求められる。ここで、ｐは透過投影データ、ｒは透過投影データの座標、θは対象物の相対回転角、Ｒはフィルタである。

【数1】

【0004】

式（１）から分かるように、一つの画素に対応する光吸収係数を算出するためには３重の積分が必要であるため、２次元の断層画像全体の光吸収係数を算出するためには５重の積分が必要となり、演算量が膨大となる。この演算量を低減するために、重畳積分法と呼ばれる方法が提案されている。式（１）を見ると、光吸収係数の演算は、透過投影データｐのフーリエ変換、フーリエ空間における乗算（フィルタリング処理）、逆フーリエ変換であることが分かる。重畳積分法では、フーリエ空間における乗算は実空間での畳み込み演算（フィルタリング処理）と等価であることに着目し、式（１）を式（２）のように変形している。式（２）では、１回の畳み込み演算と１重の積分とにより、式（１）と同等の演算を少ない演算量で行うことができる。

【数2】

【0005】

また、光源がファンビームやコーンビームの場合の高速化技術についても提案がなされている。例えば、特許文献１の技術では、ファンビームによる透過投影像から直線の投影軸に軸投影したデータを求め、その軸投影したデータから逆投影画素データを求め、逆投影画素データから断層画像の各画素の光吸収係数（画素値）を算出している。これにより、逆投影処理の簡単化および高速化を実現することができる。

【0006】

また、特許文献２の技術では、コーンビームによる透過投影像から断層画像の各画素の光吸収係数（画素値）を算出する処理において、撮影幾何条件に応じて予めコンボリューション重み係数、逆投影重み係数、逆投影データの投影画像上の座標位置、補間係数をテーブルに保持しておき、これらを用いて光吸収係数を算出している。これにより、演算処理を低減し、高速化を図っている。

【先行技術文献】

【特許文献】

【0007】

【特許文献1】特開２００３－２４５２７６号公報

【特許文献2】特開２００７－０６８７１７号公報

【発明の概要】

【発明が解決しようとする課題】

【0008】

重畳積分法や特許文献１，２の技術では、特定の条件下、具体的には断層面が光軸と直交している場合には処理を高速化ができるが、それ以外の条件には対応できない。そのため、例えば、光軸と直交していない断層面の断層画像が必要な場合には、対象物全体の光吸収係数画像を算出してから断層画像を生成する必要がある。そのため、不要な部位の光吸収係数を算出することとなり、高速化のメリットを十分に享受することができない。

【0009】

本発明は上記課題に鑑みてなされたものであり、その目的は、高速に対象物の任意の部位の光吸収係数を算出することができる光吸収係数算出技術を提供することにある。

【課題を解決するための手段】

【0010】

上記課題を解決するために、本発明に係る、光源から対象物に対して照射された光線の検出器の検出面上で得られた透過投影データに基づき、当該対象物の特定領域の光吸収係数を算出する光吸収係数算出装置は、前記透過投影データは、前記対象物を前記光線の光軸に対して相対回転させ、回転角度毎に得られたものであり、前記特定領域中の座標である特定座標をグローバル座標系における座標に変換するための作用素であるグローバル座標変換作用素を決定するグローバル座標変換作用素決定部と、前記グローバル座標系における座標を前記相対回転させた後の座標に変換するための作用素である回転変換作用素を決定する回転変換作用素決定部と、前記特定座標に前記グローバル座標変換作用素を作用させた後に、前記回転変換作用素を作用させて得られる座標に基づいて、前記検出面上での座標である検出器座標値を決定する検出器座標値決定部と、前記透過投影データに対してフィルタ関数との畳み込み演算を行う畳み込み演算部と、前記回転角度毎に前記検出器座標値の前記透過投影データに基づいて算出される前記畳み込み演算の演算値を取得し、当該演算値を積算することにより前記特定座標の前記光吸収係数を算出する積分部と、を備えている。

【0011】

この構成では、グローバル座標変換作用素と回転変換作用素との作用により、相対回転後（回転角０°を含む）の対象物中の特定座標をグローバル座標系中の座標に変換することができ、さらに、検出器座標値決定部により、その座標を検出面上の座標に変換することができる。すなわち、これらの作用により、特定座標が検出器上のどの座標に対応するかのマッピングをすることができる。そのため、直接的に式（２）と同等の計算を行うことができる。これにより、対象物の必要な領域のみの光吸収係数を算出することができ、従来の対象物全体の光吸収係数を算出した後に断層画像を生成する方法に比べ、高速に光吸収係数の算出、さらには、その画像化を行うことができる。

【0012】

本発明に係る光吸収係数算出装置で用いる透過投影データは、対象物を光軸に対して相対回転しながら取得されるものであるため、その回転軸を正確に設定する必要がある。この回転軸は一般的に対象物に対して設定される座標系の軸が使用される。しかしながら、対象物を精緻に設置することは非常に困難であるため、対象物の座標系は想定している理想的な座標系からずれている可能性がある。このような場合には、算出される光吸収係数に誤差が含まれ、画質劣化の原因となる。

【0013】

そのため、本発明に係る光吸収係数算出装置の好適な実施形態の一つでは、前記回転変換作用素決定部は、想定する座標軸に対する前記相対回転の回転軸のずれ等を補正する作用素を含む前記回転変換作用素を決定する。

【0014】

この構成では、回転変換作用素に回転軸のずれを補正するための作用素を含めているため、回転軸が想定する座標軸に対してずれている場合でも正しく光吸収係数を算出することができる。

【0015】

本発明に係る光吸収係数算出装置の好適な実施形態の一つでは、前記グローバル座標変換作用素と前記回転変換作用素とは同次座標系における正方行列であり、前記検出器座標値決定部は、前記回転変換作用素の右側に前記グローバル座標変換作用素を掛けることにより求められる合成変換行列と、前記特定座標を同次座標系で表現したベクトルと、の積に基づいて前記検出器座標値を決定する。

【0016】

この構成では、グローバル座標変換作用素と回転変換作用素とは正方行列として表現されるため、これらを連続して作用させる作用素はこれらの積である正方行列として表現することができる。そのため、検出器座標値は行列とベクトルとの積という簡易な方法で算出することができる。

【0017】

また、本発明は、光源から対象物に対して照射された光線の検出器の検出面上で得られた透過投影データに基づき、当該対象物の特定領域の光吸収係数を算出する光吸収係数算出プログラムも権利対象としており、そのような光吸収係数算出プログラムは、前記透過投影データは、前記対象物を前記光線の光軸に対して相対回転させ、回転角度毎に得られたものであり、前記特定領域中の座標である特定座標をグローバル座標系における座標に変換するための作用素であるグローバル座標変換作用素を決定するグローバル座標変換作用素決定機能と、前記グローバル座標系における座標を前記相対回転させた後の座標に変換するための作用素である回転変換作用素を決定する回転変換作用素決定機能と、前記特定座標に前記グローバル座標変換作用素を作用させた後に、前記回転変換作用素を作用させて得られる座標に基づいて、前記検出面上での座標である検出器座標値を決定する検出器座標値決定機能と、前記透過投影データに対してフィルタ関数との畳み込み演算を行う畳み込み演算機能と、前記回転角度毎に前記検出器座標値の前記透過投影データに基づいて算出される前記畳み込み演算の演算値を取得し、当該演算値を積算することにより前記特定座標の前記光吸収係数を算出する積分機能と、をコンピュータに実現させる。

【0018】

当然ながら、このような光吸収係数算出プログラムも上述した光吸収係数算出装置と同様の付加的な特徴構成を付加することができ、同様の作用効果を奏する。

【図面の簡単な説明】

【0019】

【図1】光吸収係数算出装置の機能ブロック図である。

【図2】透過投影データを取得する装置の概略図である。

【図3】表示画像座標系と対象物座標系との関係を表す図である。

【図4】回転軸のずれの例を示す図である。

【図5】グローバル座標系中の座標と検出器座標系中の座標との関係を示す図である。

【図6】光吸収係数算出装置の処理の流れを表すフローチャートである。

【発明を実施するための形態】

【0020】

以下に図面を用いて、本発明に係る光吸収係数算出装置の実施形態を説明する。本実施形態では、対象物の特定の領域の光吸収係数を画像化して表示する場合を説明する。以下、この画像化された光吸収係数を表示画像と称し、この特定の領域を表示画像生成領域Ｒと称する。なお、表示画像生成領域Ｒは平面や曲面上の閉領域、さらには３次元領域であっても構わない。表示画像生成領域が平面上の閉領域の場合には、表示画像は一般的なＣＴにおける断層画像に相当する。

【0021】

先ず、光吸収係数算出装置を説明する前に、図２を用いて光吸収係数算出装置で用いられる透過投影データを取得（撮影）するための装置および各座標系について説明する。

【0022】

図２に示すように、透過投影データの取得は、Ｘ線（本発明における光線に相当）を照射する光源Ｓと、光源Ｓに対して対向配置されたＸ線を検出する検出器Ｄとによって行われる。光源Ｓと検出器Ｄとの間に対象物Ｏが配置される。本実施形態では、光源Ｓから照射されるＸ線はコーンビームであり、検出器Ｄは平面状である。なお、本実施形態では、対象物Ｏを後述するＸｏ軸周りに回転させることによりＸ線の照射方向（光軸）に対して相対回転させ、複数方向からの透過投影データを撮影しているが、対象物Ｏを固定し、光源Ｓと検出器ＤとをＸｏ軸周りに回転させても構わない。

【0023】

光源Ｓに対しては左手系の３次元直交座標系（以下、光源座標系と称する）が設定されている。本実施形態では、この光源座標系をグローバル座標系として使用するが、別途グローバル座標系を設定しても構わない。光源座標系は光源Ｓから検出器Ｄに向かう垂線（光軸）をｚ軸（以下、Ｚｇ軸と表記する）とし、Ｚｇ軸に直交するｘ軸（以下、Ｘｇ軸と表記する）とｙ軸（以下、Ｙｇ軸と表記する）とを設定している。また、Ｘ線の照射中心を原点Ｏｇとしている。

【0024】

検出器Ｄは２次元直交座標系（以下、検出器座標系と称する）が設定されており、ｘ軸，ｙ軸をそれぞれＸｄ軸およびＹｄ軸と表記する。検出器Ｄは、Ｘｄ軸とＹｄ軸とがそれぞれＸｇ軸とＹｇ軸と平行となるように配置されている。なお、本実施形態では検出器座標系の原点Ｏｄは左下である。

【0025】

一方、対象物Ｏにも独自の左手系の３次元直交座標系（以下、対象物座標系と称する）が設定されており、ｘ軸，ｙ軸，ｚ軸をそれぞれＸｏ軸，Ｙｏ軸，Ｚｏ軸と表記し、その原点をＯｏと表記する。なお、本実施形態で使用する透過投影データは、Ｘｏ軸を回転軸として、対象物Ｏを所定角度回転させる毎にされたものである。本実施形態では、回転角度はθ１（＝０°）～θＮ（＝１８０°）のＮ段階とする。

【0026】

図１は、本実施形態に係る光吸収係数算出装置Ａの機能ブロック図である。本実施形態では、光吸収係数算出装置Ａは汎用コンピュータであり、後述する機能部は、ハードウェアまたは、読み込まれた光吸収係数算出プログラム（ソフトウェア）とハードウェアとが協働することにより構成されている。当然ながら、光吸収係数算出装置Ａは専用のハードウェアのみによって構成しても構わない。

【0027】

図１に示すように、光吸収係数算出装置Ａには、表示装置Ｂと入力装置Ｃとが接続されている。表示装置Ｂは通常のディスプレイであり、対象物や表示画像の表示等に使用される。また、入力装置Ｃは、光吸収係数算出装置Ａに対する指示等を入力するためのものであり、キーボードやマウス等を用いることができる。なお、本実施形態では、対象物の透過投影データは予め取得され、記憶されている場合を説明するが、検出器からの透過投影データを直接入力して処理することも可能である。

【0028】

本実施形態における光吸収係数算出装置Ａは、制御部１，領域特定部２，グローバル座標変換作用素決定部３，回転変換作用素決定部４，検出器座標値決定部５，畳み込み演算部６，積分部７，光吸収係数算出部８，表示画像生成部９および記録部１０を備えている。

【0029】

記録部１０はＨＤＤ（Ｈａｒｄ Ｄｉｓｋ Ｄｒｉｖｅ）やＳＳＤ（Ｓｏｌｉｄ Ｓｔａｔｅ Ｄｒｉｖｅ）等の不揮発性記録媒体であり、予め取得された対象物の透過投影データをはじめ処理に必要な静的なデータが記録されている。

【0030】

制御部１は、全体の処理の流れを制御等する機能部であり、主にＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）とソフトウェアとが協働して構成されている。

【0031】

領域特定部２は、表示装置Ｂに表示された対象物に対して、マウスやキーボード等の入力装置Ｃを用いて入力された情報から、対象物Ｏの光吸収係数を算出する領域、すなわち、表示画像生成領域Ｒを特定する機能部である。

【0032】

グローバル座標変換作用素決定部３は、表示画像の座標系（以下、表示画像座標系と称する）における座標（以下、表示画像座標値と称する）をグローバル座標系の座標（以下、グローバル座標値と称する）に変換する作用素（以下、グローバル座標変換作用素と称する）を決定する機能部である。

【0033】

回転変換作用素決定部４は、グローバル座標値をＸｏ周りに回転させた座標（以下、回転座標値と称する）に変換する作用素（以下、回転変換作用素と称する）を決定する機能部である。

【0034】

検出器座標値決定部５は、表示画像座標値に対応する検出器座標系での座標（以下、検出器座標値と称する）を決定する機能部である。

【0035】

畳み込み演算部６は、透過投影データに対してフィルタ処理としての畳み込み演算を行う機能部である。畳み込み演算部６は、主にＣＰＵとＧＰＵ（Ｇｒａｐｈｉｃｓ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）とソフトウェアとを協働させて構成することにより、高速な演算を実現することができる。

【0036】

積分部７は、畳み込み演算部６の処理結果を積分する機能部である。

【0037】

光吸収係数算出部８は、表示画像中の全ての画素に対応する光吸収係数を算出する機能部である。

【0038】

表示画像生成部９は、光吸収係数算出部８によって算出された光吸収係数に基づいて、表示画像を生成する機能部である。

【0039】

ここで、本発明における光吸収係数の算出方法の基本的な考え方を説明する。上述したように、表示画像中のある画素の光吸収係数は、その画素の座標に対応する透過投影データに対してフィルタリング（畳み込み）を施した値を、角度について積分したものである。この式（２）を計算するためには、表示画像座標値と、検出器座標値との対応をとる必要がある。しかしながら、上述したように、光源Ｓ、表示画像Ｉ、対象物Ｏ、検出器Ｄはそれぞれの座標系が設定されているため、直接的に対応をとることはできない。そこで、表示画像座標値から検出器座標値に変換する関数（作用素）をＨ（θ）として定義すると、式（２）は次の式（３）に書き換えることができる。ここで、Ｑは第１引数の座標の第２引数の回転角度の透過投影データに対するフィルタリング（畳み込み）関数であり、ベクトルｘは表示画像座標値の同次座標表現である。また、Ｑの第１引数は検出器座標値の同時座標表現である。

【数3】

【0040】

以下に、Ｈ（θ）の決定方法について説明する。図３は、表示画像座標系と対象物座標系との関係を表す図である。図には、表示画像Ｉと表示画像Ｉに対応する表示画像生成領域Ｒとを示している。ここでは説明を簡単にするために、表示画像生成領域ＲをＺｏ軸に直交する平面上の矩形領域としている。通常、表示画像Ｉは所定の大きさで生成され、表示画像生成領域Ｒは任意の大きさであるため、表示画像座標系における画素の大きさと対象物座標系における画素の大きさとは一致しない。この図では、表示画像系における画素の方が細かくなっている。また、表示画像Ｉの原点Ｏｉは左下に設定されているが、表示画像生成領域Ｒの原点Ｏｒは表示画像生成領域Ｒに依存するため、必ずしも左下とはならない。この図では、原点Ｏｒは表示画像生成領域Ｒの中心となっている。

【0041】

表示画像座標系と対象物座標系とはこのような関係性にあるため、表示画像座標値を対象物座標系における座標（以下、対象物座標値と称する）に変換するための作用素Ｖは、並進（原点移動）とスケーリング（画素の大きさの変換）とによって表すことができる。具体的には、対象物座標系における表示画像Ｉの原点Ｏｉの座標を（―ｘｖ，―ｙｖ，―ｚｖ）とし、それぞれの座標系における画素の大きさの比をａとすると、作用素Ｖは以下の式（４）の行列の形で表すことができる。なお、この説明では、表示画像生成領域Ｒは２次元ため、ｚｖ＝０である。

【数4】

【0042】

したがって、以下の式（５）に示すように、表示画像座標値ｐｉはこの作用素Ｖにより対象物座標値ｐｏに変換することができる。なお、数式中のベクトルｐｉおよびｐｏは表示画像座標値ｐｉおよび対象物座標値ｐｏの同次座標表現である。

【数5】

【0043】

次に、対象物座標値からグローバル座標値に変換するための作用素Ｍを考える。対象物座標系は対象物Ｏに対して設定された直交座標系であるため、作用素Ｍはグローバル座標系内の対象物Ｏの位置と姿勢とによって決定することができる。ここで、対象物座標系の３軸とグローバル座標系の３軸とがいずれも平行であり、対象物座標系の原点Ｏｏがグローバル座標系のＺｇ軸上にあると仮定すると、作用素ＭはＺｇ軸に沿った並進となる。グローバル座標系における対象物座標系の原点Ｏｏの座標を（０，０，ＳＯＤ）とすると作用素Ｍは以下の式（６）となり、式（７）に示すように、対象物座標値ｐｏは作用素Ｍによりグローバル座標値ｐｇに変換することができる。なお、式（７）中のベクトルｐｏおよびｐｇは対象物座標値ｐｏおよびグローバル座標値ｐｇの同次座標系表現である。

【数6】

【数7】

【0044】

式（７）に式（５）を代入すると以下の式（８）となる。すなわち、作用素Ｖおよび作用素Ｍにより表示画像座標値ｐｉをグローバル座標値ｐｇに変換することができる。

【数8】

【0045】

上述したように、透過投影データの取得は、対象物Ｏを対象物座標系のＸｏ軸周りに所定角度回転させながら行う。そのため、式（８）で求められるグローバル座標値ｐｇの回転後の座標を求める必要がある。以下、グローバル座標値ｐｇをＸｏ軸周りに角度θの回転させる回転変換作用素をＲ（θ）と表記する。

【0046】

左手系の３次元直交空間でのｘ軸周りの角度θの回転を表す作用素は

【数9】

で表すことができる。しかしながら、対象物座標系の原点Ｏｏはグローバル座標系の原点Ｏｇから並進した位置にあるため、式（９）ではグローバル座標系での回転後の座標は計算できない。この場合には、対象物座標系の原点Ｏｏをグローバル座標系の原点Ｏｇに一致させ、その後、Ｘｏ軸（Ｘｇ軸）周りの回転操作を行い、対象物座標系の原点Ｏｏを元の位置に戻すという操作を行えばよい。回転変換作用素Ｒ（θ）は以下の式（１０）となる。

【数10】

ここで、

【数11】

である。

【0047】

ただし、回転変換作用素Ｒ（θ）を式（１０）で表すことができるのは、対象物座標系の各軸がグローバル座標系の各軸と平行であり、対象物座標系の原点Ｏｏがグローバル座標系のＺｇ軸上にあるときだけである。これを理想的な対象物座標系と称する。しかしながら、実際にこのように対象物座標系、すなわち、対象物Ｏを配置することは困難であり、通常は原点Ｏｏの位置ずれや座標軸のずれが生じている。このようなずれが生じたものを現実の対象物座標系と称する。

【0048】

しかしながら、前者のずれは並進、後者のずれは回転として捉えることができるため、このようなずれが生じている場合も、上記と同様に考えることができる。すなわち、現実の対象物座標系のＸｏ軸周りの回転操作は、現実の対象物座標系を理想的な対象物座標系に一致させ、式（１０）を作用させた後に、現実の対象物座標系に戻すことにより実現することができる。現実の対象物座標系を理想的な対象物座標系に一致させる操作を行う作用素を作用素Ｚとすると、一般的な回転変換作用素Ｒ（θ）は以下の式（１２）となる。

【数12】

【0049】

例えば、図４のように、現実の対象物座標系が理想的な対象物座標系に対して、Ｙｏ軸方向にδｙ、Ｙｏ軸周りにφのずれがあった場合を考える。図中、実線が理想的な対象物座標系のＸｏ軸（本発明における想定する座標軸の例）であり、点線が実際の対象物座標系のＸｏ軸（Ｘｏ’と表記している）である。このとき、作用素Ｚは－δｙの大きさのＹｏ軸方向の並進と－φの大きさのＹｏ軸周りの回転となる。すなわち、この場合の回転変換作用素Ｒ（θ）は

【数13】

となる。ここで、

【数14】

である。

【0050】

このように、理想的な対象物座標系に対する現実の対象物座標系のずれは並進と回転とによって表現することができる。すなわち、ずれの大きさが既知であれば作用素Ｚを決定することができるため、回転変換作用素Ｒ（θ）は式（１３）で表現することができる。これにより、グローバル座標値ｐｇを対象物座標系のＸｏ軸周りに角度θ回転させた後のグローバル座標系における座標（回転座標値）ｐｇ’は

【数15】

となり、これに式（８）を代入すると、

【数16】

となる。

【0051】

次に、グローバル座標系の座標ｐを通ったＸ線が検出される検出器Ｄの検出面上での座標ｑを考える。なお、ここで、座標ｐ，ｑはグローバル座標系における座標である。図５は、点ｐと点ｑとの関係を示す図である。図に示すように、座標ｐ，ｑをそれぞれ（ｘｐ，ｙｐ，ｚｐ），（ｘｑ，ｙｑ，ｚｑ）とし、座標ｐからＺｇ軸に下ろした垂線とＺｇ軸との交点座標をｐ’、検出器Ｄの検出面とＺｇ軸との交点座標をＯｄ’とすると、三角形Ｏｇ－ｐ－ｐ’と三角形Ｏｇ－ｑ－Ｏｄ’とは相似形となる。ここで、光源Ｓから検出器Ｄまでの距離（Ｓｏｕｒｃｅ－ｔｏ－Ｉｍａｇｅ Ｄｉｓｔａｎｃｅ：ＳＩＤ）をＬとするとｚｑ＝Ｌであるため、三角形Ｏｇ－ｑ－Ｏｄ’は三角形Ｏｇ－ｐ－ｐ’のＬ／ｚｐ倍の大きさとなる。すなわち、

【数17】

【0052】

ただし、式（１７）で求められる座標ｑはグローバル座標系におけるものであるため、これを検出器座標系の座標（検出器座標値）に変換する必要がある。グローバル座標系と検出器座標系とは原点位置がずれているため、この変換は平行移動となる。具体的には、Ｏｄ’の座標を（ｘｄ，ｙｄ，Ｌ）とすると、検出器座標値ｐｄは

【数18】

となる。なお、グローバル座標系における画素の大きさと検出器座標系における画素の大きさとが異なっている場合には更にスケーリングが必要となるが、ここでは画素の大きさは同一としてスケーリングは不要とする。

【0053】

ここで、式（１８）により表されるグローバル座標値から検出器座標値に変換する作用素（関数）をＰ（）と表記すると、表示画像座標値ｐｉは以下の式（１９）により検出器座標値ｐｄに変換することができる。

【数19】

【0054】

したがって、Ｐを行列表現とするとＨ（θ）は以下の式（２０）となる。なお、ここでは便宜上Ｐを行列表現としているが、実際の演算では行列形式ではなく式（１８）で演算する方が演算量を低減することができる。

【数20】

【0055】

次に、図６のフローチャートを用いて本実施形態における光吸収係数算出装置Ａの処理の流れを説明する。先ず、制御部１は記録部１０から対象物Ｏの透過投影データを撮影した際のジオメトリ情報を読み出す（＃０１）。ジオメトリ情報は、例えば、光源Ｓから検出器Ｄまでの距離Ｌ（Ｓｏｕｒｃｅ－ｔｏ－Ｉｍａｇｅ Ｄｉｓｔａｎｃｅ：ＳＩＤ）やグローバル座標系における対象物座標系の原点Ｏｏの座標、座標軸Ｘｇに対する座標軸Ｘｏのずれ等である。

【0056】

次に、制御部１は表示装置Ｂに対象物の外形を表示し、ユーザの入力装置Ｃの操作によって指定された表示画像生成領域Ｒを特定する情報（以下、領域情報と称する）を取得する（＃０２）。領域情報としては、例えば、表示画像生成領域Ｒが平面の場合には対象物座標系における法線方向、中心座標、大きさ等を用いることができる。当然ながら、領域を特定できる情報であれば他の情報を用いることができる。ジオメトリ情報や領域情報は、ＲＡＭ（Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）等の揮発性記憶媒体に一時的に記憶される。

【0057】

また、領域情報やユーザからの入力に基づいて表示画像のサイズ等のプロパティが設定される。具体的には、表示画像は表示画像生成領域Ｒと同次元であり、アスペクト比も等しくなっている。このとき、表示画像の最長辺の画素数が、予め設定された、または、ユーザから入力された最大画素数となるよう、各軸方向の画素数が定められる。

【0058】

入力装置Ｃから表示画像生成指示が入力されると、制御部１はグローバル座標変換作用素決定部３に対して、グローバル座標変換作用素を決定するよう指示を出す。グローバル座標変換作用素決定部３は以下の処理によりグローバル座標変換作用素を決定する（＃０３）。具体的には、グローバル座標変換作用素決定部３は上述した式（４）と（６）とで表されるＭとＶとの積をグローバル座標変換作用素として決定する。

【0059】

次に制御部１は、回転変換作用素決定部４に対して回転変換作用素を決定するよう指示を出す。具体的には、回転変換作用素決定部４は、上述した式（１２）で表される回転変換作用素Ｒ（θ）を決定する（＃０４）。なお、現実の対象物座標系が理想的な対象物座標系に対してずれがない場合には、作用素Ｚは単位行列となる。

【0060】

次に、制御部１は畳み込み演算部６に対して透過投影データの畳み込み演算を行うよう指示し、畳み込み演算部６は各々の透過投影データに対して畳み込み演算を行う（＃０５）。具体的には、畳み込み演算部６は、記録部１０から回転角θにおける座標（ｘ，ｙ）の透過投影データおよび畳み込み演算に必要な透過投影データを取得して畳み込み演算を行う。なお、この畳み込み演算は通常の重畳積分法におけるものと同様であるため詳細は省略する。畳み込み演算部６は、この畳み込み演算を全ての座標位置および回転角に対して実行し、その演算結果を（座標、角度）によって定まる配列形式でＲＡＭ等に記憶させておく。この演算値は式（３）中の関数Ｑ（）に相当する。

【0061】

なお、全ての座標位置の透過投影データに対して畳み込み演算を行う場合には、表示画像の生成に寄与しないデータに対しても演算を行っている可能性がある。この不要な演算を回避するためには、領域情報に基づいて表示画像の生成に寄与するデータの座標を特定しておき、その特定された座標位置のデータに対してのみ畳み込み演算を実行すればよい。

【0062】

このようにして畳み込み演算が完了すると、制御部１は光吸収係数算出部８に対して、表示画像の全ての画素に対応する光吸収係数を算出するよう指示を出す。光吸収係数算出部８は以下の処理により、光吸収係数を算出する。

【0063】

先ず、光吸収係数算出部８は、表示画像中の一つ座標（表示画像座標値ｐｉ）を特定し（＃０６）、この表示画像座標値ｐｉを積分部７に通知する。なお、本実施形態では、この表示画像座標値ｐｉが本発明における特定座標に相当する。

【0064】

表示画像座標値ｐｉが通知された積分部７は、先ず、光吸収係数ｃを０で初期化し、回転角度のインデックスｉを１で初期化する（＃０７）。そして、表示画像座標値ｐｉに対してグローバル座標変換作用素ＭＶおよび回転変換作用素Ｒ（θｉ）を作用させ、回転座標値ｐｇ’を算出する（＃０８）。具体的には式（１６）の計算を実行することになるが、Ｒ（θｉ）ＭＶを予め計算し、一つの正方行列として記憶しておくことにより演算量を低減することができる。なお、このＲ（θｉ）ＭＶにより得られる行列が本発明における合成変換行列に相当する。

【0065】

次に、積分部７は、このようにして回転座標値ｐｇ’から検出器座標値ｐｄを算出する（＃０９）。 具体的には、式（１８）において、ｐ＝ｐｇ’とすることにより検出器座標値ｐｄが算出される。そして、光吸収係数ｃに検出器座標値ｐｄと回転角θｉにより特定される畳み込み演算値Ｑ（ｐｄ，θｉ）を加算する（＃１０）。積分部７はｉがＮに達するまでこの処理を繰り返す（＃０８～＃１２）。これにより、式（３）の計算が実現される。

【0066】

光吸収係数算出部８は、一つの表示画像座標値ｐｉについての光吸収係数の算出が完了すると、表示画像中に光吸収係数の算出が完了していない画素があるか否かを判定し（＃１３）、未算出の画素が存在する場合には（＃１３のＹｅｓ分岐）、処理を＃０６に戻し、次の未算出の表示画像座標値ｐｉを特定し（＃０６）、上述の処理を行う。

【0067】

一方、表示画像中の全ての画素の光吸収係数の算出が完了すると（＃１３のＮｏ分岐）、制御部１は表示画像生成部９に対して表示画像を生成するよう指示を出す。この指示を受けた表示画像生成部９は、各々の画素に対応する光吸収係数に基づいて表示画像を生成する（＃１４）。表示画像の生成は、例えば、光吸収係数を２５５で正規化する等により実現することができる。このようにして生成された表示画像は、表示装置Ｄに表示される。

【0068】

このように、本発明に係る光吸収係数算出装置Ａでは、表示画像の画素単位で光吸収係数を算出することができるため、必要な部位のみの光吸収係数を算出することができる。これにより、不要な部位についての演算を削減し、高速に光吸収係数の算出、表示画像の生成を行うことができる。

【0069】

〔別実施形態〕
（１）上述の実施形態では、表示画像の画素に対応する光吸収係数を算出したが、対象物座標系のボクセルの光吸収係数を算出し、そのボクセルの光吸収係数から表示画像の画素の光吸収係数を算出することもできる。

【0070】

（２）上述の実施形態では、予め畳み込み演算を行っておき、その演算結果を用いて積分を行ったが、積分時に畳み込み演算を行うこともできる。

【0071】

（３）上述の実施形態では、グローバル座標変換作用素や回転変換作用素等を行列形式としたが、これに限定されるものではなく、関数形式等他の表現であっても構わない。

【産業上の利用可能性】

【0072】

本発明は、医療や産業で用いられる、対象物にＸ線を照射して得られる透過投影データから対象物の光吸収係数を算出する技術に利用することができる。

【符号の説明】

【0073】

Ａ：光吸収係数算出装置
Ｄ：検出器
Ｉ：表示画像
Ｏ：対象物
ｐｉ：表示画像座標値
ｐｏ：対象物座標値
ｐｇ：グローバル座標値
ｐｇ’：回転座標値
ｐｄ：検出器座標値
Ｒ：表示画像生成領域
Ｓ：光源
１：制御部
２：領域特定部
３：グローバル座標変換作用素決定部
４：回転変換作用素決定部
５：検出器座標値決定部
６：畳み込み演算部
７：積分部
８：光吸収係数算出部
９：表示画像生成部
１０：記録部

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】