JP2023162713A

JP2023162713A - Image processing device, image processing method, and program

Info

Publication number: JP2023162713A
Application number: JP2022073275A
Authority: JP
Inventors: 直人高橋; Naoto Takahashi
Original assignee: Canon Inc
Current assignee: Canon Inc
Priority date: 2022-04-27
Filing date: 2022-04-27
Publication date: 2023-11-09

Abstract

To precisely reduce visually conspicuous harmonic components while suppressing the deterioration of image quality.SOLUTION: An image processing device is used to remove, from a radiation image captured using grids for scattered line removal, grid stripes caused by the grids for scattered line removal. The image processing device includes extraction means that extracts a frequency component in a frequency band including at least one of the fundamental wave and n-th harmonic of grid stripes (where n is an integer greater than 1); calculation means that calculates a contrast from the extracted frequency component; determination means that determines a processing target area on the basis of the calculated contrast; and removal means that removes the frequency component from the processing target area in the radiation image.SELECTED DRAWING: Figure 1

Description

開示の技術は、放射線画像を処理する画像処理装置、画像処理方法、およびプログラムに関する。 The disclosed technology relates to an image processing device, an image processing method, and a program that process radiation images.

従来からＸ線に代表される放射線を被写体に照射し、被写体を透過した放射線をフラットパネルディテクタ（以下、ＦＰＤと呼ぶ）にて画像化する技術が医療分野で用いられている。 2. Description of the Related Art Conventionally, a technique has been used in the medical field in which a subject is irradiated with radiation, such as X-rays, and the radiation transmitted through the subject is imaged using a flat panel detector (hereinafter referred to as FPD).

また、放射線は被写体内部で散乱線を発生させるため散乱線を除去する散乱線除去用グリッド（以下、グリッドと略す）と呼ばれる器具を被写体とＦＰＤの間に配置し撮影を行う場合がある。 Furthermore, since radiation generates scattered rays inside the subject, imaging may be performed by placing a device called a scattered ray removal grid (hereinafter abbreviated as grid) between the subject and the FPD to remove the scattered rays.

このグリッドは、鉛等の放射線遮蔽物質と、アルミニウムやカーボン等の放射線透過物質とを、所定の幅で交互に並べて構成することで散乱線を除去するが、画像上にグリッドに起因する周期的な信号（グリッド縞とも言う）を発生させ、観察者の邪魔になる場合がある。 This grid removes scattered radiation by alternately arranging radiation-shielding materials such as lead and radiation-transmitting materials such as aluminum and carbon at a predetermined width, but periodic radiation caused by the grid may appear on the image. This may generate a signal (also called grid stripes) that may disturb the observer.

そこで、従来からこのグリッド縞を除去する方法が提案されている。例えば、特許文献１では基本波（１次高調波とも呼ぶ）よりも高次の高調波成分の強度と近傍の周波数成分の強度比に基づき基本波だけでなく、視覚的に目立つ高次の高調波も除去するものである。 Therefore, methods for removing these grid stripes have been proposed. For example, in Patent Document 1, based on the intensity ratio of harmonic components higher than the fundamental wave (also called first harmonic) and the intensity of nearby frequency components, not only the fundamental wave but also visually noticeable higher harmonics are detected. It also removes waves.

国際公開第２０１７／１３８０９７号International Publication No. 2017/138097 特開２０１４－１５０８４４号公報Japanese Patent Application Publication No. 2014-150844 特開２００２－３３０３４２号公報Japanese Patent Application Publication No. 2002-330342

ところで、特許文献１は近傍の周波数成分のスペクトル強度を被写体成分と見なし、高調波成分のスペクトル強度がこの被写体成分の強度よりも有意に大きい場合に高調波成分を除去するものである。このとき、低い周波数では一般的に、被写体成分のスペクトル強度は高くなる傾向がある。その一方で、低い周波数の高調波成分は強度が低くても視覚的に目立ち易い傾向がある。そのため、高調波成分の周波数が低い場合、強度比基準では高調波成分の検出が困難となり、視覚的に目立つ低周波の高調波成分を精度良く除去できないという課題がある。 By the way, in Patent Document 1, the spectral intensity of nearby frequency components is regarded as a subject component, and when the spectral intensity of a harmonic component is significantly larger than the intensity of this subject component, the harmonic component is removed. At this time, the spectral intensity of the object component generally tends to be high at low frequencies. On the other hand, harmonic components of low frequencies tend to be visually noticeable even if their intensity is low. Therefore, when the frequency of the harmonic component is low, it becomes difficult to detect the harmonic component using the intensity ratio criterion, and there is a problem that visually noticeable low-frequency harmonic components cannot be removed with high precision.

そこで、開示の技術の例示的な目的は、上記の課題を解決するためになされたものであり、画質劣化を抑えつつ視覚的に目立つグリッドの高調波成分を精度良く低減することにある。 Therefore, an exemplary purpose of the disclosed technology is to solve the above-mentioned problems, and is to accurately reduce visually noticeable harmonic components of the grid while suppressing image quality deterioration.

開示の技術の目的を達成するために、開示の技術による例示的な画像処理装置は、
散乱線除去用グリッドを用いて撮影した放射線画像から前記散乱線除去用グリッドに起因して生じるグリッド縞を除去する画像処理装置であって、
前記グリッド縞の基本波またはｎ次高調波（ｎは１より大きい整数）の少なくとも何れか１つを含む周波数帯域の周波数成分を抽出する抽出手段と、
前記抽出した周波数成分からコントラストを算出する算出手段と、
前記算出したコントラストに基づき処理対象領域を決定する決定手段と、
前記放射線画像における前記処理対象領域から前記周波数成分を除去する除去手段と、を備える。 To achieve the objectives of the disclosed technology, an exemplary image processing device according to the disclosed technology includes:
An image processing device that removes grid stripes caused by the scattered ray removal grid from a radiation image photographed using the scattered ray removal grid, the image processing device comprising:
Extracting means for extracting a frequency component of a frequency band containing at least one of the fundamental wave or the n-th harmonic (n is an integer greater than 1) of the grid stripe;
Calculation means for calculating contrast from the extracted frequency components;
determining means for determining a processing target area based on the calculated contrast;
A removal means for removing the frequency component from the processing target area in the radiation image.

開示の技術の一つによれば、画質劣化を抑えつつ視覚的に目立つ高調波成分を精度良く低減することができる。 According to one of the disclosed techniques, it is possible to accurately reduce visually noticeable harmonic components while suppressing image quality deterioration.

実施例１による放射線撮影装置全体の構成図Configuration diagram of the entire radiation imaging apparatus according to Example 1 実施例２による放射線撮影装置全体の構成図Configuration diagram of the entire radiation imaging apparatus according to Example 2 実施例３による放射線撮影装置全体の構成図Configuration diagram of the entire radiographic apparatus according to Example 3 実施例１による処理手順を示すフローチャートFlowchart showing the processing procedure according to Example 1 フィルタ特性の例Example of filter characteristics フィルタリングを説明する図Diagram explaining filtering コントラストの算出方法を説明する図Diagram explaining how to calculate contrast 実施例２による処理手順を示すフローチャートFlowchart showing the processing procedure according to the second embodiment 実施例２の角度検出の処理手順を示すフローチャートFlowchart showing the processing procedure of angle detection in Example 2 角度の検出方法を説明する図Diagram explaining how to detect angle 実施例３による処理手順を示すフローチャートFlowchart showing the processing procedure according to the third embodiment 補間関数のフィルタ特性を説明する図Diagram explaining the filter characteristics of the interpolation function

以下、本開示を実施するための例示的な実施例を、図面を参照して詳細に説明する。ただし、以下の実施例で説明する寸法、材料、形状、及び構成要素の相対的な位置等は任意であり、本開示が適用される装置の構成又は様々な条件に応じて変更できる。また、図面において、同一であるか又は機能的に類似している要素を示すために図面間で同じ参照符号を用いる。 Hereinafter, exemplary embodiments for carrying out the present disclosure will be described in detail with reference to the drawings. However, the dimensions, materials, shapes, relative positions of components, etc. described in the following examples are arbitrary and can be changed depending on the configuration of the device to which the present disclosure is applied or various conditions. Additionally, in the drawings, the same reference numerals are used between the drawings to indicate elements that are identical or functionally similar.

（実施例１）
本実施例は、例えば図１に示すような放射線撮影装置１００に適用される。即ち、放射線撮影装置１００は、放射線撮影により得られた放射線画像からグリッド縞を除去する画像処理機能を有する放射線撮影装置であり、放射線発生部１０１、放射線検出器１０４、データ収集部１０５、前処理部１０６、ＣＰＵ１０８、記憶部１０９、操作部１１０、表示部１１１、画像処理部１１２、後処理部１１８を備えており、これらはＣＰＵバス１０７を介して互いにデータ授受が可能に接続されている。ここで、グリッド縞は、散乱線除去用グリッドに起因して生じる周期的な信号のことである。 (Example 1)
This embodiment is applied to a radiographic apparatus 100 as shown in FIG. 1, for example. That is, the radiographic apparatus 100 is a radiographic apparatus having an image processing function for removing grid stripes from a radiographic image obtained by radiographic imaging, and includes a radiation generating section 101, a radiation detector 104, a data collecting section 105, and a preprocessing section. The image forming apparatus includes a section 106, a CPU 108, a storage section 109, an operation section 110, a display section 111, an image processing section 112, and a post-processing section 118, which are connected to each other via a CPU bus 107 so as to be able to exchange data. Here, the grid stripes are periodic signals generated due to the grid for removing scattered radiation.

また、画像処理部１１２は、放射線検出器１０４で撮影した放射線画像からグリッド縞を除去するものであり、検出部１１３、抽出部１１４、コントラスト算出部１１５、対象領域決定部１１６、除去部１１７を備えている。なお、除去部１１７は、低減手段の一例である。また、「除去」は、「低減」の一例である。 The image processing unit 112 also removes grid stripes from the radiation image taken by the radiation detector 104, and includes a detection unit 113, an extraction unit 114, a contrast calculation unit 115, a target area determination unit 116, and a removal unit 117. We are prepared. Note that the removal unit 117 is an example of a reduction unit. Further, "removal" is an example of "reduction".

上述のような放射線撮影装置１００において、まず、記憶部１０９は、ＣＰＵ１０８での処理に必要な各種のデータを記憶すると共に、ＣＰＵ１０８のワーキングメモリとして機能する。ＣＰＵ１０８は、記憶部１０９を用いて、操作部１１０からの操作に従った装置全体の動作制御等を行う。これにより放射線撮影装置１００は、以下のように動作する。 In the radiographic apparatus 100 as described above, first, the storage unit 109 stores various data necessary for processing by the CPU 108, and functions as a working memory of the CPU 108. The CPU 108 uses the storage unit 109 to control the operation of the entire apparatus according to operations from the operation unit 110. As a result, the radiographic apparatus 100 operates as follows.

まず、操作部１１０を介して操作者が複数の撮影プロトコルの中から所望の１つを選択することで装置への撮影指示がなされる。ここで、撮影プロトコルとは、所望の検査を行う際に使用される一連の動作パラメータセットのことであり、複数の撮影プロトコルを予め作成しておくことで検査に応じた条件設定を簡便に行うことができる。撮影プロトコルの情報としては撮影部位や撮影条件（管電圧、管電流、照射時間など）、グリッドの有無やグリッドの仕様（集束距離、格子比、グリッド密度など）、画像処理パラメータなど様々な設定が紐付けされている。なお、本実施例では撮影プロトコルに紐付けられたグリッドの情報を利用することでグリッド縞の除去を行う。グリッド縞を除去する方法については後述する。 First, an operator selects a desired one from a plurality of imaging protocols via the operation unit 110, thereby instructing the apparatus to perform imaging. Here, an imaging protocol is a series of operating parameter sets used when performing a desired examination, and by creating multiple imaging protocols in advance, conditions can be easily set according to the examination. be able to. Information on the imaging protocol includes various settings such as the area to be imaged, imaging conditions (tube voltage, tube current, irradiation time, etc.), presence or absence of grids, grid specifications (focal distance, grid ratio, grid density, etc.), and image processing parameters. It is tied. Note that in this embodiment, grid stripes are removed by using grid information linked to the imaging protocol. A method for removing grid stripes will be described later.

上述の如く、操作者により入力された撮影指示はＣＰＵ１０８によりデータ収集部１０５に伝えられる。ＣＰＵ１０８は、撮影指示を受けると、放射線発生部１０１及び放射線検出器１０４を制御して放射線撮影を実行させる。 As described above, the imaging instruction input by the operator is transmitted to the data collection unit 105 by the CPU 108. Upon receiving the imaging instruction, the CPU 108 controls the radiation generating unit 101 and the radiation detector 104 to execute radiation imaging.

放射線撮影では、まず放射線発生部１０１が、被写体１０３に対して放射線ビーム１０２を照射する。放射線発生部１０１から照射された放射線ビーム１０２は、被写体１０３を減衰しながら透過して、放射線検出器１０４に到達する。そして、放射線検出器１０４は到達した放射線強度に応じた信号を出力する。なお、本実施例では被写体１０３を人体とする。よって、放射線検出器１０４から出力される信号は人体を撮影したデータとなる。 In radiography, first, the radiation generating unit 101 irradiates the subject 103 with a radiation beam 102 . A radiation beam 102 emitted from a radiation generating unit 101 passes through a subject 103 while being attenuated, and reaches a radiation detector 104 . Then, the radiation detector 104 outputs a signal according to the intensity of the radiation that has arrived. Note that in this embodiment, the subject 103 is a human body. Therefore, the signal output from the radiation detector 104 becomes data obtained by photographing the human body.

データ収集部１０５は、放射線検出器１０４から出力された信号を所定のデジタル信号に変換して放射線画像（以下、画像と略す）として前処理部１０６に供給する。前処理部１０６は、データ収集部１０５から供給された画像に対して、オフセット補正やゲイン補正の前処理を行う。この前処理部１０６で前処理が行われた画像は、ＣＰＵ１０８の制御により、ＣＰＵバス１０７を介して、記憶部１０９、画像処理部１１２に順次転送される。 The data collection unit 105 converts the signal output from the radiation detector 104 into a predetermined digital signal and supplies it to the preprocessing unit 106 as a radiation image (hereinafter abbreviated as image). The preprocessing unit 106 performs preprocessing such as offset correction and gain correction on the image supplied from the data collection unit 105. The image preprocessed by the preprocessing unit 106 is sequentially transferred to the storage unit 109 and the image processing unit 112 via the CPU bus 107 under the control of the CPU 108.

画像処理部１１２は、画像からグリッド縞を除去する画像処理を実行する。また、画像処理部１１２にて処理された画像は、さらに診断に適した画像にするために、後処理部１１８で階調処理や強調処理などの各種処理を実行する。後処理部１１８で処理された画像は表示部１１１にて表示され操作者の確認後、図示しないプリンタ等に出力され一連の撮影動作が終了する。 The image processing unit 112 performs image processing to remove grid stripes from the image. Further, the image processed by the image processing unit 112 is subjected to various processes such as gradation processing and emphasis processing in a post-processing unit 118 in order to make the image more suitable for diagnosis. The image processed by the post-processing unit 118 is displayed on the display unit 111 and, after confirmation by the operator, is output to a printer (not shown) or the like, and the series of photographing operations is completed.

以上のような構成を備えた放射線撮影装置１００において、本実施の形態の特徴である画像処理部１１２の動作、すなわち撮影した画像からグリッド縞を除去する動作に関して、図４のフローチャートを用いて具体的に説明する。なお、本実施例における入力画像は対数変換後のデータ、すなわち線量の対数に比例した画素値をもつデータとする。 In the radiation imaging apparatus 100 having the above configuration, the operation of the image processing unit 112, which is a feature of the present embodiment, that is, the operation of removing grid stripes from the captured image will be explained in detail using the flowchart of FIG. Explain in detail. Note that the input image in this embodiment is data after logarithmic transformation, that is, data having pixel values proportional to the logarithm of the dose.

上述のごとく前処理部１０６によって得られた画像は、ＣＰＵバス１０７を介して画像処理部１１２に転送され、検出部１１３において、グリッドの有無・方向、およびグリッドの基本波の周波数を検出する（ｓ４０１）。 The image obtained by the preprocessing unit 106 as described above is transferred to the image processing unit 112 via the CPU bus 107, and the detection unit 113 detects the presence/absence and direction of the grid and the frequency of the fundamental wave of the grid ( s401).

ここで、グリッドの基本波の検出方法は特に限定するものではないが、例えば、特許文献２などの方法を用いる。この方法は、定めた測定領域のパワースペクトルの比較によりグリッドの有無・方向、および周波数を検出するものであるが、詳細は特許文献２にて公知であるためここでは説明を省略する。 Here, the method of detecting the fundamental wave of the grid is not particularly limited, but for example, the method described in Patent Document 2 is used. This method detects the presence, direction, and frequency of a grid by comparing the power spectra of a defined measurement area, but the details are known from Patent Document 2, so a description thereof will be omitted here.

なお、本実施例では検出部１１３にてグリッドの有無・方向を画像から検出する構成としたが、グリッドの有無・方向が予め既知の場合はその情報を用いても良い。例えば、撮影プロトコルにグリッドの有無・方向を予め紐付けて設定しておき、撮影時にその情報を取得しても良い。また、グリッドの基本波の周波数についても撮影プロトコルに紐付けられたグリッド密度から算出しても良い。具体的には、グリッドの密度をＤ（本／ｃｍ）、ＦＰＤのサンプリングピッチをＳ（ｍｍ）とすれば、グリッドの基本波の周波数ｆｇ（ｒａｄ／ｓａｍｐｌｅ）は下記式にて算出できる。 In this embodiment, the detection unit 113 detects the presence/absence and direction of the grid from the image, but if the presence/absence/direction of the grid is known in advance, that information may be used. For example, the presence/absence and direction of the grid may be linked and set in the imaging protocol in advance, and the information may be acquired at the time of imaging. Furthermore, the frequency of the fundamental wave of the grid may also be calculated from the grid density linked to the imaging protocol. Specifically, if the density of the grid is D (lines/cm) and the sampling pitch of the FPD is S (mm), then the frequency fg (rad/sample) of the fundamental wave of the grid can be calculated using the following formula.

ただし、ｋは下記の条件式を満たす整数である。 However, k is an integer that satisfies the following conditional expression.

なお、グリッドは製造ばらつきやＦＰＤとの距離に応じた拡大が生じるため、上記式による算出では精度誤差が大きくなる。そのため、本実施例のように画像から解析する方がより好適である。 Note that since the grid expands depending on manufacturing variations and the distance from the FPD, calculation using the above formula increases accuracy error. Therefore, it is more suitable to analyze from images as in this embodiment.

次に、検出部１１３において、除去対象のｎ次高調波（ｎは１より大きい整数）の周波数を算出する（ｓ４０２）。具体的には、グリッドの基本波の周波数ｆｇ（ｒａｄ／ｓａｍｐｌｅ）から、ｎ次高調波の周波数ｆｇ，ｎ（ｒａｄ／ｓａｍｐｌｅ）を下記式にて算出する。 Next, the detection unit 113 calculates the frequency of the n-th harmonic (n is an integer greater than 1) to be removed (s402). Specifically, the frequency fg,n (rad/sample) of the n-th harmonic is calculated from the frequency fg (rad/sample) of the fundamental wave of the grid using the following formula.

なお、ｓ４０１において検出結果が「グリッド無し」の場合、本処理はスキップする。 Note that if the detection result in s401 is "no grid", this process is skipped.

次に、検出部１１３において、グリッドを除去するか否かの判定を行う（ｓ４０３）。具体的には、検出結果が「グリッド無し」であれば、処理をする必要がないため以降の処理をスキップし終了する。また、グリッドの有無だけでなく、グリッド密度や求めたｎ次高調波の周波数に基づき判定しても良い。例えば、８０（本／ｃｍ）などの高密度グリッドを使用している場合は、高調波の縞が目立たないため、グリッド密度に応じて処理をスキップさせても良い。また、縞が目立つのは周波数が低い場合であるため、算出した周波数が所定の周波数よりも低い場合のみ処理を実行しても良い。 Next, the detection unit 113 determines whether or not to remove the grid (s403). Specifically, if the detection result is "no grid", there is no need to perform any processing, so the subsequent processing is skipped and the process ends. Further, the determination may be made based not only on the presence or absence of a grid, but also on the grid density and the frequency of the determined n-th harmonic. For example, when a high-density grid such as 80 lines/cm is used, harmonic fringes are not noticeable, so processing may be skipped depending on the grid density. Further, since the stripes are noticeable when the frequency is low, the process may be executed only when the calculated frequency is lower than a predetermined frequency.

次に、抽出部１１４において、ｓ４０４～ｓ４０５の処理を実行することで、ｎ次高調波成分を抽出する。まず、ｓ４０４では求めたｎ次高調波の周波数を通過域とするＢＰＦ（ＢａｎｄＰａｓｓＦｉｌｔｅｒ）を設計する。フィルタの設計方法は特に限定するものではないが、本実施例ではウィンドウ法を用いてフィルタの設計を行う。具体的には、ｎ次高調波の周波数をｆｇ，ｎ（ｒａｄ／ｓａｍｐｌｅ）とし、下記式にてＢＰＦ（ｋ）を設計する。 Next, the extraction unit 114 extracts the n-th harmonic component by executing the processes of s404 to s405. First, in s404, a BPF (Band Pass Filter) whose passband is the frequency of the determined n-th harmonic is designed. Although the filter design method is not particularly limited, in this embodiment, the filter is designed using a window method. Specifically, the frequency of the n-th harmonic is set to fg,n (rad/sample), and BPF(k) is designed using the following formula.

ここで、ｉは虚数を表し、Ｎはフィルタの次数を表す。Δｆは通過域の帯域幅を規定するパラメータであり、所望の周波数帯域が抽出できるように任意に設定すれば良い。また、ｗ（ｋ）がウィンドウ関数を表しており本実施例ではカイザーウィンドウを用いる。なお、カイザーウィンドウについては公知であるため説明は省略する。 Here, i represents an imaginary number and N represents the order of the filter. Δf is a parameter that defines the bandwidth of the passband, and may be set arbitrarily so that a desired frequency band can be extracted. Further, w(k) represents a window function, and in this embodiment, a Kaiser window is used. Note that since the Kaiser window is well known, its explanation will be omitted.

図５（ａ）に設計したＢＰＦのフィルタ特性の例を示す。図５（ａ）はｎ次高調波の周波数を０．１π（ｒａｄ／ｓａｍｐｌｅ）（５０１）、Δｆを０．０５π（ｒａｄ／ｓａｍｐｌｅ）（５０２）としてフィルタを設計したものであり、このＢＰＦによるフィルタリングでｎ次高調波を中心とした周波数成分のみを抽出することが可能である。 FIG. 5(a) shows an example of filter characteristics of the designed BPF. Figure 5(a) shows a filter designed with the frequency of the n-th harmonic as 0.1π (rad/sample) (501) and Δf as 0.05π (rad/sample) (502). It is possible to extract only frequency components centered on the n-th harmonic by filtering.

なお、フィルタリングはグリッド縞と直交する方向に行う。例えば図６のように、グリッド縞が縦方向（６０１）の場合は、図５（ａ）のようにｘ方向にグリッド縞の周期信号が現れるためｘ方向に１次元のフィルタリングを行う。一方、図示しないが、グリッド縞が横方向の場合はｙ方向にグリッド縞の周期信号が現れるためｙ方向に１次元のフィルタリングを行う。 Note that filtering is performed in a direction perpendicular to the grid stripes. For example, when the grid stripes are in the vertical direction (601) as shown in FIG. 6, a periodic signal of the grid stripes appears in the x direction as shown in FIG. 5(a), so one-dimensional filtering is performed in the x direction. On the other hand, although not shown, when the grid stripes are in the horizontal direction, a periodic signal of the grid stripes appears in the y direction, so one-dimensional filtering is performed in the y direction.

次に、ｓ４０５にてグリッド縞と並行な方向にＬＰＦ（ＬｏｗＰａｓｓＦｉｌｔｅｒ）によるフィルタリングを行う。図６（ｂ）に示したようにグリッド縞と並行な方向ではグリッド縞は直流信号となる。そこで、並行な方向にＬＰＦを適用することで、ｎ次高調波以外の被写体成分を除去する。なお、用いるＬＰＦは特に限定するものではないが、例えばボックスフィルタやガウシアンフィルタなどを用いれば良い。 Next, in s405, filtering is performed using an LPF (Low Pass Filter) in a direction parallel to the grid stripes. As shown in FIG. 6(b), the grid stripes become a DC signal in the direction parallel to the grid stripes. Therefore, by applying LPF in parallel directions, object components other than the n-th harmonic are removed. Note that the LPF to be used is not particularly limited, but a box filter, a Gaussian filter, or the like may be used, for example.

次に、コントラスト算出部１１５においてｓ４０６～ｓ４０７を実行することで、抽出したｎ次高調波のコントラストを算出する。まず、ｎ次高調波のｐ－ｐ値（ピークピーク値）を求めるために包絡線を算出する（ｓ４０６）。ここで、包絡線は信号のピーク点に接する曲線を指し、例えば図７（ａ）の点線で示した信号７０１の包絡線は実線７０２で示した曲線となる。また、抽出部１１４で抽出したｎ次高調波は直流成分が除去された信号（平均が０の信号）であるため、この包絡線の値を２倍することでｎ次高調波のｐ－ｐ値を求めることができる。 Next, the contrast calculation unit 115 calculates the contrast of the extracted n-th harmonic by executing steps s406 to s407. First, an envelope is calculated to obtain the pp value (peak-to-peak value) of the n-th harmonic (s406). Here, the envelope refers to a curve that is in contact with the peak point of a signal; for example, the envelope of the signal 701 shown by the dotted line in FIG. 7A becomes the curve shown by the solid line 702. Furthermore, since the n-th harmonic extracted by the extraction unit 114 is a signal from which the DC component has been removed (a signal with an average of 0), by doubling the value of this envelope, the n-th harmonic is You can find the value.

なお、包絡線の算出方法は公知であるため、詳細は省略するがｎ次高調波の信号をｓ（ｔ）、ｓ（ｔ）をヒルベルトへ変換した信号をｈ（ｔ）とすれば、下記式にて包絡線ｒ（ｔ）を算出することができる。 Note that the method for calculating the envelope is well known, so the details will be omitted, but if the n-th harmonic signal is s(t) and the signal obtained by converting s(t) to Hilbert is h(t), then the following The envelope r(t) can be calculated using the formula.

また、具体的な処理としては抽出部１１４で抽出したｎ次高調波は２次元信号であるため、包絡線の算出はグリッド縞と直交する方向にライン単位で実行し、２次元の包絡線画像を生成する。例えば、図６のようにグリッド縞が縦方向の場合は、全てのｙ行に対しｘ方向の包絡線を算出する。一方、図示しないが、グリッド縞が横方向の場合は全てのｘ列に対し、ｙ方向の包絡線を算出する。 Further, as a specific process, since the n-th harmonic extracted by the extraction unit 114 is a two-dimensional signal, the calculation of the envelope is performed line by line in the direction orthogonal to the grid stripes, and the two-dimensional envelope image is generate. For example, if the grid stripes are vertical as shown in FIG. 6, the envelopes in the x direction are calculated for all y rows. On the other hand, although not shown, if the grid stripes are in the horizontal direction, envelopes in the y direction are calculated for all x columns.

次に求めた包絡線画像を用いてｎ次高調波のコントラストを算出する（ｓ４０７）。ここで、Ｘ線の線量に比例した画素値をもつ信号のコントラストはＭｉｃｈｅｌｓｏｎＣｏｎｔｒａｓｔの定義に従うと下記式となる。 Next, the contrast of the nth harmonic is calculated using the obtained envelope image (s407). Here, the contrast of a signal having a pixel value proportional to the X-ray dose is expressed by the following formula according to the definition of Michelson Contrast.

ただし、Ｖｍａｘ，Ｖｍｉｎ，Ｖｅは図７（ｂ）に示した通り、それぞれ、信号の最大値、最小値、平均値を表す。 However, as shown in FIG. 7(b), Vmax, Vmin, and Ve represent the maximum value, minimum value, and average value of the signal, respectively.

なお、本実施例では対数変換後のデータ、すなわち線量の対数に比例した画素値をもつデータを入力としているため、上記（７）式をそのまま適用することはできない。そこで、（７）式を対数変換後のデータに対する式に変形する。具体的には以下の変形式によってコントラストを算出する。 Note that in this embodiment, data after logarithmic transformation, that is, data having pixel values proportional to the logarithm of the dose, is input, so the above equation (7) cannot be applied as is. Therefore, equation (7) is transformed into an equation for data after logarithmic transformation. Specifically, the contrast is calculated using the following modified form.

ここで、Ｌｍａｘは対数変換後のＶｍａｘを表し、Ｌｍｉｎが対数変換後のＶｍｉｎを表している。また、ｆ（）はデータを対数変換した時に用いた対数関数であり、例えば下記式のようなものである。 Here, Lmax represents Vmax after logarithmic transformation, and Lmin represents Vmin after logarithmic transformation. Further, f( ) is a logarithmic function used when data was logarithmically transformed, and is, for example, as shown in the following equation.

ただし、ａ，ｂは任意の定数である。 However, a and b are arbitrary constants.

ここで、（９）式のｘ＝Ｖｅにおける微分係数を（８）式に代入すると、下記式となり対数変換後のデータではＬｍａｘ－Ｌｍｉｎを定数倍することでコントラストを算出することができる。 Here, by substituting the differential coefficient of equation (9) at x=Ve into equation (8), the following equation is obtained.In the data after logarithmic transformation, the contrast can be calculated by multiplying Lmax-Lmin by a constant.

また、Ｌｍａｘ－Ｌｍｉｎはｐ－ｐ値と同義であるため、ｓ４０６で求めた包絡線画像に（０．５／ａ）を乗算すればコントラストを算出することができる。 Furthermore, since Lmax-Lmin is synonymous with pp value, the contrast can be calculated by multiplying the envelope image obtained in s406 by (0.5/a).

次に、対象領域決定部１１６において、ｎ次高調波成分を除去する領域を決定する（ｓ４０９）。本実施例では求めたコントラストを基準に対象領域を決定する。具体的には、各画素におけるコントラストが所定の閾値未満であれば処理対象とするものであり、コントラスト算出部１１５で求めた各画素（ｘ，ｙ）のコントラストをＣ（ｘ，ｙ）、コントラストの閾値をＴｈとすれば、下記式にて処理対象を１とするｆｌａｇ（ｘ，ｙ）を生成する。 Next, the target region determination unit 116 determines a region from which the n-th harmonic component is to be removed (s409). In this embodiment, the target area is determined based on the obtained contrast. Specifically, if the contrast in each pixel is less than a predetermined threshold, it is processed, and the contrast of each pixel (x, y) calculated by the contrast calculation unit 115 is calculated as C(x, y), contrast If the threshold value is Th, flag (x, y) with the processing target set to 1 is generated using the following formula.

ここで、対象領域を決定する閾値Ｔｈはｎ次高調波のコントラストを基準に設定すれば良い。具体的にはｎ次高調波のコントラストは画像内で略一定であり、コントラストが大きい領域は被写体成分を多く含む領域と見なすことができる。そこで、Ｔｈをｎ次高調波のコントラストの２倍に設定すれば、被写体成分が支配的となる領域では処理を実行しない。一方で、ｎ次高調波が支配的となる領域のみに処理を実行することができる。これにより、ｎ次高調波が目立たない被写体成分が支配的な領域での不要な劣化を抑えつつ、視覚的に目立つ領域にのみ処理を適用することができる。 Here, the threshold Th for determining the target area may be set based on the contrast of the n-th harmonic. Specifically, the contrast of the n-th harmonic is approximately constant within the image, and areas with high contrast can be considered to be areas containing many subject components. Therefore, if Th is set to twice the contrast of the n-th harmonic, no processing will be performed in the region where the subject component is dominant. On the other hand, processing can be performed only in the region where the n-th harmonic is dominant. This makes it possible to apply processing only to visually noticeable areas while suppressing unnecessary deterioration in areas where the n-th harmonic is dominated by inconspicuous subject components.

なお、ｎ次高調波のコントラストは使用するグリッドのみを撮影した画像から予め算出することが可能であり、この値に基づきＴｈを設定しておけば良い。 Note that the contrast of the n-th harmonic can be calculated in advance from an image taken of only the grid to be used, and Th may be set based on this value.

次に、除去部１１７において、処理対象となった領域からｎ次高調波を除去する（ｓ４０９）。具体的には、画像処理部１１２に入力された画像をＩｎ（ｘ，ｙ）とし、抽出部１１４で抽出したｎ次高調波成分をＭ（ｘ，ｙ）、対象領域決定部１１６で算出したフラグをｆｌａｇ（ｘ，ｙ）とすれば、下記式にて処理後画像Ｏｕｔ（ｘ，ｙ）を生成する。 Next, the removal unit 117 removes the n-th harmonic from the region to be processed (s409). Specifically, the image input to the image processing unit 112 is In(x,y), and the n-th harmonic component extracted by the extraction unit 114 is M(x,y), which is calculated by the target area determining unit 116. If the flag is flag(x,y), a processed image Out(x,y) is generated using the following formula.

以上、実施例１では、ｎ次高調波のコントラストを基準に処理対象領域を決定することで、被写体成分の不要な劣化を抑えつつ、視覚的に目立つｎ次高調波成分を除去することが可能となる。 As described above, in Example 1, by determining the processing target area based on the contrast of the n-th harmonic, it is possible to remove the visually noticeable n-th harmonic component while suppressing unnecessary deterioration of the subject component. becomes.

なお、実施例１では、一つのｎ次高調波を除去する方法を説明したが、これに限定されるものではない。例えば基本波を含む複数のｎ次高調波（例えば、ｎ＝１，２，３）を除去する場合は、図４のｓ４０２～ｓ４０９をｎの値を変えて繰り返し実行すれば容易に実現可能である。また、公知のグリッド縞の除去技術を組み合わせることも可能である。例えば、ｓ４０１を実行した後に、基本波については、特許文献３などの公知の除去技術を用いて除去し、ｎ次高調波についてはｓ４０２～ｓ４０９を実行するようなことも可能である。 Note that in the first embodiment, a method of removing one n-th harmonic has been described, but the present invention is not limited to this. For example, if you want to remove multiple n-th harmonics (for example, n = 1, 2, 3) including the fundamental wave, this can be easily achieved by repeating steps s402 to s409 in Figure 4 while changing the value of n. be. It is also possible to combine known grid stripe removal techniques. For example, after executing s401, it is possible to remove the fundamental wave using a known removal technique such as Patent Document 3, and execute s402 to s409 for the n-th harmonic.

（実施例２）
本実施例は、例えば図２で示すような放射線撮影装置２００に適用される。この放射線撮影装置２００は、放射線撮影装置１００に対し、角度検出部２０１を備える構成としている。また、本実施例は、画像処理部１１２の処理手順を実施例１とは異なる図８，図９に示したフローチャートに従った動作とする。 (Example 2)
This embodiment is applied to a radiographic apparatus 200 as shown in FIG. 2, for example. This radiographic apparatus 200 is configured to include an angle detection section 201 in contrast to the radiographic apparatus 100. Further, in this embodiment, the processing procedure of the image processing unit 112 is different from that in the first embodiment, and operates according to the flowcharts shown in FIGS. 8 and 9.

ここで、図２は放射線撮影装置２００において、図１の放射線撮影装置１００と同様に動作する箇所は同じ符号を付し、その詳細は省略する。また、図８に示すフローチャートにおいて、図４に示したフローチャートと同様に処理実行するステップは同じ符号を付し、ここでは、上述した実施例１とは異なる構成についてのみ具体的に説明する。 Here, in the radiation imaging apparatus 200 shown in FIG. 2, parts that operate in the same way as the radiation imaging apparatus 100 of FIG. 1 are given the same reference numerals, and the details thereof are omitted. In addition, in the flowchart shown in FIG. 8, the same reference numerals are attached to the steps that perform the same processing as in the flowchart shown in FIG. 4, and only the configurations different from the first embodiment described above will be specifically explained here.

なお、実施例２では実施例１に対しｎ次高調波の角度を検出する手段を備える構成としている。これにより、グリッドが傾いた場合においても精度良くｎ次高調波を抽出することができる。 The second embodiment is different from the first embodiment in that it includes means for detecting the angle of the n-th harmonic. Thereby, even when the grid is tilted, the n-th harmonic can be extracted with high accuracy.

まず、ｓ４０１～ｓ４０４を実施例１と同様に実行することで、ｎ次高調波を抽出する。次に、角度検出部２０１においてｎ次高調波の角度を検出する（ｓ８０１）。ここで、実施例１では図６（ｂ）に示した通り、並行な方向ではグリッド縞は直流成分となることを前提にフィルタリングを行ったが、グリッドが傾いた場合はその仮定が成り立たない。そこで本実施例では、ｓ８０１にてｎ次高調波の角度を検出し、角度に応じてＬＰＦを設計する。以下、図９のフローチャートを用いて角度検出の詳細について説明する。 First, by executing steps s401 to s404 in the same manner as in the first embodiment, the n-th harmonic is extracted. Next, the angle of the n-th harmonic is detected in the angle detection unit 201 (s801). Here, in Example 1, as shown in FIG. 6(b), filtering was performed on the assumption that the grid stripes become DC components in parallel directions, but this assumption does not hold if the grid is tilted. Therefore, in this embodiment, the angle of the n-th harmonic is detected in s801, and the LPF is designed according to the angle. Hereinafter, details of angle detection will be explained using the flowchart of FIG.

まず、ｓ４０６～ｓ４０８を実施例１と同様に実行することで抽出したｎ次高調波から包絡線と処理対象領域を示すフラグを生成する。次に角度を算出するために角度算出画像を生成する（ｓ９０１）。具体的には、ｎ次高調波成分をＭ（ｘ，ｙ）、包絡線画像をＥ（ｘ，ｙ）、処理対象領域を示すフラグをｆｌａｇ（ｘ，ｙ）として、角度算出画像Ａ（ｘ，ｙ）を下記式にて算出する。 First, by executing steps s406 to s408 in the same manner as in the first embodiment, an envelope and a flag indicating the processing target area are generated from the extracted n-th harmonic. Next, an angle calculation image is generated to calculate the angle (s901). Specifically, the angle calculation image A(x, , y) is calculated using the following formula.

なお、ここで求めた角度算出画像Ａ（ｘ，ｙ）は、被写体成分が支配的な領域を除外し、かつ振幅を１に正規化したデータである。これにより被写体成分の影響を受けずにｎ次高調波成分の角度を精度良く検出することが可能となる。 Note that the angle calculation image A(x, y) obtained here is data in which a region where the subject component is dominant is excluded and the amplitude is normalized to 1. This makes it possible to accurately detect the angle of the n-th harmonic component without being affected by the subject component.

次に、角度算出画像を用いてｎ次高調波の角度を検出する（ｓ９０２）。角度の検出方法については特に限定するものではないが、本実施例ではフーリエ変換を用いた角度検出を行う。具体的には、図１０（ａ）に示した通り傾いたｎ次高調波Ａ（ｘ，ｙ）の角度１００１は、図１０（ｂ）の２次元フーリエ変換後の画像Ｈ（ｕ，ｖ）において、｜Ｈ（ｕ，ｖ）｜が最大となる位置１００２と原点１００３を結んだ直線とｕ軸がなす角となる。よって、第一象限または第四象限において｜Ｈ（ｕ，ｖ）｜が最大となる位置を（ｕ０，ｖ０）とすれば、角度θ（ｒａｄ）は以下の式で算出することができる。 Next, the angle of the nth harmonic is detected using the angle calculation image (s902). Although there are no particular limitations on the angle detection method, in this embodiment angle detection is performed using Fourier transform. Specifically, the angle 1001 of the tilted n-th harmonic A(x,y) as shown in FIG. 10(a) corresponds to the two-dimensional Fourier transformed image H(u,v) in FIG. 10(b). , it is the angle formed by the u-axis and the straight line connecting the position 1002 where |H(u,v)| is maximum and the origin 1003. Therefore, if the position where |H(u,v)| is maximum in the first or fourth quadrant is (u0, v0), the angle θ(rad) can be calculated using the following formula.

次に、角度検出で求めたθに基づきＬＰＦの設計を行う（ｓ８０２）。まず、グリッド縞と並行な方向のｎ次高調波の周波数ｆ（ｒａｄ／ｓａｍｐｌｅ）を算出する。ここで、グリッド縞と直交する方向のｎ次高調波の周波数をｆｇ，ｎ（ｒａｄ／ｓａｍｐｌｅ）とすれば、ｆ（ｒａｄ／ｓａｍｐｌｅ）は下記式で求めることができる。 Next, an LPF is designed based on θ obtained by angle detection (s802). First, the frequency f (rad/sample) of the n-th harmonic in the direction parallel to the grid stripes is calculated. Here, if the frequency of the n-th harmonic in the direction perpendicular to the grid stripes is fg,n (rad/sample), then f(rad/sample) can be determined by the following formula.

次に求めた周波数ｆを通過域とするＬＰＦを設計する。フィルタの設計方法は特に限定するものではないが、本実施例ではウィンドウ法を用いてフィルタの設計を行う。具体的には、下記式にてＬＰＦ（ｋ）を設計する。 Next, an LPF whose passband is the determined frequency f is designed. Although the filter design method is not particularly limited, in this embodiment, the filter is designed using a window method. Specifically, LPF(k) is designed using the following formula.

ここで、Ｎはフィルタの次数を表す。Δｆは通過域の帯域幅を規定するパラメータであり、所望の周波数帯域が抽出できるように任意に設定すれば良い。また、ｗ（ｋ）がウィンドウ関数を表しており本実施例ではカイザーウィンドウを用いる。なお、カイザーウィンドウについては公知であるため説明は省略する。 Here, N represents the order of the filter. Δf is a parameter that defines the bandwidth of the passband, and may be set arbitrarily so that a desired frequency band can be extracted. Further, w(k) represents a window function, and in this embodiment, a Kaiser window is used. Note that since the Kaiser window is well known, its explanation will be omitted.

図５（ｂ）に設計したＬＰＦのフィルタ特性の例を示す。図５（ｂ）はｎ次高調波の周波数を０．１π（ｒａｄ／ｓａｍｐｌｅ）（５０３）、Δｆを０．０５π（ｒａｄ／ｓａｍｐｌｅ）（５０４）としてフィルタを設計したものであり、このＬＰＦによるフィルタリングで、グリッドが傾いた場合でもｎ次高調波の周波数成分を劣化させずに被写体成分のみを除去することができる。 FIG. 5(b) shows an example of filter characteristics of the designed LPF. Figure 5(b) shows a filter designed with the frequency of the n-th harmonic as 0.1π (rad/sample) (503) and Δf as 0.05π (rad/sample) (504). By filtering, even if the grid is tilted, only the subject component can be removed without deteriorating the n-th harmonic frequency component.

以上、実施例２ではグリッドの角度に応じてＬＰＦを設計することで、グリッドが傾いた場合でも精度良く縞を除去できる効果がある。 As described above, in the second embodiment, by designing the LPF according to the angle of the grid, it is possible to remove fringes with high precision even when the grid is tilted.

（実施例３）
本実施例は、例えば図３で示すような放射線撮影装置３００に適用される。この放射線撮影装置３００は、放射線撮影装置２００に対し、縮小部３０１、拡大部３０２を備える構成としている。また、本実施例は、画像処理部１１２の処理手順を実施例２とは異なる図１１に示したフローチャートに従った動作とする。 (Example 3)
This embodiment is applied to a radiographic apparatus 300 as shown in FIG. 3, for example. This radiographic apparatus 300 is configured to include a reduction section 301 and an enlargement section 302 in contrast to the radiography apparatus 200. Further, in this embodiment, the processing procedure of the image processing unit 112 is operated according to the flowchart shown in FIG. 11, which is different from that in the second embodiment.

ここで、図３は放射線撮影装置３００において、図２の放射線撮影装置２００と同様に動作する箇所は同じ符号を付し、その詳細は省略する。また、図１１に示すフローチャートにおいて、図８に示したフローチャートと同様に処理実行するステップは同じ符号を付し、ここでは、上述した実施例２とは異なる構成についてのみ具体的に説明する。 Here, in the radiographic apparatus 300 shown in FIG. 3, parts that operate in the same manner as the radiographic apparatus 200 of FIG. 2 are given the same reference numerals, and the details thereof will be omitted. In addition, in the flowchart shown in FIG. 11, the same reference numerals are attached to the steps that perform the same processing as in the flowchart shown in FIG. 8, and only the configuration different from the second embodiment described above will be specifically explained here.

なお、本実施例はｎ次高調波の周波数に応じて画像を縮小する手段を備え、縞の除去にかかる処理時間を短縮することができる。 Note that this embodiment includes a means for reducing the image according to the frequency of the n-th harmonic, and can shorten the processing time required to remove stripes.

まず、ｓ４０１～ｓ４０２を実施例２と同様に実行することで、ｎ次高調波の周波数を算出する。次に、ｓ４０３の判定で処理を実行する場合は、縮小部３０１にて画像の縮小を実行する（ｓ１１０１）。 First, by executing steps s401 to s402 in the same manner as in the second embodiment, the frequency of the n-th harmonic is calculated. Next, if the process is to be executed based on the determination in s403, the reduction unit 301 executes image reduction (s1101).

ｓ１１０１では、最初に縮小率をｎ次高調波の周波数に応じて決定する。具体的にはｎ次高調波の周波数が縮小後のナイキスト周波数よりも高い場合、縮小によりｎ次高調波の成分が除去されることになる。そこで、ｎ次高調波の周波数ｆｇ，ｎ（ｒａｄ／ｓａｍｐｌｅ）がナイキスト周波数よりも低くなるように、下記式にて縮小率ｒａｔｉｏを決定する。 In s1101, the reduction rate is first determined according to the frequency of the n-th harmonic. Specifically, when the frequency of the n-th harmonic is higher than the Nyquist frequency after reduction, the n-th harmonic component is removed by reduction. Therefore, the reduction ratio ratio is determined using the following formula so that the frequency fg,n (rad/sample) of the n-th harmonic is lower than the Nyquist frequency.

ここで、ａは補間関数による信号の減衰を考慮した定数（ただし、０＜ａ＜＝１）であり、縮小に用いる補間関数のフィルタ特性に応じて任意に決定する。 Here, a is a constant (0<a<=1) considering signal attenuation due to the interpolation function, and is arbitrarily determined according to the filter characteristics of the interpolation function used for reduction.

例えば、縮小の補間関数としてＬａｎｃｚｏｓ４を用いる場合、そのフィルタ特性は図１２に示したものとなる。図１２において１１０１が縮小後のナイキスト周波数を表しており、ナイキスト周波数では少なからず信号の減衰が生じる。一方、１１０２よりも低い周波数、すなわちナイキスト周波数×０．７以下の周波数では信号の減衰が起きない特性となっている。そこで、ａを０．７の逆数（１／０．７）以上とすれば、ｎ次高調波を減衰させない縮小率の設定が可能となる。 For example, when Lanczos4 is used as the reduction interpolation function, its filter characteristics are as shown in FIG. In FIG. 12, 1101 represents the Nyquist frequency after reduction, and at the Nyquist frequency, a considerable amount of signal attenuation occurs. On the other hand, at frequencies lower than 1102, that is, frequencies below the Nyquist frequency x 0.7, signal attenuation does not occur. Therefore, by setting a to the reciprocal of 0.7 (1/0.7) or more, it becomes possible to set a reduction rate that does not attenuate the n-th harmonic.

次に、設定した縮小率で画像を縮小し、縮小画像を生成する。なお、縮小によってｎ次高調波の周波数が変わるため、縮小後のｎ次高調波ｆ′ｇ，ｎ（ｒａｄ／ｓａｍｐｌｅ）を下記式にて算出する。 Next, the image is reduced at the set reduction ratio to generate a reduced image. Note that since the frequency of the n-th harmonic changes due to reduction, the n-th harmonic f'g,n (rad/sample) after reduction is calculated using the following formula.

次に、求めた縮小画像とｎ次高調波の周波数ｆ′ｇ，ｎ（ｒａｄ／ｓａｍｐｌｅ）を入力とし、ｓ４０４～ｓ４０９の処理を実施例２と同様に実行し、縮小画像から縞を除去する。 Next, the obtained reduced image and the n-th harmonic frequency f'g,n (rad/sample) are input, and the processes from s404 to s409 are executed in the same manner as in Example 2 to remove stripes from the reduced image. .

次に、拡大部３０２において、縞を縮小した画像を拡大し縮小前の画像から縞を除去する（ｓ１１０２）。具体的には、縮小前の入力画像をＩｎ（ｘ，ｙ）、縞を除去した縮小画像を拡大した拡大画像をＭ′（ｘ，ｙ）、縞を除去する前の縮小画像を拡大した拡大画像をＭ（ｘ，ｙ）として、下記式にて処理後画像Ｏｕｔ（ｘ，ｙ）を生成する。 Next, the enlarging unit 302 enlarges the image with the stripes reduced and removes the stripes from the image before reduction (s1102). Specifically, the input image before reduction is In (x, y), the enlarged image obtained by enlarging the reduced image with stripes removed is M'(x, y), and the enlarged image obtained by enlarging the reduced image before removing stripes is In(x, y). Assuming that the image is M(x,y), a processed image Out(x,y) is generated using the following formula.

以上、実施例３では縮小画像に対し縞の除去を行うことで、処理時間を短縮する効果がある。 As described above, in the third embodiment, removing stripes from the reduced image has the effect of shortening the processing time.

以上、開示の技術の好ましい実施例について説明したが、開示の技術はこれらの実施例に限定されないことはいうまでもなく、その要旨の範囲内で種々の変形及び変更が可能である。 Although preferred embodiments of the disclosed technology have been described above, it goes without saying that the disclosed technology is not limited to these examples, and various modifications and changes can be made within the scope of the gist.

なお、本開示の技術は、以下の構成及び方法を含んでいてもよい。 Note that the technology of the present disclosure may include the following configurations and methods.

（構成１）
散乱線除去用グリッドを用いて撮影した放射線画像から前記散乱線除去用グリッドに起因して生じるグリッド縞を除去する画像処理装置であって、
前記グリッド縞の基本波またはｎ次高調波（ｎは１より大きい整数）の少なくとも何れか１つを含む周波数帯域の周波数成分を抽出する抽出手段と、
前記抽出した周波数成分からコントラストを算出する算出手段と、
前記算出したコントラストに基づき処理対象領域を決定する決定手段と、
前記放射線画像における前記処理対象領域から前記周波数成分を除去する除去手段と、
を備える画像処理装置であってもよい。 (Configuration 1)
An image processing device that removes grid stripes caused by the scattered ray removal grid from a radiation image photographed using the scattered ray removal grid, the image processing device comprising:
Extracting means for extracting a frequency component of a frequency band containing at least one of the fundamental wave or the n-th harmonic (n is an integer greater than 1) of the grid stripe;
Calculation means for calculating contrast from the extracted frequency components;
determining means for determining a processing target area based on the calculated contrast;
removal means for removing the frequency component from the processing target area in the radiation image;
The image processing apparatus may include:

（構成２）
前記算出手段は、前記抽出した周波数成分の包絡線から求まるｐ－ｐ値（ピークピーク値）から算出する構成１に記載の画像処理装置であってもよい。 (Configuration 2)
The image processing device according to configuration 1 may be one in which the calculation means calculates from a pp value (peak-to-peak value) determined from an envelope of the extracted frequency component.

（構成３）
前記決定手段は、前記コントラストが所定の閾値よりも小さい領域を前記処理対象領域として設定する構成１から２の何れか一項に記載の画像処理装置であってもよい。 (Configuration 3)
The determining means may be the image processing device according to any one of configurations 1 to 2, in which a region in which the contrast is smaller than a predetermined threshold value is set as the processing target region.

（構成４）
前記所定の閾値は、グリッドのみを撮影した画像から算出したコントラストに基づき定める構成３に記載の画像処理装置であってもよい。 (Configuration 4)
The image processing device according to configuration 3 may be configured such that the predetermined threshold value is determined based on a contrast calculated from an image of only the grid.

（構成５）
前記散乱線除去用グリッドのグリッド密度と前記撮影した放射線画像のサンプリングピッチに基づき前記グリッド縞の除去を実行するか否かが決定される構成１から４の何れか一項に記載の画像処理装置であってもよい。 (Configuration 5)
The image processing device according to any one of configurations 1 to 4, wherein it is determined whether or not to remove the grid stripes based on the grid density of the scattered radiation removal grid and the sampling pitch of the photographed radiation image. It may be.

（構成６）
前記グリッド縞の基本波を検出する検出手段を更に備え、
前記抽出手段は、前記検出した基本波の周波数に基づきグリッド縞の基本波またはｎ次高調波の少なくとも何れか１つの周波数成分を抽出する構成１から５の何れか一項に記載の画像処理装置であってもよい。 (Configuration 6)
Further comprising a detection means for detecting the fundamental wave of the grid stripes,
The image processing device according to any one of configurations 1 to 5, wherein the extraction means extracts at least one frequency component of the fundamental wave or the n-th harmonic of the grid stripe based on the frequency of the detected fundamental wave. It may be.

（構成７）
前記グリッド縞の基本波が検出できない場合には、前記グリッド縞の除去は実行されない構成６に記載の画像処理装置であってもよい。 (Configuration 7)
The image processing apparatus according to configuration 6 may not perform the removal of the grid stripes when the fundamental wave of the grid stripes cannot be detected.

（構成８）
前記検出した基本波と前記撮影した放射線画像のサンプリングピッチに基づき前記グリッド縞の除去を実行するか否かが決定される構成６又は７に記載の画像処理装置であってもよい。 (Configuration 8)
The image processing apparatus according to configuration 6 or 7, wherein it is determined whether or not to remove the grid stripes based on the detected fundamental wave and the sampling pitch of the photographed radiographic image.

（構成９）
前記抽出手段は、前記グリッド縞の基本波またはｎ次高調波の少なくとも何れか１つの周波数が通過域となるようなフィルタで周波数成分を抽出する構成１から８の何れか一項に記載の画像処理装置であってもよい。 (Configuration 9)
9. The image according to any one of configurations 1 to 8, wherein the extraction means extracts the frequency component using a filter whose passband is at least one frequency of the fundamental wave or the n-th harmonic of the grid stripes. It may also be a processing device.

（構成１０）
前記フィルタは、前記通過域の周波数に応じてウィンドウ法で設計する構成９に記載の画像処理装置であってもよい。 (Configuration 10)
In the image processing device according to configuration 9, the filter may be designed using a window method according to the frequency of the passband.

（構成１１）
前記ウィンドウ法で用いるウィンドウは、カイザーウィンドウである構成１０に記載の画像処理装置であってもよい。 (Configuration 11)
The window used in the window method may be the image processing device according to configuration 10, which is a Kaiser window.

（構成１２）
前記グリッド縞の基本波またはｎ次高調波の少なくとも何れか１つの角度を算出する角度算出手段を更に備え、
前記算出した角度に応じて前記抽出する周波数成分の周波数が変更される構成１から１１の何れか一項に記載の画像処理装置であってもよい。 (Configuration 12)
Further comprising angle calculation means for calculating at least one angle of the fundamental wave or the n-th harmonic of the grid stripes,
The image processing device according to any one of configurations 1 to 11 may be configured such that the frequency of the frequency component to be extracted is changed according to the calculated angle.

（構成１３）
前記放射線画像から縮小画像を生成する縮小手段と、
前記縮小画像を拡大する拡大手段と、を更に備え、
前記縮小画像からグリッドに起因する縞を除去し、前記グリッド縞を除去した縮小画像を拡大した拡大画像に基づき前記放射線画像からグリッド縞を除去する構成１から１２の何れか一項に記載の画像処理装置であってもよい。 (Configuration 13)
Reducing means for generating a reduced image from the radiographic image;
further comprising an enlarging means for enlarging the reduced image,
The image according to any one of configurations 1 to 12, wherein stripes caused by grids are removed from the reduced image, and grid stripes are removed from the radiation image based on an enlarged image obtained by enlarging the reduced image from which the grid stripes have been removed. It may also be a processing device.

（構成１４）
前記縮小手段は、縮小後のナイキスト周波数が前記抽出手段で抽出する周波数よりも高くなる縮小率で縮小する構成１３に記載の画像処理装置であってもよい。 (Configuration 14)
The image processing device according to configuration 13 may be such that the reduction means performs reduction at a reduction rate such that the Nyquist frequency after reduction is higher than the frequency extracted by the extraction means.

（構成１５）
前記縮小手段は、縮小に用いる補間関数と前記抽出手段で抽出する周波数に基づいて縮小率を決定する構成１３に記載の画像処理装置であってもよい。 (Configuration 15)
The image processing device according to configuration 13, wherein the reduction means determines the reduction rate based on an interpolation function used for reduction and a frequency extracted by the extraction means.

（構成１６）
前記縮小手段は、前記抽出手段で抽出する周波数が前記補間関数の通過域に含まれるように縮小率を決定する構成１５に記載の画像処理装置であってもよい。 (Configuration 16)
The image processing device according to configuration 15, wherein the reduction means determines a reduction rate such that the frequency extracted by the extraction means is included in a passband of the interpolation function.

（構成１７）
グリッドに起因して生じるグリッド縞の基本波またはｎ次高調波（ｎは１より大きい整数）の少なくとも何れか１つを含む周波数帯域の周波数成分を用いて、放射線画像における処理対象領域を決定する決定手段と、
前記放射線画像における前記処理対象領域から前記周波数成分を低減する低減手段と、
を備える画像処理装置であってもよい。 (Configuration 17)
A processing target region in a radiation image is determined using a frequency component of a frequency band that includes at least one of the fundamental wave or n-th harmonic (n is an integer greater than 1) of grid stripes caused by the grid. a means of determining;
reduction means for reducing the frequency component from the processing target area in the radiation image;
The image processing apparatus may include:

（方法１）
散乱線除去用グリッドを用いて撮影した放射線画像から前記散乱線除去用グリッドに起因して生じるグリッド縞を除去する画像処理方法であって、
前記グリッド縞の基本波またはｎ次高調波（ｎは１より大きい整数）の少なくとも何れか１つを含む周波数帯域の周波数成分を抽出する抽出工程と、
前記抽出した周波数成分からコントラストを算出する算出工程と、
前記算出したコントラストに基づき処理対象領域を決定する決定工程と、
前記放射線画像における前記処理対象領域から前記周波数成分を除去する除去工程と、
を含む画像処理方法であってもよい。 (Method 1)
An image processing method for removing grid stripes caused by the scattered ray removal grid from a radiation image photographed using the scattered ray removal grid, the method comprising:
an extraction step of extracting a frequency component of a frequency band containing at least one of the fundamental wave or n-th harmonic (n is an integer greater than 1) of the grid stripes;
a calculation step of calculating contrast from the extracted frequency components;
a determining step of determining a processing target area based on the calculated contrast;
a removal step of removing the frequency component from the processing target area in the radiation image;
An image processing method including:

（方法２）
グリッドに起因して生じるグリッド縞の基本波またはｎ次高調波（ｎは１より大きい整数）の少なくとも何れか１つを含む周波数帯域の周波数成分を用いて、放射線画像における処理対象領域を決定する決定工程と、
前記放射線画像における前記処理対象領域から前記周波数成分を低減する低減工程と、
を含む画像処理方法であってもよい。 (Method 2)
A processing target region in a radiation image is determined using a frequency component of a frequency band that includes at least one of the fundamental wave or n-th harmonic (n is an integer greater than 1) of grid stripes caused by the grid. a decision process;
a reduction step of reducing the frequency component from the processing target area in the radiation image;
An image processing method including:

（プログラムや記憶媒体）
方法１又は２に記載の画像処理方法をコンピュータに実行させるプログラムであってもよく、また、該プログラムを記憶する記憶媒体であってもよい。 (programs and storage media)
It may be a program that causes a computer to execute the image processing method described in Method 1 or 2, or it may be a storage medium that stores the program.

（その他の実施例）
また、開示の技術は、以下の処理を実行することによっても実現される。即ち、開示の技術は、上述した様々な実施例の１以上の機能を実現するソフトウェア（プログラム）を、ネットワーク又は記憶媒体を介してシステム又は装置に供給し、そのシステム又は装置のコンピュータ（またはＣＰＵやＭＰＵ等）がプログラムを読み出して実行する処理でも実現可能である。コンピュータは、１つ又は複数のプロセッサー若しくは回路を有し、コンピュータ実行可能命令を読み出し実行するために、分離した複数のコンピュータ又は分離した複数のプロセッサー若しくは回路のネットワークを含みうる。 (Other examples)
Further, the disclosed technology is also realized by performing the following processing. That is, the disclosed technology supplies software (programs) that implements one or more functions of the various embodiments described above to a system or device via a network or a storage medium, and the computer (or CPU) of the system or device. This can also be realized by a process in which a program is read out and executed by a computer, MPU, etc. A computer has one or more processors or circuits and may include separate computers or a network of separate processors or circuits for reading and executing computer-executable instructions.

このとき、プロセッサー又は回路は、中央演算処理装置（ＣＰＵ）、マイクロプロセッシングユニット（ＭＰＵ）、グラフィクスプロセッシングユニット（ＧＰＵ）、特定用途向け集積回路（ＡＳＩＣ）、又はフィールドプログラマブルゲートウェイ（ＦＰＧＡ）を含みうる。また、プロセッサー又は回路は、デジタルシグナルプロセッサ（ＤＳＰ）、データフロープロセッサ（ＤＦＰ）、又はニューラルプロセッシングユニット（ＮＰＵ）を含みうる。 The processor or circuit may include a central processing unit (CPU), microprocessing unit (MPU), graphics processing unit (GPU), application specific integrated circuit (ASIC), or field programmable gateway (FPGA). The processor or circuit may also include a digital signal processor (DSP), a data flow processor (DFP), or a neural processing unit (NPU).

１００放射線撮影装置
１０１放射線出力部
１０２放射線ビーム
１０３被写体
１０４放射線検出器
１０５データ収集部
１０６前処理部
１０７ＣＰＵバス
１０８ＣＰＵ
１０９記憶部
１１０操作部
１１１表示部
１１２画像処理部
１１３検出部
１１４抽出部
１１５コントラスト算出部
１１６対象領域決定部
１１７除去部
１１８後処理部
２０１角度検出部
３０１縮小部
３０２拡大部 100 Radiography apparatus 101 Radiation output section 102 Radiation beam 103 Subject 104 Radiation detector 105 Data collection section 106 Preprocessing section 107 CPU bus 108 CPU
109 Storage unit 110 Operation unit 111 Display unit 112 Image processing unit 113 Detection unit 114 Extraction unit 115 Contrast calculation unit 116 Target area determination unit 117 Removal unit 118 Post-processing unit 201 Angle detection unit 301 Reduction unit 302 Enlargement unit

Claims

散乱線除去用グリッドを用いて撮影した放射線画像から前記散乱線除去用グリッドに起因して生じるグリッド縞を除去する画像処理装置であって、
前記グリッド縞の基本波またはｎ次高調波（ｎは１より大きい整数）の少なくとも何れか１つを含む周波数帯域の周波数成分を抽出する抽出手段と、
前記抽出した周波数成分からコントラストを算出する算出手段と、
前記算出したコントラストに基づき処理対象領域を決定する決定手段と、
前記放射線画像における前記処理対象領域から前記周波数成分を除去する除去手段と、
を備える画像処理装置。 An image processing device that removes grid stripes caused by the scattered ray removal grid from a radiation image photographed using the scattered ray removal grid, the image processing device comprising:
Extracting means for extracting a frequency component of a frequency band containing at least one of the fundamental wave or the n-th harmonic (n is an integer greater than 1) of the grid stripe;
Calculation means for calculating contrast from the extracted frequency components;
determining means for determining a processing target area based on the calculated contrast;
removal means for removing the frequency component from the processing target area in the radiation image;
An image processing device comprising:

前記算出手段は、前記抽出した周波数成分の包絡線から求まるｐ－ｐ値（ピークピーク値）から算出する請求項１に記載の画像処理装置。 The image processing apparatus according to claim 1, wherein the calculation means calculates from a pp value (peak-to-peak value) determined from an envelope of the extracted frequency component.

前記決定手段は、前記コントラストが所定の閾値よりも小さい領域を前記処理対象領域として設定する請求項１に記載の画像処理装置。 The image processing apparatus according to claim 1, wherein the determining means sets an area where the contrast is smaller than a predetermined threshold value as the processing target area.

前記所定の閾値は、グリッドのみを撮影した画像から算出したコントラストに基づき定める請求項３に記載の画像処理装置。 The image processing device according to claim 3, wherein the predetermined threshold value is determined based on a contrast calculated from an image of only the grid.

前記散乱線除去用グリッドのグリッド密度と前記撮影した放射線画像のサンプリングピッチに基づき前記グリッド縞の除去を実行するか否かが決定される請求項１に記載の画像処理装置。 The image processing device according to claim 1, wherein whether or not to remove the grid stripes is determined based on the grid density of the scattered radiation removal grid and the sampling pitch of the photographed radiation image.

前記グリッド縞の基本波を検出する検出手段を更に備え、
前記抽出手段は、前記検出した基本波の周波数に基づきグリッド縞の基本波またはｎ次高調波の少なくとも何れか１つの周波数成分を抽出する請求項１に記載の画像処理装置。 Further comprising a detection means for detecting the fundamental wave of the grid stripes,
The image processing apparatus according to claim 1, wherein the extraction means extracts at least one frequency component of the fundamental wave or the n-th harmonic of the grid stripe based on the frequency of the detected fundamental wave.

前記グリッド縞の基本波が検出できない場合には、前記グリッド縞の除去は実行されない請求項６に記載の画像処理装置。 The image processing apparatus according to claim 6, wherein if the fundamental wave of the grid stripes cannot be detected, the removal of the grid stripes is not performed.

前記検出した基本波と前記撮影した放射線画像のサンプリングピッチに基づき前記グリッド縞の除去を実行するか否かが決定される請求項６に記載の画像処理装置。 The image processing apparatus according to claim 6, wherein it is determined whether or not to remove the grid stripes based on the detected fundamental wave and the sampling pitch of the photographed radiographic image.

前記抽出手段は、前記グリッド縞の基本波またはｎ次高調波の少なくとも何れか１つの周波数が通過域となるようなフィルタで周波数成分を抽出する請求項１に記載の画像処理装置。 The image processing apparatus according to claim 1, wherein the extraction means extracts frequency components using a filter whose passband is at least one frequency of the fundamental wave or the n-th harmonic of the grid stripes.

前記フィルタは、前記通過域の周波数に応じてウィンドウ法で設計する請求項９に記載の画像処理装置。 The image processing apparatus according to claim 9, wherein the filter is designed using a window method according to the frequency of the passband.

前記ウィンドウ法で用いるウィンドウは、カイザーウィンドウである請求項１０に記載の画像処理装置。 The image processing apparatus according to claim 10, wherein the window used in the window method is a Kaiser window.

前記グリッド縞の基本波またはｎ次高調波の少なくとも何れか１つの角度を算出する角度算出手段を更に備え、
前記算出した角度に応じて前記抽出する周波数成分の周波数が変更される請求項１に記載の画像処理装置。 Further comprising angle calculation means for calculating at least one angle of the fundamental wave or the n-th harmonic of the grid stripes,
The image processing device according to claim 1, wherein the frequency of the frequency component to be extracted is changed according to the calculated angle.

前記放射線画像から縮小画像を生成する縮小手段と、
前記縮小画像を拡大する拡大手段と、を更に備え、
前記縮小画像からグリッドに起因する縞を除去し、前記グリッド縞を除去した縮小画像を拡大した拡大画像に基づき前記放射線画像からグリッド縞を除去する請求項１に記載の画像処理装置。 Reducing means for generating a reduced image from the radiographic image;
further comprising an enlarging means for enlarging the reduced image,
2. The image processing apparatus according to claim 1, wherein the grid stripes are removed from the reduced image, and the grid stripes are removed from the radiation image based on an enlarged image obtained by enlarging the reduced image from which the grid stripes have been removed.

前記縮小手段は、縮小後のナイキスト周波数が前記抽出手段で抽出する周波数よりも高くなる縮小率で縮小する請求項１３に記載の画像処理装置。 14. The image processing apparatus according to claim 13, wherein the reduction means performs reduction at a reduction rate such that the Nyquist frequency after reduction is higher than the frequency extracted by the extraction means.

前記縮小手段は、縮小に用いる補間関数と前記抽出手段で抽出する周波数に基づいて縮小率を決定する請求項１３に記載の画像処理装置。 14. The image processing apparatus according to claim 13, wherein the reduction means determines the reduction rate based on an interpolation function used for reduction and a frequency extracted by the extraction means.

前記縮小手段は、前記抽出手段で抽出する周波数が前記補間関数の通過域に含まれるように縮小率を決定する請求項１５に記載の画像処理装置。 16. The image processing apparatus according to claim 15, wherein the reduction means determines a reduction rate such that the frequency extracted by the extraction means is included in a passband of the interpolation function.

グリッドに起因して生じるグリッド縞の基本波またはｎ次高調波（ｎは１より大きい整数）の少なくとも何れか１つを含む周波数帯域の周波数成分を用いて、放射線画像における処理対象領域を決定する決定手段と、
前記放射線画像における前記処理対象領域から前記周波数成分を低減する低減手段と、
を備える画像処理装置。 A processing target region in a radiation image is determined using a frequency component of a frequency band that includes at least one of the fundamental wave or n-th harmonic (n is an integer greater than 1) of grid stripes caused by the grid. a means of determining;
reduction means for reducing the frequency component from the processing target area in the radiation image;
An image processing device comprising:

散乱線除去用グリッドを用いて撮影した放射線画像から前記散乱線除去用グリッドに起因して生じるグリッド縞を除去する画像処理方法であって、
前記グリッド縞の基本波またはｎ次高調波（ｎは１より大きい整数）の少なくとも何れか１つを含む周波数帯域の周波数成分を抽出する抽出工程と、
前記抽出した周波数成分からコントラストを算出する算出工程と、
前記算出したコントラストに基づき処理対象領域を決定する決定工程と、
前記放射線画像における前記処理対象領域から前記周波数成分を除去する除去工程と、
を含む画像処理方法。 An image processing method for removing grid stripes caused by the scattered ray removal grid from a radiation image photographed using the scattered ray removal grid, the method comprising:
an extraction step of extracting a frequency component of a frequency band containing at least one of the fundamental wave or n-th harmonic (n is an integer greater than 1) of the grid stripes;
a calculation step of calculating contrast from the extracted frequency components;
a determining step of determining a processing target area based on the calculated contrast;
a removal step of removing the frequency component from the processing target area in the radiation image;
image processing methods including;

グリッドに起因して生じるグリッド縞の基本波またはｎ次高調波（ｎは１より大きい整数）の少なくとも何れか１つを含む周波数帯域の周波数成分を用いて、放射線画像における処理対象領域を決定する決定工程と、
前記放射線画像における前記処理対象領域から前記周波数成分を低減する低減工程と、
を含む画像処理方法。 A processing target region in a radiation image is determined using a frequency component of a frequency band that includes at least one of the fundamental wave or n-th harmonic (n is an integer greater than 1) of grid stripes caused by the grid. a decision process;
a reduction step of reducing the frequency component from the processing target area in the radiation image;
image processing methods including;

請求項１８又は１９に記載の画像処理方法をコンピュータに実行させるプログラム。 A program that causes a computer to execute the image processing method according to claim 18 or 19.