JP6465763B2 - Image processing apparatus, image processing system, image processing method, and program - Google Patents

Description

  The disclosure of the present specification relates to an image processing apparatus that performs image processing on a radiation image, and an image processing program.

The radiation image obtained by irradiating the subject with radiation includes a component due to scattered radiation, which is radiation scattered in the subject, in addition to a component due to primary radiation traveling straight from the radiation source. Such scattered radiation components may reduce the contrast of the radiation image.
As a technique for reducing the scattered radiation component of an X-ray image, Patent Document 1 discloses that a scattered radiation component is modeled using an expression obtained by convolving two Gaussian functions with respect to one expression represented by primary X-rays. Obtaining a primary X-ray image is disclosed.

JP 2010-188113 A

  Incidentally, FIG. 9A is a diagram showing an example in which a scattered radiation component included in a radiographic image is experimentally obtained by measurement using a chest phantom. According to FIG. 9A, the scattered ray component includes a so-called skin line 903 at the boundary between the region 901 and the subject region 902 where the X-rays that have passed through the periphery of the subject have entered without passing through the subject. Exists. This is because the manner in which the radiation incident on the subject is scattered differs depending on the thickness of the subject, and the scattered X-ray intensity distribution is discontinuous because it does not scatter in the blank region.

  The inventors of the present invention have found that one expression expressed by the above-described primary X-rays disclosed in Patent Document 1 represents the intensity of the primary scattered X-rays once scattered in the subject. It has been found that only one equation represented by X-rays is insufficient to accurately represent the intensity distribution of scattered X-rays as illustrated in FIG. The inventors of the present invention have intensively studied in view of the above problems, and as a result, have found that the scattered radiation component can be estimated more accurately with the following configuration.

  An image processing apparatus according to the present invention includes an image acquisition unit that acquires a radiation image obtained by irradiating a subject with radiation, and a scattered ray including multi-order scattered rays generated by the radiation being scattered a plurality of times in the subject. The radiation acquired by the image acquisition means based on the first function corresponding to the intensity of the first scattered radiation and the second function corresponding to the intensity of the primary scattered radiation generated by once scattering the radiation in the subject An estimation unit that estimates a scattered radiation component included in an image, and an output unit that outputs a corrected image in which the scattered radiation component estimated by the estimation unit is reduced from the radiation image.

  According to the present invention, the first function corresponding to the intensity of the scattered radiation including the multi-order scattered radiation generated by scattering the radiation a plurality of times in the subject, and the primary scattered radiation generated by scattering the radiation once in the subject. Based on the second function corresponding to the intensity of, the intensity distribution of scattered radiation can be estimated more accurately.

It is a figure which shows the structure of the information system containing the image processing apparatus which concerns on embodiment of this invention. It is a figure which shows the structure of the image processing apparatus which concerns on embodiment of this invention. It is a flowchart which shows the process which concerns on embodiment of this invention. It is a figure which illustrates the scattered X-ray approximated with the 1st function which concerns on embodiment of this invention. It is a figure which illustrates the scattered X-ray approximated with the 2nd function which concerns on embodiment of this invention. It is a flowchart which shows the process which concerns on another embodiment of this invention. It is a figure which shows the frequency decomposition and image reconstruction by the Laplacian pyramid concerning another embodiment of this invention. It is a figure which shows the frequency characteristic of the frequency decomposition by the Laplacian pyramid which concerns on another embodiment of this invention. It is a figure which shows an example of the scattered-radiation component contained in the radiographic image measured using the phantom.

  An image processing apparatus according to an embodiment of the present invention and an information system 120 including a radiation imaging system 100 will be described with reference to FIG. An image processing apparatus according to an embodiment corresponds to the control unit 105 included in the radiation imaging system 100. The information system 120 referred to here includes, for example, HIS (Hospital Information System) 111, RIS (Radiography Information System) 112, WS (Work Station and 113), PACS (Picture Archiving and Community, 115). Printer 116. The HIS 111 is a system that comprehensively manages medical information including patient information and examination by radiography. The RIS 112 is a system for managing the order of radiography. WS 113 is an image processing terminal that performs image processing on a radiographic image captured by the radiation imaging system 100. WS may be replaced by one or a plurality of computers installed with software having similar functions. The PACS 114 is a database system that holds images obtained by radiation imaging or other medical imaging apparatuses in the information system 120. The PACS 114 includes a storage unit (not shown) that stores medical images and incidental information such as imaging conditions of such medical images and patient information, and a controller (not shown) that manages information stored in the storage unit. The Viewer 115 is a terminal for image diagnosis, reads out an image stored in the PACS 114 and displays it for diagnosis. The printer 116 is, for example, a film printer, and outputs an image stored in the PACS 114 or the like to a film.

  In the radiation imaging system 100 in the embodiment, X-rays are used as radiation. The radiation imaging system 100 includes an X-ray source 101 that is an example of a radiation generation apparatus, an FPD (Flat Panel Detector) 102, and a control unit 105. These are connected via a cable or a communication system. The control unit 105 attaches the imaging condition of the imaging, patient information, and the like to the captured radiographic image. For example, information is attached according to the DICOM (Digital Imaging and Communications in Medicine) standard, and a DICOM image file including information such as radiographic image data, patient information, and imaging conditions is generated. The control unit 105 transmits the image to the WS 113 or the PACS 114. The imaging order is transmitted from the RIS 112 to the control unit 105, for example. The control unit 105 reads out imaging conditions from a storage unit (not shown) according to input information from the RIS 112.

  The x-ray source 101 may be an x-ray tube or any other radiation source suitable for acquiring medical or other images. When the operator depresses the exposure switch, the high-voltage generation unit 104 gives a high-pressure pulse to the X-ray source 101 so that X-rays are exposed from the X-ray source 101 to the area where the subject 103 is disposed. X-rays that have passed through the subject 103 or have passed through the surroundings enter an FPD 102 that is an X-ray detector. The FPD 102 is controlled by the control unit 105, converts incident X-rays into electrical signals, and transmits them to the control unit 105 as digital images. For example, in the FPD 102, a phosphor (not shown) converts incident X-rays into visible light, a visible light is detected by a photodiode (not shown), and converted into an electrical signal by an A / D converter (not shown). . Alternatively, the FPD 102 converts X-rays into electrical signals using amorphous selenium (not shown). The pixel value of the radiation image is obtained by the output from the radiation detection element 102a constituting the FPD 102. The radiation detection element 102a includes, for example, a phosphor (not shown) and a photodiode (not shown). In another example, it is composed of amorphous selenium (not shown).

  The digital image is subjected to image processing by the control unit 105 or WS 113 and stored in the PACS 114 or the like. Each unit included in the information system 120 only needs to be connected to each other by a bus or other communication system, and can be remotely installed.

  Next, the configuration of the image processing apparatus according to the embodiment of the present invention will be described in detail based on FIG. The image processing apparatus according to the embodiment of the present invention is a control unit 105 connected to the information system 120 and the radiation imaging system 100, and includes one or a plurality of computers. A computer constituting the control unit 105 includes a CPU (Central Processing Unit) 201 as a main control unit, a RAM (Random Access Memory) 202 as a storage unit, a ROM (Read Only Memory) 205, an SSD (Solid State Drive) 206, Graphics processing unit GPU (Graphics Processing Unit) 208, communication unit NIC (Network Interface Card) 203 and 204, connection unit USB (Universal Serial Bus) 207, HDMI (registered trademark) High Definition Multimedia Intensity Multimedia 209, which are connected to each other via an internal bus Has been. The CPU 201 is a control circuit that integrally controls the control unit 105 and each unit connected thereto. A RAM 202 is a memory for storing a program for executing processing in the control unit 105 and each unit connected thereto, and various parameters used in image processing. The CPU 201 sequentially executes instructions included in the program expanded in the RAM 202, thereby realizing image processing to be described later. For example, a first NIC 203 as a communication unit is connected to an access point of a facility that performs radiation imaging, and a second NIC 204 is connected to an access point that relays communication in the information system 120. The SSD 206 stores a program as described above, a radiographic image obtained by imaging, incidental information, and other various parameters. The USB 207 is connected to the operation unit 107. The GPU 208 is an image processing unit, and executes image processing in accordance with control from the CPU 201. An image obtained as a result of the image processing is output to the monitor 106 via the HDMI (registered trademark) 209 and displayed. The monitor 106 and the operation unit 107 may be integrated with a touch panel monitor.

  The programs stored in the SSD 206 are, for example, an imaging control module 211, a communication control module 212, an image acquisition module 213, an output module 214, a display control module 215, and an estimation module 220.

  The imaging control module 211 is a program for causing the CPU 201 to integrally control the execution of radiation imaging. The imaging control module 211 performs, for example, specification of imaging conditions according to operation input, transmission of a signal requesting the state of the FPD 102, and the like.

  The communication control module 212 controls communication by the first NIC 203 and the second NIC 204. For example, in response to an input from the operation unit 107, the communication control module 212 transmits a signal for causing the FPD 102 to transition to a state in which photographing can be performed.

  The image acquisition module 213 is executed by the CPU 201 to control a process of acquiring an image used for image processing according to the embodiment of the present invention. For example, the NIC 203 is caused to receive a radiation image captured by the FPD 102. In reception of a radiographic image, after receiving a reduced image of a radiographic image with a small amount of data in advance, reception of the radiographic image may be completed by receiving data other than the reduced image in the radiographic image. . The reduced image is obtained by using only output signals selectively read from some elements, such as reading out even columns of a plurality of radiation detection elements that constitute the FPD 102 and give the pixel values of the radiation image. Or you may read several elements collectively. The read image may be divided into a plurality of small areas, and a reduced image may be obtained using representative knowledge of the small areas. Alternatively, the NIC 203 is caused to receive a radiation image stored in the PACS 114 or other storage unit on the network. Or the radiographic image memorize | stored in SSD206 of the image processing apparatus 105 or another memory | storage part is read. In addition, before the image processing according to the embodiment of the present invention is performed, it may be controlled to perform known image processing.

  The estimation module 220 is executed by the CPU 201, and a condition that the radiation image is represented by the sum of the primary radiation component and the scattered radiation component, and a condition that the scattered radiation component can be obtained from the primary radiation component by the scattered radiation model. , The primary radiation component or scattered radiation component of the radiation image is estimated. The estimation module 220 includes, for example, a scattered radiation intensity approximation module 221, a frequency processing module 222, and an optimization module 223. These modules are executed by the CPU 201. As a result, the CPU 201 controls the GPU 208, and the GPU 208 executes a scattered radiation image estimation process.

The scattered radiation intensity approximation module 221 is executed by the CPU 201 and obtains the intensity of the scattered radiation component by inputting the primary radiation component P. For example, for the first function, the n-th order approximate solution P _n of the primary radiation component is input, and the first output data is obtained. For the second function, the n-th order approximate solution P _n of the primary radiation component is input, and second output data is obtained.

The frequency processing module 222 is executed by the CPU 201 to obtain a specific frequency component for the output of the function. For example a Gaussian function G ₁ and convolution first output data to obtain a specific frequency component for the first output data. The Gaussian function G ₂ convolution to the second output data, to obtain a specific frequency component for the second output data. By synthesized weighted by these predetermined ratios to obtain n-order approximate solution S _n of the scattered ray components.

Optimization module 223 is executed by the CPU 201, the recurrence formula based on successive approximation, the _{P n + 1} is determined based on _{P n} and _{S n.} The obtained value is input to the scattered radiation intensity approximation module 221 again. The optimization module 223 determines convergence, and the process is repeated in the same order until it is determined that convergence has occurred. Approximate solutions at the time when it is determined to have converged, for example, P _n and S _n are used as estimation results of the primary radiation component and scattered radiation component of the radiation image. Details of this processing will be described later with reference to FIG.

  The output module 214 is executed by the CPU 201 to control output of a corrected image in which scattered radiation components are reduced by performing image processing according to the embodiment of the present invention. For example, the output module 214 outputs the corrected image to the monitor 106 to display the corrected image on the monitor 106. Further, for example, the output module 214 outputs the corrected image to the PACS 114 and the printer 116 via the NIC 204. As a result, the corrected image is stored in the PACS 114, and the corrected image is output to a film or the like by the printer 116. The output module 214 may output the corrected image to other storage units (not shown) inside and outside the control unit 105 and store the corrected image. Furthermore, it is preferable to output a corrected image with various information attached according to the DICOM standard. Modality is an image generation device that images a patient and generates a medical image. In the information system 120 according to the embodiment of the present invention, for example, the radiation imaging system 100 including the X-ray source 101 and the FPD 102 corresponds to a modality. At this time, DX indicating Digital Radiography is attached as a Modality tag (0008, 0060). In the case of moving image shooting, an RF indicating Radio Fluoroscopy is attached. Further, in the case of saving in the PACS 114, an SOP Class UID (0008, 0016) tag, which is a tag for designating Service and Object Pairs, indicates a combination of the Object Digital X-ray Image and the Service Storage 1.2.840. .100088.5.4.1.1.1.1.1 is attached.

  The display control module 215 controls the contents displayed on the monitor 106. For example, control is performed to display on the monitor 106 patient information, information on imaging conditions, information indicating the state of the FPD 102, and the like. These are displayed on the monitor 106 together with the above-described corrected image.

  In another embodiment, the display control module 215 may execute display control for displaying the corrected image on the monitor 106 by the output control module 214. In this case, the display control module 215 displays a radiographic image captured on the display screen or displays a corrected image.

  The main body that executes image processing according to the embodiment of the present invention is the control unit 105. Hereinafter, the control unit 105 may be referred to as an image processing apparatus 105 from the viewpoint of describing image processing according to the embodiment of the present invention.

  Note that some or all of the components included in the control unit 105 are not fixed to the control unit 105 and may be realized as an image processing system included in the information system 120. For example, an image processing apparatus that executes an image processing program including the image acquisition module 213, the output control module 214, and the estimation module 220 may be provided separately from the control unit 105 that executes the imaging control module 211. Further, for example, the WS 113 may have some or all of these modules. The PACS 114 may have some or all of these modules. The FPD 102 may include, for example, a field-programmable gate array (FPGA) having an estimation module 220. The components included in the control unit 105 may be included in different devices, and a device that performs processing may be selected according to an instruction from the operator. Furthermore, it may be configured by a workstation, a server, and a storage device connected via a network, and may perform image processing according to the embodiment of the present invention by communicating with these devices as necessary. Good.

  Next, image processing executed by the control unit 105 will be described in detail based on FIG. If there is no notice in the following processing, the subject of processing is the CPU 201 or the GPU 208, and the estimation module 220 is executed by the CPU 201, whereby the processing for estimating the following scattered radiation component is executed.

  In an image obtained by radiography, a scattered X-ray image of scattered X-rays scattered in the subject 103 is superimposed on a primary X-ray image of primary X-rays that linearly reaches each element of the FPD 102 from the X-ray source 101. doing.

Here, M ^′ (x, y) is an image (input image) obtained by radiography, P ^′ (x, y) is a primary X-ray component included in the input image, and S ^′ (x, y) is a scattered X-ray component. y). (X, y) represents the pixel position of the image or the element position of the FPD 102 that gives the pixel value of the pixel position. The relationship between these images can be expressed by Equation 1.

Scattered X-rays are generated in the process in which the X-rays irradiated from the X-ray source pass through the subject 103, and the component that has not been scattered is the primary X-ray, so S ^′ (x, y) is P ^′ ( x, y). However, this correlation is non-linear as represented by the following formula. Therefore, the image processing apparatus 105 optimizes S ^′ (x, y) and P ^′ (x, y) by iterative processing, and uses the maximum likelihood estimation method. There are various other methods for nonlinear optimization, and methods such as the least square method, Newton method, and convex analysis method may be used.

  In step S <b> 300, the CPU 201 executes the image acquisition module 213 to acquire a radiographic image captured by irradiating the subject 103 with radiation. Such a radiation image is an input image that is a target of the scattered radiation estimation process. If a reduced image with a small amount of data is acquired as such an input image and used for estimation processing, data transmission from the FPD 102 and subsequent image processing can be performed at higher speed. Further, since the scattered radiation component is mainly a low frequency component, even if it is estimated from the reduced image, the influence on the estimation accuracy of the scattered radiation component is small.

  In step S301, the CPU 201 normalizes the input image by executing the scattered radiation intensity approximation module 221. For example, the total value of the primary X-ray dose and the scattered X-ray dose that are the radiation dose reaching the FPD 102 is normalized by the incident dose to the subject 103. The incident dose may be a radiation dose detected by the FPD 102 when it is assumed that the subject 103 does not exist. The incident dose may be estimated from the imaging conditions such as the tube voltage, the tube current, the irradiation time, and the imaging region, or may be obtained from the pixel value of the region (elementary region) where the X-ray that has not passed through the subject 103 has arrived. Good. If normalized by the incident dose Q, each can be expressed as follows.

Hereinafter, M (x, y), P (x, y), and S (x, y) are treated as normalized.

In step S302, the CPU 201 sets the initial value P ⁰ (x, y) of P (x, y) in the maximum likelihood estimation method by executing the scattered radiation intensity approximation module 221. P ⁰ (x, y) may be a positive value, but the primary X-ray component P (x, y) is considered to be close to the input image M (x, y). Therefore, in this embodiment, M (x, y) Let y) be an initial value.

  In S303, the scattered radiation intensity approximating module 221 causes the CPU 201 to calculate the first scattered X-ray intensity as the first output data obtained by substituting the input image M into the first function, for example. Here, the first function is a scattered X-ray generated by scattering of X-rays, and is a function (approximate expression) that approximates the scattered X-ray intensity of scattered X-rays that is assumed not to be attenuated after the X-rays are scattered. ). The first output data is the output of the first function. FIG. 4 is a diagram illustrating scattered X-rays approximated by the first function. Scattered X-rays generated on a path L403 passing from a focal point 406 of the X-ray source 101 to a certain position 402 in the subject 103 are not attenuated after the scattering, and are not attenuated after the scattering. Assume that y) is reached. Here, the z axis is set on the path from the focal point 406 to the pixel position 404, and the positive direction of the z axis is the traveling direction of the primary X-ray. The origin is a position where X-rays enter the subject 103. The attenuation coefficient is a constant value μ in the subject 103. At this time, scattered X-rays generated at a position z402 on the path L403 inside the subject 103 are

Is proportional to. For normalized primary X-ray P

Holds. Therefore, the total intensity of the first scattered X-ray generated on the path L403 and reaching the element position 404 (x, y) without being attenuated after scattering is

And is given as a function of the primary X-ray component P (x, y). That is, for example, 1-P is the first function.

  The scattered X-ray approximated by the first function is obtained by scattering the radiation incident on the subject at the position 402, for example, reaching the pixel position 405 in the unexposed region with little attenuation after scattering. It may be considered that it approximates a straight path from the radiation source to the pixel position 404.

  In the image processing apparatus 105, it is preferable to estimate the scattered radiation image based on a plurality of functions different from the first function. In step S304, the scattered radiation intensity approximation module 221 causes the CPU 201 to calculate the second scattered X-ray intensity as second output data obtained by substituting the input image M into the second function, for example. Here, the second function is a scattered X-ray generated by scattering of X-rays, and is a function (approximate expression) that approximates the scattered X-ray intensity of scattered X-rays assumed to be further attenuated in the subject 103 after the scattering. ). The second output data is the output of the second function. FIG. 5 is a diagram illustrating scattered X-rays approximated by the second function. Here, the scattered X-rays mainly transmitted through the subject 103 are approximated. Scattered X-rays generated on a path L503 passing through a certain element position 504 (x, y) of the radiation detection element 102a from the focal point 506 of the X-ray source 101 pass through the path L503, attenuate, and enter the FPD 102. Assume that. Here, the z axis is set on the path from the focal point 506 to the element position 504, and the positive direction of the z axis is the traveling direction of the primary radiation. The origin is a position where X-rays enter the subject 103. The attenuation coefficient is μ. At this time, the scattered X-rays generated at a position z502 on the path L503 inside the subject 103 are proportional to Equation 2. Further, assuming that the object 103 is transmitted through the distance Lz and attenuated, the total intensity of the second scattered X-ray generated in the path L is

And is given as a function of the primary X-ray component P (x, y). That is, for example, -PlnP is the second function.

Considering that scattered X-rays are scattered at a wide angle in the subject 103, it is preferable to estimate the scattered radiation component by combining specific frequency components for the output data of each function. In step S305, the frequency processing module 222 causes the CPU 201 to perform frequency processing, thereby decomposing output data obtained by the processing up to step S304 into a plurality of frequency components. First, the scattered X-ray component Since components of low frequency is mainly, the Gaussian function G ₁ to convolution. The low frequency component for the first output data is

It can be expressed. Furthermore, in consideration of the discontinuity of scattered X-rays, it is preferable to add components other than low frequencies. Components other than the low frequency of the scattered radiation image by the first scattered X-ray model are

It can be expressed. Here, F represents Fourier transform, and F ⁻¹ represents Fourier inverse transform.

The CPU 201 executes the frequency processing module 222 to synthesize a plurality of frequency components obtained by the above-described processing by weighting them with the parameters A _{1 and} C ₁ . This gives a frequency component S ₁ for the first output data.
According to the processing results so far, S ₁ is

It can be expressed. Here, A ₁ is the weight of the low-frequency component of the first scattered X-ray, B ₁ is the spread of the first scattered X-ray, and C ₁ is the parameter representing the weight of the component other than the low-frequency of the first scattered X-ray. It is.

In step S306, a specific frequency component is obtained for the second output data as in step S305. Low-frequency component is a Gaussian function G ₂ convolution,

It can be expressed. Components other than low frequency are

It can be expressed.

By executing the frequency processing module 222, the CPU 201 synthesizes the plurality of frequency components obtained by the above processing by weighting them with the parameters A ₂ and C ₂ . This gives the frequency components S ₂ for the second output data.
According to the processing result up to this point, S ₂ is

It can be expressed. Here A ₂ is the weight of the low-frequency component of the second scattered X-rays, B ₂ is the spread of the second scattered X-rays, C ₂ represents the weight of components other than the low-frequency of the second scattered X-ray parameter It is.

Here, G ₁ and G ₂ may be the same, but it is preferable to divide the first output data and the second output data into a plurality of different frequency bands, and combine them by weighting.

  The two types of components represented by Formula 9 and Formula 11 are exemplified as the frequency components of the first output data, and the two types of components represented by Formula 13 and Formula 15 are illustrated as the frequency components of the second output data. At least one of the two types of frequency components of the first output data and at least one of the two types of frequency components of the second output data may be weighted and combined to form a scattered radiation component. .

The above processing can be a frequency component S ₁ of the first output data based on the frequency components S ₂ of the second output data to obtain the scattered ray components. That is, the scattered radiation model in the image processing apparatus 105 is based on a specific frequency component for the output of the first function and a specific frequency component for the output of the second function different from the first function. This is a scattered radiation model for obtaining scattered radiation components.

A ₁ , B ₁ , C ₁ , A ₂ , B ₂ , and C ₂ may be stored in the RAM 202 in a plurality of combinations, and may be set automatically or appropriately by the operator. You may have the setting module which sets the weight of a synthesis | combination according to the tube current of the X-ray tube, tube voltage, the X-ray irradiation time, the imaging | photography site | part etc. which are the imaging conditions of the said radiographic image. For example, the setting module acquires information on the above-described shooting conditions in response to an operation input to the operation unit 107, and sets a parameter set (A ₁ , B ₁ , C ₁ , A ₂ , B ₂ , C corresponding to the shooting conditions). ₂ ) is obtained. The table information indicating the correspondence between the shooting conditions and the parameter set is stored in the SSD 208, and the CPU 201 sets the parameter set by referring to the table information. The table information is experimentally determined in advance for each shooting condition, for example.

As described above, since the scattered X-ray image mainly includes low-frequency components, the parameters A ₁ and A ₂ representing the weights of the low-frequency components are derived from the parameters C ₁ and C ₂ representing the weights of the components other than the low-frequency components. Is preferably a large value.

  Note that the frequency processing module 222 may create a low frequency component and a component other than the low frequency using Fourier transform, or may be performed in real space. The filter used for frequency processing is not limited to a Gaussian function. For example, if the weighted average filter is used, the spread around the periphery of the target pixel can be expressed, and if the moving average filter is used, the spread to the pixel far from the target pixel can be expressed.

  Moreover, although the order of S303, S304, S305, and S306 was illustrated, the processing may be performed in the order of S303, S305, S304, and S306, and the processing of S303 and S305 and the processing of S304 and S306 are performed in parallel. You may do it.

In S307, the optimization module 223 causes the CPU 201 to combine the first output frequency component calculated up to S306 and the second output frequency component to obtain the scattered radiation component, and to obtain the maximum of the primary X-ray component. Let the likelihood estimate be performed. The maximum likelihood estimation is, for example, as expressed in Expression 17, with the input image M (x, y), the primary X-ray component P ⁿ (x, y) obtained in the previous stage of the iteration, and the first output. The data is divided by the sum of the frequency component S ₁ ⁿ (x, y) of the data and the frequency component S ₂ ⁿ (x, y) of the second output data, and is multiplied by the primary X-ray component P ⁿ (x, y). Here, the superscript n indicates that it was obtained in the nth iteration.

In maximum likelihood estimation based on Equation 17, every iteration

Is known to converge to 1. Sequentially, the primary X-ray component P (x, y), the first output data frequency component S ₁ (x, y), and the second output data frequency component S ₂ so as to satisfy the relationship of Equation 1. (X, y) is determined.

  In S308, the optimization module 223 causes the CPU 201 to determine whether to continue iteratively. It is determined by a square residual or the like whether Equation 18 is sufficiently close to 1 or whether the relationship of Equation 1 is sufficiently satisfied. Alternatively, the determination may be made based on whether a fixed number of iterations have been performed.

If it is determined in S308 to continue the iteration, the process returns to S302 with the calculated P ^{n + 1} (x, y) as an initial value and the process is continued. If it is determined that the repetition is not continued, the process is terminated, and P ^{n + 1} (x, y) is set as a corrected image in which the scattered X-ray component is reduced. Alternatively, the calculated sum of S ₁ ^{n + 1} (x, y) and S ₂ ^{n + 1} (x, y) may be multiplied by a parameter to be subtracted from the input image. That is, the effect of reducing the scattered radiation image may be adjustable.

As described above, in the image processing apparatus 105, the scattered radiation component S ₁ (x, y) + S ₂ (x, y) obtained by combining specific frequency components for each of the first output data and the second output data. Is output as a scattered radiation component. Further, under the assumption that the sum of the scattered ray component and the primary X-ray component matches the input image, the scattered ray component included in the radiation image is derived by successive approximation.

In step S309, the CPU 201 executes the output control module 214 to output the corrected image obtained by the above-described processing. In this process, the corrected image is stored in the PACS 114 and displayed on the monitor 106. The estimated scattered X-ray image S ₁ ^{n + 1} (x, y) + S ₂ ^{n + 1} (x, y) may be stored in the PACS 114 as image data or an image file different from the corrected image.

  In another aspect, the first function and the second function are used to calculate the total intensity of the scattered radiation generated by the scattering of the radiation on the path until the radiation reaches the radiation detection element 102a that provides the pixel value of the radiation image. This is a function approximated for each path by a plurality of scattered ray models with different distances attenuated on the path. Furthermore, in the above-described processing, the approximated processing is uniformly performed on each pixel value corresponding to the output from the radiation detection element to estimate the scattered radiation image. Thereby, it is not necessary to switch the model for each pixel value included in the radiation image, and as a result, time and resources required for the scattered radiation image estimation process can be reduced.

  Further, as shown in Equation 19, the distance attenuated after scattering can be adjusted by a parameter. For example, assuming that after scattering, the object 103 is attenuated at distance x and not otherwise attenuated:

It can be expressed. It is assumed that x satisfies 0 <x <(Lz).
Assuming that it is attenuated in the subject 103 after scattering but not at distance y,

It can be expressed. It is assumed that y satisfies 0 <y <(Lz).

X and Y can be included in the parameters A ₁ , B ₁ , C ₁ , A ₂ , B ₂ , C ₂ . Therefore, the present invention is not limited to these two models, but includes models approximated with various attenuation distances by appropriately selecting appropriate parameters.

The parameters A ₁ , B ₁ , C ₁ , A ₂ , B ₂ , C ₂ may be multiplied by any term within the range of equivalent deformation. That is, the first function and the second function are not limited to the concept of approximating the scattered radiation intensity as described above, and those obtained by multiplying other parameters by each parameter are also included in the present invention.

  Further, the first function in the image processing apparatus 105 is not limited to 1-P. 1-P includes (1-1) the property that the scattered radiation component is 0, that is, the minimum value when the primary radiation component is 1, that is, the maximum value in the normalized radiation image, and (1-2) the primary When the radiation component is 0, that is, the minimum value, the scattered radiation component has a property of 1, that is, a positive value, and (1-3) the property that the primary radiation component converges to a positive value as it approaches the minimum value. It is a function that has. 1-P is a function having the property that (1-4) the primary radiation component monotonously decreases within the range of 0 to 1. 1-P is a function having a property that the derivative takes a negative value when the (1-5) primary radiation component is zero. That is, a function satisfying part or all of the properties of (1-1) to (1-5), and information indicating the correspondence between the primary radiation component and the value of the scattered radiation component is the first It may be a function. For example, Formula 21 is mentioned.

Furthermore, the second function in the image processing apparatus 105 is not limited to −PlnP. -PlnP has the following characteristics: (2-1) In the normalized radiation image, when the primary radiation component is 1 or the maximum value, the scattered radiation component is 0 or the minimum value. (2-2) Primary radiation When the component is 0, that is, the minimum value, the scattered radiation component is 0, that is, the minimum value, and (2-3) When the primary radiation component is neither 0 nor 1, the scattered radiation component has a positive value. And a function having the property. -PlnP is a function having the property that (2-4) the primary radiation component takes a maximum value within the range of 0 to 1. -PlnP is a function having the property that the derivative takes a positive value when the (2-5) primary radiation component is zero. That is, a function that satisfies a part or all of such properties and that indicates the correspondence between the primary radiation component and the value of the scattered radiation component may be used as the second function. For example, Formula 22 is mentioned.

The first function or the second function may not be represented by the mathematical formula as described above, and may be represented by, for example, a lookup table that represents an output value for the value of the primary X-ray component P. .

  Based on the first function, the scattered X-ray image can be accurately estimated when primary X-rays are greatly reduced and many scattered X-rays are generated, for example, when the subject is thick. For example, the primary X-ray hardly passes through the thick portion of the subject, and the generated scattered X-ray hardly passes. Specifically, the scattered X-ray intensity when P takes a value close to zero is estimated to have a scattered X-ray intensity close to zero according to the second function. However, according to the first function, since it is estimated to be a value close to 1, it is possible to estimate with high accuracy without estimating the scattered X-rays generated originally. Furthermore, in consideration of the large spread of scattered X-rays, an area where the subject does not exist, or a portion where the subject is thin, is decomposed into a plurality of frequency components and synthesized using parameters that take into account the spread of each component. It is possible to accurately estimate scattered X-rays that spread out. As a result, it is possible to more accurately estimate the portion where the subject does not exist and the scattered radiation images of the lung field region 904 and the skin line 903 shown in FIG. 9A.

  In yet another aspect, first-order scattered X-rays are generated by attenuation of incident X-rays, and the second function −PlnP is a multi-order scattering in the process in which the first-order scattered X-rays are attenuated in the object 103. It may be considered as a function based on the assumption that no X-rays are generated. According to this viewpoint, the first function 1-P is such that the primary scattered X-rays are generated by the attenuation of the incident X-rays and the primary scattered X-rays pass through the object 103 and are attenuated. It can be considered as a function based on the assumption that multi-order scattered X-rays are generated without being absorbed by 103. Actually, scattered X-rays generated in the subject 103 can generate multi-order scattered X-rays. That is, combining the output data of the first function 1-P with the output data of the second function -PlnP is combined with the multi-order scattered X-rays that can be generated from the X-rays incident on the subject 103. In other words, the ratio of primary scattered X-rays is weighted and synthesized. Therefore, it is possible to appropriately express how the scattered X-rays generated in the subject 103 are attenuated, and improve the quantitativeness of the estimated scattered radiation component. That is, the scattered radiation component can be estimated more accurately. By obtaining and synthesizing specific frequency components for the output data of the two types of functions illustrated, it is possible to take into account the difference in the spread of scattered X-rays that are attenuated differently. Component estimation accuracy is improved.

In addition to the first function and the second function, the third function may be used from the various functions described above. In this case, the scattered radiation is obtained using the output data obtained by substituting P ⁿ (x, y), which is the approximate solution calculated by the input image or Expression 17, for each of the first, second and third functions. Estimate components. The third function may be the same as the first function or the second function, but it is preferable to synthesize different frequency components for each output data. For example, to obtain a specific frequency component using a Gaussian function G ₃ in the third output data, combined by weighting the frequency component and the second output data of the frequency component of the first output data as described above May be.

As a scattered radiation model applied to a radiographic image, it is only necessary to be able to express that the distance at which the scattered radiation is attenuated in the subject 103 is various, and that the generation of multi-order scattered radiation is varied from another viewpoint. . That is, even if the function is single, at least two frequency components obtained by decomposing into a plurality of frequency components may be weighted and synthesized. The behavior of scattered radiation that can be expressed is different due to different parameters, and it is possible to estimate the change in the distribution of scattered radiation by considering at least two ways of spreading the component on the low frequency side and the component on the high frequency side. it can. For example, the respective components expressed by Equation 9 and Equation 11 may be combined by weighting.

  The phenomenon in which X-rays are attenuated in the subject 103 is a phenomenon in which the X-ray intensity reaching the FPD 102 is reduced as compared to the X-rays incident on the subject 103 due to scattering or absorption of the X-rays. The scattered X-rays generated in the object 103 and attenuated according to the transmission distance after scattering are more scattered X-rays generated once scattered in the process of transmission after scattering. Although the next scattered X-rays are generated or absorbed, the multi-order scattered rays are scattered X-rays assuming that only the primary scattered X-rays are detected by the FPD 102 without reaching the FPD 102. Such scattered X-rays are considered to be dominant in the portion where the thickness of the subject is large. Scattered X-rays generated in the subject 103 that are not attenuated after scattering are those in which radiation is not absorbed in the subject 103 and all the generated multi-dimensional scattered X-rays are detected by the FPD 102. The scattered X-rays assuming It is considered that not only the primary scattered X-rays but also the multi-order scattered X-rays are included in the portion where the thickness of the subject is small and the scattered X-rays that have scattered and reach the unexposed region. Therefore, based on the functions assuming these two cases, the manner in which the radiation is attenuated or scattered in the subject 103 can be expressed more accurately.

  Although the scattered radiation model has been described from various viewpoints, naturally, the modified examples of the functions shown so far can be used when the scattered radiation model is viewed from any viewpoint, and are included in the present invention.

  As described above, based on the first function and the second function different from the first function, the scattered radiation image can be estimated with high accuracy. In addition to the low-frequency component, weighted and synthesized at least two frequency components obtained by decomposing into a plurality of frequency components, the scattered radiation component including the boundary where the scattered radiation distribution changes is estimated. it can.

  Next, another embodiment of the present invention will be described. Similar to the embodiment described above, the scattered radiation component is determined based on the specific frequency component for the output of the first function and the specific frequency component for the output of the second function different from the first function. By applying the obtained scattered radiation model to the radiation image, the scattered radiation component included in the radiation image is estimated. Here, the first output data and the second output data are decomposed into three or more frequency bands to obtain a frequency component.

  FIG. 6 is a flowchart showing processing according to another embodiment of the present invention. Steps S601, S602, S603, S606, S609, and S610 are the same as S301, S302, S303, S305, S307, and S308 shown in FIG.

  In step S603, first output data is obtained in the same manner as in step S303.

  In step S604, the frequency processing module 222 causes the CPU 201 to perform multiple frequency decomposition processing on the first output data calculated in step S603. In this embodiment, versatile and fast Laplacian pyramid decomposition is used. Other methods for decomposing into frequency bands include fast Fourier transform and wavelet decomposition, which may be used.

FIG. 7A is a conceptual diagram of Laplacian pyramid decomposition. If the first output data calculated in S603 is L ₁

It is expressed. To L ₁ after performing aliasing processing in the two-dimensional low-pass filter AF shown in FIG. 7 (b), it is possible to obtain a low-frequency image L ₂ by downsampling the size of the vertical and horizontal 1/2.

Here, ↓ indicates 1/2 downsampling. It can be obtained sequentially low-frequency image L _n in the same manner.

Further, the high-frequency image H _n can be obtained by up-sampling L _{n + 1} twice, performing an aliasing process using a two-dimensional low-pass filter AF as shown in FIG. 7B, and subtracting from L _n. it can.

Here, ↑ indicates double upsampling. Equation 15, by performing sequentially the processing of Equation 16 is divided into n types of frequency bands, can be created L _n image and H _n images. In this way, images of a plurality of frequency bands are created.

FIG. 8 is a diagram showing frequency characteristics of frequency decomposition by the Laplacian pyramid. Image L _n of each frequency band generated by the Laplacian pyramid decomposition of _S604, an example of the frequency characteristic of H _n. For example, if a disassembled image aspect 1000 pixels to _{L 10,} are degraded to the level of vertical and horizontal numbers of pixels in the _{L 10.} Further, the high-frequency component H _n of each level and is obtained by band-limiting the input image L _1.

In step S605, the frequency processing module 222 converts the frequency-resolved images L _n and H _n that are images obtained by decomposing the first output data obtained in step S603 into a plurality of frequency bands in step S604. Is reconfigured by being executed by the CPU 201. In this process, as shown in FIG. 7C, the L _{n + 1} image at each level is up-sampled twice in the vertical and horizontal directions, aliased by a low-pass filter AF as shown in FIG. 7B, and added to H _n . It can be reconstructed _{L n} by. For example, if we break down to L ₁₀ then L ₉ ^'

It is expressed. Therefore, L _n ^′ is

It can be expressed. This calculation is performed sequentially, and in this embodiment, the reconstructed L ₁ ^′ is set as the frequency component S ₁ of the first output data. Here, α _n is a coupling constant. By adjusting α _n , it is possible to represent S ₁ by weighting and synthesizing n frequency components.

  In step S606, second output data is obtained in the same manner as in step S305.

In step S607, the second output data calculated in S606 is decomposed into a plurality of frequency-resolved images L _n and H _n as in S604.

In step S608, similarly to S605, the frequency-resolved images L _n and H _n calculated in S607 are reconstructed. That is,

And L _n ^′ is

It is expressed. This calculation is sequentially performed, and in this embodiment, the reconstructed L ₁ ^′ is set as the first scattered X-ray image S ₂ . Here, β _n is a coupling constant. By adjusting β _n , n frequency components can be weighted and combined to represent S ₂ .

  Also in this embodiment, the first method is to approximate the total intensity of the scattered radiation generated on the path where the radiation generated from the X-ray source 101 linearly reaches the radiation detection element 102a that gives the pixel value of the radiation image. Description has been made using two of the scattered X-ray model and the second scattered X-ray model. However, as shown in Equations 19 and 20, the assumption of the attenuation distance after scattering can be adjusted by parameters. Therefore, the present invention is not limited to these two models, but can include models approximated with various attenuation distances by appropriately selecting appropriate parameters. Further, for example, if a weighted average filter is used as the low-pass filter, the spread around the periphery of the target pixel can be expressed, and if the moving average filter is used, the spread to a pixel far from the target pixel can be expressed.

In the present embodiment, the scattered radiation component can be expressed by a plurality of frequency components by adjusting α _n and β _n . Further, at least one of the first output data and the second output data may be decomposed into three or more frequency band components. By decomposing into three or more frequency band components, the low-frequency components scattered and greatly spread in the body, and the structure in the skin line and lung fields where the absorption of X-rays in the subject 103 changes greatly, It can be estimated more accurately.

  Also in this embodiment, modified examples of the scattered radiation model and the function captured from various viewpoints described above can be used, and are included in the present invention.

  In the present embodiment, the scattered radiation component is estimated by performing a fitting process using a plurality of functions on a two-dimensional image obtained by radiography. The present embodiment includes an example in which the above-described scattered radiation image estimation processing is applied to a tomographic projection image. Alternatively, an example in which the above-described scattered radiation image estimation processing is applied to a radiation image obtained by so-called general imaging in which radiation is applied to the subject 103 from one direction and imaging using a two-dimensional detector is also applied to the present invention. included. Since it is not necessary to extract a portion where radiation absorption changes by an analysis process such as segmentation, the scattered radiation component can be estimated by a more robust process.

  In addition, it is preferable that on / off of the scattered radiation reduction process can be switched according to the imaging region.

  FIG. 9B is a diagram illustrating an example in which the scattered radiation component included in the radiation image acquired using the chest phantom is acquired using the image processing according to the present invention. Similar to (a), which illustrates an experimentally acquired scattered radiation component, the intensity distribution in the unexposed region 905 and the lung region 908 and the skin line 907 are represented.

  Although the image processing apparatus in the above-described embodiment is a single apparatus, a form in which the above-described processing is executed by an image processing system in which apparatuses including a plurality of information processing apparatuses are communicably combined is also included in the present invention. . Or it is good also as performing the above-mentioned process with the server apparatus or server group which is common by several modalities. In this case, the common server device corresponds to the image processing device according to the embodiment, and the server group corresponds to the image processing system according to the embodiment. The information system 120 or the plurality of devices constituting the image processing system need only be able to communicate at a predetermined communication rate, and does not need to exist in the same facility or in the same country.

  Further, the present invention provides a form in which a software program that realizes the functions of the above-described embodiments is supplied to a system or apparatus, and a computer of the system or apparatus reads and executes the code of the supplied program. Including.

  Accordingly, the program code itself installed in the computer in order to realize the processing according to the present embodiment by the computer is also one aspect of the present invention. Further, based on instructions included in the program read by the computer, an OS or the like running on the computer performs part or all of the actual processing, and the functions of the above-described embodiments are also realized by the processing. .

  A form in which the above-described embodiments are appropriately combined is also included in the present invention.

Claims (27)

Image acquisition means for acquiring a radiation image obtained by irradiating the subject with radiation;
A first function corresponding to the intensity of scattered radiation including multi-order scattered radiation generated by scattering the radiation a plurality of times in the subject, and the intensity of primary scattered radiation generated by scattering the radiation once in the subject. based on the second function corresponding to, acquired by the image acquiring unit by using the specific frequency component of the output of the second function with a particular frequency component of the output of said first function An estimation means for estimating the scattered radiation component contained in the radiographic image;
And an output unit that outputs a corrected image in which a scattered radiation component estimated by the estimation unit is reduced from the radiation image.   The estimation means estimates the scattered radiation component by weighting and synthesizing a specific frequency component for the output of the first function and a specific frequency component for the output of the second function. The image processing apparatus according to claim 1, wherein the image processing apparatus is configured as described above.   The image processing apparatus according to claim 2, wherein the specific frequency component for the output of the first function and the specific frequency component for the output of the second function have different frequency bands.   The estimation unit decomposes the output of the first function and the output of the second function into a plurality of different frequency components, and weights and synthesizes the plurality of different frequency components, thereby obtaining the radiation image. The image processing apparatus according to claim 3, wherein the image processing apparatus is configured to estimate a scattered radiation component included in the image processing apparatus.   5. The plurality of different frequency components are components decomposed into three or more frequency bands for at least one of the output of the first function and the output of the second function. An image processing apparatus according to 1. The first function is a correspondence relationship between a value of a primary radiation component that is a component of radiation that travels straight from a radiation generation device to a radiation detection element that provides a pixel value of the radiation image, and a value of a scattered radiation component of the radiation image. The scattered radiation component has a minimum value when the primary radiation component has the maximum value in the radiation image, and the scattered radiation component has a positive value when the primary radiation component has the minimum value in the radiation image. 6. The image processing according to claim 1, wherein the scattered radiation component converges to a positive value as it approaches a minimum value in the radiation image. apparatus. The first function is a correspondence relationship between a value of a primary radiation component that is a component of radiation that travels straight from a radiation generation device to a radiation detection element that provides a pixel value of the radiation image, and a value of a scattered radiation component of the radiation image. A function in which, when the primary radiation component is the maximum value in the radiation image, the scattered radiation component has a minimum value, and monotonically decreases in the range of the primary radiation component in the radiation image. The image processing apparatus according to claim 1, wherein the image processing apparatus is an image processing apparatus.   The image processing apparatus according to claim 1, wherein the first function is a function that assumes that scattered radiation generated by scattering of the radiation is not attenuated. Normalizing means for normalizing a pixel value of the radiation image with a value corresponding to an incident dose to the subject;
5. The first function is represented by an approximate expression 1-P, where P is a primary radiation component in the normalized radiation image. 6. An image processing apparatus according to 1. The second function is a correspondence relationship between a value of a primary radiation component that is a component of radiation that travels straight from a radiation generator to a radiation detection element that provides a pixel value of the radiation image, and a value of a scattered radiation component of the radiation image. The scattered radiation component is a minimum value when the primary radiation component is the maximum value in the radiation image, and the scattered radiation component is positive when the primary radiation component is not the maximum value and the minimum value. the image processing apparatus according to any one of claims 1 to 5, characterized in that a function as a value. The second function is a correspondence relationship between a value of a primary radiation component that is a component of radiation that travels straight from a radiation generator to a radiation detection element that provides a pixel value of the radiation image, and a value of a scattered radiation component of the radiation image. The primary radiation component is the minimum value when the primary radiation component is the maximum value and the minimum value in the radiation image, and the scattered radiation component is the minimum value, and when the primary radiation component is not the maximum value and the minimum value, scattered radiation component is a positive value, the image of any one of claims 1 to 5, characterized in that a function has a maximum value in a range of the primary radiation component in the radiographic image in Processing equipment. The second function, scattered radiation the radiation is generated by scattering in any one of claims 1 to 5, characterized in that a function which is assumed to be further attenuated in the subject The image processing apparatus described. Normalizing means for normalizing a pixel value of the radiation image with a value corresponding to an incident dose to the subject;
The said 2nd function is represented by approximate expression -PlnP, when P is made into the primary radiation component in the said normalized radiographic image, It is any one of Claim 1 thru | or 5 characterized by the above-mentioned. The image processing apparatus described. The plurality of different frequency components are obtained by decomposing the output of the first function and the output of the second function into two different frequency components,
The said estimation means is comprised so that the scattered-radiation component contained in the said radiographic image may be estimated by weighting and combining two said different frequency components, The Claim 1 thru | or 5 characterized by the above-mentioned. The image processing apparatus according to any one of claims.   The estimation means is configured to estimate a scattered radiation component included in the radiographic image by combining a frequency component on a low frequency side among the two different frequency components by increasing a weight. The image processing apparatus according to claim 14, wherein The first function is a correspondence relationship between a value of a primary radiation component that is a component of radiation that travels straight from a radiation generation device to a radiation detection element that provides a pixel value of the radiation image, and a value of a scattered radiation component of the radiation image. Wherein the primary radiation component has a minimum value in the radiation image, and the derivative has a negative value.
The second function is a correspondence between a value of a primary radiation component that is a component of radiation that travels straight from the radiation generation device to a radiation detection element that provides a pixel value of the radiation image, and a value of a scattered radiation component of the radiation image. a function showing a relation, the primary radiation component, any of claims 1 to 5, wherein the derivative when the minimum value in the radiographic image in is a function which is a positive value An image processing apparatus according to claim 1. Image acquisition means for acquiring a radiation image obtained by irradiating the subject with radiation;
A scattered radiation model that obtains a scattered radiation component based on a specific frequency component for the output of the first function and a specific frequency component for the output of the second function different from the first function is the radiation image. Applying to the estimation means for estimating the scattered radiation component included in the radiation image;
And an output unit that outputs a corrected image in which a scattered radiation component estimated by the estimation unit is reduced from the radiation image. Image acquisition means for acquiring a radiation image obtained by irradiating the subject with radiation;
Estimating means for estimating a scattered radiation component included in the radiographic image by weighting and synthesizing at least two frequency components obtained by decomposing an output of a specific function into the radiographic image into a plurality of frequency components; ,
Output means for outputting a corrected image in which the scattered radiation component is reduced from the radiation image. Tube current, irradiation time, using an imaging condition including at least one of imaging regions, the image processing apparatus according to claim 18, characterized in that it comprises a setting means for setting the parameters of the synthesis. Using the primary radiation component of the radiation image as an input, the first function output obtained by applying the first function to the primary radiation component and the output of the second function obtained by applying the second function Further comprising processing means for obtaining as an output a scattered radiation component of the radiation image by combining the specific frequency components for each of the
The estimation means is the radiation acquired by the image acquisition means under the assumption that the sum of the primary radiation component of the radiation image and the scattered radiation component of the radiation image obtained by the processing means matches the radiation image. 20. The apparatus according to claim 1, wherein the scattered radiation component of the radiation image acquired by the image acquisition unit is derived by successive approximation using the image as an input. An image processing apparatus according to 1.   The image processing apparatus according to claim 20, wherein a maximum likelihood estimation method is used as the successive approximation method.   The image processing apparatus according to any one of claims 1 to 21, wherein the image acquisition unit acquires a radiographic image captured by irradiating a subject with radiation from one direction. . The image acquisition means first receives a reduced image of the radiation image having a data amount smaller than that of the radiation image from a radiation detector for obtaining the radiation image, and then completes reception of the radiation image. ,
The estimating means estimates a scattered radiation component of the radiation image based on the reduced image;
The image processing according to any one of claims 1 to 22, wherein the output unit is configured to output a corrected image in which the scattered radiation component is reduced from the radiation image. apparatus.   The image processing apparatus according to any one of claims 1 to 23, wherein the output unit outputs a DICOM image file in which DX or RF is attached to the corrected image as a Modality tag. Image acquisition means for acquiring a radiation image obtained by irradiating the subject with radiation;
A first function corresponding to the intensity of scattered radiation including multi-order scattered radiation generated by scattering the radiation a plurality of times in the subject, and the intensity of primary scattered radiation generated by scattering the radiation once in the subject. based on the second function corresponding to, acquired by the image acquiring unit by using the specific frequency component of the output of the second function with a particular frequency component of the output of said first function An estimation means for estimating the scattered radiation component contained in the radiographic image;
And an output unit that outputs a corrected image in which the scattered radiation component estimated by the estimation unit is reduced from the radiation image. An image acquisition step of acquiring a radiation image obtained by irradiating the subject with radiation;
A first function corresponding to the intensity of scattered radiation including multi-order scattered radiation generated by scattering the radiation a plurality of times in the subject, and intensity of primary scattered radiation generated by scattering the radiation once in the subject. based on the second function corresponding to, acquired by the image acquiring unit by using the specific frequency component of the output of the second function with a particular frequency component of the output of said first function An estimation step for estimating the scattered radiation component contained in the radiographic image;
And an output step of outputting a corrected image in which the scattered radiation component estimated by the estimating means is reduced from the radiation image. An image acquisition process for acquiring a radiation image obtained by irradiating a subject with radiation;
A first function corresponding to the intensity of scattered radiation including multi-order scattered radiation generated by scattering the radiation a plurality of times in the subject, and the intensity of primary scattered radiation generated by scattering the radiation once in the subject. On the basis of the second function corresponding to the first function and the specific frequency component for the output of the first function and the specific frequency component for the output of the second function. An estimation process for estimating the scattered radiation component contained in the radiographic image,
An image processing program for causing a computer to execute an output process for outputting a corrected image in which a scattered radiation component estimated by the estimation unit is reduced from the radiation image.