WO1999061893A1

WO1999061893A1 - Procede et appareil permettant de mesurer la repartition de caracteristique d'un corps d'absorption/diffusion

Info

Publication number: WO1999061893A1
Application number: PCT/JP1999/002696
Authority: WO
Inventors: Yutaka Tsuchiya
Original assignee: Hamamatsu Photonics K.K.
Priority date: 1998-05-26
Filing date: 1999-05-24
Publication date: 1999-12-02
Also published as: EP1081486B1; DE69928392D1; EP1081486A4; EP1081486A1; JP3887486B2; AU3850999A; DE69928392T2; US6335792B1; JPH11337476A

Description

明糸田書

散乱吸収体の内部特性分布の計測方法及び装置技術分野

本発明は、散乱吸収体の内部特性分布を計測するための方法及びそのための装置に関する。より詳細には、本発明は、光入射位置及び光検出位置を測定対象物の表面に沿って移動させて内部の情報を得る装置に適用可能な内部特性分布の計測方法並びにそのための装置に関する。背景技術

光 C T (computer tomography)は、生体内部の光学特性あるいは吸収成分の濃度分布を計測する手法又は装置を意味し、光を生体組織に入射したときにその中を移動して検出される光（信号光）を利用する。

現在までに直進光、準直進光、及び散乱光を利用する 3種の光 C T技術が報告されている。このうち、直進光を利用する方法は光の利用率が極めて悪く、極小さい媒体にしか適用できない。例えば、近赤外線を用いる場合、標準的な生体組織では l cm の厚さに対して透過光は約 1/10 ⁵倍に減衰する。これに比べて、散乱光を利用する方法は、媒体から出射する全ての光を利用するため、信号対雑音比が向上するので、より大きい媒体への適用が期待される。準直進光を利用する方法は、前記二つの方法の中間にある。

そして、実用的な光 CT は、生体に対する最大許容光入射パワー、計測感度、所要計測時間などの制限から散乱光を利用する必要があるが、以下の技術上の問題に起因して上記光 CTは未だ実用化に至っていない。

第 1の問題は、散乱吸収体内部の光あるいはフォトンの振る舞いを十分な精度で記述する方法が開発されていないことである。散乱吸収体内部を移動するフォトンの振る舞いの解析には、従来から、輸送方程式を近似したもの、あるいは輸送理論に拡散近似を適用した光拡散方程式が用いられてきた。しかし、拡散近似は媒体中のフォトンの平均自由光路長に対して十分に大きな媒体においてのみ成立するため、比較的小さい媒体、内部の複雑な形状の組織、及び複雑な形状の媒体を取り扱うことができない。また、拡散近似は等方散乱を前提としているため、非等方な散乱特性をもつ実際の生体組織の計測に適用すると、非等方散乱に起因する無視できない誤差が発生する。さらに、拡散方程式は、解析的あるいは数値的（有限要素法など）手法のいずれを用いても、前もって境界条件を設定しないとその解が求められない。つまり、計測に先立って、光入射，検出する各部の境界条件、つまり媒体形状や界面での反射特性を設定する必要がある。また、これらが個体差などによって変化すれば、それに応じて境界条件を変えて、計算をやりなおす必要がある。従って、輸送方程式の近似式や光拡散方程式から導出される信号光と散乱吸収体の光学特性との関係を利用する光 C Tは、精度及び使い勝手に大きな問題がある。

また、上記とは別に、輸送方程式の近似式や光拡散方程式に摂動理論を適用して信号光と散乱吸収体の光学特性との関係を導出し、この関係を用いて光 C T画像を再構成する方法がある。しかし、この方法は、非線形効果（2次以上の項）の取り扱いが極めて複雑になる。この際、 2次以上の項を含めた計算は、理論的にはコンピュータで計算することが可能である力 s、現在の最速のコンピュータを用いても演算時間が膨大なものになり、実用化は不可能である。このため、常套的には 2次以上の項を無視している。そのため、この方法においては、複数個の比較的強い吸収領域が存在する媒体の光 C T画像を再構成する際などに吸収領域の相互作用を無視することができなくなり、これに起因する大きな誤差が発生する。

第 2の問題は、従来の光 C Tが、狭義の重み関数、すなわち平均光路長あるいはそれに相当する位相遅れを利用していることである。このため、吸収係数に依存して変化する検出光の平均光路長の取り扱いが極めて複雑になる。そこで通常は近似を用いるが、近似を用いると誤差が増大するという重大な問題がある。このような狭義の重み関数を利用した方法に関しては、例えば以下の文献に記載されてレヽる。 (1) S. Arridge: SPIE Institutes for Advanced Optical Technologies,

Vol. IS11, Medical Optical Tomograph ： Functional Imaging and Monitoring, 35-64 (1993); (2) R. L. Barbour and H. L. Graber: ibid. 87-120 (1993); (3)

H. L. Graber, J. Chang, R. Aronson and R. L. Barbour: ibid. 121-143 (1993);

(4) J. C. Schotland, J. C. Hasel rove and J. S. Leigh: Applied Optics, 32，

448-5453 (1993); (5) Chang, R. Aronson, H. L. Graber andR. L. Barbour: Proc.

SPIE, 2389， 448-464 (1995); (6) B. W. Pogue, M. S. Patterson, H. Jiang and K. D. Paulsen: Phys. Med. Biol. 40， 1709-1729 (1995); (7) S. R. Arridge:

Applied Optics, 34, 7395-7409 (1995); (8) H. L. Graber, J. Chang, and R.

L. Barbour: Proc. SPIE, 2570, 219-234 (1995); (9) A. Maki and H. Koizumi:

0SA TOPS, Vol.2, 299-304 (1996); (10) H. Jiang, K. D. Paulsen and Ulf L.

Osterberg: J. Opt. Soc. Am. A13, 253-266 (1996); (11) S. R. Arridge and J. C. Hebden: Phys. Med. Biol. 42， 841-853 (1997); (12) S. B. Colak, D. G.

Papaioannou, G. W. ' t Hooft, M. B. van der Mark, H. Schomberg, J. C. J.

Paasschens, J. B. M. Melissen and N. A. A. J. van Astten: Applied Optics,

36， 180-213 (1997)。

なお、光 CTの画像再構成に際しては、信号光が散乱媒体のどの部分を通過してきたか、つまり媒体内の信号光の光路分布を知ることが最も重要である。この観点から、信号光の光路分布をランダムウォーク理論などで導出する方法もあるが、上記の問題は解決されていない。

以上のように、従来の光 CTにおいては十分な精度の再構成画像が得られず、空間解像度、画像ひずみ、定量性、測定感度、所要計測時間などに大きな問題がある。

本発明者は、以上のような状況を打破するため、次のことが重要であると考えて一連の研究を推進してきた。すなわち、光 C Tを実現させるために重要なことは、強い散乱媒体である生体組織の中を移動する光の振る舞いを明らかにし、検出される信号光と吸収成分を含む散乱媒体（散乱吸収体）の光学特性との関係を明らかにし、前記信号光と前記関係とを利用して光 C T画像を再構成するァルゴリズムを開発するということである。

そして本発明者は、①マイク口'ベア'ランバート則（Microscopic Beer-Lambert Law, 以下「M B L」という）に基づくモデルを提案すると共に、 ②散乱吸収体の光学特性と信号光との関係を表す解析式を導出し、例えば以下の文献に掲載している。 ( 13) Y. Tsuchiya and T. Urakami： Jpn. J. Appl . Phys. 34， L79-81 ( 1995)

( 14) Y. Tsuchiya and T. Urakami： Jpn. J. Appl . Phys. 35, 4848-4851 ( 1996)

( 15) Y. Tsuchiya and T. Urakami： Optics Communications, 144, 269-280 ( 1997)

( 16) H. Zhang, M. Miwa, Y. Yamashita and Y. Tsuchiya: Ext. Abstr. Optics Japan ' 97, 30aA08 ( 1997)。この M B Lに基づく方法は原理的に、媒体形状、境界条件及び散乱の影響を受けないという大きな特長があり、非等方散乱媒体や小さい媒体への適用も可能である。この結果、散乱吸収体内の光子移動が明らかになり、この理論に基づいて散乱媒体の吸収係数あるレ、は吸収成分濃度を計測することが可能になった。発明の開示

しかしながら、本発明者は、上記従来の方法はいずれも上記従来の狭義の重み関数を利用していたため、不均一系への適用が困難であり、生体などの散乱吸収体内の吸収成分についての計測精度が未だ充分なものではないことを見出した。本発明は、上記従来の課題に鑑みてなされたものであり、生体などの不均一系の散乱吸収体について適用した場合であっても吸収成分についてより高い精度の計測が可能な方法及び装置を提供することを目的とする。

本発明者は、上記目的を達成すべく鋭意研究した結果、上記本発明者による知見をさらに発展させて、 ③上記マイクロ ·ベア ·ランバート則を不均一系に適用し、 ④不均一散乱媒体の光学特性と信号光との関係を表す解析式を導出し、 ⑤前記信号光と不均一散乱吸収体の光学特性との関係から従来とは異なる定義の重み関数を導出し、 ⑥この重み関数を用いて光 C T画像を再構成するためのアルゴリズム及び装置を提供することによって上記目的が達成されることを見出し、本発明に到達した。

すなわち、本発明の散乱吸収体の内部特性分布の計測方法は、

散乱吸収体である被計測媒体中に一以上の光入射位置から順次光線を入射する光入射ステップと、

前記被計測媒体内部を透過した光線を複数の光検出位置で検出する光検出ステップと、

検出された各光線に基づいて、前記光入射位置と前記光検出位置との組み合わせ毎に前記光線の所定パラメ一夕の測定値を取得する測定値取得ステップと、前記被計測媒体の吸収係数の基準値を設定する基準値設定ステツプと、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記光入射位置と前記光検出位置との組み合わせ毎に前記パラメ一夕の推定値を求める推定値算出ステップと、

前記吸収係数の基準値を用レ、て、複数のポクセルに分割された前記被計測媒体の各ポクセルにおける重み関数をマイクロ ·ベア ·ランバート則に基づいて求める重み関数演算ステップと、

前記パラメータの測定値、前記パラメ一夕の推定値及び前記重み関数に基づいて、前記各ポクセルにおける前記吸収係数の基準値に対する吸収係数の偏差を算出する吸収係数偏差算出ステツブと、

前記吸収係数の基準値及び前記吸収係数の偏差に基づいて前記各ポクセルにおける吸収係数の絶対値を算出して前記被計測媒体における吸収係数の絶対値分布を求める吸収係数絶対値算出ステップと、を含むことを特徴とする方法である。

また、本発明の散乱吸収体の内部特性分布の計測装置は、

散乱吸収体である被計測媒体中に一以上の光入射位置から順次光線を入射する光入射手段と、

前記被計測媒体内部を透過した光線を複数の光検出位置で検出する光検出手段と、

検出された各光線に基づいて、前記光入射位置と前記光検出位置との組み合わせ毎に前記光線の所定パラメ一夕の測定値を取得する測定値取得手段と、前記被計測媒体の吸収係数の基準値を設定する基準値設定手段と、

前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記光入射位置と前記光検出位置との組み合わせ毎に前記パラメ一夕の推定値を求める推定値算出手段と、

前記吸収係数の基準値を用いて、複数のポクセルに分割された前記被計測媒体の各ポクセルにおける重み関数をマイクロ ·ベア ·ランバート則に基づいて求める重み関数演算手段と、

前記パラメ一夕の測定値、前記パラメータの推定値及び前記重み関数に基づいて、前記各ポクセルにおける前記吸収係数の基準値に対する吸収係数の偏差を算出する吸収係数偏差算出手段と、

前記吸収係数の基準値及び前記吸収係数の偏差に基づいて前記各ポクセルにおける吸収係数の絶対値を算出して前記被計測媒体における吸収係数の絶対値分布を求める吸収係数絶対値算出手段と、

を備えることを特徴とする装置である。

本発明の方法及び装置においては、個々の計測に際して本発明によって初めて開示される各ポクセルにおける重み関数をマイクロ ·ベア ·ランバート則に基づいて直接的に求め、その重み関数、所定パラメ一夕の測定値及び前記パラメータの推定値に基づいて各ポクセルにおける吸収係数の偏差が算出される。この際、重み関数は計測状況とは無関係な一種の式で表わされる。このように本発明においては、個々の計測状況が変化しても適切な重み関数に基づいて吸収係数の偏差が算出されるため、近似法の採用に伴う誤差の発生が防止される。そのため、本発明の方法及び装置においては、生体などの不均一系の散乱吸収体について適用した場合であっても、各ポクセルにおける吸収係数の偏差が正確に求まることとなり、このような吸収係数偏差に基づいて得られる内部特性分布（光 C T画像）の測定精度が向上する。

上記本発明の方法及び装置において採用される好適な重み関数として、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなしたときの、前記各ボクセルにおける平均光路長とその分布の分散との関数が挙げられる。この場合、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記各ポクセルにおける平均光路長を求める平均光路長取得ステップ（平均光路長取得手段）を更に含んでいることが好ましく、前記重み関数演算ステップ（重み関数演算手段）において、下記の式：

[式中、 W iは重み関数、 Z a v )は平均光路長、 A〃_{a i}は吸収係数の偏差、び _{i v} ²は平均光路長の分布の分散をそれぞれ示す]

に基づいて前記重み関数を求めることが好ましい。

また、上記本発明の方法及び装置において採用される他の好適な重み関数として、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなしたときの、前記各ポクセルにおける所定時間領域内の平均光路長とその分布の分散との関数が挙げられる。

この場合、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記各ポクセルにおける所定時間領域内の平均光路長を求める平均光路長取得ステップ（平均光路長取得手段）を更に含んでいることが好ましく、前記重み関数演算ステップ（重み関数演算手段）において、下記の式：

[式中、は重み関数、 Z i ( _{a v} )は平均光路長、 A〃_{a i}は吸収係数の偏差、び _{i v} ²は平均光路長の分布の分散をそれぞれ示す]

に基づいて前記重み関数を求めることが好ましい。

更に、上記本発明の方法及び装置において採用される更に他の好適な重み関数として、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなしたときの、前記各ポクセルにおける群遅延とその分布の分散との関数が挙げられる。

この場合、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記各ポクセルにおける群遅延を求める群遅延取得ステップ（群遅延取得手段）を更に含んでおり、前記重み関数演算ステップ（重み関数演算手段）において、下記の式：

[式中、 W iは重み関数、

は c (媒体内の光速）倍の群遅延、 A〃_{a i}は吸収係数の偏差、び _f ²は c倍の群遅延の分布の分散をそれぞれ示す]

に基づいて前記重み関数を求めることが好ましい。

上記本発明の方法及び装置に好適な重み関数を採用した場合、吸収係数に依存して変化する平均光路長又は群遅延のみならずその分布の分散をも考慮されて重み関数が演算されるため、このような重み関数を用いて得られる内部特性分布の測定精度がより向上する傾向にある。

本発明の方法及び装置は、更に、前記吸収係数の絶対値を用いて前記各ボクセルにおける吸収成分の濃度を算出して前記被計測媒体における吸収成分の濃度分布を求める濃度算出ステップ（濃度算出手段）を含んでもよい。このような方法及び装置によれば、前述のように正確に求められた吸収係数偏差に基づいて吸収成分の濃度が求められるため、高精度の濃度分布が得られる。

上記本発明の方法及び装置を少なくとも 2つの吸収成分を含有している被計測媒体に適用する場合には、前記光入射ステップ（光入射手段）において前記被計測媒体中に入射される光線が、前記吸収成分に対する吸収係数が互いに相違する少なくとも 2つの波長を有していることが好ましい。この場合、前記光検出ステヅプ（光検出手段）において前記少なくとも 2つの波長を有する光線をそれそれ検出し、前記測定値取得ステップ（測定値取得手段）において前記少なくとも 2 つの波長を有する光線に関してそれそれ前記測定値を取得し、前記基準値設定ステツプ（基準値設定手段）において前記少なくとも 2つの波長を有する光線に関してそれそれ前記基準値を設定し、前記推定値算出ステップ（推定値算出手段）において前記少なくとも 2つの波長を有する光線に関してそれぞれ前記推定値を求め、前記重み関数演算ステップ（重み関数演算手段）において前記少なくとも 2つの波長を有する光線に関してそれぞれ前記重み関数を求め、前記吸収係数偏差算出ステップ（吸収係数偏差算出手段）において前記少なくとも 2つの波長を有する光線に関してそれぞれ前記吸収係数の偏差を算出し、前記吸収係数絶対値算出ステップ（吸収係数絶対値算出手段）において前記少なくとも 2つの波長を有する光線に関してそれそれ前記吸収係数の絶対値を算出し、前記濃度算出ステップ（濃度算出手段）において前記少なくとも 2つの波長を有する光線に関してそれぞれ前記吸収成分の濃度を算出することによって、前記被計測媒体における前記各吸収成分の濃度分布が高精度で求められる。上述の本発明の方法及び装置は、前記の求められた分布に基づいて、前記被計測媒体における該分布を示す光 C T画像を表示する画像表示ステップ（画像表示手段）を更に含んでもよい。このような本発明の方法及び装置によれば、前述のように精度良く求められた内部特性分布を画像化して高精度な光 C T画像として表示することが可能である。図面の簡単な説明

図 1は、不均一媒体内の光子移動に関するモデルを示す模式図である。

図 2は、本発明の散乱吸収体の内部特性分布の計測装置の一実施形態を示す模式図である。

図 3は、本発明の散乱吸収体の内部特性分布の計測方法の一実施形態を示すフ口一チヤ一卜である。

図 4は、本発明にかかる重み関数の演算方法の一実施形態を示すフローチヤ一卜である。

図 5は、各ポクセルにおける吸収係数分布（吸収係数の基準値と吸収係数の偏差との関係）を示すグラフである。

図 6は、ヘモグロビン及びミオグロビンの吸収スぺクトルを示すグラフである。発明を実施するための最良の形態

以下、図面を参照しつつ本発明の原理及び好適な実施形態について詳細に説明する。尚、図面中、同一又は相当部分には同一符号を付することとする。なお、散乱されつつ進行する光線に関しては 3次元座標を用いて考える必要があるが、以下では説明を簡単にするために場合により 2次元座標を用いて説明する。先ず、本発明の原理について説明する。

[本発明の原理]

本発明は、マイクロ 'ベア 'ランバート則に基づいて導出される新規な重み関数を用いる。すなわち、本発明にかかわる重み関数は、従来の重み関数、つまり各ポクセル（voxel、体積素）内の光路長や平均光路長、光子滞在時間、あるいはそれらに対応する位相遅延量などとは異なる。本発明で採用する重み関数は、媒体形状、境界条件、光入射 ·検出位置及びそれらの間の距離、媒体の散乱係数などを繰り込んだ形で表される。また、この重み関数は散乱媒体に対して一意的に定まる。本発明では上記のような重み関数を用いて光 C Tの画像再構成を行うため、後述するように、従来の問題点である誤差が大きく低減され、高精度の光 C T画像を得ることができる。

1 - 光子移動モデル

先ず、本発明の基礎である散乱媒体内の光子移動のモデルを説明する。

A ) 散乱媒体は再入射不可能な任意の 3次元形状をもつ。従って、媒体から出た光子が再度その媒体に入射することはない。この制限内で、散乱媒体と外部媒体との境界条件、例えば形状や屈折率のマッチング条件などは任意である。

B ) 散乱媒体内で、光子と媒体との非線形作用や、光子と光子との相互作用はない。従って、媒体の光学定数は光子密度に依存しない。

C ) マクロに見たとき、散乱特性及び屈折特性は一様（均一）である。し力 ^し、吸収係数分布は一様でない（不均一である）。以下、このような媒体を便宜的に不均一媒体あるいは不均一系とよぶ。従って、以下で議論する散乱媒体内の光速度は一定値 cになる。

D ) 上記不均一散乱媒体に対して、光入射 ·検出位置 p ( u , V ) (表面位置

Uに光子を入射し、他の表面位置 Vで散乱光子を検出する）を定めたとき、光路分布（ / _a =0 )がー意的に決まる。この光路分布（〃_a =0 )は、散乱媒体に吸収がないと仮定したときに光子が取り得る光路の分布あるいは光路長分布である。なお、このように吸収がないと仮定した散乱媒体を「仮想媒体」という。

E ) 上記の光路分布（〃_a =0 )は、散乱媒体に吸収がないと仮定した仮想媒体に対して定義したから、散乱媒体の吸収あるいは吸収分布には影響されない。また、この光路分布（ _a=0)は、散乱媒体の散乱特性、媒体形状、境界条件、及び光入射'検出位置 p(u，v)に依存する。そして、これらの影響は光路分布（〃_a=0)の中に繰り込まれている。つまり、散乱等の影響は、吸収がない仮想媒体に対する光路分布（ z _a=0)に反映されるため、以下では光路分布（〃_a=0)や z _a=0 である仮想媒体に対する応答を、単に散乱等の影響とよぶこともある。

F ) また、前記 B ) 項と関連するが、上記の光路分布（ _a=0 )は、複数あるいは多数光子（光子集合）入射や多数光子の移動に対して不変である。

また、媒体内の光子移動と減衰に関しては、以下の関係：

G ) 吸収係数が^ _aである均一媒体内を光路長 1(エル） = c tだけ移動する単一光子に対する減衰量（attenuation)は〃 _a l となるという関係（以下、この関係を「MBL第 1基本則」とよぶ）；及び

H ) 散乱媒体内の多数の光子（光子集合）に対する減衰量は、それぞれの光子に対する減衰量の重ね合せ（superposition)で表されるという関係；

が成立する。

2 . 不均一系内の光子移動に関する基礎事項

最初に不均一系内の光子移動に関する基礎事項を説明する。

2.1 不均一媒体内の光子移動

吸収（又は吸収成分）が不均一に分布する 3次元の不均一媒体全体を N個のボクセルに分割して番号 iを付け、それぞれのポクセル iについて非等方散乱係数 s、散乱角の余弦の平均値 g、吸収係数 _{a i}を定義する（図 1参照）。光子入射'検出位置は p(u，v)である。ここで、〃_sと _{a i}は不均一媒体固有の値で、かつ _sは定数、また^ _{a i}は位置の関数である。この際、ポクセル i は任意の大きさと形状にすることができ、大きい体積をもつポクセルではそのポクセル内の平均吸収係数や平均光路長を考えればよい。また、空間的に離れた 2個以上のポクセルを合わせて、 1個のポクセルと考えてもよい。いま、単一光子を繰り返して入射するとして、 m番目の入射光子（m- th光子）に対して時間 tに検出される出力]! _m(t)を考える。そして、この m- th光子のボクセル iにおける滞在（飛行）時間を t _im、また対応する光路長を

_m)とする。この場合、検出される光子の光路長 1 と l_imとの関係は、 l=ct= ∑l_imである。ここで、添字 mは m-th光子に関する量を表す。上記の系で m- th 光子の光路、つまり光路分布（〃_a=0)が決まると、 l_imが一意的に決まる。また、ポクセル番号 iを決めると、そのポクセル iの吸収係数/ z_aiがー意的に決まる。この場合、吸収係数は任意に選べるから、ポクセルの光路長と吸収係数は独立である。

時間 tに検出される m-th光子の減衰量（減光量ともよばれる） B_mは、

と表される。但し、光子が通過しなかったポクセル iに対して、 l_im=0 すなわち t_im=0を定義した（光路長分布がゼロ）。従って、（2.1.1)式の中の B_mには光子が通過しなかったポクセルの影響は含まれていない。また、 _{a i}と _mはともに非負の値をもっから、 B_m=∑(/z_{a i}l_im)=0 となるのは l_im≠0 であるポクセル i の吸収係数が// _{a i}=0のときに限られる。この際、 l_im=0 となるポクセルの影響は（2.1.1)式の中に含まれていないから、 ∑ (〃_{a i}l_im)=0 なる条件は∑ _{a i}=0と置き換えてよい。

この（2.1.1)式は、不均一系における MB L第 1基本則、つまり「吸収係数が〃_{a i}である散乱体の中をジグザグ光路に沿って長さ l=∑l_imだけ移動した光子の減衰量は、媒体の散乱特性とは無関係な値∑(〃_{a i}l_im)になる」ことを表す。この関係は、光子移動の際に吸収係数が変化する位置で MBL を接続したものである。また、この関係は、ジグザグ光路が交差（光子が 2 回以上同一場所を通過）する場合にも成立する。ここで、 ∑ /_{a i}=0のときの m-th光子に対する応答を s_m( z_s，t)とおけば、 m-th光子に対する応答 h_m(t)は、

In h_m (t) = In s_m (μ₅,ή-∑f_=l (μ_αϊί _im ) (2.1.2) となる。各ボクセルの光路長と吸収係数は独立であるから、 (2.1.3)

が得られる。但し、上式中の B_{i m}は、ポクセル i における減衰量である。この（2.1.3)式は、不均一系内の m-th光子に関する MB Lを表わす。

すると、（2.1.3)式から、不均一系内の m-th光子に対して、

In ^( = In ， ή-∑^ ^ I _imd _a = In s_w ，り-∑'"₌₁ {μ_αί£ _im ) (2.1.4)

が得られる。

不均一媒体に多数の光子（光子集合）を入射したときの応答を構成する光子で、光路長が 1となる光子を考える。この場合、光路長が 1となる光路は媒体内に無限といえるほど多数存在する。そして、多数の光子（光子集合）に対する減衰量は、単一光子に対する減衰量の重ね合せ（superposition)で表される。また、通常の応答、例えばインパルス応答には、光路長の異なる光子が含まれている。この場合にも、減衰量は単一光子に対する減衰量の重ね合せで表される。また、先に定義した光路分布（ _a=0)は、多数の光子入射に対して不変である。従つて、多数光子入射に対する不均一系の各種応答の減衰量は、単一光子に対する減衰量の重ね合せ、つまり線形問題として解析することができる。 2.2 不均一媒体の光路分布（ _a=0)

再入射不可能な任意の不均一媒体に対して光入射 ·検出位置 p (u, V) を定めたとき、光子が取り得る光路分布は一意的に決まる。この光路分布（〃_a=0) は、その不均一媒体に吸収がない（ _a=0)と仮定したとき（仮想媒体、 „=0) の光路分布で定義されるから、不均一媒体内の吸収やその分布には影響されない。しかし、光路分布（〃_a二 0)は、不均一媒体の光学特性、媒体形状や境界条件、光入射 ·検出位置 P(u，v)、及び光路長に依存する。

不均一系に対するインパルス応答、インパルス応答の時間分解ゲート積分信号、及び周波数領域の応答を構成する光子集合に対して、それぞれの光路分布（ /_a =0)を定義することができる。また、光路分布（ _a=0)が決まれば、各ボクセルの平均光路長（ z_a=0)がー意的に決まる。

以上のような光路分布（ /_a=0)は、モンテカルロシミュレーションなどによつて計算することができる。さらに、モンテカルロシミュレーションでは、光路分布（〃_a=0)を基準にして、任意の一様な吸収（例えば〃 _a = /i_{a v})がある場合の検出光の光路分布や平均光路長を計算することもできる。

2.3 検出光の光路分布

インパルス応答には、種々の光路長をもつ光子が含まれる。媒体に吸収があるとき、個々の光子の減衰量は B_m=∑(〃_{a i} l_{i m})になるが、光路長が異なる光子集合の平均光路長は、吸収係数の関数になる。従って、応答の時間波形は、吸収に依存して変化することになる。つまり、吸収があると光子移動の途中で光子が消滅し、その消滅の程度が光路長あるいは光路分布（〃_a=0) によって異なるため、その分だけ検出光に対する平均光路長や光路分布が変化する。この結果、検出光の光路分布は、通常、吸収に依存することになり、前出の光路分布（ _a=0)に等しくならない。

2.4 重み関数

不均一媒体内を 1だけ移動した単一光子の減衰量は、前述したように∑ ( _a i l_{i ra})で表される。そして、不均一媒体における多数光子の減衰量は、単一光子の減衰量の重ね合せで表される。ここで、多光子で構成される不均一媒体の各種応答の減衰量をポクセル i における寄与係数 \^とポクセル i の吸収係数//„ iとの積で表すことを考える。そして、前記の寄与係数を「重み関数（Weight Function) j とよぶ。

前出の（2.1.4)式と（2.1.5)式の場合には、光路長 1 =c t =∑ 1 _imが吸収係数 _{a i}の関数にならないから、光路分布（ _a二 0)と重み関数は等しい。しかし、ィンパルス応答の時間積分値のように平均光路長が吸収係数の関数になる場合、後述するように、重み関数は光路分布（〃_a=0)や検出光の光路分布に等しくない。また、一般に重み関数は吸収に依存する。

3. 不均一系における各種の応答関数

前述のことを踏まえて、不均一系における各種の応答を求めて、不均一系に対する MBLを導出する。この際、不均一媒体に対する光入射 ·検出位置を p (u, V) とする。すると、時間領域や周波数領域の各種計測で得られる各種応答に対して、それそれの光路分布（ _a=0)がー意的に決まる。そして、各種応答の減衰量は、各種応答を構成する単一光子に対する減衰量の重ね合せになる。なお、光路分布（ _a=0)は吸収分布に対して不変であるから、計測系の設定、つまり光入射 ·検出位置 p(u，v)や計測方式が同一であれば、被計測媒体が不均一媒体であっても均一媒体であっても、同一の光路分布（ _a=0)になる。 MBL に基づく記述'解析方法は、不均一系の応答の場合でも、 _a=0 としたときの応答（つまり散乱や境界条件の影響）と吸収による減衰量とを分離した形で記述することができるという点で優れている。

3.1 不均一系のインパルス応答

不均一系のインパルス応答 h(t)は、種々の光路長をもつ多数の検出光子で構成される。この場合、観測時間 tを決めると l=c tとして光路長 1が決まる。いま、光路長が 1である検出光子の数 Mが十分に大きいとする。この場合、全ての光子の光路長は l_m二 1であるから、 M個の光子に対する平均光路長も 1、つまり、

^ である。このような光路長が 1となる光子集合の光路分布（ _a=0)は、媒体の形状や散乱特性及び光入射 ·検出位置 p(u，v)に対して一意的に決まり、不均一媒体の吸収係数や吸収分布には関係しない。そこで、 z_a=0のときのインパルス応答に対して光路長が 1となる十分に大きい数 M個の光子（光子集合）を考え、この光子集合の光路分布（〃_a=0)を各ポクセル i(i=l〜N)の平均光路長（〃 0)で表す。ここで（3.1.1)式を参照すると、この平均光路長 l i に関して、

(³¹·²)

が成立する。つまり、多数光子に対する各ポクセルの平均光路長 1 iの全ボクセルに対する総和は多数光子に対する平均光路長 1に等しい。

また、インパルス応答を構成する光子（複数）の中で光路長が 1 (エル）となる光子集合に対する減衰量 B_hとポクセル iにおける減衰量 B_ihは、

B_ih = μ_αϊίι (3.1.4)

となる。従って、 ∑〃_{a i}=0のときのィンパルス応答を s(〃_s，t)とおけば、不均一系のィンパルス応答 h(t)は、

（3丄 5)

と表される。すると、不均一系のインパルス応答に対する MB L、つまり din (t) = _ δ¾_ ₌ _^ih_ = __g , (_{3 1 6})

dMai ^ ai ^d ai '

が得られる。

以上から、不均一系のインパルス応答 h(t)は、

（3丄7)

⁽³丄⁸⁾ となる。これらの式では、先に述べたように、 _a=0 としたときのインパルス応答、つまり媒体形状や散乱の影響を表す項 s( z _s，t )と、吸収に依存する減衰量を表す項 B_hとが分離されている。

なお、上記 (3.1.7)式から、不均一系のインパルス応答 h (t)に対する重み関数は、先に定義した光路分布（〃_a=0)で決まるポクセル i ( i=l〜N)の光路長 1 _a=0)に等しいことがわかる。従って、時間分解計測による光 CT画像再構成アルゴリズムは、比較的簡単になる。

以下では、上記で求めた不均一系のィンパルス応答 h (t)を構成する光子集合の光路分布（〃_a=0)、つまり前記で定義した各ポクセルの光路長

あるいはそれらの時間分布（異なる値をもつそれぞれの l(〃_a=0)に対する 1 i (； _a=0)であり、時間分解型の光路分布（〃_a=0 )とよぶ）を用いて各種の応答とそれらに対する重み関数を導出する。

3.2 不均一系のィンパルス応答の時間分解ゲート積分信号

ここでは、先ず不均一系のィンパルス応答の時間分解ゲー卜積分信号に対する解析式を導出し、減衰量や平均光路長の吸収依存性を考慮した重み関数を導出する。なお、時間域の積分範囲 [ t を [0，∞]にすると、定常（CW)光に対する応答を表すことになる。

3.2.1 時間分解ゲート積分信号の解析式

ィンパルス応答の時間分解ゲ一ト積分信号 I _τは、（3. 1.8)式を積分して

= I = I ² expl-

= t)_exp[—∑ ₌₁(^ ) （3.2.1) と表される。

この式から、

ニー _Γ〉 = - ）（3.2.2)

が得られる。ここで〈1 _{i T}〉=Li(〃_a)はインパルス応答の時間分解ゲート積分信号に対するポクセル iの平均光路長である。この Li( _a)は媒体の散乱特性や境界条件、及び，に依存する。また、後出の（3.2.9)式に示すように、 Li( /_a)は _{a i}に対する単調非増加関数である。さらに、上記の Li( /_a) はボクセル iの吸収係数/ _{a i}の関数であるが、ポクセル i を含めた全ボクセル

セルの吸収に対する依存性については後述する。

すると、（3.2.2)式から

(₃₂₃

が得られる。つまり、インパルス応答の時間分解ゲート積分信号 I _τに対する各ボクセルの平均光路長くパ ^の全ボクセルに対する総和はくこ等しい。

次に、（3.2.2)式を〃 _{a i}で積分して、時間分解ゲート積分信号 Ι_τに対する基本式を求めると、後で詳述するように、

1η/_Γ (μ_α)α_μα (3.2.4)

が得られる。右辺第 1項は積分定数であり、（3.2.1)式に∑〃_&1=0 を代入して得られる。この際、右辺第 1 項は吸収とは無関係、つまり散乱や境界の影響を表し、右辺第 2項は減衰量 Β_τを表す。ここで、上記までの結果をまとめると、

が得られる。

ここで、先に述べた L i ( _a )の全ボクセルの吸収係数に対する依存性が問題になる。つまり、この依存性があると、（3.2.5)式をこれ以上簡単にすることができない。この事実は重要であり、光 C T画像再構成アルゴリズムを複雑にしている。

ところで一般には、各ポクセル間の吸収係数 / _{a i}の差、つまり平均吸収係数からの偏差はそれほど大きくない場合が多い。この場合、注目しているポクセル i の平均光路長

に対する他部分（又は他のポクセル）の吸収係数偏差の影響は小さい。換言すれば、各ポクセル間の吸収係数の差が小さい場合、注目しているポクセル i の平均光路長 Li (〃_a )は、全ポクセルに対する吸収係数の平均値で大まかに決まり、他部分の吸収係数の平均吸収係数に対する偏差に依存するの変化は無視できる程度に小さく、 L i (〃_a )はボクセル i の吸収係数 z _{a i}のみの関数であると見なしてよい。この場合、ポクセル iの吸収係数の平均吸収係数との差（吸収係数偏差）に対応する（ _a )の差分は、線形近似して記述することができる。従って、以下では、 ( _a )はポクセル iの吸収係数〃 _{a i}のみの関数であると仮定する。つまり、ポクセル iの平均光路長（ _a )を考える際に、ポクセル i の吸収係数〃 _{a i}以外のパラメ一夕は一定（不変）であると仮定する。なお、 2 つ以上のポクセルの吸収係数が同時に変化する場合、これらの 2つのポクセルを合わせたボクセルを 1つの新しいポクセルと考える。このことは、「1つの測定値から複数ポクセルの吸収係数の変化を定量することはできない」ことと等価であり、理論上の矛盾はない。

以上の線形近似を適用すると、ポクセル iにおける減衰量 B _{i T}を定義すると d\nI_T 一 dB_T d

dM_ai μ_αί δμ_αί

0 CUai )φ_α = -^ = ,丁〉=一） (3-2-6) d _ai ^J0 一" μ_αί

が得られる。この式は、不均一系の時間分解ゲート積分信号 Ι_τを構成する光子集合に対する MBL を表す。但し、先に述べた線形近似を用いた。また、前述したように、積分範囲 ]を [0，∞]にすると、通常の時間分解積分信号、つまり定常（CW)光に対する MBLを表す。

次に、時間分解ゲート積分信号 Ι_τに対する重み関数について考える。先ず、平均値定理を用いて（3.2.4)式と（3.2.6)式の減衰量 Β_τと B_iTを変形すると、

Β_τ 。】 ⁷¹)

B_iT (3.2.7.2)

となる。但し、 _χ。は 0≤〃_x。≤Ad_{a i}なる適宜の値である。従って、 2.4節に記載した重み関数の定義から、時間分解ゲート積分信号 I _τに対するポクセル iの重み関数は。=Li( _x。）であることがわかる。この重み関数 ^。は、時間分解ゲート積分信号に対するポクセル i の平均光路長 Li(〃_ai)=<l_{i T}> とは異なる。

なお、前記 Li( i_a )の全ポクセルの吸収に対する依存性が問題になる場合には、（3.2.6)式の代わりに（3.2.5)式を使用する必要がある。

以上のことも、積分範囲 [tい ]を [0，∞]にすることによって、通常の時間分解積分信号、つまり定常（CW)光に対する応答に適用することができる。

[(3.2.4)式の導出]

前出の（3.2.2)式から、

1η/_Γ =一 ^(μα^μα +1η/_Γ|^ ₌₀ が得 ηられる。従って、 i = 1〜Nについての漸化式は、

- ln/_r =- ²L₂^ ?/_a+ln/_T

+ (Λ2)

to

となる。そして、上式の両辺の和を求めて整理すると、

^1η

3) となる。但し、右辺第 2項、第 3項は積分定数であり、前出の（3.2.1)式に∑ / _ai=0を代入して得られる。すると結局、 ln/_r )^ (3-2.4)

が得られる。

3.2.2 平均光路長の吸収依存性

各ポクセルの平均光路長は、そのポクセルの吸収に依存する。また、重み関数もそのポクセルの吸収に依存する。本節ではポクセル iの吸収係数〃 _{a i}が（ _a 1)に変化（厶 ^=11)したときの減衰量や平均光路長の変化を解析して、それらの変化量を表す解析式を導出する。この結果、光 CTの画像再構成に際して問題となる平均光路長と重み関数の差、及びそれらの吸収依存性が明らかになる。

先ず、（3.2.7.2)式に示したポクセル iの減衰量 B_iTを〃 _{a i}の周りにテイラ-展開すると、

となる。但し、 B_{i T} ^{( n}) ( _a)は 8^( ^)の ^に対する n階導関数である。また、（3.2.6)式の関係から、

）（3.2.9) が得られる。従って、（3.2.2)式を用ぃて6^ ( ² ) (〃_{a i})を求めると、

_D (2)₍ 、 _T {X)₍ 、

expト—∑₌₁ 。 )

=—〈½²〉 +〈丁〉² =—ひ, ( T) (3-2.10) となる。つまり、一 B_{i T} ( ² ) ( _{a i})二一 ( ¹ ) (〃_{a i})は l_iTの平均値まわりの分散を示す。従って、減衰量 B_{i T} ( _{a i})は、

Β_ίΤ{μ_αί +h)=B_iT{ _ai)+hL_{i ai})-—a_i ²[t_iT)+ (3.2.11) となる。

ここで、吸収係数が〃 _a iから（〃_a i+h)に変化したときの減衰量の差 ΔΒ_{ί Τ} を（3.2.7.2)式から求めると、

ΑΒ_ιΤ

(3.2.12) となる。この〃 _xは前出の _x。と異なり、 min i/ia a

C z _{a i},〃_{a i}+h)なる適宜の値である。また、この（〃は減衰量の差 ΔΒ_{ί Τ} に対する重み関数 ^である。ここで、（3.2.11)式と（3.2.12)式を比べると、

となる。他方、平均光路長 Li(〃_{a i})を〃 _{a i}の周りにティラー展開すると、 ,)+ ⁽²⁾ ,)+ ."）となる。但し、 (ⁿ) ( _a)は Li(//_a)の〃 &に対する n階導関数である。従って、（3.2.13)式及び (3.2.14)式から、

W, (/½)- ²( ) = ' + + , )】（32.15)

が得られる。なお、（3.2, 15)式の最右辺（第 4式）に示す近似式は、本発明者によって既に均一散乱媒体内の吸収成分の濃度定量に用いられ、その有効性が確認されている。

以上を踏まえてポクセル iの吸収係数〃 _{a i}が（〃_{a i}+h)に変化したとき、 (3.2.4)式を用いて変化前後の時間分解ゲート積分信号 1ηΙ_τの差 ΔΙ を求めると、

M = \ ！ = \^^{ai +h} L_t μ_α )άμ_α =hW_f =/( 。,)ーび,² _ίΓ) (3.2.16) となる。つまり、インパルス応答の時間分解ゲート積分信号 lnl _τの差 ΔΙ と、吸収係数が〃 _{a i}であるときの重み関数 Wiから、ポクセル iの吸収係数変化分 h=A〃_{a i}を求めることができる。この場合、重み関数は平均光路長 L パ /^)=<1₁ >及びその分散び₁ ²の関数である。

また、複数ポクセル i ( 1〜Nの中の複数の値）の吸収係数〃 _{a i}が（〃_{a i}+h Jに変化したとき、（3.2.4)式を用いて変化前後の時間分解ゲート積分信号 In

I _τの差 ΔΙを求めると、

となる。つまり、インパルス応答の時間分解ゲート積分信号 1ηΙ_τの差 ΔΙ と、吸収係数が〃 _{a i}であるときの重み関数 Wiから、複数のポクセル i(i=l〜N) の吸収係数変化分

〃_{a i}を求めることができる。但し、この場合には N個の連立方程式を解く必要があり、 N個の独立な計測値 ΔΙが必要である。また、前記の重み関数は平均光路長 Li(〃_{a i})=<l_iT>及びその分散び丄の関数である。

以上の場合、吸収係数 x _{a i}を被計測媒体の平均吸収係数あるいはそれに近い値に選んでおくと、精度が向上する。つまり、さきに述べたの全ポクセルの吸収に対する依存性の問題が大きく緩和される。

なお、各ポクセルの吸収係数 _{a i}の平均吸収係数に対する偏差が極端に大きい場合で、他のポクセルの吸収係数偏差に対する Li(〃_a)の依存性が問題になる場合には、（3.2.5)式に戻って再計算をするような繰り返し計算が必要になる。

以上のことは、積分範囲 [t を [0，∞]にすると、通常の時間分解積分信号、つまり定常（CW)光に対する応答に適用することができる。

3.2.3 アベレージ ·バリュ一 'メッツドへの応用

光 C Tは、不均一媒体中の吸収成分の濃度分布を計測しょうとするものである。この場合、先に述べた 2種類の平均光路長の吸収依存性、つまりの全ポクセルの吸収に対する依存性、及び Li(〃_a)のポクセル i 自身の吸収に対する依存性が問題を複雑にしている。本発明者は、これらの吸収依存性に基づく誤差を低減する一つの方法として、ァへ、、レ-シ、、'ハ、、リュ-づソ、ソに（Average Value Method, 以下「AVM」という）を開発した。この AVMは、計測すベき不均一媒体に対して近似的な平均吸収係数を推定あるいは計測して、この値を基準にして、媒体各部の吸収係数の偏差を定量する。この場合、上記で推定あるいは計測する平均吸収係数は、真値に等しい必要はない。この方法は、平均光路長の吸収係数依存性の問題を大きく緩和するもので、前記で推定あるいは計測した平均吸収係数が真値に近いほど、再構成光 CT画像の精度が向上する。

いま、被計測媒体の平均吸収係数に近い値の吸収係数^ _{a v}をもつ仮想媒体（ a= _av)を仮定する。すると計測対象である不均一媒体の吸収係数は〃 _ai二 a_V+A z_{a i}(i=l〜N)となる。時間分解ゲート積分計測を考えると、仮想媒体に対する計測値 I(〃_av)と実際の媒体に対する計測値 I (〃_{a i}= _av + 厶〃^)( 1〜? との関係は（3.2.4)式及び（3.2.15)式、あるいは（3.2.17)式から、

となる。つまり、仮想媒体と実際の媒体に対するインパルス応答の時間分解ゲート積分信号 1ηΙ_τの差 ΔΙ と、吸収係数が〃 _avであるときの重み関数から、複数のポクセル i(i=l〜N)の吸収係数偏差 A /_aiを求めることができる。但し、この場合には N個の連立方程式を解く必要があり、 N個の独立な計測値 ΔΙ が必要である。また、前記重み関数 Wiは、吸収係数が// _avであるときの平均光路長 L z aJ li とその分散び _{i v} ²との関数である。なお、上記のような一様な吸収がある媒体（仮想媒体）のポクセル i の平均光路長及び分散は、モンテカルロシミュレーションなどで計算することができる。前記 (3.2.18)式を行列式で表せば、

となる。但し、重み関数 Wiは、

W_{t +}Δμ_α,)₊Ι_ί{μ_αν)] (3.2.20)

である。

3.3 不均一系の周波数領域の応答

周波数領域の応答は、前出の（3.1.8)式に示すィンパルス応答 h(t)をフーリエ変換して、

=J₀ ^sK s » ^exp (- *)exp (- ]ωί)ώ = (δ,ω) + jX(b,a)

= 4， ) exp [—ゾ^^)】 (3.3.1.1)

Σ'^ ） (3.3.1.2) となる。ここで、 R, X, Α, øはそれぞれ、

R = /₀ ^a> Cu_J,/)exp[-∑¾_{1 a/}c/_/)] co mtdt (3.3.2.1)

A = {R²+X² (3.3.2.3)

-iX

φ = -tan — (3.3.2.4)

^Ψ R

である。すなわち、 R及び Xはそれぞれ実部及び虚部であり、 A及び øはそれぞれ振幅及び位相遅れであり、これらは口ックインアンプなどで容易に計測することができる。

すると、後で詳述するように、 dR dXi

- c (3.3.3.1)

8 Mai δω dx dRi

=一 c (3.3.3.2)

dm

なる関係が得られる。ここで、各式の右辺の添字 iが付いたパラメ一夕は、ポクセル iにおける量であることを示し、前節に述べた平均光路長に相当する。また、ここでは、後で詳述するように、

- i₀" (^ exp[— Σ'^Ο^ )] smwtdt (3.3.4.1)

δΐη^' 1 ,dRi , _vdXi

R—± + X- (3-3.4.3)

~ A¹ δω οω

を定義した。以上では、

d\ A vN 31η 4 (3.3.51)

(3-3.5.2)

の関係が成り立つ。

従って、（3.3.3.1)式ないし（3.3.3.4)式の関係式から、前出の（3.2.4)式の導出と同様にして、

in 1η _{ι ½ )= 0})__{c∑ i}^ (³·³.⁶,³)

δω ^ φ。 (3.3.6.4)

が得られる。右辺第 1項は積分定数であり、それそれ、左辺に示すパラメータに ∑ _{a i}=0を代入して得られる。

以上の（3.3.6.1)式ないし（3.3.6.4)式は、前節に記載した時間分解ゲート積分信号に対する式と相似であるから、前節と同様にして、これらの式を変形することができる。そこで以下では、（3.3.6.3)式を例として取り上げる。先ず、前節と同様に、周波数領域の応答に対する重み関数を考える、先ず、平均値定理を用いると（3.3.6.3)式の減衰項（右辺第 2項）は、

となる。こで、

(3.3.7.3)

d ω

は検出光に対するポクセル iの重み関数（cx群遅延）、つまり周波数領域の応答の重み関数を表す。但し、〃_xは 0≤〃_x≤ _{a i}なる適宜の値である。この重み関数（3.3.7.3)は、前節に述べた重み関数 Wi (

〃_x)に相当し、ポクセル iの平均群遅延 3 _ί

(3.3.7.4)

d ω

の c倍の値

3

c- (3.3.7.5)

όω

とは異なる値をもつ。

さらに、ポクセル iの減衰 ≠量 B i _f及び平均群遅延（3.3.7.4)の吸収係数依存性は、前節に述べた減衰量 B_{i T}及び平均光路長 Ι^(〃_χ)と同様になる。つまり、群遅延（3.3.7.4)の吸収依存性を考えると、（3.3.3.4)式及び（3.3.5.2)式よ、

In Ai_ d δφ d δ d _N δφ_ί = d δφ_ί (₃ , ₈) ^_ω ² ' δωδμ_αί μ_αί δω δμ_αί ^ί=1 δω d μ_αί δω となる。但し、異なるポクセル間において（i≠ jのとき）、

=ο (3.3.9)

^Mai

を仮定した。つまり、平均群遅延（3.3.7.4)はポクセル i以外のポクセルの減群 ; に依存しないとした。これらは、前節における仮定と同様である。なお、上式中の d InA

όω

は、 3種の変調周波数を利用して計測することができる <

以上から、前述の説明（3.2.3節）と同様にして、 h ^c

c (3.3.10)

i —び/ =一

δω δω Ma 2 Mai

σ f (3.3.11)

が得られる。この際、一様な吸収がある媒体（仮想媒体）（^a ^ Z av )のポクセル i の平均群遅延及び分散は、モンテカルロシミュレーションなどで計算することができる。従って、 _{a i}が既知であるとき、吸収係数が/ Ζ _χであるときの平均群遅延 d ω β、

を計算で求めることができる。

[(3.3.3.1)式〜（3.3.3.4)式の導出]

先ず、 t =∑ t iであることを考慮して、（3.3.2.1)式及び（3.3.2.2)式から tdt

( 1)

dX = ~V^tsi , exp -∑_/=1 Maj tj ) cosotdt

COSG) di (ΒΛ.2)

が得られる。そこで、ョ- ( 5' ^exP -Hi- aiCti) (5-2.1)

を定義するの場合、

ex一 _VN dxj

(β.3.2)

dm "ⁱ⁼⁼i d

の関係が成り立つ。また、（3.3.2.1)式及び (3.3.2.2)式から、 dR

= -c ij ί (^，り exp| -∑;: a^) ]∞_&ωίώ .4.1) αι

が得られる。以上から、 dR dX

- c (3.3.3.1)

δω dX dR

= -c- (3.3.3.2)

dm

が得られる。

さらに、（3.3.2.3)式及び上式から、

が得られる。但し、こでは、

を定義した。この場合、上式に（Β.3.1)式及び (Β.3.2)式を代入すると.

A^A δω όω

2 d X 2 . d . δφ

-cos φ = cos φ tan^ = ( .7) δω R οω δω

の関係が得られ、 δ 一 N

(5-8) όω― / =1

となる。以上により、（3.3.3.3)式が得られる。

また、同様にして、（3.3.2.4)式、 (3.3.3.1)式及び（3.3.3.2)式から

が得られる。但し、ここでは、 dlnAj dRi _+xdX_f

R ( 10) δω A δω δω

を定義した。従って、 δώ cln

—^J— = c ^L (3.3.3.4) δω

が得られる。

なお、上記のとき、

ΛΓ OUlAj ΣΝ dRi dXi

R + X

¹ δω δω δω

+ 、

A' 厶 ¹ δω ム' ⁼¹ δω

dR ^δΧλ din A

R— + X

A' δω σω が成立する。従って、

の関係が成り立つ。

4. 不均一系に関する一般式

前章までの説明から明らかなように、不均一系における MB Lは、微分形で表すと有効である。つまり、不均一系に対するインパルス応答 h(t)の減衰量 B_h を表す（3.1.6)式、インパルス応答の時間分解ゲート積分信号 I_T( ^_S、 //_a)の減衰量 Β_τを表す（3.2.6)式、及び周波数領域の計測で検出される振幅 Aの減衰量 B_f を表す（3.3.3.3)式等は、微分形で表され、かつ、それらの式は相似である。

これらの応答（検出光の所定パラメータ）を Y、減衰量を Β、各ポクセルの平均光路長を Ziとおけば、上記の各種応答は一般式、 d\nY dB

= -^i\ ai) (4-1)

δμ αι ομ αι

で統一して表すことができる。

上記で、対象とする応答が光路長の異なる複数光子からなる場合、平均光路長 Ζは〃 _{a i}の関数となる。また、対象とする応答が等しい光路長 1 をもつ光子からなる場合、平均光路長は光路分布（〃_a=0)で決まる

となり、この liは〃 _{a i}に依存しない。以上から明らかなように、（4.1)式は不均一系に対する MB Lの一般式を表す。すなわち、この（4.1)式は、不均一媒体における減衰量が、各ポクセルの吸収係数と平均光路長のみによって決まることを表わしている。従って、不均一系に対する MB Lはこの事実を記述したものと考えてよい。

なお、 ₃=0のときの応答を S( /_s)とおけば、上記の（4.1)式から、

が得られる。但し、 Wiは各ポクセルの重み関数である。

ここで、不均一媒体の平均吸収係数に近い値を基準にして、媒体各部の吸収係数分布を定量する A VMによる光 CT画像再構成との関係を説明する。この場合、不均一媒体の平均吸収係数に近い値 ^_avが既知であるとし、ポクセル i の吸収係数を〃 _av+A _{a i}とする。また、 3.2.3節と同様にして、吸収係数の偏差 △〃_aiの範囲内で、 Zi ( _{a i} )が〃 _{a i}のみの関数であると仮定すれば、 (3.2.18)式と同様にして、

ΔΚ = Inノ _Λ) 、 =∑^ Z_t {μ_α μ_α

—∑ jA "。' ∑ J Α"_σ/ί - A "。' v²〕] (4.4) が得られる。但し、び i _v ²は（〃 )の分散を表す。また、（4.4)式を行列で表示すれば、

となる。但し、重み関数 Wiは、

と表される, なお、一様な吸収がある媒体（仮想媒体）のポクセル iの平均光路長及び分散は、モンテカルロシミュレーションなどで計算することができる。また、先に述ベたような Z i ( _a i )の全ポクセルの吸収に対する依存性が問題になる場合には、次に示す（4.7)式に戻って再計算するような繰り返し計算が必要である。

3 lnF dB _{ίΑ Ί},

~" (⁴·⁷)

(第 1実施形態）

図 2〜図 6を参照して本発明の好適な一実施形態である第 1実施形態について説明する。図 2には、散乱吸収体である被計測媒体 1 0の内部に分布する吸収成分の濃度分布を計測する内部特性分布計測装置（光 C T画像計測装置） 1が示してある。先ず、内部特性分布計測装置 1の構成について説明する。

図 2に示す装置は、光入射ファイバー 2 1を備えており、光入射ファイバー 2

1の先端が被計測媒体 1 0の表面の異なる位置 u j ( j=l〜p)に移動可能となつている。そして、光入射フアイ Γ- 21 には波長選択器 22を介して光源 23が光学的に接続されており、光源 23 から発せられた光線は波長選択器 22 で波長選択され、光入射フアイハ"- 21 を通して異なる位置 U j (j=l~p)から順次被計測媒体 10中に入射される。光源 23には、発光タ"ィォ -に、レ- Γ-タ"ィォ -に、 He-Ne レ- -等種々のものが使用できる。また、光源 23は、パルス光線や方形波光線、又はそれらの変調光線を発生するものでもよい。本実施形態において使用する光源 23 は単波長の光線を発生するものであってもよいが、 2波長以上の光線を発生可能なものであってもよい。

また、図 2に示す装置は、複数の光検出ファイバー 2 5を備えており、光検出ファイバー 2 5の先端が被計測媒体 1 0の表面の複数の位置 v _k (k=l〜q)にそれそれ設置されている。そして、光検出フアイ; - 25には光検出器 26がそれぞれ光学的に接続されており、被計測媒体 10中で散乱されつつ透過した光線は光検出フアイ; - 25 を介して光検出器 26に導かれ、光検出器 26で受光信号を検出信号（電気信号）に変換して増幅し、それぞれに対応する検出信号が出力される。光検出器 26 は光電子増倍管のほか、光電管、フォトイオ-ド、ァ Γ、ランシェフオト夕、、ィオ-ド、 PINフォトタ、'ィオ-ド等、あらゆる種類の光検出器を使用することができる。光検出器 26の選択に際しては、使用される測定光（測定光線）の波長の光が検出できる分光感度特性をもっていれば良い。また、光信号が微弱であるときは高感度あるいは高利得の光検出器を使用することが好ましい。

なお、被計測媒体 1 0の表面における光入射ファイバー 2 1による光入射面及び光検出ファイバー 2 5による光検出面以外の場所は、光を吸収あるいは遮光する構造にしておくことが望ましい。また、被計測媒体 1 0の内部を拡散伝搬した光線が複数の波長の光線を含む場合には、光検出器 2 6と光検出ファイバー 2 5 との間に波長選択フィルタ（図示せず）を適宜配置してもよい。

光源 2 3及び光検出器 2 6には C P U (制御及び演算処理部） 3 0が電気的に接続されており、光入射に同期した光検出のタイミング等が C P U 3 0によって制御されると共に、光検出器 2 6から出力された検出信号は C P U 3 0に導かれる。また、複数の波長を有する測定光線を使用する場合は、入射光線の波長も C P U 3 0によって制御される。具体的な手法としては、異なる波長の光線を時分割で入射させて使用する手法と、異なる波長の光線を同時に含む光を使用する手法とがある。具体的な波長選択手段としては、ミラ一を用いた光ビーム切り換え器、フィルターを用いた波長切り換え器、光スィッチを用いた光切り換え器等がある。

上記の光入射ファイバ一 2 1、波長選択器 2 2、光源 2 3及び C P U 3 0が本発明にかかる光入射手段を構成し、上記の光検出ファイバ一 2 5、光検出器 2 6 及び C P U 3 0が本発明にかかる光検出手段を構成する。

図 2に示す内部特性分布計測装置 1は、更に、オペレーティングシステム（0

S ) 4 1及び後で詳述する内部特性分布計測プログラム 4 2が記憶された第 1記憶装置 4 0と、各種ファイルが記憶された第 2記憶装置 5 0と、得られた内部特性分布を示す光 C T画像デ一夕を記憶する画像メモリ 6 1と、作業用データを一時的に記憶する作業用メモリ 6 2と、データの入力を受け付けるキーボード 7 1 及びマウス 7 2を備える入力装置 7 0と、得られたデータを出力するディスプレィ 8 1及びプリン夕 8 2を備える出力装置 8 0とを備えており、これらも電気的に接続されている C P U 3 0によって制御される。なお、上記の記憶装置及びメモリは、コンピュータの内部メモリ（ハードディスク）であっても、フレキシブルディスクであってもよい。

第 2記憶装置 5 0は、基本光路長分布データファイル 5 1と、測定データファィル 5 2と、入力データファイル 5 3と、吸収係数データファイル 5 4と、吸収成分濃度データファイル 5 5とを備えている。

この内、基本光路長分布データファイル 5 1には、前もって用意された基本光路長分布（種々の計測に使う種々の重み関数演算の元になる共通の光路長分布）が格納されている。このような基本光路長分布は、既知の手法、例えばモンテ力ルロシミュレーションや光拡散方程式を利用して計算することができ、例えば (17) J. Haselgrove, J. Leigh, C. Yee, N-G, Wang, M. Maris and B. Chance： Proc. SPIE, Vol . 1431 , 30-41 (1991)； (18) J. C. Schotland, J. C. Haselgrove and J. S. Leigh: Appl . Opt. 32, 448-453 (1993) ； ( 19) Y. Tsuchiya, K. Ohta and T. Urakami： Jpn. J. Appl . Opt. 34, 2495-2501 (1995)； (20) S. R. Arrige : Appl . Opt. 34, 7395-7409 (1995) に記載されている。

また、入力装置 7 0を用いて予め設定された計測条件や既知値が入力され、このような入力デ一夕が入力データファイル 5 3に格納されている。このような入力データとしては、予め任意に設定されたポクセル（体積素）の数、形状及び大きさ、被計測媒体の形状、光入射位置、光検出位置、散乱係数、平均吸収係数、計測に用いる光線の波長、計測の種類（時間分解積分計測、時間分解ゲート積分計測、位相変調計測など）、計測対象となる吸収成分の所定の波長における消光比などがある。

更に、測定データファイル 5 2には、内部特性分布計測プログラム 4 2の実行過程において光検出器 2 6から発せられた検出信号に基づいて求められる所定パラメ一夕の測定値が光入射位置と光検出位置との組み合わせ毎に記憶される。また、吸収係数データファイル 5 4及び吸収成分濃度データファイル 5 5には、内部特性分布計測プログラム 4 2の実行によって得られた吸収係数データ及び吸収成分濃度データが記憶される。

上記の C P U 3 0、第 1記憶装置 4 0及び第 2記憶装置 5 0が、本発明にかかる測定値取得手段と、基準値設定手段と、推定値算出手段と、重み関数演算手段と、吸収係数偏差算出手段と、吸収係数絶対値算出手段と、濃度算出手段とを構成し、上記の出力装置 8 0が画像表示手段を構成する。これらの本発明にかかる諸手段については、図 3及び図 4に示す本発明の方法の一実施形態のフローチヤ —ト（図 2に示す内部特性分布計測プログラム 4 2の処理を示すフローチヤ一ト）に基づいて以下に詳細に説明する。

1 )図 3に示すフローチャートにおいては、先ず、光源 2 3で生成した光線（光ビーム）を光入射ファイバー 2 1を介して被計測媒体 1 0の異なる光入射位置 u j ^：^〜パこ順次入射し ^ 、被計測媒体内部で散乱されつつ透過した各光線（各光入射位置 Ujから入射された光線の一部）をそれぞれ複数の光検出位置 v_k (k二 l〜q)に設置した光検出フアイ Γ'- 25 を介して光検出器 26で検出する（S102)。

そして、各光入射位置 U jから入射された光線に対して、複数の光検出位置 v_k (k=l〜q)で検出された光線に対応する光検出信号が光検出器 26 から発せられ、 CPU30において検出光線の所定パラメ一夕の測定値 Y_nに変換される。このようにして、光入射位置 U j ( j=l〜p)と光検出位置 v_k (k=l〜q)との組み合わせ P_n毎に検出光線の所定パラメ一夕の測定値 Y_nが得られ、測定デ一夕ファイル 52 に η セヅト（n=p x q)の測定デ一夕セヅトとして記憶される (S103 )。この際、 nセットのデータに対する光入射'検出位置の組み合わせ P _nはそれそれ異なる。

本発明にかかる所定パラメータの測定値 Y _nとしては、測定光線の測定対象物内部での散乱及び吸収に関係する所定パラメ一夕の測定値が好ましく、例えば検出光線の光量、位相差（又は位相遅れ）、振幅、時間分解波形等のパラメータの測定値が好適に用いられる。また、 CPU30において、光源 23からの光線の発生に同期した信号を利用して光検出信号に対する所定時間域での積分演算を行えば時間分解ゲ一ト積分計測が可能となり、他方、積分範囲を 0→∞とすれば通常の時間分解積分計測が可能となる。なお、パルス光線等を利用する場合には、この同期信号を省略することができる。また、 CPU30 において、測定値を平均化フィルタリンク"や最小二乗フィ、ソテインク等を利用して修正してもよい。

2 ) 次に、被計測媒体 1 0の吸収係数の基準値 z _{a v}及び輸送散乱係数の基準値〃' _sを設定する（S104)。なお、吸収係数の基準値 _{a v}は、仮想媒体に対する真の平均吸収係数に等しくする必要はなく、 z _{a v}=0 と近似してもよいが、一般的には^ _{a v}として真の平均吸収係数に近い値を用いる方が高い計測精度が得られる傾向にある。従って、吸収係数の基準値 _{a v}として、被計測媒体 10 を巨視的に見たとき、すなわち被計測媒体 10 が全体的に一様な吸収係数を有するとみなした場合の平均的な吸収係数（近似値）を採用することが好ましい。また同様に、輸送散乱係数の基準値〃' _sとしても、被計測媒体 10が全体的に一様な吸収係数を有するとみなした場合の平均的な輸送散乱係数（近似値）を採用することが好ましい。

このような吸収係数〃 _{a v}と輸送散乱係数/ ' _sは、前記の複数の光検出位置 v_kで検出された光検出信号（又は所定パラメ一夕の測定値 Y_n )から求めることができる。例えば、前記光検出信号がインパルスレスポンス（光検出信号の時間波形に対して十分に短いと見なせるパルス光入射に対する光検出信号）である場合には、その時間波形を光拡散方程式にフイツティンク"することで、所定の光入射'検出位置の組み合わせ P_nに対応する媒体内部の吸収係数 / _a _vと輸送散乱係数 ' _sを求めることができる。従って、各々の光入射'検出位置の組み合わせ P_nに対して求めた前記の吸収係数〃 _avの平均値及び輸送散乱係数/ _sの平均値を採用すればよい。さらに、これらの値は、時間分解ゲート積分計測や位相変調計測に利用することができる。また、このような吸収係数〃 _av及び輸送散乱係数〃' _sの取得方法に関しては、例えば

(21) R. Berg, S. Andersson-Engels, 0. Jar 1 man and S. Svanberg: Proc. SPIE, Vol. 1431， 110-119 (1991)； (22) M. Miwa, Y. Ueda andB. Chance: Proc. SPIE, Vol. 2389, 142-149 (1995)； (23) R. Cubeddu, A. Pifferi, P. Taroni, A. Torricelli and G. Valentini, Phy. Med. Biol. 42, 1971-1979 (1997) に記載されている。

また、上記基準値を設定する際に、被計測媒体 10を巨視的に見たときの平均的な他の光学定数を更に求めてもよく、このような光学定数としては、被計測媒体の屈折率 n_e、散乱係数 _s、散乱角の余弦の平均値 g が挙げられる。なお、被計測媒体の屈折率 n_eは、通常は水の屈折率で近似してよい。

3) 次に、被計測媒体 10が全体的に一様な吸収係数の基準値 _avを有するとみなして、吸収係数の基準値 _{a v}に基づいて、光入射位置と光検出位置との組み合わせ P_n毎に前記の所定パラメータの推定値 Y_avnを求める（S105)。すなわち、例えば、上記の平均的な光学定数（吸収係数の基準値〃 _av、輸送散乱係数の基準値 ' _s等）を与えて光拡散方程式を解けば、各光入射 ·検出位置の組み合わせ P_nにおける前記の所定パラメータ（光検出信号）の推定値 Y_avnが求められる。また、任意の光学定数（例えば、吸収係数）を与えたモンテカルロ計算の結果から、光学定数を変えたときの結果を推定する方法が開発されている（例えば、 A. Kienle and M. S. Patterson: Phys. Med. Biol. 41, 2221-2227 (1996)に記載されている）。従って、このようなモンテカルロ計算の結果から、前記吸収係数の基準値^ _{a v}を用いて各光入射'検出位置の組み合わせ P_nにおける前記の所定パラメ一夕（光検出信号）の推定値 Y_avnを計算することができる。

4) 他方、被計測媒体 10が全体的に一様な吸収係数の基準値〃 _avを有するとみなして、吸収係数の基準値// _avに基づいて、複数のポクセルに分割された被計測媒体 10の各ポクセル i における平均光路長

を求める (S106)。

なお、本発明では、個々の計測に際してその状況に合致した重み関数を後述の方法で求める必要がある。この場合、計測のリアルタイム性を考慮して、重み関数をその場で（計測時に）迅速に求める必要があるため、種々の計測に使う種々の重み関数演算の元になる共通の基本光路長分布を前もって用意しておき、計測時に状況に応じてその基本光路長分布を変換することが好ましく、それによつて必要な重み関数を迅速に求めることが可能となる傾向にある。

従って、基本光路長分布データファイル 51からモンテカルロ計算等によって予め求められている前記基本光路長分布を読み出し、前記の吸収係数の基準値 _av及び輸送散乱係数の基準値〃' _sに基づいて、これらの基準値を有する仮想媒体に対する検出光線の光路長分布を求めることができる。これが基本光路長分布の変換であり、例えば（24) A. Kienle andM. Patterson: Phys. Med. Biol. 41, 2221-2227 (1996) に示されている方法を利用することができる。なお、ここで求められた仮想媒体に対する検出光線の光路長分布 Zi( _av) に基づいて、前記の各光入射'検出位置の組み合わせ P_nにおける所定パラメ一夕（光検出信号）の推定値 Y_avnを計算することもできる。このような場合における処理の流れを図 3中に鎖線 Aで示す。

5) 次に、上記のように吸収係数の基準値〃 _avを用いて得られた各ボクセル i における平均光路長

に基づいて、各ポクセル i における重み関数 Wiをマイク Π·へ、、ァ'ラン; Γ-ト則に則して（4.6)式にしたがって求める（S 107) と共に、前記の各光入射 ·検出位置の組み合わせ P_nにおける所定パラメータの測定値 Y_n及び推定値 Y_avnと前記重み関数 Wiとに基づいて、各ボタセル iにおける吸収係数の基準値〃 _avに対する吸収係数の偏差 A _aiを（4.5)式にしたがって算出する（S 108)。

なお、求めるべき重み関数 Wiは（4.6)式に示したもので、

である。ここで、み関数 Wiはそのポクセルのもつ吸収係数偏差 A〃_ai=〃 _{a i}- /_av(i=l〜！の関数になっていることに注意する必要がある。従って、このような重み関数 Wiは、以下のような繰り返し演算によって吸収係数偏差 A〃_{a i}と共に求めることが好ましい。すなわち、 S107及び S108 の演算処理（図 3 中に鎖線 B で示す）は、通常、繰り返し演算されることとなる。図 4は、このような繰り返し演算によって重み関数 \^及び吸収係数偏差 Δ a_iを求める処理の好適な一実施形態を示すフローチャートである。以下、図 4に示すフローチャートについて詳述する。

先ず、次数 Kの初期値を K= l、吸収係数偏差 A^_{a i}の初期値を Δ _aiK=0 と設定する（S201)。次いで、各ポクセル i における平均光路長 Zi ( _av)及び吸収係数偏差 A〃_{a i K}に基づいて（4.6)式にしたがって重み関数 W_{i K}を求める（S202)。なお、（4.6)式中のび^ ²は

の分散を表し、例えば (3.2.10) 式や（3.3.8) 式を計算して求められる。

次に、前記の所定パラメータの測定値 Y_n及び推定値 Y_avnと前記重み関数 W_iKとに基づいて、各ポクセル i における吸収係数の偏差 A /_ai (_{κ+ 1}) を (4.5)式にしたがって算出する（S 203)。なお、前記 Y_nと Y_avnとの比の自然対数値 [In (Y_avn/Y_n)]が（4.5)式中の厶¥ に相当する。

そして、 S 203によって求められた吸収係数の偏差 (_κ+ と S202 で用いた吸収係数偏差 A〃_{a i K}との残差 E を算出し（S204)、その値が所定値（err) 以下であるか否かを判断する（S 205)。

その結果、上記残差 Eが所定値以下である場合はその時の重み関数 W_{i K}を重み関数 Wiとして、吸収係数偏差 A〃_{a i (κ+}ΐ)を吸収係数偏差 A _{a i}としてそれぞれ出力する（S206)。他方、上記残差 Eが所定値以下でない場合は、 S202に用いる吸収係数偏差 A _{a i K}を前記 S203によって求められた吸収係数偏差 A〃_{a i (κ+ 1 )}に置き換えて（次数 Κに 1を付加して： S207)、重み関数 W_iKの演算（S202)、吸収係数偏差 A〃_{a i} (_{κ+ 1})の演算（S203)及び前記残差 Ε の算出（S204)を再度行なう。そして、得られた残差 Εが所定値（err) 以下であるか否かを判断し（S 205 )、残差 Eが所定値以下となるまで上記の S 20 7→S 202→S 203→S 204→S 205を繰り返す。なお、上記の S 20 1を初期値設定手段、 S 202を重み関数演算手段、 S 203を吸収係数偏差演算手段、 S 204、 S 205及び S 206を評価手段、 S 207を吸収係数偏差修正手段でそれぞれ行なってもよい。

6) 以上によって、吸収係数偏差の分布（各ポクセル iにおける吸収係数の偏差） ^=〃^ -〃 _av(i=l〜N)が得られる。そして、 A〃_{a i}は先に述べた吸収係数の基準値（好ましくは平均吸収係数）〃 _avからの偏差を示す。従つて、各ポクセル iにおける吸収係数の偏差 Δ _{a i}及び吸収係数の基準値// _a _vを用いて、次式：

a_i ^{= A}〃_{a i}+ _av(ⁱ⁼¹〜^N)

にしたがって各ポクセル iにおける吸収係数の絶対値 _{a i}を算出することができ（S109)、吸収係数データファイル 54に記憶される。図 5に、吸収係数の基準値〃 _{a v}と偏差 A〃_{a i}と絶対値〃 _{a i}との関係を示す。そして、このようにして求められた各ポクセル i における吸収係数の絶対値〃 _{a i}に基づいて被計測媒体内部における吸収係数絶対値に関する分布が求められ、その分布を示す光 CT画像が出力装置 80に表示される（S110)。

7) 次に、吸収成分の濃度演算について説明する。ここでは、被計測媒体が一種類の吸収成分を含む場合を考える。すると、その吸収成分のポクセル iにおける濃度偏差 A C iは、ベア 'ランパート則から導かれる次式：

で与えられる。但し、 £は吸収成分の単位濃度当たりの吸収係数（又は消光係数）であり、分光光度計で測定することができる。さらに、前記吸収成分のボタセル iにおける濃度 C iは、ベア ·ランバート則から導かれる次式：

£

_{a i} (C.3)

で与えられる。従って、吸収成分の既知の吸収係数を用いて上記の各ポクセル i における吸収係数の絶対値^ _{a i}から各ポクセル i における吸収成分の濃度を算出することができ（Sl l l )、吸収成分濃度データファイル 55 に記憶される。そして、このようにして求められた各ポクセル iにおける吸収成分の濃度に基づいて被計測媒体内部における吸収成分濃度に関する分布が求められ、その分布を示す光 CT画像が出力装置 80に表示される（S112)。

(第 2実施形態）

上記の第 1実施形態の方法において 2種類以上の波長の光線を用いれば、 2種類以上の吸収成分の濃度偏差や濃度絶対値が測定できる。本実施形態では、波長え iとえ ₂の 2つの光線で計測する 2波長分光計測法について説明する。すなわち、波長人ェ及び人 ₂の光線に対応する各吸収成分の単位濃度当たりの吸収係数（又は消光係数） £ i及び ε ₂が既知であれば、前記（C .2)式及び (C. 3 ) 式に対して 2本の式が成立する。従って、これらの式を連立させることによって、 2種類の吸収成分の濃度偏差や濃度絶対値を測定することができる。

より具体的には、散乱吸収体が少なくとも 2つの吸収成分（例えば酸素化型及び脱酸素化型ヘモグロビン）を含有する場合は、それらの吸収成分に対する吸収係数が互いに相違する少なくとも 2つの波長を有する測定光線を使用し、各波長を有する測定光線のそれぞれに関して測定値 Υ _η、推定値 Y_{a v n}、基準値 / _{a v} 及び重み関数を求め、それらに基づいて各波長を有する測定光線のそれそれに関して吸収係数偏差 A^_{a i}、吸収係数絶対値 _{a i}を求めることによつて、各吸収成分の濃度の分布が求められる。

以下に、上記 2波長分光法を利用したヘモグロビンの濃度の計測について説明する。哺乳類の脳における吸収成分の主なものは、水、チトクローム（cytochrom)、酸素化型及び脱酸素化型ヘモグロビンである。近赤外線領域での水及びチトクロームの吸収は、酸素化型及び脱酸素化型ヘモグロビンに対して、ほぼ無視することができる程度に少ない。また、酸素化型及び脱酸素化型ヘモグロビンは、図 6 に示すように、吸収スぺクトルが異なる。さらに、頭蓋骨は、近赤外線に対して、散乱体と考えてよい。

いま前節までに述べた方法で波長え iと人 ₂の 2種の波長の光線に対して吸収係数// _a iと〃 _{a 2}が求められたとすれば、ランバー卜 'ベール（Lambert- Beer) 則によって、次式が成立する。

^ a i ^^ H b , i[Hb]+£_Hb0; , 〔HbO〕

〃_{a2 Hb}, ₂[Hb]+£_{Hb0> 2} 〔HbO〕

但し、

£_Hb, 脱酸素化型ヘモグロビンの波長え iに対するモル吸収係数〔mm — 1._M— 1 ]

_{H b0}, ！；酸素化型ヘモグロビンの波長え iに対するモル吸収係数〔mm 1 ·Μ— ¹ ]

£_Hb, ₂；脱酸素化型ヘモグロビンの波長え ₂に対するモル吸収係数〔mm

H h 0 2 ；酸素化型へモグロビ >~ン zの波長え 2₉に対するモル吸収係数〔mm

- ¹.M- ¹ ]

〔Hb〕：脱酸素化型ヘモグロビンのモル濃度〔M〕

〔HbO〕：酸素化型ヘモグロビンのモル濃度〔M〕

である。従って、既知のパラメ一夕 £_{H b}, い £ _{H b0}, い _Hb, ₂、 _Hb0, ₂及び計測値から演算された _{a l}と〃 _{a 2}から、脱酸素化型ヘモグロビンのモル濃度 [Hb]、及び酸素化型ヘモグロビンのモル濃度 [HbO]を求めることができる。

また、上記に対してチトクロームを考慮する場合のように、吸収スペクトルが既知である 3成分のそれぞれの濃度の定量は、 3波長以上の光線を使用すればよレ、。一般的には、吸収スペクトルが既知である n個の成分の濃度の定量計測は、 n個又は（n+ 1) 個の波長に対する吸収係数の計測値から、上記と同様にして求めることができる。

さらに、飽和度 Yは、

Y= 〔Hb〇〕 / (〔Hb〕 + 〔HbO〕）

であるから、

〃a l //" a 2 ⁼ [ ^£ H b , 1 °' H b O , 1一 ^£ H b , i )]丁

[ ^£ H b , 2 +Y( ^£ H b 0 , 2 ~ ^£ H b , 2 ^

を用いて、既知のパラメ一タ £_{H b}, い e_{H b0}, い £_Hb, ₂、 £ _Hb0, ₂と計測値から演算された _{a l}及び〃 _{a 2}とから、飽和度 Yが容易に算出される。

以上の方法では、本発明によって各波長の光線に対する吸収係数/ i_{a l}と _a2 が高精度に求められるため、各濃度も高精度に求められる。なお、酸素化型及び脱酸素化型ヘモグロビンに対して吸収が同一値になる波長（800皿、 isosbestic wavelength)を使用すれば上記の式はさらに簡単になる。

(第 3実施形態）

本実施形態は本発明を位相変調計測に応用する例を示す。この場合、前記の第 1実施形態における入射光線を振幅変調光線、所定パラメータを検出光線の振幅 Aと位相遅れ ø、検出光線の光路長分布 Zi ( _av)を c (媒体内の光速）倍の群遅延 3 _;

C の分布、び _{i v} ²を上記 c 倍の群遅延の分散にそれぞれ変更する以外は第 1 実施形態と同様にして測定値 Y_n、推定値 Y_{a v n}、基準値 _{a v}及び重み関数 Wiを求め、それらに基づいて吸収係数偏差 Δ _{a i}、吸収係数絶対値/ z _{a i}を求めることによって、吸収成分濃度の分布が求められる。

なお、この場合は、求めるべき重み関数 Wiは次式：

に示したものであり、び _f ²は上記 C倍の群遅延の分布の分散を表す。また、第 1 実施形態における基本光路長分布に代えて、 C 倍の群遅延分布が用いられる。

以上、本発明の好適な実施形態について説明したが、本発明は勿論上記実施形態に限定されるものではない。

すなわち、上記実施形態においては吸収係数の基準値及び輸送散乱係数の基準値を光 C T装置自身で得られたデ一夕から求めているが、これらの基準値を別の装置で求めるようにしてもよい。この場合の利点は、例えば光 C T装置で得られるデ一夕は C W (連続光線）で測定し、上記基準値を求める装置でのみパルス光線や変調光線を用いればよくなるため、光 C T装置のシステムの構成が簡単になる。なお、別の装置で上記基準値を求める手法は位相変調法や時間分解分光法であってもよい。

また、上記実施形態においては光入射位置を走査させていたが、光入射位置と同期して光検出位置を走査させてもよい。また、被計測媒体の周囲に複数の光入射ファイバ一及び複数の光検出ファイバーを配置しておき、光入射に使用するフアイバ—をそれぞれ順次変えていくことによって光入射位置を移動させてもよい _c 更に、以上説明した本発明の装置及び方法は、一般的な断層像を得る光 C Tだけでなく、 3次元媒体内部の吸収成分を 3次元計測する装置及び方法に利用できることは明らかである。また、本発明の装置及び方法は、頭部の皮膚、頭蓋骨、灰白質、白質のように、多数の層状の構造をもつ媒体の計測に利用できることも明白である。なお、この場合は、それぞれの層をそれぞれのポクセルに対応させればよい。産業上の利用可能性

以上説明したように、本発明によれば、新規な重み関数を使用して個々の計測に際して各ポクセルにおける重み関数がマイクロ ·ベア ·ランバート則に基づいて直接的に求められるため、計測状況に応じた適切な重み関数に基づいて吸収係数の偏差が算出され、従って近似法の採用に伴う誤差の発生が防止される。そのため、本発明によれば、生体などの不均一系の散乱吸収体について適用した場合であっても各ボクセルにおける吸収係数の偏差を正確に求めることが可能となり、このような吸収係数偏差に基づいて高精度の内部特性分布を得て光 C T画像化することが可能となる。

このように、本発明によって再入射不可能な様々な外形をもつ種々の散乱吸収体の内部吸収成分の分布を高精度かつ迅速に計測することが可能となり、例えば、光マンモグラフィ一、頭部光 C T、全身用光 C T、臨床モニター、診断や解析、手術や治療への応用も可能となる。

Claims

請求の範囲

1 . 散乱吸収体である被計測媒体中に一以上の光入射位置から順次光線を入射する光入射ステップと、

前記被計測媒体内部を透過した光線を複数の光検出位置で検出する光検出ステヅプと、

検出された各光線に基づいて、前記光入射位置と前記光検出位置との組み合わせ毎に前記光線の所定パラメータの測定値を取得する測定値取得ステップと、前記被計測媒体の吸収係数の基準値を設定する基準値設定ステ、ソプと、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記光入射位置と前記光検出位置との組み合わせ毎に前記パラメ一夕の推定値を求める推定値算出ステップと、

前記吸収係数の基準値を用いて、複数のポクセルに分割された前記被計測媒体の各ポクセルにおける重み関数をマイクロ ·ベア ' ランバート則に基づいて求める重み関数演算ステップと、

前記パラメータの測定値、前記パラメ一夕の推定値及び前記重み関数に基づいて、前記各ポクセルにおける前記吸収係数の基準値に対する吸収係数の偏差を算出する吸収係数偏差算出ステップと、

前記吸収係数の基準値及び前記吸収係数の偏差に基づいて前記各ポクセルにおける吸収係数の絶対値を算出して前記被計測媒体における吸収係数の絶対値分布を求める吸収係数絶対値算出ステップと、

を含む、散乱吸収体の内部特性分布の計測方法。

2 . 前記重み関数は、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなしたときの、前記各ポクセルにおける平均光路長とその分布の分散との関数である、請求項 1記載の方法。

3 . 前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記各ポクセルにおける平均光路長を求める平均光路長取得ステツプを更に含んでおり、

前記重み関数演算ステップにおいて、下記の式：

^·=^·(^ν)-^^σ ₍ ² (1)

[式中、 Wiは重み関数、 Zi(〃_{a v})は平均光路長、 A〃_{a i}は吸収係数の偏差、び _{i v} ²は平均光路長の分布の分散をそれぞれ示す]

に基づいて前記重み関数を求める、請求項 2記載の方法。

4. 前記重み関数は、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなしたときの、前記各ポクセルにおける所定時間領域内の平均光路長とその分布の分散との関数である、請求項 1記載の方法。

5. 前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記各ポクセルにおける所定時間領域内の平均光路長を求める平均光路長取得ステツプを更に含んでおり、

前記重み関数演算ステップにおいて、下記の式：

[式中、は重み関数、 Zi(〃_{a v})は平均光路長、 A〃_{a i}は吸収係数の偏差、び _{i v} ²は平均光路長の分布の分散をそれそれ示す]

に基づいて前記重み関数を求める、請求項 4記載の方法。

6. 前記重み関数は、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなしたときの、前記各ポクセルにおける群遅延とその分布の分散との関数である、請求項 1記載の方法。

7. 前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記各ポクセルにおける群遅延を求める群遅延取得ステツプを更に含んでおり、

前記重み関数演算ステップにおいて、下記の式：

[式中、 W iは重み関数、

θ φ,· は c (媒体内の光速）倍の群遅延、 Δ _{a i}は吸収係数の偏差、び _f ²は c倍の群遅延の分布の分散をそれぞれ示す]

に基づいて前記重み関数を求める、請求項 6記載の方法。

8 . 前記吸収係数の絶対値を用いて、前記各ポクセルにおける吸収成分の濃度を算出して前記被計測媒体における吸収成分の濃度分布を求める濃度算出ステップを更に含む、請求項 1記載の方法。

9 . 前記被計測媒体が少なくとも 2つの吸収成分を含有しており、前記光入射ステップにおいて前記被計測媒体中に入射される光線が、前記吸収成分に対する吸収係数が互いに相違する少なくとも 2つの波長を有しており、前記光検出ステップにおいて前記少なくとも 2つの波長を有する光線をそれぞれ検出し、

前記測定値取得ステップにおいて前記少なくとも 2つの波長を有する光線に関してそれそれ前記測定値を取得し、

前記基準値設定ステップにおいて前記少なくとも 2つの波長を有する光線に関してそれそれ前記基準値を設定し、

前記推定値算出ステップにおいて前記少なくとも 2つの波長を有する光線に関してそれぞれ前記推定値を求め、

前記重み関数演算ステツプにおいて前記少なくとも 2つの波長を有する光線に関してそれぞれ前記重み関数を求め、

前記吸収係数偏差算出ステップにおいて前記少なくとも 2つの波長を有する光線に関してそれぞれ前記吸収係数の偏差を算出し、

前記吸収係数絶対値算出ステップにおいて前記少なくとも 2つの波長を有する光線に関してそれぞれ前記吸収係数の絶対値を算出し、前記濃度算出ステップにおいて前記少なくとも 2つの波長を有する光線に関してそれぞれ前記吸収成分の濃度を算出し、前記被計測媒体における前記各吸収成分の濃度分布を求める、請求項 8記載の方法。

1 0 . 前記の求められた分布に基づいて、前記被計測媒体における該分布を示す光 C T画像を表示する画像表示ステップを更に含む、請求項 1記載の方法。

1 1 . 散乱吸収体である被計測媒体中に一以上の光入射位置から順次光線を入射する光入射手段と、

前記吸収係数の基準値を用いて、複数のボクセルに分割された前記被計測媒体の各ポクセルにおける重み関数をマイクロ 'ベア ' ランバート則に基づいて求める重み関数演算手段と、

前記パラメータの測定値、前記パラメ一夕の推定値及び前記重み関数に基づいて、前記各ポクセルにおける前記吸収係数の基準値に対する吸収係数の偏差を算出する吸収係数偏差算出手段と、

を備える、散乱吸収体の内部特性分布の計測装置。

1 2 . 前記重み関数は、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなしたときの、前記各ポクセルにおける平均光路長とその分布の分散との関数である、請求項 1 1記載の装置。

1 3 . 前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記各ポクセルにおける平均光路長を求める平均光路長取得手段を更に含んでおり、

前記重み関数演算手段において、下記の式：

[式中、 W iは重み関数、

_iは吸収係数の偏差、び _{i v} ²は平均光路長の分布の分散をそれぞれ示す]

に基づいて前記重み関数を求める、請求項 1 2記載の装置。

1 4 . 前記重み関数は、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなしたときの、前記各ポクセルにおける所定時間領域内の平均光路長とその分布の分散との関数である、請求項 1 1記載の装置。

1 5 . 前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記各ポクセルにおける所定時間領域内の平均光路長を求める平均光路長取得手段を更に含んでおり、

前記重み関数演算手段において、下記の式：

[式中、 W iは重み関数、

は平均光路長、厶〃 _{a i}は吸収係数の偏差、び _{i v} ²は平均光路長の分布の分散をそれそれ示す]

に基づいて前記重み関数を求める、請求項 1 4記載の装置。

1 6 . 前記重み関数は、前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなしたときの、前記各ポクセルにおける群遅延とその分布の分散との関数である、請求項 1 1記載の装置。

1 7 . 前記被計測媒体が全体的に一様な前記吸収係数の基準値を有するとみなして、前記吸収係数の基準値に基づいて、前記各ポクセルにおける群遅延を求める群遅延取得手段を更に含んでおり、

前記重み関数演算手段において、下記の式：

[式中、 W iは重み関数、

は c (媒体内の光速）倍の群遅延、 A / _{a i}は吸収係数の偏差、び _f ²は c倍の群遅延の分布の分散をそれそれ示す]

に基づいて前記重み関数を求める、請求項 1 6記載の装置。

1 8 . 前記吸収係数の絶対値を用いて、前記各ポクセルにおける吸収成分の濃度を算出して前記被計測媒体における吸収成分の濃度分布を求める濃度算出手段を更に備える、請求項 1 1記載の装置。

1 9 . 前記被計測媒体が少なくとも 2つの吸収成分を含有しており、前記光入射手段において前記被計測媒体中に入射される光線が、前記吸収成分に対する吸収係数が互いに相違する少なくとも 2つの波長を有しており、前記光検出手段において前記少なくとも 2つの波長を有する光線をそれぞれ検出し、

前記測定値取得手段において前記少なくとも 2つの波長を有する光線に関してそれぞれ前記測定値を取得し、

前記基準値設定手段において前記少なくとも 2つの波長を有する光線に関してそれぞれ前記基準値を設定し、

前記推定値算出手段において前記少なくとも 2つの波長を有する光線に関してそれぞれ前記推定値を求め、

前記重み関数演算手段において前記少なくとも 2つの波長を有する光線に関してそれそれ前記重み関数を求め、前記吸収係数偏差算出手段において前記少なくとも 2つの波長を有する光線に関してそれぞれ前記吸収係数の偏差を算出し、

前記吸収係数絶対値算出手段において前記少なくとも 2つの波長を有する光線に関してそれぞれ前記吸収係数の絶対値を算出し、

前記濃度算出手段において前記少なくとも 2つの波長を有する光線に関してそれぞれ前記吸収成分の濃度を算出し、前記被計測媒体における前記各吸収成分の濃度分布を求める、請求項 1 8記載の装置。

2 0 . 前記の求められた分布に基づいて、前記被計測媒体における該分布を示す光 C T画像を表示する画像表示手段を更に備える、請求項 1 1記載の装置。