[go: up one dir, main page]

WO2020241030A1 - 画像処理装置、画像処理方法及びプログラム - Google Patents

画像処理装置、画像処理方法及びプログラム Download PDF

Info

Publication number
WO2020241030A1
WO2020241030A1 PCT/JP2020/014095 JP2020014095W WO2020241030A1 WO 2020241030 A1 WO2020241030 A1 WO 2020241030A1 JP 2020014095 W JP2020014095 W JP 2020014095W WO 2020241030 A1 WO2020241030 A1 WO 2020241030A1
Authority
WO
WIPO (PCT)
Prior art keywords
image
plane distribution
radiation
image processing
acquisition means
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.)
Ceased
Application number
PCT/JP2020/014095
Other languages
English (en)
French (fr)
Inventor
明 佃
貴司 岩下
聡太 鳥居
野田　剛司
晃介 照井
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.)
Canon Inc
Original Assignee
Canon Inc
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 Canon Inc filed Critical Canon Inc
Publication of WO2020241030A1 publication Critical patent/WO2020241030A1/ja
Priority to US17/526,697 priority Critical patent/US12112473B2/en
Anticipated expiration legal-status Critical
Ceased legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/70—Denoising; Smoothing
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/0002—Inspection of images, e.g. flaw detection
    • G06T7/0012—Biomedical image inspection
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10116—X-ray image
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20076—Probabilistic image processing
    • G—PHYSICS
    • G06—COMPUTING OR CALCULATING; COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30004—Biomedical image processing
    • G06T2207/30008—Bone

Definitions

  • the present invention relates to an image processing apparatus, an image processing method and a program.
  • spectral imaging technology which is an imaging technology that uses radiation energy information, obtains thickness images of multiple substances from multiple different radiation energy images, or obtains images of surface density and effective atomic number images. can do.
  • These thickness images, areal density images, and effective atomic number images (an example of "plane distribution related to a substance") are pixels of radiation after transmission, assuming that radiation after material penetration is attenuated as an exponential function such as thickness. It can be obtained by calculating back from the value.
  • Patent Document 1 discloses a technique for improving the image quality of a bone image by smoothing an image of soft tissue and subtracting the image from the accumulated image.
  • an object of the present invention is to provide an image processing technique capable of reducing stochastically distributed noise and acquiring an image with better image quality.
  • the image processing apparatus uses a plurality of radiation images taken with different radiation energies as input images, and outputs an output image showing a planar distribution of substances contained in the input images.
  • the plane distribution acquisition means is provided with the plane distribution acquisition means for acquiring the It is characterized in that the output image is determined so that the probability that the output becomes the value of the plane distribution is maximized.
  • radiation includes not only X-rays but also ⁇ -rays, ⁇ -rays, ⁇ -rays, and various particle beams.
  • FIG. 1 is a diagram showing a configuration example of the radiography apparatus 100 according to the embodiment of the present invention.
  • the radiation imaging device 100 includes a radiation tube 101, an FPD 103 (Flat Panel Detector (radiation detection unit)), a computer 109 (image processing device), an operation panel 104, a synchronization device 118, and an irradiation instruction unit 108.
  • the radiation tube 101 generates radiation based on the radiation irradiation instruction signal transmitted from the synchronization device 118, and irradiates the subject 102 with radiation.
  • the radiation tube 101 may be provided with a filter for removing low-energy radiation. When radiation is emitted from the radiation tube 101 to the subject 102, the radiation passes through the subject 102 and enters the FPD 103.
  • the FPD 103 uses a phosphor that converts radiation into visible light, an image sensor that converts visible light into electric charge and voltage (there are a plurality of elements, and one element corresponds to one pixel) and the converted voltage. Further, it has an image forming unit that converts a digital value into an image, and an image correction unit that corrects the image (offset correction, gain correction, defect correction).
  • image sensor pixels that output a signal corresponding to the incident light are arranged in an array (two-dimensional region).
  • the image forming unit of each pixel converts the radiation converted into visible light by the phosphor into an electric signal and outputs it as an image signal.
  • the present embodiment can also be applied to an FPD in which radiation is directly converted into a voltage, for example, an FPD using cadmium telluride (CdTe) or amorphous selenium (a-Se).
  • the FPD 103 accumulates electric charges based on the image signal, acquires a radiographic image, and outputs the radiographic image to the computer 109.
  • the operation panel 104 can specify the irradiation conditions of the radiation emitted from the radiation tube 101, and the operation panel 104 has a condition designation unit 105, an imaging site designation unit 106, and a condition display unit 107.
  • the condition designation unit 105 can specify the tube voltage, tube current, and irradiation time based on the operation input from the operator via the operation unit such as a button.
  • the operation panel 104 holds an image quality adjustment unit (dial, etc.) for adjusting the image quality, and the operator operates the dial to set the image quality so that an image confirmer such as a doctor can easily confirm the radiation image. Can be done.
  • the imaging region designation unit 106 is associated with the imaging region and the irradiation conditions of radiation in advance, and when the operator specifies the imaging region of the subject, the imaging region designation unit 106 is associated with the designated imaging region. It is possible to call and set the irradiation conditions of the radiation.
  • the condition display unit 107 is configured to display the irradiation conditions of the designated radiation.
  • the irradiation instruction unit 108 is a foot pedal, and when the operator steps on the pedal, an irradiation instruction signal is input to the synchronization device 118, and radiation can be emitted from the radiation tube 101 under the control of the synchronization device 118. Become.
  • the image (radiation image) output from the FPD 103 is transferred to the computer 109.
  • the computer 109 can function as an image processing device that processes a radiation image acquired from the FPD 103, and as a functional configuration, an image acquisition unit 110 that acquires an image (radiation image) from the FPD 103, a prior information acquisition unit 111, It has a plane distribution acquisition unit 112 and a prior information storage unit 113.
  • These functional configurations include, for example, a central processing unit (CPU) of a computer 109 (not shown), a main storage device, a storage device such as a hard disk, a graphic processing unit for high-speed calculation, and a LAN (Local Area Network). It can be configured by hardware such as an adapter.
  • each function can be implemented as software.
  • the central processing unit and the program read from the storage device may be composed of integrated circuits (equivalent electric circuits) or the like as long as they perform similar functions.
  • the image correction unit mounted on the FPD 103 may be provided on the computer 109.
  • the operation panel 104 it is also possible to store the radiation irradiation conditions in the computer 109 in advance. In this case, when the radiation irradiation conditions are specified on the operation panel 104, the irradiation conditions stored in the operation panel 104 may be read from the computer 109 and set.
  • the image acquisition unit 110 acquires an image (radiation image) from the FPD 103.
  • the prior information acquisition unit 111 may include, for example, information on the probability distribution related to radiation (for example, information such as Poisson distribution, normal distribution function, mean value which is a parameter of Poisson distribution or normal distribution, standard deviation, etc.).
  • Tube information tube voltage, tube current, irradiation time, information including presence / absence of radiation filter), information about FPD103 (accumulation time, etc.), imaging site information (for example, site such as chest, head, etc.), subject information (Identification number, name, gender, age, height, weight, etc. of the subject input by the input device 114), contrast agent information (type, concentration, etc. of contrast agent), information on C-arm 116 (material, position, angle, etc.) ), Information on the gantry 117 (position of the gantry, etc.) is included.
  • information on the probability distribution related to radiation for example, information such as Poisson distribution, normal distribution function, mean value which is a parameter of
  • the plane distribution acquisition unit 112 uses a plurality of radiation images taken with different radiation energies as input images, and acquires an output image showing the plane distribution of the substances contained in the input images. The calculation process by the plane distribution acquisition unit 112 will be described later.
  • the input device 114 connected to the computer 109 functions as a user interface for inputting instructions to the computer 109 and the FPD 103.
  • the input device 114 can be configured by, for example, the keyboard of the computer 109.
  • the operator can input the imaged portion information and the subject information by using the input device 114.
  • the imaging site information and the subject information may be input from an external device via a network, for example, through a LAN adapter mounted on the computer 109.
  • the display device 115 connected to the computer 109 can display a radiation image (digital image) received by the computer 109 from the FPD 103 and an image (image processed image) output from the computer 109.
  • a radiation image digital image
  • image image processed image
  • the display device 115 it is possible to express the planar distribution related to the substance as a color, and it is possible to improve the diagnostic imaging performance.
  • the input device 114 and the display device 115 may include a part of the functions of the operation panel 104.
  • the input device 114 may be used to input the radiation irradiation conditions performed through the condition designation unit 105 or the imaging site designation unit 106 of the operation panel 104, or the radiation irradiation conditions such as the input tube voltage may be input.
  • the display may be performed by the display device 115 instead of the condition display unit 107.
  • the synchronization device 118 is connected to the radiation tube 101, the FPD 103, the operation panel 104, the irradiation instruction unit 108, and the computer 109, and irradiates radiation from the state of the FPD 103, the pressed state of the irradiation instruction unit 108, and the processing state of the computer 109. Decide whether or not to accept.
  • the synchronization device 118 functions as an irradiation control unit, and controls irradiation and imaging of radiation based on a determined timing.
  • step S201 the operator of the radiography apparatus 100 inputs the radiation irradiation conditions.
  • the irradiation condition input of radiation can be performed, for example, through the condition designation unit 105 of the operation panel 104 or the imaging site designation unit 106.
  • information such as a subject identification number is also input through the input device 114.
  • step S202 the operator instructs the irradiation of radiation through the irradiation instruction unit 108.
  • step S203 the synchronization device 118 determines whether or not irradiation is possible from the states of the radiation tube 101, the FPD 103, the operation panel 104, the irradiation instruction unit 108, and the computer 109.
  • step S204 the synchronization device 118 causes the display device 115 to display a warning.
  • the process proceeds to step S205.
  • step S205 the advance information acquisition unit 111 acquires information on the probability distribution regarding radiation, tube information, information regarding FPD 103, imaging site information, subject information, and the like as advance information.
  • the number of images to be acquired (for example, two images when separating substances by two energies (high energy and low energy) having different energy levels) is determined according to information such as the number of identification substances.
  • the advance information acquisition unit 111 notifies the synchronization device 118 of the determined number of image acquisitions, information (accumulation time, etc.) in the FPD 103, and the like.
  • step S206 the synchronization device 118 generates a radiation irradiation instruction signal based on the information acquired from the prior information acquisition unit 111, and transmits the radiation irradiation instruction signal to the radiation tube 101 and the FPD 103.
  • the irradiation of radiation by the radiation tube 101 and the image acquisition by the FPD 103 are controlled.
  • the FPD 103 captures a plurality of radiation images (eg, low energy radiation images, high energy radiation images) taken at different radiation energy levels.
  • step S207 when the required number of images to be acquired is acquired by the FPD 103, the shooting is ended in step S207.
  • the image acquisition unit 110 acquires an image (radiation image) from the FPD 103. That is, the image acquisition unit 110 acquires a plurality of radiation images (low-energy radiation image, high-energy radiation image) taken at different radiation energy levels.
  • step S208 the plane distribution acquisition unit 112 calculates the plane distribution related to the substance based on the prior information such as the information of the probability distribution regarding radiation obtained in step S205 and the image acquired in step S206. .. The detailed calculation process of the plane distribution will be described later.
  • step S209 the CPU of the computer 109 stores the calculation result of the plane distribution in the computer 109, and the CPU causes the display device 115 to display the calculation result of the plane distribution.
  • the CPU displays, for example, an image (substance identification image) in which the calculation result of the plane distribution such as the thickness, effective atomic number, and surface density of the separated substances is associated with the color of visible light on the display device 115 in color. It is possible to perform such display control.
  • Image acquisition unit 110 by subtracting the radiation distribution information (X L) from the radiation distribution information (X L + X H), can be obtained radiation distribution information (X H).
  • the low-energy radiation distribution information ( XL ) becomes an image by low-energy radiation (low-energy radiation image)
  • the high-energy radiation distribution information (X H ) becomes an image by high-energy radiation (high-energy radiation).
  • the plane distribution acquisition unit 112 acquires substance property images of substances contained in a plurality of radiation images taken with different radiation energies.
  • the plane distribution acquisition unit 112 can acquire a substance characteristic image such as a substance separation image and a substance identification image by using a plurality of radiation images acquired by one irradiation.
  • the substance-separated image is an image in which a plurality of substances contained in a plurality of radiation images (low-energy radiation image, high-energy radiation image) are separated, and when the subject is represented by two or more specific substances, An image separated into two or more images composed of the thickness of the substance.
  • the substance identification image is an image showing the distribution of effective atomic numbers and surface densities of substances contained in a plurality of radiation images (low-energy radiation image, high-energy radiation image), and represents a subject with one specific substance. When this is done, it means an image decomposed into the distribution of the effective atomic number and surface density of the substance.
  • the plane distribution related to substances includes an image showing the distribution of thickness of separated substances (thickness image), an image showing the distribution of effective atomic numbers for a specific substance (effective atomic number image), and the distribution of areal density.
  • the plane distribution calculation process can be applied when dividing into two types of images (area density) to be shown, or when dividing into three or more types of an effective atomic number image, an area density image, and a thickness image. ..
  • equation 1 may not always hold.
  • the noise is n
  • the relationship between the input image (I), the output image (f (d)) and the noise is as shown in the following equation (2), and the noise n is the input image (I) and the output image (I). It is expressed as a difference from f (d)).
  • noise n The probability distribution followed by noise n can be considered as the probability that the image will be I when the thickness of the substance is d. Therefore, the probability distribution that noise n follows can be expressed as p (I
  • the required thickness can be considered to be the thickness d of the substance having the maximum probability when the input image is I. That is, the thickness of the substance satisfying the following equation 3 equation.
  • the plane distribution acquisition unit 112 determines the output image so that the probability (posterior probability) that the output becomes the value (d) of the plane distribution is maximized. To do.
  • the thickness of the substance to be obtained may be considered as the thickness d of the substance having the maximum logarithm of the probability distribution.
  • the form of the equation can be expressed by the following equations 4 or 5.
  • the plane distribution acquisition unit 112 outputs a plane distribution when the input is the input image (I) from the likelihood and the prior probability that the pixel value of the pixel of interest or the pixel value of the surrounding pixels is input.
  • the output image is determined so that the probability (posterior probability) of the value (d) is maximized.
  • the likelihood indicates the probability that the image corresponding to the input becomes the input image when the input value is the value of the plane distribution (thickness d of the substance).
  • d) is a term called likelihood
  • log prior (d) is a term called prior probability.
  • the extremum when it is the minimum value, it can be obtained by, for example, the steepest descent method.
  • Newton's method and other generally known methods can be used.
  • the steepest descent method When the steepest descent method is applied, the successive approximation solution of d n + 1 can be obtained in the form of the following equation 7 and the solution of equation 6 can be obtained from the obtained successive approximation solution.
  • the plane distribution acquisition unit 112 can acquire an output image by a successive approximation solution by iterative calculation.
  • Equation 7 d n is the thickness of the substance obtained in the nth iteration.
  • ⁇ is a coefficient representing the speed of repetition.
  • Grad is a symbol that represents the derivative (gradient) of a vector. That is, when the input image is I, the thickness d n + 1 of the substance having the maximum probability can be sequentially obtained by repeating the repetition as in the equation (7).
  • FIG. 3 is a flowchart illustrating the flow of the detailed process (process of acquiring the planar distribution of the substance) in step S208 of FIG.
  • the image acquisition unit 110 acquires an image (radiation image) from the FPD 103.
  • the image acquisition unit 110 acquires a plurality of radiation images (low-energy radiation image, high-energy radiation image) captured at different radiation energy levels.
  • Two radiographic images, a low-energy radiographic image and a high-energy radiographic image are the input image I in step S301.
  • step S302 the prior information acquisition unit 111 acquires prior information such as information on the probability distribution regarding radiation obtained in step S205.
  • the plane distribution acquisition unit 112 sets the initial value (d 1 ) before performing the sequential calculation.
  • the initial value to be set can be, for example, the plane distribution (thickness of the substance) assuming that there is no noise in the system. Since the relationship between the input image (I) in the absence of noise and the output image (f (d)) relating to the thickness d of the substance can be expressed as in Equation 1, the plane distribution (thickness of the substance).
  • the initial value d 1 of can be obtained from the following equation (8). Based on the equation (8), the plane distribution acquisition unit 112 performs an iterative operation with the value of the plane distribution assuming that there is no noise as an initial value.
  • the energy subtraction method smoothes an output image (for example, a soft tissue thickness image or a bone thickness image) separated into a plurality of substances contained in a plurality of radiation images taken with different radiation energies.
  • step S304 is the iterative calculation (successive approximation), and the plane distribution acquisition unit 112 performs the successive approximation in the equation 7 in step S305 until the preset number of iterations n is reached (S304-Yes).
  • step S305 the process of the differential term (grad) of the equation 7 may be obtained numerically or analytically. The details of the method of processing the differential term (grad) will be described later.
  • the plane distribution acquisition unit 112 updates the coefficient ⁇ n representing the speed of repetition.
  • the solution in the first few iterations should have a larger difference than the final solution, so the plane distribution acquisition unit 112 sets the coefficient ⁇ n , which represents the speed of the iteration, to a large value (first coefficient value).
  • the plane distribution acquisition unit 112 sets the coefficient ⁇ n representing the speed of the iteration to a small value, that is, a second coefficient value smaller than the first coefficient value.
  • the coefficient ⁇ n representing the speed of repetition is in the form of, for example, k ⁇ ⁇ n ⁇ ⁇ n + 1 .
  • the notation of " ⁇ " means that k ⁇ ⁇ n is substituted into ⁇ n + 1 , and the coefficient k is a real number constant greater than 0 and less than 1 (0 ⁇ k (real number constant) ⁇ 1. ).
  • the first coefficient value is ⁇ n
  • the second coefficient value ( ⁇ n + 1 ) can be obtained by k ⁇ the first coefficient value.
  • Flat distribution obtaining unit 112 depending on the number of iterations of the iterative calculation (successive approximation), the coefficient alpha n representing the rate of repeating updates the setting of the coefficients alpha n to be smaller.
  • step S304-No the plane distribution acquisition unit 112 ends this process.
  • Equation 9 I'/ I 0 is the attenuation rate
  • ⁇ is the linear attenuation coefficient
  • d is the thickness of the substance
  • E is the energy of radiation
  • ⁇ (E) is the spectrum of radiation.
  • I is an input image, which is an image (low energy radiation image, high energy radiation image) taken with two different radiation energies.
  • the plane distribution acquisition unit 112 uses either a normal distribution or a Poisson distribution as the probability distribution of likelihood.
  • Equation 6 The likelihood (logp (I
  • Equation 12 approximated by the Poisson distribution can be used.
  • ⁇ indicates the variance.
  • the plane distribution acquisition unit 112 takes a plurality of radiation images taken with different radiation energies as input images, acquires an output image showing the plane distribution of the substance contained in the input image, and acquires the pixel value of the pixel of interest in the output image and the surroundings. The continuity of the plane distribution is evaluated by the difference from the pixel value of the pixel.
  • Prior probabilities for continuity in the output image are obtained based on the evaluation of the continuity of the plane distribution.
  • the prior probability term (logp prior (d)) shown in Equation 6 for example, a term for comparing (difference) between the pixel of interest and the surrounding pixels can be introduced.
  • the following equation 13 can be used as a term of prior probability.
  • the plane distribution acquisition unit 112 acquires the prior probability based on the difference between the pixel value of the pixel of interest in the output image and the pixel value of the surrounding pixels with respect to the pixel of interest.
  • " indicates an absolute value.
  • the plane distribution acquisition unit 112 can use the pixel value of the pixel of interest and the pixel value of the surrounding pixels in the output image (substance-separated image).
  • the spatial continuity of the plane distribution in the output image can be evaluated by comparison (difference) with.
  • the plane distribution acquisition unit 112 performs noise reduction processing on the output image based on the evaluation result of the spatial continuity of the plane distribution. As a result, the noise that is stochastically distributed is reduced, and it becomes possible to acquire an image with better image quality.
  • FIG. 4 is a diagram schematically showing the distribution of the difference between the pixel values and the prior probability.
  • the difference between the pixel values is small, such as the pixel of interest 401 and the surrounding pixels 402
  • the prior probability tends to increase.
  • the difference between the pixel values is large, such as the pixel of interest 401 and the surrounding pixels 403, the prior probability tends to be small.
  • the prior probability term (regularization term) of the equation 13 becomes larger as the evaluation value of continuity (difference between the pixel value of the pixel of interest and the pixel value of the surrounding pixels) is smaller. Means.
  • the output of the prior probability term can be increased.
  • the assumption introduced with prior probability (preliminary assumption that the output image is smooth due to the continuity of pixel values) has a strong influence.
  • Regularization coefficient mu i, eta i is more certain when the large value, the difference is reduced in the pixel value between the target pixel and the surrounding pixels by the acquisition processing of the output image (successive approximation operation), a smoother output image ( Material separation image) can be obtained. Since smooth output images tend to be obtained, the noise of the plane distribution (substance separation image showing the distribution of the thickness of the substance) related to the substance should be reduced by using the term of prior probability (regularization term). Can be done.
  • the plane distribution acquisition unit 112 has a regular coefficient ⁇ i , which represents the degree (strength) of the influence of the prior probability term (regularization term) based on the shooting conditions included in the prior information acquired by the prior information acquisition unit 111.
  • ⁇ i can be changed. For example, when the irradiation condition (dose) of radiation is smaller than the predetermined threshold value, the plane distribution acquisition unit 112 changes the regular coefficients ⁇ i and ⁇ i to values larger than the predetermined reference value to obtain the prior probability.
  • noise removal noise reduction
  • the plane distribution acquisition unit 112 can acquire information on the radiation irradiation condition (dose) based on the prior information acquired by the prior information acquisition unit 111.
  • the tube information (tube voltage, tube current, irradiation time) included in the prior information is information related to the dose
  • the plane distribution acquisition unit 112 is a term of the prior probability based on the radiation irradiation condition (dose).
  • the regular coefficients ⁇ i and ⁇ i which represent the degree of influence of (regularization term), can be changed.
  • the plane distribution acquisition unit 112 has a prior probability term (regularization term) based on the value indicated by the adjustment device. )
  • Regular coefficients ⁇ i and ⁇ i can be adjusted.
  • the subject information (height, weight, etc.) included in the prior information is information related to the body shape of the subject and information related to the dose.
  • the plane distribution acquisition unit 112 has a regular coefficient ⁇ i , which indicates the degree of influence of the prior probability term (regularization term) based on the subject information (information on the body shape of the subject) as information on the radiation irradiation condition (dose).
  • ⁇ i can be changed. For example, the heavier the body weight of the subject, the lower the transmitted dose. Therefore, the plane distribution acquisition unit 112 changes the regular coefficients ⁇ i and ⁇ i to values larger than the predetermined reference values, and the term of prior probability (regular). By increasing the coefficient), noise removal (noise reduction) can be set with higher priority.
  • the information on the position setting of the C arm 116 and the gantry 117 is also information related to the dose. For example, when the distance through which radiation passes through the subject 102 becomes longer depending on the position and angle of the C arm 116, that is, the relative position or angle between the radiation generator (radiation tube 101) and the radiation detection unit (FPD103). It is expected that the permeated dose will be relatively lower than when the distance is short. Further, when the distance through which the radiation passes through the subject is short, it is expected that the transmitted dose will be relatively higher than when the distance is long.
  • the plane distribution acquisition unit 112 has a prior probability term (pre-probability term) based on the relative position or angle between the radiation generator (radiation tube 101) and the radiation detection unit (FPD103) as information on the radiation irradiation condition (dose).
  • pre-probability term based on the relative position or angle between the radiation generator (radiation tube 101) and the radiation detection unit (FPD103) as information on the radiation irradiation condition (dose).
  • the regular coefficients ⁇ i and ⁇ i which represent the degree of influence of the regularization term), can be changed.
  • the gantry 117 is one of the elements that block radiation
  • the presence or absence of the gantry 117 is also information related to the dose
  • the plane distribution acquisition unit 112 provides information on the gantry 117 as information on the radiation irradiation condition (dose).
  • the regular coefficients ⁇ i and ⁇ i which represent the degree of influence of the prior probability term (regularization term), can be changed.
  • the inventor may introduce a prior probability term (regularization term) into a substance separation image that shows the distribution of the thickness of a substance particularly related to soft tissues such as fat (thickness image).
  • a prior probability term for the thickness image related to the soft tissue
  • the regular coefficient of the prior probability term for the thickness image related to the soft tissue
  • the prior coefficient term larger than the regular coefficient of the prior probability term (regularization term)
  • noise is removed (noise reduction) and the image quality is better.
  • the plane distribution acquisition unit 112 evaluates the continuity in the soft tissue thickness image based on the difference between the pixel value of the pixel of interest and the pixel values of the surrounding pixels in the soft tissue thickness image. Then, the plane distribution acquisition unit 112 performs noise reduction processing on the output image based on the evaluation result of the spatial continuity of the plane distribution. As a result, the noise that is stochastically distributed is reduced, and it becomes possible to acquire an image with better image quality.
  • the regular coefficient of the prior probability term (regularization term) for the thickness image of the hard tissue such as bone is used as the soft tissue such as fat.
  • the thickness may be set larger than the regular coefficient of the prior probability term (regularization term) for the image.
  • equations 14 and 15 can be expressed by the following equations 16.
  • the integral of the function of the following equation 17 including the sigmoid function for example, the function in which the derivative is continuous instead of the term containing the absolute value, is used. Discontinuity can be avoided by using it instead of the term shown in the equation.
  • the prior probability term has described an example in which the difference between the pixel of interest and the surrounding pixels (for example, the equation 13), but the prior probability term (regularization term) is not necessarily the focus. It is not always the difference between the pixel and the surrounding pixels.
  • the plane distribution acquisition unit 112 can also acquire the prior probability based on the pixel value of the pixel of interest in the output image (thickness image of the substance). For example, since the material thickness image (d 1 or d 2 ) has a finite value of zero or more, it is possible to introduce a regularization term as shown in the following equation 18 with the material thickness as an input. it can.
  • Equation 18 is a positive constant and is an example of a regular coefficient.
  • the prior probability term (regularization term) shown in Eq. 18 also becomes smaller than zero, and therefore the likelihood and the prior (d 1i, j , d 2i, j ).
  • the total posterior probability including the probability term (regularization term) becomes smaller.
  • Material thickness images (d 1 , d 2 ) are more likely to be images with pixel values greater than zero. Since this prior probability also has the effect of limiting the thickness image that should have varied between positive and negative values to only positive values, it is similar to the prior probability term (regularization term) that has a term of difference from surrounding pixels. It is an example of a term for evaluating the continuity of a plane distribution, and also gives a result of reducing noise.
  • the plane distribution acquisition unit 112 evaluates the continuity of the plane distribution based on the pixel values of the pixels of interest in the output image (thickness image of the substance), and reduces the noise in the output image based on the evaluation result. As a result, the noise that is stochastically distributed is reduced, and it becomes possible to acquire an image with better image quality.
  • Equation 11 is an equation that approximates the probability distribution of likelihood indicating the probability that an image becomes I when the thickness of the substance is d with a normal distribution. However, when the logarithm of Equation 11 is taken, Equation 11 is It can be expressed as the following equation.
  • logp prior (d) can be expressed by the following equation 20.
  • Equation 12 is an equation that approximates the probability distribution of likelihood with a Poisson distribution, but if the logarithm of Equation 12 is taken, Equation 12 can be expressed as Equation 22.
  • logp prior (d) can be expressed by the following equation 23.
  • the plane distribution acquisition unit 112 acquires an output image by approximating the radiation spectra corresponding to different radiation energies as a monochromatic radiation spectrum.
  • the number 9 or number 10 becomes the following number 26, and the number 26 can be expressed as the number 27 when the coordinates of the thickness d of the substance are taken into consideration. ..
  • each component of the equation 21 can be expressed as the following equations 30 and 31.
  • each component of the equation 24 relating to the Poisson distribution can be obtained when the monochromatic approximation is introduced.
  • the monochromatic approximation of the radiation spectrum has been described here, the reproducibility of the normal distribution and Poisson distribution, that is, the independent random variables X 1 and X 2 follow the normal distribution and Poisson distribution (parameters ⁇ 1 and ⁇ 2 ), respectively.
  • the sum of random variables X 1 + X 2 also follows a normal distribution or Poisson distribution
  • the plane distribution acquisition unit 112 acquires an output image by adding the pixel values of the acquired images based on the radiation spectra corresponding to different radiation energies, and obtains the output image with respect to the acquired output image of the present embodiment.
  • the process can be applied.
  • the present invention supplies a program that realizes one or more functions of the above-described embodiment to a system or device via a network or storage medium, and one or more processors in the computer of the system or device reads and executes the program. It can also be realized by the processing to be performed. It can also be realized by a circuit (for example, ASIC) that realizes one or more functions.
  • a circuit for example, ASIC
  • 100 Radiation imaging device, 101: Radiation tube (radiation generator), 103: FPD (radiation detection unit): 104: Operation panel, 108: Irradiation instruction unit, 109: Computer (image processing device), 110: Image acquisition Unit, 111: Advance information acquisition unit, 112: Plane distribution acquisition unit, 113: Advance information storage unit, 114: Input device, 115: Display device, 116: C arm, 117: Stand, 118: Synchronous device

Landscapes

  • Engineering & Computer Science (AREA)
  • Health & Medical Sciences (AREA)
  • Medical Informatics (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Radiology & Medical Imaging (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Quality & Reliability (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Biophysics (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Optics & Photonics (AREA)
  • Pathology (AREA)
  • Biomedical Technology (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • Animal Behavior & Ethology (AREA)
  • Public Health (AREA)
  • Veterinary Medicine (AREA)
  • Apparatus For Radiation Diagnosis (AREA)
  • Image Processing (AREA)

Abstract

画像処理装置は、異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、入力画像に含まれる物質について平面分布を示す出力画像を取得する平面分布取得部を備え、平面分布取得部は、尤度と、注目画素の画素値または周囲の画素の画素値を入力とする事前確率から、入力が入力画像であるときに、出力が平面分布の値となる確率が最大となるように出力画像を決定する。

Description

画像処理装置、画像処理方法及びプログラム
　本発明は、画像処理装置、画像処理方法及びプログラムに関する。
　近年、放射線のエネルギー情報を利用した画像化技術であるスペクトラル・イメージング技術では複数の異なる放射線エネルギー画像から複数の物質の厚さ画像を求めたり、あるいは面密度と実効原子番号画像の画像を求めたりすることができる。こうした厚さ画像、面密度画像、実効原子番号画像(「物質にかかわる平面分布」の一例)は物質透過後の放射線が厚さなどの指数関数として減衰すると仮定して、透過後の放射線の画素値から逆算することで求めることができる。
　特許文献1では軟部組織の画像を平滑化し、その画像を蓄積画像から減算することで、骨部画像の画質を改善する技術が開示されている。
特開平３－２８５４７５号公報
　しかしながら、物質にかかわる平面分布の推定演算においては、システムノイズや放射線に由来する量子ノイズが計算過程で伝播し、出力された画像はノイズを含んだものとなり得る。特許文献１の技術においても、単純な四則演算を用いて平面分布の推定演算をしており、放射線の線量が低くなる場合などにおいては、ノイズを除去しきれず、また、平滑化により画像が不鮮明になる場合が生じ得る。
　本発明は、上記の課題に鑑み、確率的に分布するノイズを低減し、より優れた画質の画像を取得することが可能な画像処理技術を提供することを目的とする。
　本発明の目的を達成するために、本発明の一態様による画像処理装置は、異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、前記入力画像に含まれる物質について平面分布を示す出力画像を取得する平面分布取得手段を備え、前記平面分布取得手段は、尤度と、注目画素の画素値または周囲の画素の画素値を入力とする事前確率から、入力が前記入力画像であるときに、出力が前記平面分布の値となる確率が最大となるように前記出力画像を決定することを特徴とする。
　本発明によれば、ノイズを低減し、より優れた画質の画像を取得することが可能になる。
 　本発明のその他の特徴及び利点は、添付図面を参照とした以下の説明により明らかになるであろう。なお、添付図面においては、同じ若しくは同様の構成には、同じ参照番号を付す。
　添付図面は明細書に含まれ、その一部を構成し、本発明の実施形態を示し、その記述と共に本発明の原理を説明するために用いられる。
実施形態に係る放射線撮影装置の構成例を示す図。 実施形態に係る撮影処理の流れを説明する図。 物質の平面分布を取得する処理の流れを説明する図。 画素値の差分と事前確率の分布を模式的に示す図。
　以下、添付図面を参照して実施形態を詳しく説明する。なお、以下の実施形態は特許請求の範囲に係る発明を限定するものではない。実施形態には複数の特徴が記載されているが、これらの複数の特徴の全てが発明に必須のものとは限らず、また、複数の特徴は任意に組み合わせられてもよい。さらに、添付図面においては、同一若しくは同様の構成に同一の参照番号を付し、重複した説明は省略する。以下の実施形態及び特許請求の範囲において、放射線は、Ｘ線の他、α線、β線、γ線、及び各種粒子線なども含む。
　図１は、本発明の実施形態に係る放射線撮影装置１００の構成例を示す図である。放射線撮影装置１００は、放射線管球１０１、ＦＰＤ１０３（Flat　Panel　Detector（放射線検知部））、コンピュータ１０９（画像処理装置）、操作盤１０４及び同期装置１１８、照射指示部１０８を有する。
　放射線管球１０１は、同期装置１１８から送信された放射線照射の指示信号に基づいて、放射線を発生させ、被写体１０２に放射線を照射する。放射線管球１０１には、低エネルギーの放射線を除去するためのフィルタを設けることがある。放射線管球１０１から放射線が被写体１０２に照射されると、放射線は被写体１０２を透過してＦＰＤ１０３に入射する。
　ＦＰＤ１０３は、放射線を可視光に変換する蛍光体と、可視光を電荷、電圧に変換する撮像素子（素子は複数個あり、1つの素子と1つの画素が対応する）と、変換された電圧をさらにデジタル値の画像へ変換する画像形成部と、画像に補正（オフセット補正、ゲイン補正、欠陥補正）を行う画像補正部を有する。撮像素子は、入射光に応じた信号を出力する画素がアレイ状（二次元の領域）に配置されている。各画素の画像形成部は蛍光体により可視光に変換された放射線を電気信号に変換し、画像信号として出力する。尚、放射線を電圧に直接変換する方式のＦＰＤ、例えば、テルル化カドミウム(ＣｄＴｅ)や非晶質セレン(ａ-Ｓｅ)を用いたＦＰＤにおいても本実施形態を適用することが可能である。ＦＰＤ１０３は画像信号に基づく電荷の蓄積を行って放射線画像を取得し、放射線画像をコンピュータ１０９に出力する。
　操作盤１０４は、放射線管球１０１から照射する放射線の照射条件を指定することが可能であり、操作盤１０４は、条件指定部１０５、撮影部位指定部１０６、条件表示部１０７を有する。条件指定部１０５は、ボタン等の操作部を介した操作者からの操作入力に基づいて、管電圧、管電流、照射時間を指定可能である。操作盤１０４は画質調整を行うための画質調整部(ダイヤルなど)を保持しており、操作者はダイヤルを操作することで医師などの画像確認者が放射線画像を確認しやすい画質に設定することができる。
　撮影部位指定部１０６には、撮影部位と放射線の照射条件が予め関連付けられていて、操作者は、被写体の撮影部位を指定すると、撮影部位指定部１０６は、指定されたさ撮影部位に関連付けられている放射線の照射条件を呼び出し、設定することができるようになっている。条件表示部１０７には、指定されている放射線の照射条件が表示されるように構成されている。
　照射指示部１０８はフットペダルで、操作者はペダルを踏むことで、照射指示信号が同期装置１１８に入力され、同期装置１１８の制御に基づいて放射線管球１０１から放射線を照射することが可能になる。
　ＦＰＤ１０３から出力された画像（放射線画像）は、コンピュータ１０９へ転送される。コンピュータ１０９は、ＦＰＤ１０３から取得した放射線画像を処理する画像処理装置として機能することが可能であり、機能構成として、ＦＰＤ１０３から画像（放射線画像）を取得する画像取得部１１０、事前情報取得部１１１、平面分布取得部１１２、事前情報保存部１１３を有している。これらの機能構成は、例えば、図示されていないコンピュータ１０９の中央演算処理装置（ＣＰＵ）、主記憶装置、ハードディスクなどの記憶装置、高速計算のためのグラフィック・プロセッシング・ユニット、ＬＡＮ(Local Area Network)アダプタなどのハードウェアにより構成することが可能である。
　また、各機能はソフトウェアとして実装することが可能である。例えば、中央演算処理装置と、記憶装置から読み込んだプログラムを用いて、各部の機能を構成することが可能である。同様の機能を果たすのであれば、それらは集積回路（同等の電気回路）などで構成してもよい。尚、ＦＰＤ１０３に搭載されている画像補正部を、コンピュータ１０９に設けてもよい。また、操作盤１０４の代わりに、コンピュータ１０９に予め放射線の照射条件を保存しておくことも可能である。この場合、操作盤１０４で放射線の照射条件を指定する際に、操作盤１０４が保存されている照射条件をコンピュータ１０９から読み出して、設定すればよい。
　画像取得部１１０は、ＦＰＤ１０３から画像（放射線画像）を取得する。また、事前情報取得部１１１は、事前情報として、例えば、放射線に関する確率分布の情報（例えばポアソン分布、正規分布の関数やポアソン分布や正規分布のパラメタである平均値、標準偏差などの情報）、管球情報(管電圧、管電流、照射時間、放射線フィルタの有無を含む情報)、ＦＰＤ１０３に関する情報(蓄積時間など)、撮影部位情報(例えば、被写体の胸部、頭部などの部位)、被写体情報(入力装置１１４で入力される被写体の識別番号、氏名、性別、年齢、身長、体重など)、造影剤情報(造影剤の種類、濃度など)、Ｃアーム１１６の情報(材質、位置、角度など)、架台１１７の情報(架台の位置など)が含まれる。
　平面分布取得部１１２は、異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、前記入力画像に含まれる物質について平面分布を示す出力画像を取得する。平面分布取得部１１２による計算プロセスは後に説明する。
　コンピュータ１０９に接続する入力装置１１４は、コンピュータ１０９やＦＰＤ１０３に対する指示を入力するユーザーインターフェイスとして機能する。入力装置１１４は、例えば、コンピュータ１０９のキーボードなどにより構成することが可能である。操作者は、入力装置１１４を使用して、撮影部位情報や被写体情報を入力することが可能である。撮影部位情報や被写体情報は、入力装置１１４による入力の他に、例えば、コンピュータ１０９に実装されるＬＡＮアダプタを通じ、ネットワークを介して外部の装置から入力してもよい。
　また、コンピュータ１０９に接続する表示装置１１５は、コンピュータ１０９がＦＰＤ１０３から受信した放射線画像（デジタル画像）やコンピュータ１０９から出力された画像（画像処理された画像）を表示することが可能である。表示装置１１５として、特にカラー・ディスプレイが用いられた場合は、物質にかかわる平面分布を色として表現することが可能になり、画像診断性能の向上を図ることが可能になる。
　また、操作盤１０４の機能の一部を入力装置１１４、表示装置１１５が備えてもよい。例えば、操作盤１０４の条件指定部１０５あるいは撮影部位指定部１０６を通じて行った放射線の照射条件の入力操作を入力装置１１４にて行ってもよいし、入力された管電圧などの放射線の照射条件の表示を、条件表示部１０７の代わりに表示装置１１５で行ってもよい。
　同期装置１１８は、放射線管球１０１、ＦＰＤ１０３、操作盤１０４、照射指示部１０８及びコンピュータ１０９に接続され、ＦＰＤ１０３の状態と、照射指示部１０８の押下状態と、コンピュータ１０９の処理状態から放射線の照射可否を決定する。同期装置１１８は照射制御部として機能し、決定したタイミングに基づいて放射線の照射及び撮影を制御する。
　本発明の実施形態に係る撮影処理の流れを図２のフローチャートを参照して説明する。放射線撮影装置１００の操作者は、ステップＳ２０１において、放射線の照射条件を入力する。放射線の照射条件入力は、例えば、操作盤１０４の条件指定部１０５あるいは撮影部位指定部１０６を通じて行うことが可能である。また、入力装置１１４を通じて、被写体の識別番号などの情報も併せて入力する。
　ステップＳ２０２において、操作者は、照射指示部１０８を通じて放射線の照射を指示する。その後、ステップＳ２０３において、同期装置１１８は、放射線管球１０１、ＦＰＤ１０３、操作盤１０４、照射指示部１０８、コンピュータ１０９の状態から照射可否を判断する。照射不可の場合（Ｓ２０３－Ｎｏ）、ステップＳ２０４において、同期装置１１８は表示装置１１５に警告を表示させる。一方、ステップＳ２０３において、同期装置１１８が照射可と判定した場合（Ｓ２０３－Ｙｅｓ）、処理をステップＳ２０５に進める。
　ステップＳ２０５において、事前情報取得部１１１は、事前情報として、放射線に関する確率分布の情報、管球情報、ＦＰＤ１０３に関する情報、撮影部位情報、被写体情報等を取得する。この時、識別物質数などの情報に応じて、画像取得枚数(例えば、エネルギーレベルの異なる２つのエネルギー（高エネルギー、低エネルギー）による物質分離を行うなら2枚)を決定する。そして、事前情報取得部１１１は、決定した画像取得枚数やＦＰＤ１０３における情報(蓄積時間など)等を同期装置１１８に通知する。
　ステップＳ２０６において、同期装置１１８は、事前情報取得部１１１から取得した情報に基づいて、放射線照射の指示信号を生成し、放射線管球１０１とＦＰＤ１０３に放射線照射の指示信号を送信する。これにより、放射線管球１０１による放射線の照射とＦＰＤ１０３による画像取得が制御される。ＦＰＤ１０３は、異なる放射線エネルギーレベルで撮影された複数の放射線画像（例えば、低エネルギー放射線画像、高エネルギー放射線画像）を撮影する。
　ＦＰＤ１０３により必要な画像取得枚数が取得されたら、ステップＳ２０７において、撮影を終了する。ステップＳ２０７において、撮影が終了すると、画像取得部１１０は、ＦＰＤ１０３から画像（放射線画像）を取得する。すなわち、画像取得部１１０は、異なる放射線エネルギーレベルで撮影された複数の放射線画像（低エネルギー放射線画像、高エネルギー放射線画像）を取得する。
　ステップＳ２０８において、平面分布取得部１１２は、ステップＳ２０５で求められた放射線に関する確率分布の情報等の事前情報と、ステップＳ２０６にて取得された画像とに基づいて、物質にかかわる平面分布を算出する。平面分布の詳細な算出処理については後述する。
　ステップＳ２０９において、コンピュータ１０９のＣＰＵは平面分布の算出結果をコンピュータ１０９内の記憶装置し、また、ＣＰＵは、平面分布の算出結果を表示装置１１５に表示させる。ＣＰＵは、例えば、分離した物質の厚さ、実効原子番号、面密度などの平面分布の算出結果と可視光の色とを対応させた画像（物質識別画像）を表示装置１１５にカラーで表示するような表示制御を行うことが可能である。
　（平面分布算出処理）
　次に、ステップＳ２０８の物質にかかわる平面分布の算出処理方法を詳細に説明する。ここでは、物質にかかわる平面分布算出の一例として、２色(高エネルギー放射線画像と低エネルギー放射線画像の２種類)の画像から、物質ごとの厚さの分布を示す画像（物質分離画像）を出力するプロセスについて説明するが、画像が３色以上の場合においても同様に平面分布の算出処理を適用することができる。
　ＦＰＤ１０３内に一時保存された放射線分布情報は、サンプリングホールド実施後、読み出し可能となり、画像取得部１１０は、ＦＰＤ１０３から異なるタイミングで、放射線分布情報（ＸＬ）と放射線分布情報（ＸＬ＋ＸＨ）の読み出しを実施する。画像取得部１１０は、放射線分布情報（ＸＬ＋ＸＨ）から放射線分布情報（ＸＬ）を差し引くことで、放射線分布情報（ＸＨ）を得ることができる。ここで、低エネルギーの放射線分布情報（ＸＬ）が低エネルギーの放射線による画像（低エネルギー放射線画像）になり、高エネルギーの放射線分布情報（ＸＨ）が高エネルギーの放射線による画像（高エネルギー放射線画像）になる。
　平面分布取得部１１２は、異なる放射線エネルギーで撮影した複数の放射線画像に含まれる物質についての物質特性画像を取得する。平面分布取得部１１２は、１回の放射線照射で取得された複数の放射線画像を用いて、物質分離画像、物質識別画像といった物質特性画像を取得することが可能である。
　ここで物質分離画像とは、複数の放射線画像（低エネルギー放射線画像、高エネルギー放射線画像）に含まれる複数の物質を分離した画像であり、被写体を特定の２以上の物質で表した場合に、その物質の厚さで構成された２以上の画像に分離した画像をいう。
　また、物質識別画像とは、複数の放射線画像（低エネルギー放射線画像、高エネルギー放射線画像）に含まれる物質の実効原子番号と面密度の分布を示す画像であり、被写体を特定の１物質で表した場合に、その物質の実効原子番号と面密度の分布に分解した画像をいう。
　また、物質にかかわる平面分布は、分離した物質の厚さの分布を示す画像（厚さ画像）、特定の物質に関する実効原子番号の分布を示す画像（実効原子番号画像）と面密度の分布を示す画像（面密度）の２種類に分ける場合や、実効原子番号画像と面密度画像と厚さ画像の３種類以上の画像に分ける場合においても同様に平面分布の算出処理を適用することができる。
　入力画像として、異なる放射線エネルギーの画像として、２つの放射線画像（低エネルギー放射線画像、高エネルギー放射線画像）をＩ=(Ｉ(Ｅ1、Ｉ(Ｅ2））とする。ここで、Ｅ1、Ｅ2は放射線のエネルギーを表す。また、出力画像として、物質の厚さ画像をｆ（ｄ）=(ｄ1、ｄ2）とする。ノイズのない場合の入力画像（Ｉ）と出力画像（ｆ（ｄ））との関係は、以下の数１式により表すことができる。
Figure JPOXMLDOC01-appb-M000001
　実際に得られる画像において、放射線の統計的性質に由来する量子ノイズや、検出器の電気的な雑音に由来するシステムノイズが存在する場合には、入力画像（Ｉ）と出力画像（ｆ（ｄ））との関係において、数１式の関係は必ずしも成立しない場合が生じ得る。一般にノイズをｎとすれば、入力画像（Ｉ）、出力画像（ｆ（ｄ））とノイズの関係は、以下の数２式のように、ノイズｎは、入力画像（Ｉ）と出力画像（ｆ（ｄ））との差分として表される。
Figure JPOXMLDOC01-appb-M000002
　ノイズｎの従う確率分布は、物質の厚さがｄの時、画像がIとなる確率と考えることができる。よって、ノイズｎの従う確率分布はｐ（Ｉ｜ｄ）と表すことができる。
　求めようとしている物質にかかわる平面分布は、入力画像がI（=(Ｉ(Ｅ1、Ｉ(Ｅ2））)の時、入力画像に含まれ得る各物質の厚さｄの平面分布である。各物質の厚さｄについて、求める厚さは入力画像がIであるとき、確率が最大となる物質の厚さｄであると考えることができる。つまり以下の数３式を満たす物質の厚さｄとなる。平面分布取得部１１２は、入力が入力画像（Ｉ）であるときに、出力が平面分布の値（ｄ）となる確率（事後確率）が最大となるように出力画像を決定する。
Figure JPOXMLDOC01-appb-M000003
　また、求める物質の厚さは確率分布の対数が最大となる物質の厚さｄと考えてもよく、その場合、式の形は以下の数４式、あるいは数５式で表すことができる。
Figure JPOXMLDOC01-appb-M000004
Figure JPOXMLDOC01-appb-M000005
　数５式をベイズの定理を用いて変形すると、以下の数６式のように変形することができる。確率（事後確率）は、尤度と事前確率に基づいて取得することができる。すなわち、平面分布取得部１１２は、尤度と、注目画素の画素値または周囲の画素の画素値を入力とする事前確率とから、入力が入力画像（Ｉ）であるときに、出力が平面分布の値（ｄ）となる確率（事後確率）が最大となるように出力画像を決定する。ここで、尤度は、入力される値が平面分布の値（物質の厚さｄ）のときに、入力に対応する画像が入力画像となる確率を示すものである。
Figure JPOXMLDOC01-appb-M000006
　ここで、数６式において、ｌｏｇｐ（Ｉ｜ｄ）は、尤度と呼ばれる項であり、ｌｏｇｐprior（ｄ）は事前確率と呼ばれる項である。事前確率の項の形式を変更することで、多様な事前情報を反映した確率を計算することが可能になる。この事前確率は、数６式の解をより頑健にする効果を持ち、正則化項と呼ばれることもある。
　数６式を求める方法はいくつか存在するが、極値が最小値の場合、例えば、最急降下法で求めることが可能である。この他、ニュートン法やその他一般的に知られている方法を用いることが可能である。最急降下法を適用する場合、以下の数７式の形でｄｎ＋１の逐次近似解を求め、求めた逐次近似解から数６式の解を求めることができる。平面分布取得部１１２は、出力画像を反復演算による逐次近似解により取得することができる。
Figure JPOXMLDOC01-appb-M000007
　数７式において、ｄｎはn回目の反復で求められた物質の厚さである。αは反復の速度を表す係数である。ｇｒａｄはベクトルの微分（勾配）を表す記号である。つまり、入力画像がIであるとき、確率が最大となる物質の厚さｄｎ＋１は、数７式のように反復を繰り返すことで逐次的に求めることができる。
　図３は、図２のステップＳ２０８の詳細な処理（物質の平面分布を取得する処理）の流れを説明するフローチャートである。ステップＳ３０１において、画像取得部１１０は、ＦＰＤ１０３から画像（放射線画像）を取得する。図２のステップＳ２０７において撮影が終了すると、画像取得部１１０は、異なる放射線エネルギーレベルで撮影された複数の放射線画像（低エネルギー放射線画像、高エネルギー放射線画像）を取得する。低エネルギー放射線画像、高エネルギー放射線画像の２つの放射線画像が、ステップＳ３０１における入力画像Ｉとなる。
　ステップＳ３０２において、事前情報取得部１１１は、ステップＳ２０５で求められた放射線に関する確率分布の情報等の事前情報を取得する。
　次にステップＳ３０３において、平面分布取得部１１２は、逐次的演算を行う前の初期値（ｄ１）の設定を行う。設定する初期値は、例えば、系にノイズがないと仮定したときの平面分布（物質の厚さ）とすることができる。ノイズのない場合の入力画像（Ｉ）と物質の厚さｄに関する出力画像（ｆ（ｄ））との関係は、数１式のように表すことができるので、平面分布（物質の厚さ）の初期値ｄ１は、以下の数８式から求めることができる。数８式に基づいて、平面分布取得部１１２は、ノイズがないと仮定したときの平面分布の値を初期値として反復演算を行う。
Figure JPOXMLDOC01-appb-M000008
　別の初期値の設定方法として、数８式で得られた画像に更に画像処理を加えた画像を物質の厚さの初期値として入力することも可能である。あるいは、エネルギーサブトラクション法により、異なる放射線エネルギーにより撮影された複数の放射線画像に含まれる複数の物質に分離した出力画像（例えば、軟部組織の厚さ画像や骨部の厚さ画像）について平滑化する手法で物質にかかわる平面分布を計算し、その出力を初期値として入力してもよい。すなわち、平面分布取得部１１２は、出力画像に対する平滑化処理により取得した平面分布の値を初期値として反復演算を行うことも可能である。
　ステップＳ３０４以降の処理が反復演算（逐次近似）の処理となり、平面分布取得部１１２は、予め設定された繰り返し回数ｎに到達するまで（Ｓ３０４－Ｙｅｓ）、ステップＳ３０５において、数７式における逐次近似解の更新を繰り返す。逐次近似解の更新処理において、数７式の微分項(ｇｒａｄ)の処理は数値的に求めてもよいし、解析的に求めてもよい。微分項(ｇｒａｄ)の処理方法の詳細は後述する。
　そして、ステップＳ３０６において、平面分布取得部１１２は、反復の速度を表す係数αｎを更新する。一般に繰り返しの最初の数回における解は最終的な解と比較して差分が大きいはずであるから、平面分布取得部１１２は、反復の速度を表す係数αｎを大きい値（第１係数値）に設定する。反復演算を行うにつれて誤差が小さくなるはずだから、平面分布取得部１１２は、反復の速度を表す係数αｎを小さい値、すなわち、第１係数値に比べて小さい第２係数値に設定する。反復の速度を表す係数αｎは、例えば、ｋ×αｎ→αｎ＋１という形式となる。ここで、「→」の表記は、ｋ×αｎをαｎ＋１に代入するという意味であり、係数ｋは0より大きく1より小さい実数の定数である（０＜ｋ（実数の定数）＜１）。第１係数値をαｎとすると、第２係数値（αｎ＋１）は、ｋ×第１係数値で求めることができる。平面分布取得部１１２は、反復演算（逐次近似）の繰り返し回数に応じて、反復の速度を表す係数αｎが小さくなるように係数αｎの設定を更新する。
　ステップＳ３０４の判定処理で、設定された繰り返し回数ｎに到達した場合（Ｓ３０４－Ｎｏ）、平面分布取得部１１２は、本処理を終了する。
　（具体的な数式を用いた処理の例）
　（物質の厚さと画像（出力画像）の関係ｆ（ｄ）について）
　次に、具体的な数式を用いた処理の例について説明する。まず、数１式、数２式などに現れる物質の厚さと画像（出力画像）の関係ｆ（ｄ）については、以下の数９式、あるいは数１０式に示した透過に関するランベルト・ベールの法則の式を用いることができる。数９式、数１０式において、Ｉ'／Ｉ０は減弱率、μは線減弱係数、ｄは物質の厚さであり、Ｅは放射線のエネルギー、Φ（Ｅ）は放射線のスペクトルを示す。Ｉは入力画像であり、異なる２つの放射線エネルギーで撮影した画像（低エネルギー放射線画像、高エネルギー放射線画像）である。独立した２つの連立方程式を解くことで、未知の変数である物質の厚さｄ１、ｄ２の値を求めることができる。
Figure JPOXMLDOC01-appb-M000009
Figure JPOXMLDOC01-appb-M000010
　（数６式の尤度の項について）
　入力される値が平面分布の値（物質の厚さｄ）のときに、入力に対応する画像が入力画像となる確率を示す尤度の項について説明する。平面分布取得部１１２は、尤度の確率分布として、正規分布またはポアソン分布のいずれかを用いる。
　物質の厚さｄの時、画像がIとなる確率を示す尤度（ｌｏｇｐ（Ｉ｜ｄ）：数６式）は、確率分布を正規分布で近似した以下の数１１式や、確率分布をポアソン分布で近似した以下の数１２式を用いることができる。ここで、数１１式において、σは分散を示している。また、数１２式において、所定の時間中に平均でλ１回発生する事象がｋ１回（ｋは正の整数）発生する確率、λ2回発生する事象がｋ2回（k は正の整数）発生する確率を示し、Ｉ１＝ｋ１となる確率と、Ｉ２＝ｋ２となる確率は独立である。
Figure JPOXMLDOC01-appb-M000011
Figure JPOXMLDOC01-appb-M000012
　（数６式の事前確率の項について）
　平面分布取得部１１２は、異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、入力画像に含まれる物質について平面分布を示す出力画像を取得し、出力画像における注目画素の画素値と周囲の画素の画素値との差分により平面分布の連続性を評価する。
　出力画像における連続性に関する事前確率は、平面分布の連続性の評価に基づいて取得される。数６式で示した事前確率の項（ｌｏｇｐprior（ｄ））については、例えば、注目画素と周囲の画素との比較（差分）を行う項を導入することができる。一例として、例えば、以下の数１３式を事前確率の項として用いることが可能である。ここで、平面分布取得部１１２は、事前確率を、出力画像における注目画素の画素値と注目画素に対する周囲の画素の画素値との差分に基づいて取得する。
Figure JPOXMLDOC01-appb-M000013
　数１３式において、「｜　｜」の表記は絶対値を示す。数１３式における、ｋ、ｌは、注目画素に対して上下、左右、斜めの隣接画素の項を表現するための添え字であり、一例として、ｋ＝±１、ｌ＝±１の組み合わせが挙げられる。比較（差分）の範囲は、隣接画素に限られず、ｋ＝±１、ｌ＝±１の組み合わせの隣接画素を含む、より広い範囲にわたった周囲の画素との比較（差分）の項を、数１３式の事前確率の項に追加することも可能である。
　数１３式で示したような事前確率の項(正則化項)を導入することで、平面分布取得部１１２は、出力画像（物質分離画像）における注目画素の画素値と周囲の画素の画素値との比較（差分）により、出力画像における平面分布の空間的な連続性を評価することができる。そして、平面分布取得部１１２は、平面分布の空間的な連続性の評価結果に基づいて、出力画像に対するノイズ低減処理を行う。これにより、確率的に分布するノイズが低減され、より優れた画質の画像を取得することが可能になる。
　数１３式で示したような事前確率の項(正則化項)における正則係数μi、ηiは事前確率の項(正則化項)の影響の度合い（強さ）を表す係数である。正則係数μi、ηiが大きい値をとるほど事前確率で導入した仮定が強く影響することになる。図４は、画素値の差分と事前確率の分布を模式的に示す図である。注目画素４０１及び周囲の画素４０２のように画素値の差分が小さい場合、事前確率は大きくなる傾向を示す。一方、注目画素４０１及び周囲の画素４０３のように画素値の差分が大きい場合、事前確率は小さくなる傾向を示す。
　すなわち、数１３式の事前確率の項(正則化項)は、連続性の評価値（注目画素の画素値と周囲の画素の画素値との差分）が小さければ小さいほど事前確率が大きくなることを意味する。
　また、正則係数μi、ηiに大きい値を設定しても事前確率の項(正則化項)の出力を大きくすることができる。これにより、事前確率で導入した仮定（画素値の連続性により出力画像は滑らかであるという事前の仮定）が強く影響する。正則係数μi、ηiが大きい値であればあるほど、出力画像の取得演算処理（逐次近似演算）により注目画素と周囲の画素との画素値の差分が小さくなり、より滑らかな出力画像（物質分離画像）が得られるようになる。滑らかな出力画像が得られる傾向にあるから、事前確率の項(正則化項)を用いることで、物質にかかわる平面分布（物質の厚さの分布を示す物質分離画像）のノイズを低減することができる。
　平面分布取得部１１２は、事前情報取得部１１１が取得した事前情報に含まれる撮影条件に基づいて、事前確率の項(正則化項)の影響の度合い（強さ）を表す正則係数μi、ηiを変更することができる。例えば、放射線の照射条件(線量)が所定の閾値に比べて小さい場合、平面分布取得部１１２は、正則係数μi、ηiを所定の基準値に比べて大きい値に変更して事前確率の項(正則化項)を大きくすることで、ノイズ除去（ノイズ低減）をより優先した設定にすることができる。
　平面分布取得部１１２は、放射線の照射条件(線量)に関する情報を、事前情報取得部１１１が取得した事前情報に基づいて取得することができる。例えば、事前情報に含まれる管球情報(管電圧、管電流、照射時間)は線量に関わる情報であり、平面分布取得部１１２は、放射線の照射条件(線量)に基づいて、事前確率の項(正則化項)の影響の度合いを表す正則係数μi、ηiを変更することができる。尚、放射線撮影装置１００の入力装置１１４に画質調整を行うための調整装置が付属している場合は、平面分布取得部１１２は、調整装置が示す値に基づいて事前確率の項(正則化項)の正則係数μi、ηiの設定を調整することができる。
　また、事前情報に含まれる被写体情報(身長、体重など)は被写体の体形に関わる情報であり、線量に関わる情報である。平面分布取得部１１２は、放射線の照射条件(線量)に関する情報として被写体情報（被写体の体形に関する情報）に基づいて、事前確率の項(正則化項)の影響の度合いを表す正則係数μi、ηiを変更することができる。例えば、被写体の体重が大きいほど透過する線量は低くなるので、平面分布取得部１１２は、正則係数μi、ηiを所定の基準値に比べて大きい値に変更して事前確率の項(正則化項)を大きくすることで、ノイズ除去（ノイズ低減）をより優先した設定にすることができる。
　また、Ｃアーム１１６や、架台１１７の位置設定の情報も線量に関わる情報である。例えば、Ｃアーム１１６の位置や角度、すなわち、放射線発生装置（放射線管球１０１）と放射線検知部（ＦＰＤ１０３）との相対的な位置または角度によって、放射線が被写体１０２を通過する距離が長くなる場合には、距離が短い場合に比べて透過する線量は相対的に低くなることが予想される。また、放射線が被写体を通過する距離が短くなる場合は、距離が長い場合に比べて透過する線量は相対的に高くなることが予想される。
　平面分布取得部１１２は、放射線の照射条件(線量)に関する情報として放射線発生装置（放射線管球１０１）と放射線検知部（ＦＰＤ１０３）との相対的な位置または角度に基づいて、事前確率の項(正則化項)の影響の度合いを表す正則係数μi、ηiを変更することができる。
　また、架台１１７は放射線を遮る要素の一つであるため、架台１１７の有無も線量に関わる情報であり、平面分布取得部１１２は、放射線の照射条件(線量)に関する情報として架台１１７の情報に基づいて、事前確率の項(正則化項)の影響の度合いを表す正則係数μi、ηiを変更することができる。
　発明者は、物質分離画像のうち、特に脂肪などの軟部組織に関する物質の厚さの分布を示す画像（厚さ画像）に対して、事前確率の項(正則化項)を導入することが、ノイズ低減等の画質向上に効果的であることを見出した。そのため、着目する部位が脂肪などの軟部組織である場合には、軟部組織に係る厚さ画像に対する事前確率の項(正則化項)の正則係数を、骨などの硬部組織に係る厚さ画像に対する事前確率の項(正則化項)の正則係数に比べて大きく設定して事前確率の項(正則化項)を大きくすることで、ノイズが除去（ノイズ低減）され、より優れた画質の画像（物質分離画像）を提供することが可能になる。
　この場合、平面分布取得部１１２は、軟部組織の厚さ画像における、注目画素の画素値と周囲の画素の画素値との差分により軟部組織の厚さ画像における連続性の評価を行う。そして、平面分布取得部１１２は、平面分布の空間的な連続性の評価結果に基づいて、出力画像に対するノイズ低減処理を行う。これにより、確率的に分布するノイズが低減され、より優れた画質の画像を取得することが可能になる。
　尚、着目する部位が、骨などの硬部組織である場合には、骨などの硬部組織に係る厚さ画像に対する事前確率の項(正則化項)の正則係数を、脂肪などの軟部組織に係る厚さ画像に対する事前確率の項(正則化項)の正則係数に比べて大きく設定してもよい。
　数１３式を解く場合に最急降下法においては、事前確率の項(正則化項)の厚さに関する微分を計算することになるが、数１３式を微分した場合、数１３式のうち、以下の数１４式、数１５式に示す２つの項は不連続な解となる。
Figure JPOXMLDOC01-appb-M000014
Figure JPOXMLDOC01-appb-M000015
　数１４式及び数１５式に示す２つの項は、具体的には以下の数１６式でそれぞれ表すことができる。
Figure JPOXMLDOC01-appb-M000016
　数１４式及び数１５式における不連続性を回避するためには、絶対値を含む項の代わりに微分が連続になる関数、例えば、シグモイド関数を含む以下の数１７式の関数の積分を、数１４式で示す項の代わりに用いることで不連続性を回避することができる。
Figure JPOXMLDOC01-appb-M000017
　以上、事前確率の項(正則化項)は注目画素と周囲の画素との差分の項（例えば、数１３式）を持つ例を説明したが、事前確率の項(正則化項)は必ずしも注目画素と周囲の画素の差分であるとは限らない。平面分布取得部１１２は、事前確率を、出力画像（物質の厚さ画像）における注目画素の画素値に基づいて取得することも可能である。例えば、物質の厚さ画像(ｄ1やｄ2)はゼロ以上の有限の値となるから、物質の厚さを入力として、以下の数１８式に示すような正則化項を導入することができる。数１８式において、κは正の定数で、正則係数の一例である。
Figure JPOXMLDOC01-appb-M000018
　物質の厚さがゼロより小さくなると、数１８式に示す事前確率の項(正則化項)であるｌｏｇprior(ｄ1i,j、ｄ2i,j)もゼロより小さくなり、ゆえに尤度と事前確率の項(正則化項)とを合わせた全体の事後確率は小さくなる。物質の厚さ画像(ｄ1，ｄ2)はゼロより大きい画素値を持つ画像となる可能性が向上する。この事前確率も、正負の値にばらつくはずだった厚さ画像を正の値だけにする効果を持つから、周囲の画素との差分の項を持つ事前確率の項(正則化項)と同様に平面分布の連続性を評価する項の一例であり、またノイズを低減する結果をもたらす。すなわち、平面分布取得部１１２は、出力画像（物質の厚さ画像）における注目画素の画素値により平面分布の連続性を評価し、評価結果に基づいて出力画像におけるノイズを低減する。これにより、確率的に分布するノイズが低減され、より優れた画質の画像を取得することが可能になる。
　（処理の高速化）
　数１式、数８式から数１３式を、数７式に代入することで、数７式の逐次近似解（更新解）を求めることができるが、式の変形により解析解をもとめることが可能となり、結果演算をより高速化することができる。
　数１１式は、物質の厚さｄの時、画像がIとなる確率を示す尤度の確率分布を正規分布で近似した式であるが、数１１式の対数をとると、数１１式は以下の数１９式のように表すことができる。
Figure JPOXMLDOC01-appb-M000019
　尤度に関しては数１９式のように展開することができ、事前確率の項(正則化項)に関しては数１３式のように展開できるので、数７式中の－logp(I｜ｄ)－logｐprior（ｄ）という数式部分は、以下の数２０式で表すことができる。
Figure JPOXMLDOC01-appb-M000020
　数２０式を数７式に代入すると、ベクトルの微分（勾配）を表すｇｒａｄの項は、以下の数２１式のように表すことができる。
Figure JPOXMLDOC01-appb-M000021
　数２１式の解を解析的に求めることにより、演算精度を維持しつつ演算処理をより高速化することができる。
　数１２式は尤度の確率分布をポアソン分布で近似した式であるが、数１２式の対数をとると、数１２式は数２２式のように表すことができる。
Figure JPOXMLDOC01-appb-M000022
　尤度に関しては数２２式のように展開することができ、事前確率の項(正則化項)に関しては数１３式のように展開できるので、数７式中の－logp(I｜ｄ)－logｐprior（ｄ）という数式部分は、以下の数２３式で表すことができる。
Figure JPOXMLDOC01-appb-M000023
　数２３式を数７式に代入すると、ベクトルの微分（勾配）を表すｇｒａｄの項は、以下の数２４式のように表すことができる。
Figure JPOXMLDOC01-appb-M000024
　数２４式の解を解析的に求めることにより、演算精度を維持しつつ演算処理をより高速化することができる。
　（単色近似）
　尤度の確率分布を正規分布で近似した数１１式に現れるｆi（ｄ）、あるいは尤度の確率分布をポアソン分布で近似した数１２式に現れるλi（ｄ）（i＝１、２）は平均画素値となる。これはノイズがないときの平均画素値であり、数９式に基づいて決定することができる。数９式あるいは数１０式を、数２１式、数２４式に代入することで最急降下法による逐次近似解を数７式から求めることができる。数９式及び数１０式における放射線のスペクトル（Φ（Ｅ））の幅を考慮することで、放射線が透過したときのエネルギー変化の影響(ビームハードニング)を抑制することが可能になる。
　一方、以下の数２５式に示すように放射線スペクトルを、デルタ関数（ｅ＝１、２）を用いて単色近似することで、より簡潔に解を表現することができる。計算精度が必要な場合は数９式や数１０式を用いて計算を行い、計算速度が必要である場合、数９式及び数１０式における放射線のスペクトル（Φ（Ｅ））に単色近似を導入にして計算を行えばよい。平面分布取得部１１２は、異なる放射線エネルギーに対応した放射線スペクトルを単色の放射線スペクトルとして近似して出力画像を取得する。
Figure JPOXMLDOC01-appb-M000025
　単色近似を導入したとき、数９式、あるいは、数１０式は、以下の数２６式となり、物質の厚さｄの座標も考慮すると、数２６式は数２７式のように表すことができる。
Figure JPOXMLDOC01-appb-M000026
Figure JPOXMLDOC01-appb-M000027
　数２７式を微分すると数２８式、数２９式のように表すことができる。
Figure JPOXMLDOC01-appb-M000028
Figure JPOXMLDOC01-appb-M000029
　数２７式、数２８式、数２９式を正規分布に関する数２１式に代入すると、数２１式の各成分は、以下の数３０式及び数３１式のように表すことができる。
Figure JPOXMLDOC01-appb-M000030
Figure JPOXMLDOC01-appb-M000031
　同様に数２７式、数２８式、数２９式をポアソン分布に関する数２４式に代入することで、単色近似を導入した場合におけるポアソン分布に関する数２４式の各成分を求めることができる。尚、ここでは、放射線スペクトルの単色近似について説明したが、正規分布やポアソン分布の再生性、つまり独立な確率変数Ｘ１とＸ２がそれぞれ正規分布やポアソン分布(パラメタλ１とλ２)に従うとき、確率変数の和Ｘ１＋Ｘ２も正規分布やポアソン分布に従うことを考慮すると、放射線スペクトルの幅を考慮したときの画素値Ｉ＝Σｉi（Ｅi）についても同様の議論を展開することができる。この場合、平面分布取得部１１２は、異なる放射線エネルギーに対応した放射線スペクトルに基づいて取得された各画像の画素値を加算して出力画像を取得し、取得した出力画像に対して本実施形態の処理を適用することができる。
　（その他の実施形態）
　本発明は、上述の実施形態の１以上の機能を実現するプログラムを、ネットワーク又は記憶媒体を介してシステム又は装置に供給し、そのシステム又は装置のコンピュータにおける１つ以上のプロセッサーがプログラムを読出し実行する処理でも実現可能である。また、１以上の機能を実現する回路（例えば、ＡＳＩＣ）によっても実現可能である。
　発明は上記実施形態に制限されるものではなく、発明の精神及び範囲から離脱することなく、様々な変更及び変形が可能である。従って、発明の範囲を公にするために請求項を添付する。
　本願は、２０１９年５月３０日提出の日本国特許出願、特願２０１９－１０１４８４を基礎として優先権を主張するものであり、その記載内容の全てを、ここに援用する。
　１００：放射線撮影装置、１０１：放射線管球（放射線発生装置）、１０３：ＦＰＤ（放射線検知部）：１０４：操作盤、１０８：照射指示部、１０９：コンピュータ（画像処理装置）、１１０：画像取得部、１１１：事前情報取得部、１１２：平面分布取得部、１１３：事前情報保存部、１１４：入力装置、１１５：表示装置、１１６：Ｃアーム、１１７：架台、１１８：同期装置

Claims (19)

  1. 　異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、前記入力画像に含まれる物質について平面分布を示す出力画像を取得する平面分布取得手段を備え、
    　前記平面分布取得手段は、尤度と、注目画素の画素値または周囲の画素の画素値を入力とする事前確率とから、入力が前記入力画像であるときに、出力が前記平面分布の値となる確率が最大となるように前記出力画像を決定する
    　ことを特徴とする画像処理装置。
  2. 　前記平面分布取得手段は、前記事前確率を、前記出力画像における前記注目画素の画素値と前記周囲の画素の画素値との差分に基づいて取得することを特徴とする請求項１に記載の画像処理装置。
  3. 　前記平面分布取得手段は、前記事前確率を、前記出力画像における前記注目画素の画素値に基づいて取得することを特徴とする請求項１または２に記載の画像処理装置。
  4. 　前記平面分布取得手段は、前記差分により前記平面分布の連続性を評価し、前記評価に基づいて前記出力画像におけるノイズを低減することを特徴とする請求項２に記載の画像処理装置。
  5. 　前記平面分布取得手段は、前記出力画像を反復演算による逐次近似解により取得することを特徴とする請求項１乃至４のいずれか１項に記載の画像処理装置。
  6. 　前記平面分布取得手段は、前記反復演算の繰り返し回数に応じて、反復の速度を表す係数が小さくなるように前記係数の設定を更新することを特徴とする請求項５に記載の画像処理装置。
  7. 　前記平面分布取得手段は、ノイズがないと仮定したときの平面分布の値を初期値として前記反復演算を行うことを特徴とする請求項５または６に記載の画像処理装置。
  8. 　前記平面分布取得手段は、前記出力画像に対する平滑化処理により取得した平面分布の値を初期値として前記反復演算を行うことを特徴とする請求項５または６に記載の画像処理装置。
  9. 　前記平面分布取得手段は、前記異なる放射線エネルギーに対応した放射線スペクトルに基づいて取得された各画像の画素値を加算して前記出力画像を取得することを特徴とする請求項１乃至８のいずれか１項に記載の画像処理装置。
  10. 　前記平面分布取得手段は、前記異なる放射線エネルギーに対応した放射線スペクトルを単色の放射線スペクトルとして近似して前記出力画像を取得することを特徴とする請求項１乃至８のいずれか１項に記載の画像処理装置。
  11. 　前記平面分布取得手段は、前記尤度の確率分布として、正規分布またはポアソン分布のいずれかを用いることを特徴とする請求項１乃至１０のいずれか１項に記載の画像処理装置。
  12. 　前記平面分布取得手段は、前記出力画像として、前記複数の放射線画像に含まれる複数の物質に分離した画像を取得することを特徴とする請求項１乃至１１のいずれか１項に記載の画像処理装置。
  13. 　前記平面分布取得手段は、前記出力画像として、被写体を構成する軟部組織の厚さ画像と、前記被写体を構成する硬部組織の厚さ画像を取得することを特徴とする請求項１２に記載の画像処理装置。
  14. 　前記平面分布取得手段は、前記軟部組織の厚さ画像における、注目画素の画素値と周囲の画素の画素値との差分により前記軟部組織の厚さ画像における連続性の評価を行うことを特徴とする請求項１３に記載の画像処理装置。
  15. 　前記平面分布取得手段は、前記事前確率の影響の度合いを表す正則係数を撮影条件に基づいて変更することを特徴とする請求項１乃至１４のいずれか１項に記載の画像処理装置。
  16. 　前記平面分布取得手段は、前記事前確率の影響の度合いを表す正則係数を放射線の照射条件に基づいて変更することを特徴とする請求項１乃至１５のいずれか１項に記載の画像処理装置。
  17. 　前記平面分布取得手段は、前記事前確率の影響の度合いを表す正則係数を被写体の体形に関する情報に基づいて変更することを特徴とする請求項１乃至１６のいずれか１項に記載の画像処理装置。
  18. 　画像処理装置の画像処理方法であって、
    　異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、前記入力画像に含まれる物質について平面分布を示す出力画像を取得する工程を有し、
    　前記工程では、尤度と、注目画素の画素値または周囲の画素の画素値を入力とする事前確率とから、入力が前記入力画像であるときに、出力が前記平面分布の値となる確率が最大となるように前記出力画像を決定する
    　ことを特徴とする画像処理方法。
  19. 　コンピュータを、請求項１乃至１７のいずれか１項に記載の画像処理装置の手段として機能させるプログラム。
PCT/JP2020/014095 2019-05-30 2020-03-27 画像処理装置、画像処理方法及びプログラム Ceased WO2020241030A1 (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
US17/526,697 US12112473B2 (en) 2019-05-30 2021-11-15 Image processing apparatus, image processing method, and non-transitory computer-readable storage medium that use a plurality of captured radiation images to obtain, by sequential approximation solution by iterative calculation, an output image

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2019101484A JP7352382B2 (ja) 2019-05-30 2019-05-30 画像処理装置、画像処理方法及びプログラム
JP2019-101484 2019-05-30

Related Child Applications (1)

Application Number Title Priority Date Filing Date
US17/526,697 Continuation US12112473B2 (en) 2019-05-30 2021-11-15 Image processing apparatus, image processing method, and non-transitory computer-readable storage medium that use a plurality of captured radiation images to obtain, by sequential approximation solution by iterative calculation, an output image

Publications (1)

Publication Number Publication Date
WO2020241030A1 true WO2020241030A1 (ja) 2020-12-03

Family

ID=73545442

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2020/014095 Ceased WO2020241030A1 (ja) 2019-05-30 2020-03-27 画像処理装置、画像処理方法及びプログラム

Country Status (3)

Country Link
US (1) US12112473B2 (ja)
JP (1) JP7352382B2 (ja)
WO (1) WO2020241030A1 (ja)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2022118536A (ja) * 2021-02-02 2022-08-15 富士フイルムヘルスケア株式会社 医用画像処理装置および医用画像処理方法

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2021162026A1 (ja) 2020-02-10 2021-08-19 キヤノン株式会社 画像処理装置及び画像処理方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2011156302A (ja) * 2010-02-03 2011-08-18 Kyoto Univ Ｘ線ｃｔ画像処理方法，ｘ線ｃｔプログラムおよび該プログラムが搭載されたｘ線ｃｔ装置
WO2014185078A1 (ja) * 2013-05-15 2014-11-20 国立大学法人京都大学 Ｘ線ｃｔ画像処理方法，ｘ線ｃｔ画像処理プログラム及びｘ線ｃｔ画像装置
WO2018105493A1 (ja) * 2016-12-05 2018-06-14 キヤノン株式会社 放射線撮影装置、放射線撮影システム、放射線撮影方法、及びプログラム
JP2019068953A (ja) * 2017-10-06 2019-05-09 キヤノン株式会社 画像処理装置、画像処理方法及びプログラム

Family Cites Families (23)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2651625B2 (ja) 1990-02-14 1997-09-10 富士写真フイルム株式会社 エネルギーサブトラクション画像生成方法
ATE209344T1 (de) * 1995-04-28 2001-12-15 Gyula Domjan Verfahren und vorrichtung zur analyse eines gegenstandes
US7028899B2 (en) * 1999-06-07 2006-04-18 Metrologic Instruments, Inc. Method of speckle-noise pattern reduction and apparatus therefore based on reducing the temporal-coherence of the planar laser illumination beam before it illuminates the target object by applying temporal phase modulation techniques during the transmission of the plib towards the target
US6988660B2 (en) * 1999-06-07 2006-01-24 Metrologic Instruments, Inc. Planar laser illumination and imaging (PLIIM) based camera system for producing high-resolution 3-D images of moving 3-D objects
US6959870B2 (en) * 1999-06-07 2005-11-01 Metrologic Instruments, Inc. Planar LED-based illumination array (PLIA) chips
US6754298B2 (en) * 2002-02-20 2004-06-22 The Regents Of The University Of Michigan Method for statistically reconstructing images from a plurality of transmission measurements having energy diversity and image reconstructor apparatus utilizing the method
US20100312072A1 (en) * 2008-01-02 2010-12-09 Yeda Research And Development Co. Ltd. Method and system for detecting and grading prostate cancer
US8416986B2 (en) * 2009-10-29 2013-04-09 Raytheon Company Methods and systems for processing data using non-linear slope compensation
JP5921180B2 (ja) * 2011-12-15 2016-05-24 キヤノン株式会社 画像処理装置、画像処理方法及びプログラム
CN103744124B (zh) * 2013-12-31 2015-07-22 华中科技大学 一种平面地形中地下管状设施红外成像探测定位方法
US9873669B2 (en) 2014-01-09 2018-01-23 Eisai R&D Management Co., Ltd. Cationic lipid
JP6284898B2 (ja) * 2015-03-31 2018-02-28 富士フイルム株式会社 ノイズ抑制処理装置および方法並びにプログラム
US10149958B1 (en) * 2015-07-17 2018-12-11 Bao Tran Systems and methods for computer assisted operation
US10176642B2 (en) * 2015-07-17 2019-01-08 Bao Tran Systems and methods for computer assisted operation
JP6676338B2 (ja) * 2015-10-30 2020-04-08 キヤノン株式会社 放射線撮像システム、放射線画像の情報処理装置、放射線画像の情報処理方法、及び、そのプログラム
US9911208B2 (en) * 2016-04-11 2018-03-06 Toshiba Medical Systems Corporation Apparatus and method of iterative image reconstruction using regularization-parameter control
US20200151880A1 (en) * 2017-06-02 2020-05-14 The Johns Hopkins University Reconstruction of difference images using prior structural information
WO2019221778A2 (en) * 2017-10-25 2019-11-21 Northwestern University High speed/low dose multi-objective autonomous scanning materials imaging
US10497612B2 (en) 2017-12-11 2019-12-03 Globalfoundries Inc. Methods of forming contact structures on integrated circuit products
JP7208723B2 (ja) * 2018-07-09 2023-01-19 キヤノン株式会社 画像処理装置およびその制御方法、放射線撮影システム、プログラム
JP2020061080A (ja) * 2018-10-12 2020-04-16 キヤノン株式会社 画像処理装置、撮像装置、画像処理方法
FI3825737T3 (fi) * 2019-11-22 2025-05-22 Hellma Mat Gmbh Menetelmä ja laite gammasäteilyn moniulotteiseen suuntamittaukseen kaukokentässä
JPWO2021171394A1 (ja) * 2020-02-26 2021-09-02

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2011156302A (ja) * 2010-02-03 2011-08-18 Kyoto Univ Ｘ線ｃｔ画像処理方法，ｘ線ｃｔプログラムおよび該プログラムが搭載されたｘ線ｃｔ装置
WO2014185078A1 (ja) * 2013-05-15 2014-11-20 国立大学法人京都大学 Ｘ線ｃｔ画像処理方法，ｘ線ｃｔ画像処理プログラム及びｘ線ｃｔ画像装置
WO2018105493A1 (ja) * 2016-12-05 2018-06-14 キヤノン株式会社 放射線撮影装置、放射線撮影システム、放射線撮影方法、及びプログラム
JP2019068953A (ja) * 2017-10-06 2019-05-09 キヤノン株式会社 画像処理装置、画像処理方法及びプログラム

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
LEMMENS, CATHERINE ET AL.: "Suppression of Metal Artifacts in CT Using a Reconstruction Procedure That Combines MAP and Projection Completion", IEEE TRANSACTIONS ON MEDICAL IMAGING, vol. 28, no. 2, 30 January 2009 (2009-01-30), pages 250 - 260, XP011232599, DOI: 10.1109/TMI.2008.929103 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2022118536A (ja) * 2021-02-02 2022-08-15 富士フイルムヘルスケア株式会社 医用画像処理装置および医用画像処理方法
JP7395523B2 (ja) 2021-02-02 2023-12-11 富士フイルムヘルスケア株式会社 医用画像処理装置および医用画像処理方法

Also Published As

Publication number Publication date
JP7352382B2 (ja) 2023-09-28
US12112473B2 (en) 2024-10-08
US20220076415A1 (en) 2022-03-10
JP2020192269A (ja) 2020-12-03

Similar Documents

Publication Publication Date Title
US10292672B2 (en) Radiographic image processing device, method, and recording medium
US10631808B2 (en) Generating a lung condition map
JP5124363B2 (ja) エネルギーサブトラクション処理装置、方法、およびプログラム
JP7093233B2 (ja) 放射線撮影装置、放射線撮影方法およびプログラム
US10258305B2 (en) Radiographic image processing device, method, and program
JP7054329B2 (ja) 画像処理装置、画像処理方法及びプログラム
JP6156847B2 (ja) 放射線画像処理装置および方法並びにプログラム
EP2567659A1 (en) Calibration free dual energy radiography method.
JP2015043959A5 (ja)
JP6746676B2 (ja) 画像処理装置、画像処理方法、及びプログラム
WO2020095578A1 (ja) 情報処理装置および方法、放射線撮影システム
WO2021024715A1 (ja) 画像処理装置およびその制御方法、放射線撮影装置、プログラム
EP3301650B1 (en) Image processing apparatus, image processing system, image processing method, and program
JP2017080343A (ja) 放射線撮像システム、放射線画像の情報処理装置、放射線画像の情報処理方法、及び、そのプログラム
WO2020241030A1 (ja) 画像処理装置、画像処理方法及びプログラム
CN110113996A (zh) 用于暗场成像的系统和方法
JP6379114B2 (ja) Ｘ線ｃｔ装置、および、ｘ線ｃｔ装置用画像演算装置
JPH0448453B2 (ja)
WO2020241110A1 (ja) 画像処理装置、画像処理方法及びプログラム
US20230102862A1 (en) Fat mass derivation device, fat mass derivation method, and fat mass derivation program
Kim et al. Model-based noise reduction in scatter correction using a deep convolutional neural network for radiography
CN108780573B (zh) 图像重建
JP2021122674A (ja) 画像処理装置および方法、プログラム
Ohashi et al. Segmentation ability of pulmonary nodules using deep learning in dual-energy subtraction images
WO2018105493A1 (ja) 放射線撮影装置、放射線撮影システム、放射線撮影方法、及びプログラム

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 20813595

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 20813595

Country of ref document: EP

Kind code of ref document: A1