ＭＲＩ画像作成法 その１［位置情報： 磁界操作の妙法 をくり出して］、その２［パルス・シーケンス： 磁界操作と測定の手順 ］において、 体内の水素の原子核に与えた共鳴に、その存在状態を担った信号を付加して、
身体の外に取り出し、画像データに変換して表示する仕組みを 勉強してきた。
このページでは、共鳴信号を取り扱う方法を、少し掘り下げて解説する。
共鳴信号の受信
信号は、電波として放射される。難しく言えば「電磁波」、電界と磁界が織物のように交互に細かく絡み合いながら空間に拡がって、遠方にまで「光速」で走っていく。その正体を理解するにはかなり高等な数学を必要とするので、ここでは下図に示すように、水素原子核の「スピン」から放射され、アンテナコイルに高周波電流としてキャッチされる信号として、概念的に理解しておこう。
原子核の磁極がアンテナに近づくと、コイルに交差する磁力線が増えて、コイルに誘起する電流が多くなり、遠ざかると小さくなる。反対側の磁極が近づくと電流の向きは逆になる。１回転で１波ができる。回転が続く限り、繰り返される。このような波を「正弦波」（サイン波、Sine
wave）、および、「余弦波」（コサイン波、Cosine wave）という。
コサイン・アンテナコイルとサインアンテナコイルは、角度を９０度変えてあるので、それぞれに誘起する信号の波形は、９０度ほどずれていて、余弦波と正弦波が誘起される。
相互の波の位置のずれを「位相」と言い、余弦波は正弦波より位相が９０度ほど進んでいる。
（上の図には略したが、通常、アンテナコイルは同形を 反対側にも設けて直列に接続し、広い空間から信号を検知できるようにする。）
アンテナを図のように直角に２個設けると、回転方向を捉えることができる。１個だけでも、回転数すなわち周波数を検知することはできるが、回転方向は検知されないので、後述する検波回路で補わなければならない。
また、２個のアンテナで検出すると、信号が２倍得られることにになり、それに含まれるノイズ分は、合わせても２倍にはならない（√２倍になる）ので、Ｓ／Ｎ比（Signal/Noize
Ratio）が√２倍ほど向上して、きれいな画像になる。
ただし、両アンテナ・コイルは完全に直交していて相互の干渉がないように、機械的また電気的に 正確に作らないと、共振効果が妨げられて、感度を上げることができない。ふつう、アンテナは、インダクタンスとキャパシタンスを組み合わせて特定の周波数で共振させるようにして、信号を大きく検知できるようにするのだが、２個の共振回路が干渉する場合には、共振を消し合う電流が流れて、異常な信号になってしまう。人体が入るとコイルのバランスが崩れるので、干渉を無くすることは実際上は困難である。スピンの回転方向は決まっているので、一方の信号は回路的に発生させて検波することも多い。
アンテナコイルの一例
頭部用。初期の試作品。
送受信兼用。
左右対向１コイル式
（直交２コイルではない）。
右上の物体は
共振周波数調節（同調）用の
バリコン（可変キャパシター）。
周波数変換と検波
得られた信号は、例えば 0.5 Ｔ（ﾃｽﾗｰ）静磁界のＭＲＩ装置では、中心周波数が 21.289 MHz（ﾒｶﾞﾍﾙﾂ）の高周波である。共鳴の発生位置を知るために、視野に傾斜磁界が加えられて、共鳴周波数が場所によってシフトしているが、その量は、50cm視野で±0.1MHz＝±100kHz
である。すなわち、最大0.5%の周波数変調が加えられている。
高周波信号から、直接に、わずか0.5%の周波数の違いを検出するのはたいへんなので、周波数変換という処理をほどこして、共鳴信号の周波数を低周波領域に移してから、被変調信号のみを取り出す。
よく用いられる周波数変換回路は、「乗算回路」である。ふたつの入力信号を入れると、それらの積になる信号をアナログ信号のままで算出して出力する。乗算回路の一方の入力端子に受信信号を入れ、もう一方に安定な周波数の基準信号を加えておくと、出力には、両周波数の和の周波数の信号と、両周波数の差になる周波数の信号とが混ざって現れる。
式で示すと、
sin (fa･t) × sin (fb･t) = sin（ (fa+fb)･ｔ）＋ sin（ (fa-fb)･ｔ）
fbとして、中心周波数＝21.289 MHzを入れると、(fa+fb) は42.5MHz ほどの高周波になり、(fa-fb) は、±0.1MHz＝±100kHz以内の低周波になる。前者は低周波域通過型のフィルター回路を通せば完全に除去でき、残った後者は被変調信号そのものになる。この過程を「検波」という。
これをＡＤ（アナログ−ディジタル）変換器に渡して、さらにコンピューターに入力し、フーリエ変換すれば、周波数分析ができ、めでたく画像を得ることができる。
ところで、視野の両端で、一方は＋100kHzの周波数偏移を 与えられ、他方では−100kHzの偏移となる共鳴信号が重なったとき、ふたつは区別できるか？ どちらも同じ周波数だから、中心を折り目として、両側の画像が重なってしまわないか？
実は左右を区別するために、すなわち、周波数の正と負を区別するために、コサイン検出信号とサイン検出信号の２系列が用意されている。それぞれを周波数変換した後において、左右の対称位置では周波数偏移は同じ値だが、
回転方向が逆になる、サイン検出信号とコサイン検出信号の位相関係が逆になるので、フーリエ変換した後では、左右がきちんと分かれる。
以上の回路を図にまとめると、下図のようになる。
以上の回路は複雑なものではないが、アンテナコイルの検知信号は相当に低レベルであって、良質な画像を得るためには、かなり性能の高いものが必要である。
例えば、ノイズ発生が充分に少なく、周波数特性が広く、位相特性、過渡応答特性が 良いことが求められる。そして、直線性が悪いと混変調歪みが現れて高調波が発生することになり、ノイズが加わって画質の劣化をもたらすので、広いレンジにわたって直線性の良いことが求められる。これは、特に乗算回路の設計において厳しい課題となる。
信号発生回路の時間精度、位相精度の良いことも重要である。また、ＡＤ変換器の変換レートはエイリアスが発生するのを避けるために、信号上限周波数の数倍が必要であり、その前の低域通過フィルターが余分な高域信号を抑えるようにしっかりと作動することも
重要なことである。
信号処理過程をシミュレーションで理解しよう
以上の過程を、コンピューターの上で模擬的に計算させて追跡してみよう。実機のような波形を目で眺めることができ、計算方法を知ることができる。（このような計算を「シミュレーション
simulation 」（模擬実験）という。）
計算ソフトとして、「Mathematica」を用いる。これは、プログラム文が分かりやすいので、なじみやすいソフトである。
ここで、波を扱う数式について、ごくわずか、説明する。後日もっと詳しく勉強することがあれば、ここの説明が具体的な応用例として思い出されて、興味が深まると思う。
余弦波、正弦波を表現する関数として、三角関数 cos(t), sin(t) がある。t は時間の経過を示す。このページの冒頭の図に示したような波形である。
cos(t), sin(t) は、しかし、別々の関数だが、この例のように、一つの回転から生ずる互いに密接な関係にある波は、一つの関数で現した方が便利だ。そこで、「虚数」と、それを取り入れた「複素数」を利用する算法が考え出された。（前々世紀ころの話である。）
虚数はｉで現して、ｉ＝ √−１、すなわちｉを自乗したものは＝-１となる、と定義される。 普通の数は自乗すると必ずプラスの数になるが、ｉは自乗するとマイナスになる、けったいな数である。（ｉらしくないね
^o^）
ところが、これを使うと、
Exp（ｉ･ωｔ）＝ cos（ωｔ）＋ｉ・sin（ωｔ）
ｴｸｽﾎﾟｰﾈﾝｼｬﾙ（ｱｲ ｵﾒｶﾞ ﾃｨｰ） ｲｸｫｰﾙ ｺｻｲﾝ（ｵﾒｶﾞ ﾃｨｰ） ﾌﾟﾗｽ ｱｲ ｶｹﾙ ｻｲﾝ（ｵﾒｶﾞ
ﾃｨｰ） と読む。
Exp（ｉ･ωｔ）とは、定数ｅの（ｉ･ωｔ）乗を示す。
ω は「角速度」を現す係数、ω＝２πｆと置き直して周波数ｆとしてもよい
定数ｅは、自然対数の底(ﾃｲ)と言い、その数値は、2.71828...と、無限小数になる。
という数式によって、見事に 一つの式で現される。下の計算で使ってみよう。
（Mathematicaソフトの中では、 ｉはＩ（大文字のアイ）で書く。ついでながら、電気工学の世界の習慣では、ｉは電流を示す記号に使っているので、虚数にはｊを用いる。「ｊω(ｼﾞｪｲ ｵﾒｶﾞ) の学徒」などとしゃれて言ったりするくらいに、電気回路の状態解析計算などにｊωを便利に使っている。）
次のプログラムは、共鳴信号を模擬的に発生させる。
左図のように、４個の発生源が視野の左端から右端まで
４個あり、強さと周波数を違えて回転している。
中心周波数は、仮に 2000 kHz にしている。
上の実例では 21 MHz余だが、計算量を減らすために
下げた。
（ Mathematica では、
関数の中身は、角カッコ ［ ］でくくって示し、
掛け算の記号 × は、「スペース（空白）」で現す。
例えば、
Exp[ I (2 Pi (fa+fcenter)/r t/n)] は、
Exp（Ｉ×２π（fa+fcenter）／r ×ｔ／ｎ）
である。Pi は π(ﾊﾟｲ)を現す。
(* *) で囲んだ文は注釈のために書き込んだもので、
計算には関係しない。
(* MRI signal simulation. 2MHz carrier, +-100kHz modulation. *)
fcenter=2000; (* =2000 kHz *)
(* 1024points, 20.48uS time-length. 50MHz-sampling (0.02uS/sample). *)
n=1024;
fa=-100; fb=-50; fc=20; fd=100; (* in kHz *) r=1/0.02048;
wave=Table[ N[ ←表waveに
( Exp[ I (2 Pi (fa+fcenter)/r t/n)] 0.8
+Exp[ I (2 Pi (fb+fcenter)/r t/n)] 0.7 (...) の値を、
+Exp[ I (2 Pi (fc+fcenter)/r t/n)] 0.6 順に計算して並べる、
+Exp[ I (2 Pi (fd+fcenter)/r t/n)] 0.5
) ], {t, 0, n-1}]; ←変数tが 0 からn-1 まで。
ListPlot[Re[wave], PlotJoined -> True, ←表の実数部（real）をグラフにする。
AspectRatio -> 0.2, AxesOrigin -> {0, 0}]
ListPlot[Im[wave], PlotJoined -> True, ←表の虚数部（Imagenary）をグラフにする。
AspectRatio -> 0.2, AxesOrigin -> {0, 0}]
このように、Mathematica というプログラム言語は、数値の固まりを集団として取り扱い、その作図を簡単に指示できるので、とても使いよい。
計算された波形を上に示す。このように、複素数を利用した指数関数は、回転から現れる二つの波形を一つの式で取り扱う。
この波は、４個の発生源から出てくるわずかに周波数の違った波が合わさっているので、唸りを含んだような、やや複雑な波形になっている（お寺の鐘をゴーンと突くと、ワーン、ワーン、、と余韻がうねりながら響く、そんな感じか）。
この波形のスペクトルを計算してみても、1000点くらいのデータ量では分解能が足りなくて、４個の周波数はとても分解されない。
そこで、次のプログラムでは、周波数変換をして信号を低周波に落とし、これを256点サンプリングしてスペクトルを求める。
(* MRI signal simulation. 2MHz carrier, +-100kHz modulation. *)
(* Frequency conversion, hi-freq suppression, Fourier analysis. *)
fcenter=2000; (* =2000 kHz *)
n=256;
fa=-100; fb=-50; fc=20; fd=100; (* in kHz *) r=1/0.02048;
pa= .8; pb=.6; pc=.4; pd=.2;
lowwave=Table[Sum[ N[ ←Sum は、合計を求める関数、Ｎは、数値扱いにする関数。
Cos[ (2 Pi fcenter/r t/n)] ←中心周波数のコサイン波を、
( Exp[ -I (2 Pi (fa+fcenter)/r t/n)] pa
+Exp[ -I (2 Pi (fb+fcenter)/r t/n)] pb 信号波に掛け算する。
+Exp[ -I (2 Pi (fc+fcenter)/r t/n)] pc
+Exp[ -I (2 Pi (fd+fcenter)/r t/n)] pd
) ], { t, 50 i, 50 i+49 } ], {i, 0, n-1} ]; ←50点ごとに合計、ｎ回続ける。
（高周波成分を除去し低周波成分だけにするために移動平均処理を加えている。）
factors=Fourier[lowwave]; ←フーリエ変換の関数。
factors=Chop[factors]; ←半端な数（計算誤差）を丸める関数。
spectrum=Join[Take[factors, -n/2]],Take[factors,n/2]];
↑周波数スペクトルを表示する順序を整える（電気工学の慣習に合わせる）。
ListPlot[Re[factors], PlotRange -> {{0,n}, {0, 35}}, ←波形を描画。
AspectRatio -> 0.2, PlotJoined -> True]
ListPlot[Abs[spectrum], PlotRange -> {{0,n}, {0, 200}}, ←スペクトルを描画。
AspectRatio -> 0.2, PlotJoined -> True]
計算して得られた低周波信号の、波形とスペクトルを、上図に示す。 ４個の物体から発信した共鳴信号が、それぞれの位置に算出されている。
位相変調が加わった信号の処理過程
次に、実際の共鳴波のように、スピンの存在位置を示す周波数変調と位相変調を加えた信号を模擬して発生させ、２回のフーリエ変換を直交してほどこして画像に変換する過程を、プログラムで調べてみよう。
計算結果を先に示す。下図は、５個の共鳴スピンを視野に散在させ、画像を求めている。左は、第１回目のフーリエ変換結果、すなわち周波数変調情報をスペクトルに変換して、縦方向に並べた図。スピン物体が存在する列に波が現れ、その波の周波数が物体のＹ位置を示す。先の説明で、「Ｙはリズムで表現」と述べたリズムである。 右はこれを横方向に列ごとにフーリエ変換して、物体の位置を定め、表示している。
計算したプログラムを以下に示す。
n=32; r=1/0.02048; fcenter=2000; (* =2000 kHz *)
各共鳴点の、
pa= 1; pb=.8; pc=.6; pd=.4; pe=.2; ←大きさ（振幅）
fa=-10; fb=-5; fc=0; fd=3; fe=10; (* in kHz *) ←周波数（Ｘ位置）
ya= -1; yb=.2; yc=0; yd=-.4; ye=1; ←位相 （Ｙ位置）
preImage=Table[ i+I j, {i,1,n},{j,1,n}]; データエリアを用意する
finaleImage=Table[100 j+i, {i,1,n},{j,1,n}]; 画像エリアを用意する
Do[
lowwave=Table[Sum[
N[ Cos[(2 Pi fcenter/r t/n)] ←復調のためにコサイン波を掛ける
( Exp[-I (2 Pi (fa+fcenter)/r (t-4 ya (1-2 j/n))/n)] pa ）
+Exp[-I (2 Pi (fb+fcenter)/r (t-4 yb (1-2 j/n))/n)] pb ）
+Exp[-I (2 Pi (fc+fcenter)/r (t-4 yc (1-2 j/n))/n)] pc ）←５点の共鳴波を算出
+Exp[-I (2 Pi (fd+fcenter)/r (t-4 yd (1-2 j/n))/n)] pd ）
+Exp[-I (2 Pi (fe+fcenter)/r (t-4 ye (1-2 j/n))/n)] pe ）
)],
{ t, 50 i-50, 50 i }], ←50点ごとに積算平滑化
{i,1,n}];
factors=Fourier[lowwave]; ←第１回フーリエ変換
factors=Chop[factors];
spectrum=Join[Take[factors,-n/2],Take[factors,n/2]];
Do[ preImage[[i,j]]=spectrum[[i]], {i,1,n}], {j,1,n}];
ListDensityPlot[Re[preImage]]; ←第１データ画像を表示
Do[
Do[lowwave[[i]]=preImage[[j,i]],{i,1,n}]; ←縦列データを取り出し
factors=Fourier[lowwave]; ←第２回フーリエ変換
factors=Chop[factors];
spectrum=Join[Take[factors,-n/2],Take[factors,n/2]];
Do[finaleImage[[j,i]]=Abs[spectrum[[i]] ],{i,1,n}], ←縦列データとして収納
{j,1,n}];
ListDensityPlot[finaleImage, PlotRange -> {0,Max[finaleImage]}]; ←できた画像を濃度表示
List3DPlot [finaleImage, PlotRange -> {0,Max[finaleImage]},
ViewPoint -> { -.3, -2.4, 2}]; ←できた画像を３次元表示 （下図）
このように、コンピューターの中で、模擬的に信号を発生させたり信号処理を加えたりして、実際の状態を実験することができる。ただし、たくさんの点に対して、ていねいな計算を克明に実施するので、大きなコンピューター容量を使い、長い計算時間を必要とする。上の計算では、画像を荒くしてデータ量を減らし、大まかな現象しか模擬していないが、２時間ほどの計算時間を費やしている（手持ちの古い方のパソコン（マキントッシュCentris
650 1993年製）を使ったのでこんな時間だが、最新のマックでは1/10以下の所要時間であろう）。
さらに、画像を精細にしたり、データの拡がりを含めたり、緩和時間の影響を含めたりすれば、数百倍以上の計算時間になってしまう。現象の大事な点だけに集中して検討することが肝要。だが、予期しないところをあぶり出すのが実験の効果でもあるから、下手な集約は危険でもある。
シミュレーションは、諸種の現象の探求にはきわめて便利な手段であり、現在、多くの領域で利用が進んでいる。スーパーコンピューターが珍重される由縁である。
そして、いろんな分野に向けて、たくさんの計算機能を織り込んだシミュレーション・ソフトウェアが、開発されている。
あらまし ＭＲＩの元信号 「フーリエ変換」の効用 ＭＲＩ画像 作成法-1 ＭＲＩ画像 作成法-2 ＭＲＩの黎明期 もくじに戻る