[go: up one dir, main page]

JP2009000164A - 放射線画像処理装置 - Google Patents

放射線画像処理装置 Download PDF

Info

Publication number
JP2009000164A
JP2009000164A JP2007161526A JP2007161526A JP2009000164A JP 2009000164 A JP2009000164 A JP 2009000164A JP 2007161526 A JP2007161526 A JP 2007161526A JP 2007161526 A JP2007161526 A JP 2007161526A JP 2009000164 A JP2009000164 A JP 2009000164A
Authority
JP
Japan
Prior art keywords
filter
small section
edge
scanning
image processing
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
JP2007161526A
Other languages
English (en)
Inventor
Daisuke Notohara
大介 能登原
Kazuyoshi Nishino
和義 西野
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shimadzu Corp
Original Assignee
Shimadzu Corp
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Shimadzu Corp filed Critical Shimadzu Corp
Priority to JP2007161526A priority Critical patent/JP2009000164A/ja
Publication of JP2009000164A publication Critical patent/JP2009000164A/ja
Pending legal-status Critical Current

Links

Images

Landscapes

  • Apparatus For Radiation Diagnosis (AREA)

Abstract

【課題】エッジ部分をボケさせずにノイズ成分を低減させることができる放射線画像処理装置を提供することを目的とする。
【解決手段】フィルタを、フィルタを等方低域通過フィルタ(等方LPF)31a、等方高域通過フィルタ(等方HPF)31b、非等方高域通過フィルタ(非等方HPF)31cの3つに区分して、エッジの強度Vおよび方向θに基づいて、区分された各フィルタをそれぞれ変えることで各々のフィルタをそれぞれ決定するので、エッジ部分をボケさせずにノイズ成分を低減させることができる。
【選択図】図2

Description

この発明は、放射線画像に対して所定の処理を行う放射線画像処理装置に係り、特に、放射線画像に対してフィルタ処理を行う技術に関する。
放射線画像として、X線画像を例に採って説明する。従来、X線画像のノイズを減らすノイズリダクションとして、スムージング(画像の平滑化処理)等による空間フィルタ演算、またはリカーシブフィルタ(再帰的演算処理)等の時間フィルタ演算がある。空間フィルタ演算はノイズの振幅を減弱させ、画像が粒状にざらついて出力される画像の粒状性を向上させることができる。その一方で、空間フィルタ演算は、画像の信号強度自体も劣化させ、エッジがボケる欠点がある。時間フィルタ演算は、連続的な画像を重み付け加算することによって、ノイズ成分を低減させることができる。その一方で、時間フィルタ演算では、動きボケを生じるといった大きな欠点がある。ビデオレートは30フレーム/秒であるが、このビデオレートよりもフレームレートを15フレーム/秒、7フレーム/秒と下げると、前の画像との時間差が大きくなってしまい、前の画像による残像が残ってしまう。したがって、低線量化を狙って低フレームレートで透視画像(動画)を行った場合には、その透視画像では動きボケの弊害が顕著になる。
この発明では、空間フィルタ演算について説明する。上述したように空間フィルタ演算ではエッジがボケるので、構造物などを撮像した場合に画像に映った構造物の輪郭部分に相当するエッジまでボケて、構造物が正確に撮像されない場合が生じる。そこで、エッジの強度および方向を検出して、その検出結果に基づいて2つのフィルタのうちのいずれかに切り換えて、空間フィルタ演算、すなわちフィルタ処理を行うことで、エッジ部分をボケさせずにノイズ成分を低減させる技術がある(例えば、特許文献1参照)。
特開2001−283215号公報(第1−18頁、図1−27)
かかる特許文献1以外の手法で、エッジ部分をボケさせずにノイズ成分を低減させる手法も望まれる。
この発明は、このような事情に鑑みてなされたものであって、エッジ部分をボケさせずにノイズ成分を低減させることができる放射線画像処理装置を提供することを目的とする。
この発明は、このような目的を達成するために、次のような構成をとる。
すなわち、請求項1に記載の発明は、放射線画像に対して所定の処理を行う放射線画像処理装置であって、前記放射線画像を構成する各々の画素に基づいて画素値の統計量を算出する統計量算出手段と、その統計量算出手段によって求められた画素値の統計量に基づいて、周囲の画素に対して特異となる箇所の強度をエッジの強度とするとともに、特異となる箇所の方向をエッジの方向として算出するエッジ算出手段と、放射線画像に対してフィルタ処理を行うためのフィルタを複数に区分して、前記エッジ算出手段によって求められたエッジの強度および方向に基づいて、前記区分された各フィルタをそれぞれ変えることで各々のフィルタをそれぞれ決定するフィルタ決定手段と、そのフィルタ決定手段によって決定された各フィルタを複合して放射線画像に対してフィルタ処理を行うフィルタ処理手段とを備えることを特徴とするものである。
[作用・効果]請求項1に記載の発明によれば、統計量算出手段は、放射線画像を構成する各々の画素に基づいて画素値の統計量を算出する。その統計量算出手段によって求められた画素値の統計量に基づいて、エッジ算出手段は、周囲の画素に対して特異となる箇所の強度をエッジの強度とするとともに、特異となる箇所の方向をエッジの方向として算出する。構造物の輪郭部分などのエッジは周囲の画素に対して特異となるので、統計学的な知見から画素値の統計量に基づいてエッジの強度および方向を算出することが可能である。ここで、放射線画像に対してフィルタ処理を行うためのフィルタを複数に区分する。そして、上述したエッジ算出手段によって求められたエッジの強度および方向に基づいて、区分された各フィルタをフィルタ決定手段がそれぞれ変えることで各々のフィルタをそれぞれ決定する。そのフィルタ決定手段によって決定された各フィルタをフィルタ処理手段が複合して放射線画像に対してフィルタ処理を行うことで、上述した特許文献1のような2つのフィルタのうちのいずれかに切り換えてフィルタ処理を行う手法とは別の手法で、エッジの強度および方向に基づくフィルタ処理が可能になる。このように、フィルタを複数に区分して、エッジの強度および方向に基づいて、区分された各フィルタをそれぞれ変えることで各々のフィルタをそれぞれ決定するので、エッジの強度および方向に応じてフィルタの決定がより柔軟に細かく行われる。その結果、エッジの強度および方向をより一層考慮して、エッジ部分をボケさせずにノイズ成分を低減させることができる。
上述した発明では、放射線画像全体に対してフィルタ処理を一度で行ってもよいし、放射線画像をそれよりも小さな小区画に区分して走査するたびにフィルタ処理をその都度行ってもよい。
後者の場合には、この発明に係る放射線画像処理装置は、放射線画像をそれよりも小さな小区画に区分する小区画区分手段と、その小区画区分手段によって区分された小区画の位置を変えることで放射線画像内を走査する小区画走査手段とを備える。そして、統計量算出手段は、小区画区分手段によって区分された小区画を構成する各々の画素に基づいて画素値の統計量を算出し、フィルタ処理手段は小区画に対してフィルタ処理を行う。上述した小区画走査手段は、以下のような走査を行う。すなわち、同一の小区画に対して、(1)統計量算出手段による統計量の算出、(2) エッジ算出手段によるエッジの強度および方向の算出、(3) フィルタ決定手段による各々のフィルタの決定、(4)フィルタ処理手段によるフィルタ処理を行う(1)〜(4)までの工程を行った後に、(5)次なる小区画の位置に変えて、(1)の工程に戻って繰り返すことで、小区画走査手段は走査する(請求項2に記載の発明)。このように、小区画に区分して走査するたびにフィルタ処理が行われるので、小区画に応じたフィルタ処理を適正に行うことができる。
後者の場合における一例は、小区画走査手段は、対象となる小区画の位置から次なる小区画の位置に変えて走査するときに、上述した対象となる小区画と次なる小区画との間に走査されない領域が含まれるように走査することである(請求項3に記載の発明)。
後者の場合における別の一例は、小区画走査手段は、対象となる小区画の位置から次なる小区画の位置に変えて走査するときに、上述した対象となる小区画と次なる小区画とが接するように走査することである(請求項4に記載の発明)。
後者の場合におけるさらなる別の一例は、小区画走査手段は、対象となる小区画の位置から次なる小区画の位置に変えて走査するときに、上述した対象となる小区画の一部を上述した次なる小区画が含んで、対象となる小区画と次なる小区画とが重なるように走査することである(請求項5に記載の発明)。
このように、フィルタ処理の演算負荷は、請求項5、請求項4、請求項3の順に軽くなる。逆にフィルタ処理による演算としては、請求項3、請求項4、請求項5の順に精密になる。
また、後者の場合には、小区画区分手段によって区分される小区画は方形であるのが好ましい(請求項6に記載の発明)。放射線画像全体が方形であるので、小区画も方形にすることで各種の演算処理が簡潔になる。
上述したこれらの発明の一例として、フィルタ決定手段は、エッジ算出手段によって求められたエッジの強度および方向に基づいた複数の係数を、区分された各フィルタにそれぞれ乗じることで、各フィルタをそれぞれ変えて、各々のフィルタをそれぞれ決定することである(請求項7に記載の発明)。各係数に応じて、係数を乗じたフィルタがそれぞれ決定される。
この発明に係る放射線画像処理装置によれば、フィルタを複数に区分して、エッジの強度および方向に基づいて、区分された各フィルタをそれぞれ変えることで各々のフィルタをそれぞれ決定するので、エッジ部分をボケさせずにノイズ成分を低減させることができる。
以下、図面を参照してこの発明の実施例を説明する。
図1は、実施例に係る画像処理部を用いたX線透視撮影装置のブロック図であり、図2は、実施例に係る画像処理部のブロック図である。なお、本実施例では、X線画像処理装置として画像処理部を例に採るとともに、X線画像処理装置(画像処理部)はX線透視撮影装置に用いられるものとして、以下を説明する。
X線透視撮影装置は、図1に示すように、被検体Mを載置する天板1と、その被検体Mに向けてX線を照射するX線管2と、被検体Mを透過したX線を検出するフラットパネル型X線検出器(以下、「FPD」と略記する)3とを備えている。
X線透視撮影装置は、他に、天板1の昇降および水平移動を制御する天板制御部4や、FPD3の走査を制御するFPD制御部5や、X線管2の管電圧や管電流を発生させる高電圧発生部6を有するX線管制御部7や、FPD3から電荷信号であるX線検出信号をディジタル化して取り出すA/D変換器8や、A/D変換器8から出力されたX線検出信号に基づいて種々の処理を行う画像処理部9や、これらの各構成部を統括するコントローラ10や、処理された画像などを記憶するメモリ部11や、オペレータが入力設定を行う入力部12や、処理された画像などを表示するモニタ13などを備えている。画像処理部9は、この発明におけるX線画像処理装置に相当する。
天板制御部4は、天板1を水平移動させて被検体Mを撮像位置にまで収容したり、昇降、回転および水平移動させて被検体Mを所望の位置に設定したり、水平移動させながら撮像を行ったり、撮像終了後に水平移動させて撮像位置から退避させる制御などを行う。FPD制御部5は、FPD3を水平移動させたり、被検体Mの体軸の軸心周りに回転移動させることによる走査に関する制御などを行う。高電圧発生部6は、X線を照射させるための管電圧や管電流を発生してX線管2に与え、X線管制御部7は、X線管2を水平移動させたり、被検体Mの体軸の軸心周りに回転移動させることによる走査に関する制御や、X線管2側のコリメータ(図示省略)の照視野の設定の制御などを行う。なお、X線管2やFPD3の走査の際には、X線管2から照射されたX線をFPD3が検出できるようにX線管2およびFPD3が互いに対向しながらそれぞれの移動を行う。
コントローラ10は、中央演算処理装置(CPU)などで構成されており、メモリ部11は、ROM(Read-only Memory)やRAM(Random-Access Memory)などに代表される記憶媒体などで構成されている。また、入力部12は、マウスやキーボードやジョイスティックやトラックボールやタッチパネルなどに代表されるポインティングデバイスで構成されている。X線透視撮影装置では、被検体Mを透過したX線をFPD3が検出して、検出されたX線に基づいて画像処理部9で画像処理を行うことで被検体Mの撮像を行う。メモリ部11は、画像処理部9で処理された各々の画像を書き込んで記憶するように構成されている。
次に、画像処理部9は、図2に示すように、大別すると検出部20と処理部30とを備えている。画像処理部9は、FPD3から出力されてA/D変換器8でディジタル化されたX線検出信号をX線検出信号の信号値に応じて画素値に変換して、各々の画素から構成されたX線画像として取り込み、そのX線画像を検出部20および処理部30へ送り込む。検出部20は、後述するエッジを検出し、処理部30は、検出されたエッジの強度および方向に基づいてフィルタ処理を行う。
検出部20は、エッジを検出する前の処理を行う前処理部21と、エッジ検出部22とを備えている。エッジ検出部22は、画素値の平均値や揺らぎ(偏差)などに代表される画素値の統計量を算出する統計量算出部22aと、その統計量算出部22aによって求められた平均値や揺らぎに基づいて、周囲の画素に対して特異となる箇所の強度をエッジの強度とするとともに、特異となる箇所の方向をエッジの方向として算出するエッジ算出部22bとを備えている。前処理部21で行われる前処理としては、例えば被検体Mの体厚による輝度差(画素値の差分に相当)から生じるX線検出信号の強度検出偏りを軽減させるLog変換や、強度検出可能な画像を効率よく作成するために高域通過フィルタ(HPF: High Pass Filter)によって画像の高周波成分のみを抽出する高周波強調処理などがある。統計量算出部22aは、この発明における統計量算出手段に相当し、エッジ算出部22bは、この発明におけるエッジ算出手段に相当する。
処理部30は、大別するとフィルタ決定部31とフィルタ処理部32とを備えている。フィルタ決定部31は、X線画像に対してフィルタ処理を行うためのフィルタを複数に区分して、上述したエッジ算出部22bによって求められたエッジの強度および方向に基づいて、区分された各フィルタをそれぞれ変えることで各々のフィルタをそれぞれ決定する。フィルタ処理部32は、そのフィルタ決定部31によって決定された各フィルタを複合してX線画像に対してフィルタ処理を行う。フィルタ決定部31は、この発明におけるフィルタ決定手段に相当し、フィルタ処理部32は、この発明におけるフィルタ処理手段に相当する。
具体的には、フィルタ決定部31は、フィルタを等方低域通過フィルタ(LPF: Low Pass Filter)(以下、低域通過フィルタを「LPF」と略記する)31aと等方高域通過フィルタ(HPF: High Pass Filter) (以下、高域通過フィルタを「HPF」と略記する)31bと非等方高域通過フィルタ(非等方HPF)31cと3つに区分する。等方HPF31bの下流部分には乗算器31dを備えるとともに、非等方HPF31cの下流部分には乗算器31eを備えている。エッジ算出部22bによって求められたエッジの強度に基づいた係数(後述するブレンド率limitまたは(1−limit))を乗算器31dおよび乗算器31eに送り込む。乗算器31dおよび乗算器31eによって、区分された等方HPF31b、非等方HPF31cにエッジの強度に基づいた係数(limit,1−limit)をそれぞれ乗じる。このように乗じることで、各フィルタをそれぞれ変えて、各々のフィルタをそれぞれ決定する。
具体的には、フィルタ処理部32は、等方LPF31aの下流部分、等方HPF31bの乗算器31dの下流部分、非等方HPF31cの乗算器31eの下流部分には加算器32aを備えている。このような加算器32aを備えることで、フィルタ決定部31によって決定された各フィルタを複合する。
次に、X線画像の具体的な処理について、図3〜図13を参照して説明する。図3は、一連のX線画像処理の流れを示すフローチャートであり、図4は、X線画像内の小区画の一態様を示す模式図であり、図5は、方位角が0°や90°のときの5×5局所領域の小区画の模式図であり、図6は、方位角の定義を示す模式図であり、図7は、方位角が0°のときの小区画の模式図であり、図8は、方位角が45°のときの小区画の模式図であり、図9は、方位角が90°のときの小区画の模式図であり、図10は、方位角が135°のときの小区画の模式図であり、図11は、スレッシュを設定するための模式図であり、図12は、ブレンド率を求めるための模式図であり、図13(a)〜図13(c)は、小区画の走査の各一例を示す模式図である。
メモリ部11(図1を参照)のROMには、図3のフローを行うプログラムが予め記憶されており、そのプログラムをコントローラ10(図1を参照)と同様のCPUで構成された画像処理部9が実行する。具体的には、X線画像をそれよりも小さな小区画に区分し(ステップS1)、その区分された小区画の位置を変えることでX線画像内を走査する(ステップS2〜S7)。より具体的には、対象となる小区画で統計量を算出し(ステップS2)、対象となる小区画でエッジの強度/方向を算出し(ステップS3)、対象となる小区画で各々のフィルタをそれぞれ決定し(ステップS4)、対象となる小区画でフィルタ処理を行い(ステップS5)、次なる小区画があれば(ステップS6)、次なる小区画の位置に変えて走査し(ステップS7)、ステップS2に戻って繰り返すことで、走査する。したがって、画像処理部9は、この発明における小区画区分手段の機能、およびこの発明における小区画走査手段の機能を備えている。
(ステップS1)小区画に区分
図4に示すような方形のX線画像Pを、そのX線画像Pよりも小さな小区画P(I=1,2,…,I,…,M−1,M)に区分する。区分される小区画Pは、図4に示すように方形であるのが好ましい。X線画像P全体が方形であるので、小区画Pも方形にすることでステップS2以降の各種の演算処理が簡潔になる。もちろん、X線画像を方形以外(例えば三角形や多角形)で小区画に区分してもよい。
本実施例では、図5に示すように、後述する方位角が0°や90°のときを例に採って、縦に5つの画素、横に5つの画素からなる正方形の小区画P(5×5局所領域)で説明する。もちろん、5×5局所領域以外のn×n局所領域(nは任意の自然数)であってもよいし、長方形のm×n局所領域(m、nは任意の自然数)であってもよい。また、後述する方位角が0°や90°以外の45°や135°の場合には、図8や図10のような小区画になる。また、5×5局所領域の小区画Pの基準座標を(x,y)とし、横軸をi、縦軸をjとして、左から順にi=−2,−1,0,1,2とし、上から順にj=−2,−1,0,1,2とする。したがって、5×5局所領域の小区画Pの中心座標は(0,0)となる。本実施例では、この中心座標を基準座標とする。したがって、基準座標(x,y)は中心座標(0,0)となる。なお、基準座標は中心座標に限定されず、5×5局所領域の小区画P内の任意の座標であってもよい。
(ステップS2)対象となる小区画で統計量を算出
対象となる小区画Pに対して、その小区画Pを構成する各々の画素(図4では縦横に5×5の画素)に基づいて、統計量算出部22a(図2を参照)は画素値の統計量として画素値の平均値avg(i)や揺らぎ(偏差)を算出する。平均値avg(i)や揺らぎを求めるには、例えば画素の並びに対する方位角を、図6に示すように定義する。図6では0°,45°,90°,135°と方位角を定義する。もちろん、それ以外の角度で方位角を定義してもよいし、細かく方位角を定義してもよい。
方位角が0°のときには、図7に示すように下から上へ0°の方向に沿って平均値avg(i)および揺らぎを求める。方位角が45°のときには、図8に示すように左下から右上へ45°の方向に沿って平均値avg(i)および揺らぎを求める。方位角が90°のときには、図9に示すように左から右へ90°の方向に沿って平均値avg(i)および揺らぎを求める。方位角が135°のときには、図10に示すように左上から右下へ135°の方向に沿って平均値avg(i)および揺らぎを求める。なお、基準となる方位角0°は必ずしも図6のような方向である必要はなく、例えば図6に示す45°,90°,135°あるいはそれ以外の角度を基準となる方位角0°として定義してもよい。
各々の方位角ごとの平均値avg(i)および揺らぎに基づいて平均値との差の絶対値Sを設定する。方位角が0°のときをS0とし、方位角が45°のときをS45とし、方位角が90°のときをS90とし、方位角が135°のときをS135とする。また、平均値avg(i)は方位角ごとの1×5の局所領域であり、後述するd(x,y)は基準座標での画素値、d(x+i,y+j)は基準座標から(i,j)分移動した座標での画素値である。また、Wiは重み係数を表す。
先ず、方位角が0°のときのS0を設定する場合を例に採って説明する。図7に示すように、方位角が0°のときの平均値avg(i)および揺らぎを求めて、その平均値avg(i)および揺らぎに基づいてS0を設定する。すると、S0は下記(1)式のように表される。
Figure 2009000164
他の方位角45°,90°,135°のときも同様の手順でS45,S90,S135をそれぞれ設定する。図8に示すように、S45は下記(2)式のように表される。
Figure 2009000164
また、図9に示すように、S90は下記(3)式のように表される。
Figure 2009000164
また、図10に示すように、S135は下記(4)式のように表される。
Figure 2009000164
なお、上記(1)〜(4)式はnが5のときであるが、5以外の自然数にnを拡張したときで、nが奇数の場合には、(1)〜(4)式中の±2で計算しているΣを±(n−1)/2に置換すればよい。これらの設定されたS0,S45,S90,S135は、周囲の画素に対して特異となる、すなわち画素値が異なる箇所をエッジとして抽出したときに、そのエッジの度合いの目安となる物理量(パラメータ)となる。
例えば、エッジが、方位角が0°の方向(すなわち下から上への方向)に沿って延びている場合には、S0は極端に小さな値となるが、逆にS90はエッジを垂直に横切るので、S90は極端に大きな値となる。なお、S45,S135は、その間の値となる。逆にエッジが、方位角が90°の方向(すなわち左から右への方向)に沿って延びている場合には、S90は極端に小さな値となるが、逆にS0はエッジを垂直に横切るので、S0は極端に大きな値となる。通常は、エッジは必ずしも下から上へ、あるいは左から右への方向に沿って延びているわけでなく、方位角が45°や135°のときのように斜め方向や途中で折れ曲がって延びていることに留意されたい。
(ステップS3)対象となる小区画でエッジの強度/方向を算出
対象となる小区画Pで、ステップS2と同一の小区画Pに対して、ステップS2で求められた画素値の平均値avg(i)や揺らぎに基づいて、さらにはS0,S45,S90,S135に基づいて、エッジ算出部22b(図2を参照)は周囲の画素に対して特異となる箇所の強度をエッジの強度とするとともに、特異となる箇所の方向をエッジの方向として算出する。強度をV(x,y)とするとともに、エッジの方向、すなわち角度をθ(x,y)とすると、強度V(x,y)および角度θ(x,y)は下記(5)、(6)式のように表される。
Figure 2009000164
また、小区画P全体の平均値をavgとすると、平均値avgは下記(7)式のように表される。
Figure 2009000164
ところで、スレッシュは、ランダムノイズレベルをどこまで許容するかの係数である。スレッシュをshとしたときに、図11に示すようにスレッシュshの値によって、sh√(avg)のしきい値が変わり、しきい値sh√(avg)(図11中では「sh*(avg)1/2」)よりも上の領域ではエッジ(Edge)となり、しきい値sh√(avg)よりも下の領域ではノイズ(Noise)となる。このしきい値sh√(avg)を下記(8)式のような条件式で設定する。
Figure 2009000164
上記(8)式の条件式の意味は、下記に示す通りである。すなわち、もし、上記(5)式で求められたエッジの強度Vがしきい値sh√(avg)よりも大きければ(上記(8)式中の「if(V≧sh√avg)」を参照、エッジ(上記(8)式中の「{E}」を参照)と認定して、エッジの強度Vがしきい値sh√(avg)未満であれば、ノイズ(上記(8)式中の「{N}」を参照)と認定する。
ここで、図12に示すように、横軸としてV/{sh√(avg)}(図12中では「V/{sh*(avg)1/2}」)をとるとともに、縦軸としてブレンド率limitをとる。すなわち、ブレンド率limitはV/{sh√(avg)}を変数とした関数となる。このブレンド率limitの関数については任意の関数で予め適宜設定すればよいが、V/{sh√(avg)}の値が大きくなるのにしたがって、ブレンド率limitの値も大きくなるように設定するのが好ましい。また、ブレンド率limitの値が0〜1の値をとるように設定するのがより好ましい。さらに、V/{sh√(avg)}=1のときにブレンド率limitの値が0.5の値をとるように設定するのがより一層好ましい。
図11や上記(8)式の条件式でも述べたように、エッジの強度Vがしきい値sh√(avg)よりも大きければエッジとするとともに、エッジの強度Vがしきい値sh√(avg)未満であればノイズとしているので、エッジの強度Vをしきい値sh√(avg)で除算した値であるV/{sh√(avg)}については、図12に示すように、“1”よりも大きければエッジとなり、“1”未満であればノイズとなる。また、V/{sh√(avg)}=1のときにブレンド率limitの値が0.5の値をとるように設定した場合には、ブレンド率limit が“0.5”よりも大きければエッジとなり、ブレンド率limitが“0.5”未満であればノイズとなる。このようにして求められたブレンド率limitや、“1”からブレンド率limitを減算した値である(1−limit)は、この発明におけるエッジの強度および方向に基づいた複数の係数に相当する。これらブレンド率limitや(1−limit)を乗算器31dおよび乗算器31e(図2を参照)に送り込む。
(ステップS4)対象となる小区画で各々のフィルタをそれぞれ決定
フィルタ決定部31(図2を参照)は、X線画像に対してフィルタ処理を行うためのフィルタを複数に区分する。そして、対象となる小区画Pで、ステップS2、S3と同一の小区画Pに対して、ステップS3で求められたエッジの強度および方向に基づいて、区分された各フィルタをフィルタ決定部31がそれぞれ変えることで各々のフィルタをそれぞれ決定する。本実施例では、フィルタを等方低域通過フィルタ(等方LPF)31aと等方高域通過フィルタ(等方HPF)31bと非等方高域通過フィルタ(非等方HPF)31cと3つに区分する。ステップS3で求められたエッジの強度に基づいた係数に相当するブレンド率limitや(1−limit)を、区分された等方HPF31b、非等方HPF31cに、乗算器31dおよび乗算器31e(図2を参照)を介してそれぞれ乗じることで、各フィルタをそれぞれ変えて、各々のフィルタをそれぞれ決定する。なお、非等方HPF31cには、エッジの方向(角度)の成分が含まれているので、上記(6)式で求められたエッジの角度θを非等方HPF31cに送り込む。
(ステップS5)対象となる小区画でフィルタ処理
対象となる小区画Pで、ステップS2〜S4と同一の小区画Pに対して、ステップS4で決定された各フィルタをフィルタ処理部32(図2を参照)が複合してX線画像に対してフィルタ処理を行う。本実施例では、等方LPF31aの下流部分、等方HPF31bの乗算器31dの下流部分、非等方HPF31cの乗算器31eの下流部分に備えられた加算器32a(図2を参照)によって、ステップS4で決定された各フィルタを複合する。複合されたフィルタでフィルタ処理されたX線画像をg(x,y)とすると、g(x,y)は下記(9)式のように表される。
Figure 2009000164
なお、上記(9)式中のg(x,y)は(フィルタ処理されたX線画像の)出力であり、上記(9)式中のf(x,y)は(フィルタ処理される前のX線画像の)入力であり、上記(9)式中のh(x,y)は等方LPF(上記(9)式では「等方型ローパスフィルタ」)であり、上記(9)式中のm(x,y)は等方HPF(上記(9)式では「等方型ハイパスフィルタ」)であり、上記(9)式中のn(x,y,θ)は非等方HPF(上記(9)式では「非等方型ハイパスフィルタ」)である。また、n(x,y,θ)中のθはエッジの方向(角度)の成分であり、上記(6)式で求められたエッジの角度θに基づいてフィルタ処理される。
上記(9)式からも明らかなように非等方HPFを示すn(x,y,θ)にはブレンド率limitを乗じ、等方HPFを示すm(x,y)には(1−limit)を乗じている。また、等方LPFを示すh(x,y)には、ブレンド率limitや(1−limit)などのエッジの強度に基づいた係数を乗じないでフィルタを決定している。上記(9)式を行うことでX線画像に対してフィルタ処理を行う。
(ステップS6)次なる小区画があるか?
ステップS5で小区画P毎のフィルタ処理が終了すれば、次なる小区画PI+1があるか否かを判定する。次なる小区画PI+1があれば次のステップS7に進み、次なる小区画PI+1がなければ走査が終了したものとして一連のフローを終了する。
(ステップS7)次なる小区画の位置に変えて走査
ステップS6で次なる小区画PI+1があると判定された場合には、対象となる小区画Pの位置から次なる小区画PI+1の位置に変えて走査する。そして、次なる小区画PI+1の位置が、今度は対象となる小区画となって、ステップS2に戻って、同様のステップS2〜S7を繰り返す。すなわち、同一の小区画Pに対して、(1)統計量算出部22a(図2を参照)による統計量の算出(ステップS2)、(2)エッジ算出部22b(図2を参照)によるエッジの強度および方向の算出(ステップS3)、(3)フィルタ決定部31(図2を参照)による各々のフィルタの決定(ステップS4)、(4)フィルタ処理部32によるフィルタ処理(ステップS5)の後に、ステップS6の次なる小区画PI+1の存在の判定を経て、(5)次なる小区画PI+1の位置に変えて(ステップS7)、(1)の工程であるステップS2に戻って繰り返すことで、走査することになる。
ステップS7での小区画の走査については、特に限定されないが、図13に示す小区画の走査の各一例のように行うのが好ましい。例えば、図13(a)に示すように、対象となる小区画Pの位置から次なる小区画PI+1の位置に変えて走査するときに、対象となる小区画Pと次なる小区画PI+1との間に走査されない領域が含まれないように走査する。
または、図13(b)に示すように、対象となる小区画Pの位置から次なる小区画PI+1の位置に変えて走査するときに、対象となる小区画Pと次なる小区画PI+1とが接するように走査する。
あるいは、図13(c)に示すように、対象となる小区画Pの位置から次なる小区画PI+1の位置に変えて走査するときに、対象となる小区画Pの一部を次なる小区画PI+1が含んで、対象となる小区画Pと次なる小区画PI+1とが重なるように走査する。
このように、フィルタ処理の演算負荷は、図13(c)、図13(b)、図13(a)の順に軽くなる。逆にフィルタ処理による演算としては、図13(a)、図13(b)、図13(c)の順に精密になる。
本実施例に係るX線透視撮影装置に用いられる画像処理部9によれば、統計量算出部22aは、X線画像を構成する各々の画素に基づいて画素値の統計量(ここでは画素値の平均値avg(i)や揺らぎ)を算出する。その統計量算出部22aによって求められた画素値の統計量に基づいて、エッジ算出部22bは、周囲の画素に対して特異となる箇所の強度をエッジの強度とするとともに、特異となる箇所の方向をエッジの方向として算出する。構造物の輪郭部分などのエッジは周囲の画素に対して特異となるので、統計学的な知見から画素値の統計量に基づいてエッジの強度および方向を算出することが可能である。
ここで、X線画像に対してフィルタ処理を行うためのフィルタを複数(本実施例では3つ)に区分する。そして、上述したエッジ算出部22bによって求められたエッジの強度および方向に基づいて、区分された各フィルタをフィルタ決定部31がそれぞれ変えることで各々のフィルタをそれぞれ決定する。そのフィルタ決定部31によって決定された各フィルタをフィルタ処理部32が複合してX線画像に対してフィルタ処理を行うことで、上述した特許文献1のような2つのフィルタのうちのいずれかに切り換えてフィルタ処理を行う手法とは別の手法で、エッジの強度および方向に基づくフィルタ処理が可能になる。このように、フィルタを複数に区分して、エッジの強度および方向に基づいて、区分された各フィルタをそれぞれ変えることで各々のフィルタをそれぞれ決定するので、エッジの強度および方向に応じてフィルタの決定がより柔軟に細かく行われる。その結果、エッジの強度および方向をより一層考慮して、エッジ部分をボケさせずにノイズ成分を低減させることができる。
本実施例では、X線画像をそれよりも小さな小区画に区分する小区画区分手段の機能と、その区分された小区画の位置を変えることでX線画像内を走査する小区画走査手段の機能とを画像処理部9が備えている。そして、統計量算出部22aは、区分された小区画を構成する各々の画素に基づいて画素値の統計量(ここでは画素値の平均値avg(i)や揺らぎ)を算出し、フィルタ処理部32は小区画に対してフィルタ処理を行う。上述した小区画走査手段の機能は、以下のような走査を行う。
すなわち、同一の小区画Pに対して、(1)統計量算出部22a(図2を参照)による統計量の算出(ステップS2)、(2)エッジ算出部22b(図2を参照)によるエッジの強度および方向の算出(ステップS3)、(3)フィルタ決定部31(図2を参照)による各々のフィルタの決定(ステップS4)、(4)フィルタ処理部32によるフィルタ処理(ステップS5)の後に、ステップS6の次なる小区画PI+1の存在の判定を経て、(5)次なる小区画PI+1の位置に変えて(ステップS7)、(1)の工程であるステップS2に戻って繰り返すことで、走査することになる。このように、小区画に区分して走査するたびにフィルタ処理が行われるので、小区画に応じたフィルタ処理を適正に行うことができる。
また、フィルタ決定部32は、エッジ算出部22bによって求められたエッジの強度および方向に基づいた複数の係数(ここではブレンド率limit,1−limit)を、区分された各フィルタにそれぞれ乗じることで、各フィルタをそれぞれ変えて、各々のフィルタをそれぞれ決定している。各係数に応じて、係数を乗じたフィルタがそれぞれ決定される。
この発明は、上記実施形態に限られることはなく、下記のように変形実施することができる。
(1)上述した実施例では、X線画像処理装置(画像処理部)を用いた装置としてX線透視撮影装置を例に採って説明したが、PET(Positron Emission Tomography)装置やSPECT(Single Photon Emission CT)装置などに代表されるECT(Emission Computed Tomography)装置のように、X線以外の放射線(PET装置の場合にはγ線)を検出して、検出された放射線に基づいて放射線画像を得ることで放射線撮像を行う放射線撮像装置に適用してもよい。
(2)上述した実施例では、図1に示すようなX線透視撮影を例に採って説明したが、この発明は、例えばC型アームに配設されたX線撮像装置にも適用してもよい。また、この発明は、X線CT装置にも適用してもよい。
(3)上述した実施例では、放射線画像(実施例ではX線画像)をそれよりも小さな小区画に区分して走査するたびにフィルタ処理をその都度行ったが、放射線画像全体に対してフィルタ処理を一度で行ってもよい。
(4)上述した実施例では、上記(1)〜(9)式を利用して一連のX線画像処理を行ったが、上記(1)〜(9)式に限定されない。
(5)上述した実施例では、上記(1)〜(4)式のように重み付け(「W」を参照)を行ったが、必ずしも重み付けを行う必要はない。
実施例に係る画像処理部を用いたX線透視撮影装置のブロック図である。 実施例に係る画像処理部のブロック図である。 一連のX線画像処理の流れを示すフローチャートである。 X線画像内の小区画の一態様を示す模式図である。 方位角が0°や90°のときの5×5局所領域の小区画の模式図である。 方位角の定義を示す模式図である。 方位角が0°のときの小区画の模式図である。 方位角が45°のときの小区画の模式図である。 方位角が90°のときの小区画の模式図である。 方位角が135°のときの小区画の模式図である。 スレッシュを設定するための模式図である。 ブレンド率を求めるための模式図である。 (a)〜(c)は、小区画の走査の各一例を示す模式図である。
符号の説明
9 … 画像処理部
22a … 統計量算出部
22b … エッジ算出部
31 … フィルタ決定部
31a … フィルタを等方低域通過フィルタ(等方LPF)
31b … 等方高域通過フィルタ(等方HPF)
31c … 非等方高域通過フィルタ(非等方HPF)
31d,31e … 乗算器
32 … フィルタ処理部
32a … 加算器
P … X線画像
… 小区画
V … エッジの強度
θ … エッジの方向(角度)
limit … ブレンド率

Claims (7)

  1. 放射線画像に対して所定の処理を行う放射線画像処理装置であって、前記放射線画像を構成する各々の画素に基づいて画素値の統計量を算出する統計量算出手段と、その統計量算出手段によって求められた画素値の統計量に基づいて、周囲の画素に対して特異となる箇所の強度をエッジの強度とするとともに、特異となる箇所の方向をエッジの方向として算出するエッジ算出手段と、放射線画像に対してフィルタ処理を行うためのフィルタを複数に区分して、前記エッジ算出手段によって求められたエッジの強度および方向に基づいて、前記区分された各フィルタをそれぞれ変えることで各々のフィルタをそれぞれ決定するフィルタ決定手段と、そのフィルタ決定手段によって決定された各フィルタを複合して放射線画像に対してフィルタ処理を行うフィルタ処理手段とを備えることを特徴とする放射線画像処理装置。
  2. 請求項1に記載の放射線画像処理装置において、前記放射線画像をそれよりも小さな小区画に区分する小区画区分手段と、その小区画区分手段によって区分された小区画の位置を変えることで放射線画像内を走査する小区画走査手段とを備え、前記統計量算出手段は、小区画区分手段によって区分された小区画を構成する各々の画素に基づいて画素値の統計量を算出し、前記フィルタ処理手段は小区画に対してフィルタ処理を行い、前記小区画走査手段は、同一の小区画に対して、(1)統計量算出手段による統計量の算出、(2)前記エッジ算出手段によるエッジの強度および方向の算出、(3)前記フィルタ決定手段による各々のフィルタの決定、(4)フィルタ処理手段によるフィルタ処理を行う(1)〜(4)までの工程を行った後に、(5)次なる小区画の位置に変えて、(1)の工程に戻って繰り返すことで、走査することを特徴とする放射線画像処理装置。
  3. 請求項2に記載の放射線画像処理装置において、前記小区画走査手段は、対象となる小区画の位置から次なる小区画の位置に変えて走査するときに、前記対象となる小区画と次なる小区画との間に走査されない領域が含まれるように走査することを特徴とする放射線画像処理装置。
  4. 請求項2に記載の放射線画像処理装置において、前記小区画走査手段は、対象となる小区画の位置から次なる小区画の位置に変えて走査するときに、前記対象となる小区画と次なる小区画とが接するように走査することを特徴とする放射線画像処理装置。
  5. 請求項2に記載の放射線画像処理装置において、前記小区画走査手段は、対象となる小区画の位置から次なる小区画の位置に変えて走査するときに、前記対象となる小区画の一部を前記次なる小区画が含んで、対象となる小区画と次なる小区画とが重なるように走査することを特徴とする放射線画像処理装置。
  6. 請求項2から請求項5のいずれかに記載の放射線画像処理装置において、前記小区画区分手段によって区分される小区画は方形であることを特徴とする放射線画像処理装置。
  7. 請求項1から請求項6のいずれかに記載の放射線画像処理装置において、前記フィルタ決定手段は、前記エッジ算出手段によって求められたエッジの強度および方向に基づいた複数の係数を、前記区分された各フィルタにそれぞれ乗じることで、各フィルタをそれぞれ変えて、各々のフィルタをそれぞれ決定することを特徴とする放射線画像処理装置。
JP2007161526A 2007-06-19 2007-06-19 放射線画像処理装置 Pending JP2009000164A (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP2007161526A JP2009000164A (ja) 2007-06-19 2007-06-19 放射線画像処理装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2007161526A JP2009000164A (ja) 2007-06-19 2007-06-19 放射線画像処理装置

Publications (1)

Publication Number Publication Date
JP2009000164A true JP2009000164A (ja) 2009-01-08

Family

ID=40317241

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2007161526A Pending JP2009000164A (ja) 2007-06-19 2007-06-19 放射線画像処理装置

Country Status (1)

Country Link
JP (1) JP2009000164A (ja)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2012090978A (ja) * 2010-10-28 2012-05-17 Toshiba Corp 医用画像処理装置および医用画像処理方法
WO2013125325A1 (ja) * 2012-02-23 2013-08-29 富士フイルム株式会社 放射線画像撮影装置、放射線画像撮影システム、放射線画像撮影装置の制御方法、及び放射線画像撮影装置の制御プログラム
CN110298900A (zh) * 2019-06-25 2019-10-01 中山大学 一种基于各向异性基函数建立spect重构的方法

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPS62125481A (ja) * 1985-11-26 1987-06-06 インターナショナル・ビジネス・マシーンズ・コーポレーション パタ−ン認識装置
JPH09138849A (ja) * 1995-07-21 1997-05-27 Philips Electron Nv ディジタル画像のノイズの空間フィルタリング方法及びその装置
JP2001111835A (ja) * 1999-10-05 2001-04-20 Toshiba Corp 画像処理装置およびx線診断装置
JP2001148787A (ja) * 1999-11-19 2001-05-29 Konica Corp 画像処理装置
JP2002092588A (ja) * 2000-09-13 2002-03-29 Canon Inc 画像処理装置及び画像処理方法並びに記憶媒体
JP2002133410A (ja) * 2000-10-25 2002-05-10 Fuji Photo Film Co Ltd ノイズ抑制処理装置並びに記録媒体
JP2005168726A (ja) * 2003-12-10 2005-06-30 Ge Medical Systems Global Technology Co Llc 放射線断層撮影装置、放射線断層撮影方法、画像処理装置、および、エッジ抽出装置

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPS62125481A (ja) * 1985-11-26 1987-06-06 インターナショナル・ビジネス・マシーンズ・コーポレーション パタ−ン認識装置
JPH09138849A (ja) * 1995-07-21 1997-05-27 Philips Electron Nv ディジタル画像のノイズの空間フィルタリング方法及びその装置
JP2001111835A (ja) * 1999-10-05 2001-04-20 Toshiba Corp 画像処理装置およびx線診断装置
JP2001148787A (ja) * 1999-11-19 2001-05-29 Konica Corp 画像処理装置
JP2002092588A (ja) * 2000-09-13 2002-03-29 Canon Inc 画像処理装置及び画像処理方法並びに記憶媒体
JP2002133410A (ja) * 2000-10-25 2002-05-10 Fuji Photo Film Co Ltd ノイズ抑制処理装置並びに記録媒体
JP2005168726A (ja) * 2003-12-10 2005-06-30 Ge Medical Systems Global Technology Co Llc 放射線断層撮影装置、放射線断層撮影方法、画像処理装置、および、エッジ抽出装置

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2012090978A (ja) * 2010-10-28 2012-05-17 Toshiba Corp 医用画像処理装置および医用画像処理方法
US10896486B2 (en) 2010-10-28 2021-01-19 Toshiba Medical Systems Corporation Denoising method and system for preserving clinically significant structures in reconstructed images using adaptively weighted anisotropic diffusion filter
WO2013125325A1 (ja) * 2012-02-23 2013-08-29 富士フイルム株式会社 放射線画像撮影装置、放射線画像撮影システム、放射線画像撮影装置の制御方法、及び放射線画像撮影装置の制御プログラム
CN110298900A (zh) * 2019-06-25 2019-10-01 中山大学 一种基于各向异性基函数建立spect重构的方法

Similar Documents

Publication Publication Date Title
US8938110B2 (en) Enhanced image data/dose reduction
US6879715B2 (en) Iterative X-ray scatter correction method and apparatus
CN103156629B (zh) 图像处理设备和图像处理方法
Colijn et al. Accelerated simulation of cone beam X-ray scatter projections
RU2480821C2 (ru) Устройство обработки информации, способ обработки и считываемый компьютером носитель информации
US9265475B2 (en) Methods and apparatus for scatter correction for CBCT system and cone-beam image reconstruction
US8433154B2 (en) Enhanced contrast for scatter compensation in X-ray imaging
Van Aarle et al. Super-resolution for computed tomography based on discrete tomography
JP6124868B2 (ja) 画像処理装置及び画像処理方法
JP7575511B2 (ja) 放射線撮影により獲得された映像に含まれているノイズを低減するための映像処理方法および装置
US9922409B2 (en) Edge emphasis in processing images based on radiation images
WO2005110232A1 (ja) 画像処理装置及びその方法
JP2010054356A (ja) 画像処理装置及びそれを備えたx線異物検出装置並びに画像処理方法
Aufrichtig et al. X-ray fluoroscopy spatio-temporal filtering with object detection
US8724879B2 (en) Medical image processing apparatus, a medical imaging apparatus, and a method of processing medical images
US20220198651A1 (en) Multi-phase filter
US10891764B2 (en) Image processing apparatus, image processing method, and image processing program
US11010937B2 (en) Image processing apparatus, image processing method, and image processing program
CN103310471A (zh) Ct 图像生成装置及方法、ct 图像生成系统
JP2009000164A (ja) 放射線画像処理装置
JP6002324B2 (ja) 放射線画像生成装置及び画像処理方法
JP6778158B2 (ja) X線ct装置、画像生成方法および画像生成プログラム
JP2010509652A (ja) 画像信号のノイズ低減
EP3773214B1 (en) Cross directional bilateral filter for ct radiation dose reduction
US10950013B2 (en) Image processing apparatus, image processing method, and image processing program

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20091005

A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20110915

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20110927

A02 Decision of refusal

Free format text: JAPANESE INTERMEDIATE CODE: A02

Effective date: 20120306