WO2013089155A1 - Ｘ線ｃｔ装置および散乱ｘ線補正方法 - Google Patents

Ｘ線ｃｔ装置および散乱ｘ線補正方法 Download PDF

Info

Publication number
WO2013089155A1
WO2013089155A1 PCT/JP2012/082251 JP2012082251W WO2013089155A1 WO 2013089155 A1 WO2013089155 A1 WO 2013089155A1 JP 2012082251 W JP2012082251 W JP 2012082251W WO 2013089155 A1 WO2013089155 A1 WO 2013089155A1
Authority
WO
WIPO (PCT)
Prior art keywords
ray
subject
distribution
rays
scattered
Prior art date
Application number
PCT/JP2012/082251
Other languages
English (en)
French (fr)
Other versions
WO2013089155A8 (ja
Inventor
植木　広則
史人 渡辺
康隆 昆野
小嶋　進一
Original Assignee
株式会社 日立メディコ
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 株式会社 日立メディコ filed Critical 株式会社 日立メディコ
Priority to CN201280060843.XA priority Critical patent/CN103987320B/zh
Priority to US14/364,278 priority patent/US9307949B2/en
Priority to JP2013549295A priority patent/JP5815048B2/ja
Publication of WO2013089155A1 publication Critical patent/WO2013089155A1/ja
Publication of WO2013089155A8 publication Critical patent/WO2013089155A8/ja

Links

Images

Classifications

    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
    • A61B6/52—Devices using data or image processing specially adapted for radiation diagnosis
    • A61B6/5258—Devices using data or image processing specially adapted for radiation diagnosis involving detection or reduction of artifacts or noise
    • A61B6/5282—Devices using data or image processing specially adapted for radiation diagnosis involving detection or reduction of artifacts or noise due to scatter
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
    • A61B6/02—Devices for diagnosis sequentially in different planes; Stereoscopic radiation diagnosis
    • A61B6/03—Computerised tomographs
    • A61B6/032—Transmission computed tomography [CT]
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N23/00—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00
    • G01N23/02—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material
    • G01N23/04—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material and forming images of the material
    • G01N23/046—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material and forming images of the material using tomography, e.g. computed tomography [CT]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/00—2D [Two Dimensional] image generation
    • G06T11/003—Reconstruction from projections, e.g. tomography
    • G06T11/005—Specific pre-processing for tomographic reconstruction, e.g. calibration, source positioning, rebinning, scatter correction, retrospective gating
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
    • A61B6/42—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment with arrangements for detecting radiation specially adapted for radiation diagnosis
    • A61B6/4291—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment with arrangements for detecting radiation specially adapted for radiation diagnosis the detector being combined with a grid or grating
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
    • A61B6/52—Devices using data or image processing specially adapted for radiation diagnosis
    • A61B6/5205—Devices using data or image processing specially adapted for radiation diagnosis involving processing of raw data to produce diagnostic data
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
    • A61B6/58—Testing, adjusting or calibrating apparatus or devices for radiation diagnosis
    • A61B6/582—Calibration
    • A61B6/583—Calibration using calibration phantoms
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2211/00—Image generation
    • G06T2211/40—Computed tomography
    • G06T2211/416—Exact reconstruction

Definitions

  • the present invention relates to an X-ray CT apparatus and a scattered X-ray correction method for data photographed by an X-ray CT apparatus or the like.
  • the X-ray CT (Computed Tomography) apparatus is an apparatus for reconstructing the difference in the X-ray attenuation rate (X-ray absorption coefficient) inside the subject as an image using a data processing system.
  • the X-ray CT apparatus includes an X-ray source that irradiates an object with X-rays and an X-ray detector that detects X-rays transmitted through the object.
  • the X-ray source and the X-ray detector are arranged so as to face each other with the subject interposed therebetween.
  • the X-ray source and the X-ray detector rotate around the subject while maintaining the relationship of being placed oppositely across the subject. Take line transmission image data.
  • an X-ray tube that generates X-rays by irradiating an electron accelerated at a high voltage to an anode is usually used.
  • the X-ray detector has a structure in which X-ray detection elements are arranged two-dimensionally so as to image a wide area at a high speed at a time.
  • Projection data obtained by imaging with an X-ray CT apparatus includes incident intensity of X-rays scattered by a subject (scattered X-rays) in addition to the intensity of X-rays (direct X-rays) transmitted through the subject without scattering. Information is also included.
  • a scattered radiation prevention grid for removing scattered X-rays generated in the subject is disposed on the X-ray source side of the X-ray detector.
  • scattered X-ray correction by software is also performed (see, for example, Patent Documents 1 to 3 below).
  • X-ray CT apparatuses are capable of imaging a wide range of subjects with one imaging as the number of X-ray detectors increases.
  • the number of X-ray detectors increases the width of the X-ray irradiation on the subject, and the scattered X-ray dose increases accordingly, resulting in a false image (artifact) on the reconstructed image.
  • the ratio of the detection signal derived from scattered X-rays to the detection signal derived from direct X-rays increases relatively, and the X-rays of the subject The absorption coefficient may be underestimated.
  • FIG. 7 is an explanatory diagram schematically showing artifacts generated in the reconstructed image.
  • FIG. 7A shows a human body simulated structure (hereinafter referred to as a phantom) in which two high absorber rods H are present inside the low absorber L.
  • FIG. 7A When this phantom is photographed, a dark band artifact D with a CT value measured smaller than the actual value is observed between the two superabsorbent rods H, as in the reconstructed image shown in FIG. .
  • the magnitude of the CT value is illustrated by hatching interval (density), and the smaller the CT value, the shorter (dense) the hatching interval is.
  • the present invention provides an X-ray CT apparatus and a scattered X-ray correction method that can prevent deterioration in image quality due to scattered X-rays by estimating and correcting the scattered X-ray dose with high accuracy. Let it be an issue.
  • the present invention provides an X-ray detection in which an X-ray source that generates X-rays from an X-ray focal point and an X-ray detection element for detecting the X-rays are two-dimensionally arranged.
  • An imaging unit that rotates around the subject while maintaining a relationship of being opposed to each other across the subject, and photographing the X-ray transmission image data of the subject from a plurality of projection directions, and photographing by the photographing unit
  • An internal distribution estimation unit for estimating an X-ray absorption coefficient distribution inside the subject based on the X-ray transmission image data, and a simulated subject having the X-ray absorption coefficient distribution estimated by the internal distribution estimation unit.
  • a Monte Carlo simulation that simulates the physical interaction of the X-rays is performed, and a point spread function estimation unit that estimates a point spread function of scattering derived from the subject and a point spread function estimation unit.
  • a correction unit that corrects the X-ray transmission image data by deconvolution and integrating the estimated point spread function and the X-ray transmission image data, and the X-ray transmission image data corrected by the correction unit.
  • an imaging unit configured to form an X-ray absorption coefficient distribution image of the subject.
  • an X-ray source that generates X-rays from an X-ray focal point and an X-ray detector in which X-ray detection elements for detecting the X-rays are two-dimensionally arranged sandwich a subject.
  • the X-ray transmission image data imaged by the imaging unit that rotates around the subject while capturing the X-ray transmission image data of the subject from a plurality of projection directions, while maintaining the relationship of being opposed to each other
  • an X-ray distribution estimator Based on the distribution estimated by the X-ray distribution estimator, an X-ray distribution estimator that performs a Monte Carlo simulation that simulates an interaction, and estimates the distribution of scattered X-rays derived from the subject.
  • An X-ray absorption coefficient distribution image of the subject using a correction unit that removes the scattered X-ray component from the transmission image data and the X-ray transmission image data from which the scattered X-ray component has been removed by the correction unit An X-ray CT apparatus.
  • an X-ray CT apparatus and a scattered X-ray correction method that can prevent deterioration in image quality due to scattered X-rays by estimating and correcting the scattered X-ray dose with high accuracy. it can.
  • FIG. 1 is an explanatory diagram showing a configuration of an X-ray CT apparatus 100 according to the first embodiment.
  • the X-axis direction in FIG. 1 is the channel direction
  • the Y-axis direction is the X-ray focal direction
  • the Z-axis direction is the slice direction.
  • the X-ray CT apparatus 100 is illustrated as viewed from the body axis direction (slice direction: Z-axis direction) of the subject 3.
  • An opening 2 through which the subject 3 can enter is provided at the center of the gantry (not shown) of the X-ray CT apparatus 100.
  • the scanner device of the X-ray CT apparatus 100 includes an X-ray tube 1 that is an X-ray source and an X-ray detector 4, and the gantry can be rotated about the center of the opening 2 as a rotation center axis. It is supported by. With such a configuration, the subject 3 in the opening 2 can be rotated and photographed.
  • the X-ray tube 1 as an X-ray source generates X-rays from an X-ray focal point 9 having a finite size in the X-ray tube 1.
  • An X-ray detector 4 is disposed at a position facing the X-ray tube 1 across the subject 3.
  • the X-ray detector 4 is divided into a plurality of detector modules 8, and each detector module 8 is arranged in an arc shape or a flat panel shape around the X-ray focal point 9.
  • FIG. 2 is an explanatory diagram showing the structure of the detector module 8.
  • the X-ray detection elements 6 are two-dimensionally arranged with respect to the channel direction (scanner rotation direction) and the slice direction (body axis direction) via the separator 7. This is for obtaining a wide range of X-ray transmission image data (projection data) of the subject 3 by one X-ray irradiation.
  • the X-ray detection element 6 is composed of, for example, a combination of a scintillator and a photodiode, or a semiconductor that converts radiation into an electrical signal, and measures the X-ray incident intensity to the X-ray detection element 6. is there.
  • a scattered radiation prevention grid 5 is disposed on the X-ray tube 1 side of the detector module 8 in order to remove scattered X-rays generated in the subject 3 and the like.
  • the imaging control of the X-ray CT apparatus 100 is performed by the control device 103 through the recording device 101 and the arithmetic device 102 based on the scanning conditions set by the user through the input device 104.
  • a large number of projection data obtained by rotational imaging is recorded in the recording device 101, and an image processing calculation is executed by the calculation device 102.
  • the projection data after the image processing calculation is displayed on the output device 105 as information such as a tomographic image of the subject 3.
  • the projection data (X-ray transmission image data) photographed by the X-ray CT apparatus 100 includes scattered X-rays resulting from the subject 3. For this reason, the X-ray CT apparatus 100 performs a point spread function (PSF) that takes into account the distribution of scattered X-rays caused by the subject 3 that changes with each imaging, in addition to image correction processing in general image processing. An image is formed after estimation and correction.
  • PSF point spread function
  • FIG. 3 is an explanatory diagram showing an outline of image forming processing by the X-ray CT apparatus.
  • the X-ray CT apparatus 100 performs a correction / reconstruction process (F2) on the projection data (Raw Data) obtained by imaging (F1), and generates a reconstructed image (Image).
  • the correction / reconstruction process in F2 is a general image correction process such as sensitivity correction.
  • the computing device 102 of the X-ray CT apparatus 100 virtually subjects the subject on the computing device 102 (see FIG. 1) based on the reconstructed image (Image) generated by the correction / reconstruction processing (F2). 3 is simulated (F3). Then, a Monte Carlo simulation capable of calculating a detailed physical interaction of X-rays is performed on the virtual subject 3A simulated in F3 (F4). This makes it possible to accurately estimate a point spread function (PSF) reflecting the structure of the subject 3 that changes with each shooting. Note that there is a possibility that a false image due to scattered X-rays exists in the reconstructed image (Image) before correction of scattered X-rays. However, if the influence of the false image is not excessively serious, it is possible to read the internal structure information of the subject 3 from the reconstructed image.
  • PSF point spread function
  • the arithmetic unit 102 corrects the projection data (Raw Data) obtained by photographing (F1) based on the point spread function (PSF) estimated in F4 (F5). Specifically, the scattered X-ray component is removed from the projection data. Then, the arithmetic unit 102 corrects and reconstructs the projection data after the scattered X-ray component removal (F6). Thereby, in the arithmetic unit 102, a high-quality reconstructed image (New Image) in which the influence of scattered X-rays is further reduced as compared with the reconstructed image (Image) generated by the correction / reconstruction processing (F2). can get.
  • the correction / reconstruction process in F6 is a general image correction process such as sensitivity correction as in F2.
  • imaging (F1), correction / reconstruction processing (F2), and correction / reconstruction processing (F6) are imaging and image correction processing of a general X-ray CT apparatus 100, and will not be described.
  • the arithmetic unit 102 of the X-ray CT apparatus 100 virtually simulates the subject 3 on the computer based on the reconstructed image (Image) and generates a virtual subject 3A.
  • the reconstructed image (Image) represents the difference in the X-ray absorption coefficient inside the subject 3 as a CT value.
  • CT values are standardized so that water is 0 HU and air is -1000 HU.
  • what is required in the Monte Carlo simulation is the elemental composition and density of the substance, but the composition and density cannot be determined from the CT value alone.
  • a constituent material is defined as a predefined material or a mixture thereof for each pixel of the reconstructed image.
  • FIG. 4 is an explanatory diagram showing an example of a method for determining a substance that simulates the subject 3.
  • air, water, and Teflon registered trademark, hereinafter the same
  • these mixtures are configured as follows.
  • the horizontal axis represents the CT value of the pixel
  • the vertical axis represents the mixing ratio of the representative substance.
  • the ratio is a volume ratio.
  • a and B are expressed as in the following (formula 1) and (formula 2).
  • A ⁇ CT value [HU] / 1000 ⁇ 100 (Formula 1)
  • B CT value [HU] / 1000 ⁇ 100 (Expression 2)
  • composition of the subject 3 can be predicted to some extent, a substance that is present in a large amount in the subject 3 is selected as a representative substance, and a substance that is hardly present in the subject 3 is not selected as a representative substance, thereby improving the precision and accuracy. It is possible. However, in order to be able to interpolate an arbitrary CT value, it is desirable to widen the CT value range of the representative substance.
  • a region where the CT value is in a certain range may be replaced with a single substance.
  • the following settings are conceivable. Air (when CT value ⁇ -500HU) Water (when -500HU ⁇ CT value ⁇ 500HU) Teflon (500HU ⁇ CT value)
  • the reconstructed image (Image) is expressed as a matrix of, for example, 512 ⁇ 512 pixels.
  • the object 3 is represented by 512 ⁇ 512 voxels per slice on the arithmetic unit 102, and the constituent substance information determined by the above method is added to each voxel. That is, the virtual subject 3A is configured as a set of rectangular parallelepiped regions including surrounding air.
  • the area captured by the X-ray CT apparatus 100 is limited to a part of the subject 3 except when scanning the entire subject 3.
  • information about the subject 3 outside the reconstruction area cannot be obtained, but scattered X-rays may interact outside the reconstruction area and recurse into the reconstruction area. Therefore, it is necessary to simulate the structure of the subject 3 outside the reconstruction area.
  • a method can be used in which the structure indicated by the surrounding pixels of the reconstruction area is uniformly present outside the field of view.
  • the computing device 102 of the X-ray CT apparatus 100 obtains projection data when an X-ray projection similar to that for actual imaging is performed on the subject 3 reproduced on the computing device 102 by Monte Carlo simulation.
  • the device structure such as the X-ray source (X-ray tube 1 and X-ray focal point 9) and the X-ray detector 4 are also placed on the arithmetic device 102 according to the actual structure. Simulate it.
  • the X-ray CT apparatus 100 irradiates a simulated subject 3 with a pencil-like X-ray beam by Monte Carlo simulation, and generates an X-ray intensity distribution (point spread function p) on the surface of the X-ray detector 4 with respect to the X-ray incident direction.
  • the measured projection data g (ch, sl), ideal projection data t (ch, sl) without scattered X-rays, and point spread function p (ch, sl) are used, and these Fourier transforms are respectively G (CH, SL). , T (CH, SL), and P (CH, SL), it is known that the following relationships (Equation 3) and (Equation 4) hold.
  • CH and SL represent frequency components in the ch direction and sl direction, respectively.
  • the point spread function changes relatively slowly, and several representative projection angles / detection element 6 positions (ch ′, sl ′, ⁇ ′) are used for speeding up the calculation. Only find the point spread function.
  • the X-ray detector 4 is divided into ten in the channel direction, and a point spread function is calculated for the detection element 6 position (ch ′, sl ′) in the center of each region. From the obtained point spread function and (Equation 5), ideal projection data t in each region is obtained (F5).
  • the point spread function of the intermediate projection angle / position can be estimated by interpolation from the point spread function p of the representative projection angle / position.
  • the point spread function p is obtained by, for example, fitting with C ⁇ EXP ( ⁇ D ⁇
  • EXP represents an exponential function.
  • smoothing processing may be performed on the distribution of the point spread function or the distribution of the fitting parameters in order to remove high-frequency components due to insufficient statistics in the Monte Carlo simulation.
  • the smoothing process is performed, for example, by taking a moving average of data in the channel direction, slice direction, and projection angle direction.
  • the X-ray CT apparatus 100 obtains ideal projection t data in all projection angle directions, and performs image formation again. Thereby, a high-quality image without the influence of scattered X-rays can be obtained.
  • the X-ray CT apparatus 100 estimates the point spread function of the scattering derived from the subject 3 based on the reconstructed image of the subject 3, the structure corresponding to the subject 3 is used. Can be reproduced. Thereby, the X-ray CT apparatus 100 can estimate the point spread function with high accuracy, and can reduce deterioration in image quality caused by scattered X-rays. Further, according to the X-ray CT apparatus 100 according to the first embodiment, since the Monte Carlo simulation that can strictly simulate the physical interaction is used, the point spread function can be estimated with higher accuracy.
  • the X-ray CT apparatus 100 limits the projection angle only to the representative angle direction, performs the simulation only to the position representing the X-ray detection element 6, and obtains the point spread function, By calculating point spread functions at other projection angles and other X-ray detection element 6 positions by interpolation, the calculation time required for correcting scattered X-rays can be greatly reduced. Furthermore, since the X-ray CT apparatus 100 according to the first embodiment smoothes the distribution of the point spread function or the fitting parameter distribution obtained by the calculation, it can remove high-frequency noise caused by statistical fluctuations. Thus, the calculation time can be shortened and overcorrection and erroneous correction can be suppressed. That is, according to the X-ray CT apparatus 100 according to the second embodiment, a high-quality CT image that is not affected by scattered X-rays can be obtained within a realistic calculation time.
  • the scattered X-ray correction is performed by estimating the point spread function at the representative point.
  • the scattered X-ray correction is performed by estimating the scattered X-ray distribution by Monte Carlo simulation.
  • the configuration of the X-ray CT apparatus 100 (see FIGS. 1 and 2) is the same as that of the first embodiment, and a detailed description thereof will be omitted.
  • the projection data (X-ray transmission image data) photographed by the X-ray CT apparatus 100 includes scattered X-rays resulting from the subject 3. For this reason, the X-ray CT apparatus 100 estimates the distribution of scattered X-rays caused by the subject 3 that changes at each imaging as well as image correction processing in general image processing, and forms an image after correcting the distribution.
  • FIG. 5 is an explanatory diagram showing an outline of image forming processing by the X-ray CT apparatus according to the second embodiment.
  • the imaging (F1), the correction / reconstruction process (F2), and the simulation of the virtual subject 3 (F3) are the same as the image forming process by the X-ray CT apparatus 100 according to the first embodiment, and the description is omitted. .
  • Monte Carlo simulation capable of calculating detailed physical interaction of X-rays is performed on the virtual subject 3A simulated in F3 (F7).
  • F7 F3
  • the influence of the false image is not excessively serious, it is possible to read the internal structure information of the subject 3 from the reconstructed image.
  • the computing device 102 corrects the projection data (Raw Data) obtained by imaging (F1) based on the scattered X-ray distribution estimated in F7 (F8). Specifically, the scattered X-ray component is removed from the projection data. Then, the arithmetic unit 102 again corrects and reconstructs the projection data after the scattered X-ray component removal (F9). Thereby, in the arithmetic unit 102, a high-quality reconstructed image (New Image) in which the influence of scattered X-rays is further reduced as compared with the reconstructed image (Image) generated by the correction / reconstruction processing (F2). can get.
  • the correction / reconstruction processing in F9 is a general image correction processing such as sensitivity correction as in F2.
  • imaging (F1), correction / reconstruction processing (F2), and correction / reconstruction processing (F9) are imaging and image correction processing of a general X-ray CT apparatus 100, and will not be described.
  • the simulation (F3) of the subject 3 is the same as that in the first embodiment, and a description thereof will be omitted.
  • the computing device 102 of the X-ray CT apparatus 100 obtains projection data when an X-ray projection similar to that for actual imaging is performed on the subject 3 reproduced on the computing device 102 by Monte Carlo simulation.
  • the device structure such as the X-ray source (X-ray tube 1 and X-ray focal point 9) and the X-ray detector 4 are also placed on the arithmetic device 102 according to the actual structure. Simulate it.
  • X-rays are emitted from the X-ray focal point 9 as a large number of photons.
  • X-ray energy is distributed to each photon based on an energy spectrum measured in advance under each scanning condition.
  • the interaction received by each photon is described by probability, and the presence or absence of the interaction is determined using a pseudo-random number.
  • Typical physical processes involving X-rays include, for example, Compton scattering, Rayleigh scattering, the photoelectric effect, and characteristic X-ray emission processes.
  • the X-ray source (X-ray tube 1 and X-ray focal point 9) and the X-ray detector 4 rotate around the subject 3 from a plurality of angular directions in the virtual space, as in actual imaging. Take a picture and obtain projection data by calculation.
  • the Monte Carlo simulation since it is possible to know whether or not there is an interaction between the X-ray and the subject 3, the incident intensity of the direct X-ray and the incident intensity of the scattered X-ray to the X-ray detector 4 are separately set Can be requested.
  • the Monte Carlo simulation on the virtual subject 3A, it is possible to accurately estimate the scattered X-ray distribution (Scat. Dist.) In consideration of the subject 3.
  • the representative angle direction may be selected at random or equiangular intervals, but may be selected as follows, for example, in order to reproduce the scattered X-ray distribution with high accuracy within a limited calculation time.
  • FIG. 6 is a flowchart illustrating an example of a procedure of Monte Carlo simulation.
  • the computing device 102 performs the above-described Monte Carlo simulation for a predetermined initial projection angle (for example, 0 degrees, 90 degrees, 180 degrees, and 270 degrees), and calculates the scattered X-ray distribution.
  • a predetermined initial projection angle for example, 0 degrees, 90 degrees, 180 degrees, and 270 degrees
  • step S ⁇ b> 2 the arithmetic unit 102 calculates the magnitude of the change in the scattered X-ray dose according to the following (Equation 7).
  • data scat (ch, sl, ⁇ ) represents the scattered X-ray distribution obtained per detector rotation, and the position (ch) of the channel direction X-ray detection element 6 and the X in the slice direction It is expressed as a function of the position (sl) of the line detection element 6 and the projection angle ( ⁇ ).
  • ch, sl, and ⁇ are discrete variables.
  • step S ⁇ b> 3 the computing device 102 obtains a section in which the magnitude of the change in the scattered X-ray dose is maximized.
  • step S4 the arithmetic unit 102 sets the projection angle in the section obtained in step S3, performs the above-described Monte Carlo simulation, and calculates the scattered X-ray distribution. For example, when the change in (Equation 7) shows the maximum value between the projection angles ⁇ 1 and ⁇ 2 (S3), the arithmetic unit 102 takes, for example, ( ⁇ 1 + ⁇ 2) / 2 as the next projection angle, and again A Monte Carlo simulation is performed to calculate the scattered X-ray distribution (S4).
  • step S5 the arithmetic unit 102 returns to step S2 until the required statistical accuracy is satisfied or the allowable calculation time is exceeded (S5, No), and the subsequent processing is repeated.
  • the process proceeds to step S6.
  • the time required for the Monte Carlo simulation varies depending on the required accuracy.
  • the user can arbitrarily set an X-ray statistic or a calculation time used for the Monte Carlo simulation through the input device 104. Further, the user can interrupt the Monte Carlo simulation at an arbitrary timing.
  • the X-ray CT apparatus 100 performs scattered X-ray correction using data that has been calculated by the time of the interruption.
  • step S6 the arithmetic unit 102 performs an interpolation process and a smoothing process on the scattered X-ray distribution obtained in steps S1 and S4, and ends the process according to this flowchart.
  • the interpolation processing in step S6 is performed to interpolate intermediate angle data because the scattered X-ray distribution obtained by the above processing has a smaller number of data than actually measured in the projection angle direction.
  • data interpolation methods for example, a spline interpolation method, a Lagrange interpolation method, and the like are known.
  • the smoothing process in step S6 is performed in order to remove a high-frequency component caused by a shortage of statistics in the Monte Carlo simulation.
  • the smoothing process is performed, for example, by taking a moving average of data in the channel direction, slice direction, and projection angle direction.
  • the arithmetic unit 102 performs the scattered X-ray correction by the calculation shown in the following (Formula 8).
  • DATA is projection data obtained by imaging
  • DATA ref is reference data of projection obtained by imaging (for example, average data of the output of the X-ray detection element 6 located at the end)
  • data scat is reference data obtained scattered X-ray distribution data calculation
  • the display of the variables have been omitted .
  • ⁇ is a constant representing the correction intensity. Note that the measurement data is corrected so that the output (offset value) when X-rays are not irradiated is zero.
  • the arithmetic unit 102 also performs scattered X-ray correction of air data (output value when the subject 3 is air) by the calculation shown in the following (formula 9).
  • AIR is projection data when the subject 3 is photographed as air.
  • the meaning of each subscript is the same as in the above (Formula 8).
  • the arithmetic unit 102 uses the values calculated in the above (Expression 8) and (Expression 9), the arithmetic unit 102 performs sensitivity correction (referred to as air correction) of the X-ray detection element 6 by the following (Expression 10).
  • Air correction sensitivity correction
  • NewData / NewAir Form 10
  • the computing device 102 performs image reconstruction again based on the data after the scattered X-ray correction and air correction shown in the above (Equation 10). Thereby, the X-ray CT apparatus 100 can provide a user with a high-quality image that is not affected by scattered X-rays.
  • the scattered X-ray distribution (Scat. Dist.) Is estimated based on the reconstructed image (Image) of the subject 3.
  • the change of the structure can be reproduced.
  • the X-ray CT apparatus 100 can estimate the scattered X-ray distribution with high accuracy, and can reduce deterioration in image quality caused by the scattered X-rays.
  • the Monte Carlo simulation that can strictly simulate the physical interaction is used, the scattered X-ray distribution can be estimated with higher accuracy.
  • the X-ray CT apparatus 100 obtains a scattered X-ray distribution by performing a simulation by limiting the projection angle only to the representative angle direction, and interpolates the scattered X-ray distribution in other angular directions. Therefore, the calculation time required for estimating the scattered X-ray distribution can be greatly shortened. Furthermore, since the X-ray CT apparatus 100 according to the second embodiment smoothes the scattered X-ray distribution obtained by calculation, high-frequency noise due to statistical fluctuations can be removed, and calculation time can be reduced. It is possible to suppress overcorrection and erroneous correction. That is, according to the X-ray CT apparatus 100 according to the second embodiment, a high-quality CT image that is not affected by scattered X-rays can be obtained within a realistic calculation time.
  • the X-ray CT apparatus 100 employs elliptical columnar water as a phantom, for example.
  • the variables are only the major axis a and the minor axis b.
  • the X-ray CT apparatus 100 obtains the dependence of the scattered X-ray distribution or the point spread function on the major axis a for a plurality of elliptical cylinders having various flatness ratios f. Then, the X-ray CT apparatus 100 obtains the scattered X-ray distribution and the point spread function distribution for the arbitrary variables a and f by interpolation and function fitting, and creates a database.
  • the X-ray CT apparatus 100 approximates the imaging region shape of the object 3 to an ellipse based on the reconstructed image obtained by imaging, and sets the major axis a and the flatness ratio f. taking measurement.
  • the X-ray CT apparatus 100 selects a scattered X-ray distribution or a distribution of a point spread function that approximates the measured major axis a and flatness f from the above-described database, and by the method described in the first embodiment or 2,
  • the scattered X-ray correction is performed on the measurement data.
  • the X-ray CT apparatus 100 performs image reconstruction using the measurement data after the scattered X-ray correction. Thereby, a high-quality image without the influence of scattered X-rays can be obtained.
  • the X-ray CT apparatus 100 according to the third embodiment has the following effects in addition to the same effects as those of the X-ray CT apparatus 100 according to the first embodiment or the second embodiment. can get. That is, according to the X-ray CT apparatus 100 according to the third embodiment, since the Monte Carlo simulation that requires a long time for the calculation is executed in advance before photographing the subject 3, the time required for correcting the scattered X-ray is greatly reduced. can do.
  • X-ray tube (X-ray source) 2 Opening 3 Subject 3A Virtual subject (simulated subject, X-ray absorption coefficient distribution) 4 X-ray detector 5 Scattered ray prevention grid 6 Separator 7 X-ray detection element 8 Detector module 9 X-ray focus 100 X-ray CT apparatus 101 Recording apparatus 102 Arithmetic apparatus 103 Control apparatus 104 Input apparatus 105 Output apparatus F1 Imaging (imaging part) ) F2 Correction / reconstruction processing F3 Simulation of the subject F4 Monte Carlo simulation for virtual subjects (point spread function estimation unit) F5 Projection data correction (correction unit) F6 correction / reconstruction processing (imaging part) F7 Monte Carlo simulation for virtual subjects (X-ray distribution estimation unit) F8 Projection data correction (correction unit) F9 Correction / reconstruction processing (imaging part)

Abstract

　Ｘ線ＣＴ画像において、散乱Ｘ線による画像低下を防止する。　Ｘ線ＣＴ装置（１００）は、Ｘ線検出器（４）などによって被写体（３）のＸ線透過像データを撮影し（Ｆ１）、撮影されたＸ線透過像データに基づいて、被写体３内部のＸ線吸収係数分布を推定する（Ｆ３）。Ｘ線ＣＴ装置（１００）は、推定されたＸ線吸収係数分布を有する模擬被写体（３）に対して、モンテカルロシミュレーションを実施して、被写体に由来する点拡がり関数または散乱Ｘ線分布を推定する（Ｆ４，Ｆ７）。そして、Ｘ線ＣＴ装置（１００）は、推定された点拡がり関数または散乱Ｘ線分布に基づいて、Ｘ線透過像データを補正して（Ｆ５，Ｆ８）、被写体（３）のＸ線吸収係数分布画像を形成する（Ｆ６，Ｆ９）。

Description

Ｘ線ＣＴ装置および散乱Ｘ線補正方法
　本発明は、Ｘ線ＣＴ装置、およびＸ線ＣＴ装置などで撮影したデータに対する散乱Ｘ線補正方法に関する。
　Ｘ線ＣＴ（Ｃｏｍｐｕｔｅｄ　Ｔｏｍｏｇｒａｐｈｙ）装置とは、被写体内部のＸ線減衰率（Ｘ線吸収係数）の違いを、データ処理系を用いて画像として再構成するための装置である。Ｘ線ＣＴ装置は、被写体にＸ線を照射するＸ線源と、被写体を透過したＸ線を検出するＸ線検出器とを備える。Ｘ線源とＸ線検出器とは、被写体をはさんで対向配置されており、被写体をはさんで対向配置された関係を保ちながら被写体の周りを回転し、複数の投影方向から被写体のＸ線透過像データを撮影する。Ｘ線ＣＴ装置のＸ線源には、通常、高電圧で加速された電子を陽極に照射してＸ線を発生させるＸ線管球が用いられる。また、Ｘ線検出器は、一度に広範囲を高速に撮影するためにＸ線検出素子を２次元状に並べた構造を有している。
　Ｘ線ＣＴ装置の撮影によって得られる投影データには、被写体を非散乱で透過してきたＸ線（直達Ｘ線）強度以外に、被写体などで散乱されたＸ線（散乱Ｘ線）の入射強度の情報も含まれている。Ｘ線ＣＴ装置では、散乱Ｘ線を除去するために、Ｘ線検出器のＸ線源側に、被写体で発生した散乱Ｘ線を除去するための散乱線防止グリッドが配置されている。しかし、この方法によってもすべての散乱Ｘ線を除去することはできない。このため、ソフトウェアによる散乱Ｘ線補正も併せておこなわれている（例えば、下記特許文献１～３参照）。
特許第４２１８９０８号公報 特開２００９－８２６１５号公報 特許第３５６６７６２号公報
　近年、Ｘ線ＣＴ装置は、Ｘ線検出器の多列化にともなって、１回の撮影で被写体を広範囲に撮影することが可能になっている。一方で、Ｘ線検出器の多列化によって、被写体へのＸ線照射幅が広がり、それにともなって散乱Ｘ線量が増加することで、再構成画像上に偽像（アーチファクト）が生じて、画質が低下してしまうという問題点がある。とりわけ、骨などのＸ線の高吸収体が存在する領域を撮影する場合、直達Ｘ線由来の検出シグナルに対して散乱Ｘ線由来の検出シグナルの割合が相対的に増加し、被写体のＸ線吸収係数を過小評価してしまう場合がある。
　図７は、再構成画像に生じるアーチファクトを模式的に示す説明図である。図７（ａ）には、低吸収体Ｌの内部に２つの高吸収体ロッドＨが存在する人体模擬構造物（以下、ファントム）を示している。このファントムを撮影した場合、図７（ｂ）に示す再構成画像のように、２つの高吸収体ロッドＨの間に、ＣＴ値が実際よりも小さく計測されたダークバンドアーチファクトＤが観測される。なお、図７において、ＣＴ値の大きさをハッチングの間隔（密度）で図示し、ＣＴ値が小さいほどハッチングの間隔が短く（密に）なることで図示している。
　そこで、本発明は、高精度に散乱Ｘ線量を推定して補正することによって、散乱Ｘ線に起因する画質低下を防止することができるＸ線ＣＴ装置および散乱Ｘ線補正方法を提供することを課題とする。
　このような課題を解決するために、本発明は、Ｘ線焦点からＸ線を発生するＸ線源と、前記Ｘ線を検出するためのＸ線検出素子が２次元配列されているＸ線検出器とが、被写体をはさんで対向配置された関係を保ちながら前記被写体の周りを回転し、複数の投影方向から前記被写体のＸ線透過像データを撮影する撮影部と、前記撮影部によって撮影された前記Ｘ線透過像データに基づいて、前記被写体内部のＸ線吸収係数分布を推定する内部分布推定部と、前記内部分布推定部によって推定された前記Ｘ線吸収係数分布を有する模擬被写体に対して、前記Ｘ線の物理相互作用を模擬するモンテカルロシミュレーションを実施して、前記被写体に由来する散乱の点拡がり関数を推定する点拡がり関数推定部と、前記点拡がり関数推定部によって推定された前記点拡がり関数と、前記Ｘ線透過像データとを逆畳み込み積分して、前記Ｘ線透過像データを補正する補正部と、前記補正部によって補正された前記Ｘ線透過像データを用いて、前記被写体のＸ線吸収係数分布画像を形成する画像化部と、を備えることを特徴とするＸ線ＣＴ装置である。
　また、本発明は、Ｘ線焦点からＸ線を発生するＸ線源と、前記Ｘ線を検出するためのＸ線検出素子が２次元配列されているＸ線検出器とが、被写体をはさんで対向配置された関係を保ちながら前記被写体の周りを回転し、複数の投影方向から前記被写体のＸ線透過像データを撮影する撮影部と、前記撮影部によって撮影された前記Ｘ線透過像データに基づいて、前記被写体内部のＸ線吸収係数分布を推定する内部分布推定部と、前記内部分布推定部によって推定された前記Ｘ線吸収係数分布を有する模擬被写体に対して、前記Ｘ線の物理相互作用を模擬するモンテカルロシミュレーションを実施して、前記被写体に由来する散乱Ｘ線の分布を推定するＸ線分布推定部と、前記Ｘ線分布推定部によって推定された前記分布に基づいて、前記Ｘ線透過像データから前記散乱Ｘ線の成分を除去する補正部と、前記補正部によって前記散乱Ｘ線の成分が除去された前記Ｘ線透過像データを用いて、前記被写体のＸ線吸収係数分布画像を形成する画像化部と、を備えることを特徴とするＸ線ＣＴ装置である。
　本発明によれば、高精度に散乱Ｘ線量を推定して補正することによって、散乱Ｘ線に起因する画質低下を防止することができるＸ線ＣＴ装置および散乱Ｘ線補正方法を提供することができる。
第１実施形態に係るＸ線ＣＴ装置の構成を示す説明図である。 検出器モジュールの構造を示す説明図である。 第１実施形態に係るＸ線ＣＴ装置による画像形成処理の概要を示す説明図である。 被写体を模擬する物質の決定方法の一例を示す説明図である。 第２実施形態に係るＸ線ＣＴ装置による画像形成処理の概要を示す説明図である。 モンテカルロシミュレーションの手順の一例を示すフローチャートである。 再構成画像に生じるアーチファクトを模式的に示す説明図である。
　以下、本発明を実施するための形態（以下「実施形態」という）について、適宜図面を参照しながら詳細に説明する。なお、各図において、共通する部分には同一の符号を付し重複した説明を省略する。
≪第１実施形態≫
＜Ｘ線ＣＴ装置１００の構成＞
　図１は、第１実施形態に係るＸ線ＣＴ装置１００の構成を示す説明図である。以下の説明において、図１中のＸ軸方向をチャネル方向、Ｙ軸方向をＸ線焦点方向、Ｚ軸方向をスライス方向とする。図１では、Ｘ線ＣＴ装置１００を、被写体３の体軸方向（スライス方向：Ｚ軸方向）から見たものとして図示している。
　Ｘ線ＣＴ装置１００の図示しないガントリの中央部には、被写体３が進入できる開口部２が設けられている。また、Ｘ線ＣＴ装置１００のスキャナ装置には、Ｘ線源であるＸ線管球１と、Ｘ線検出器４とが備えられ、開口部２の中心を回転中心軸として、ガントリに回転可能に支持されている。このような構成により、開口部２内の被写体３を回転撮影することが可能となる。
　Ｘ線源であるＸ線管球１は、Ｘ線管球１内にある有限の大きさを持つＸ線焦点９からＸ線を発生する。被写体３をはさんでＸ線管球１と対向する位置には、Ｘ線検出器４が配置される。Ｘ線検出器４は、複数の検出器モジュール８に分割されており、それぞれの検出器モジュール８は、Ｘ線焦点９を中心として円弧状またはフラットパネル状に配置されるものとする。
　図２は、検出器モジュール８の構造を示す説明図である。検出器モジュール８は、Ｘ線検出素子６が、セパレータ７を介して、チャネル方向（スキャナ回転方向）およびスライス方向（体軸方向）に対して２次元配列される。これは、１回のＸ線照射で、被写体３の広範囲なＸ線透過像データ（投影データ）を得るためである。ここで、Ｘ線検出素子６は、例えば、シンチレータとフォトダイオードを組み合わせたものや、放射線を電気信号に変換する半導体で構成され、Ｘ線検出素子６へのＸ線入射強度を計測するものである。また、検出器モジュール８のＸ線管球１側には、被写体３などで発生した散乱Ｘ線を除去するために、散乱線防止グリッド５が配置される。
　図１の説明に戻り、Ｘ線ＣＴ装置１００の撮影制御は、ユーザが入力装置１０４を通して設定したスキャン条件に基づいて、記録装置１０１、演算装置１０２を通して、制御装置１０３によって行われる。回転撮影によって得られた多数の投影データは、記録装置１０１に記録され、演算装置１０２で画像処理演算が実行される。そして、画像処理演算後の投影データは、被写体３の断層画像などの情報として出力装置１０５に表示される。
＜Ｘ線ＣＴ装置１００による画像形成処理＞
　つづいて、Ｘ線ＣＴ装置１００による画像形成処理について説明する。Ｘ線ＣＴ装置１００によって撮影された投影データ（Ｘ線透過像データ）には、被写体３に起因する散乱Ｘ線が含まれている。このため、Ｘ線ＣＴ装置１００は、一般的な画像処理における画像補正処理の他、撮影ごとに変わる被写体３に起因する散乱Ｘ線の分布を考慮した点拡がり関数（ＰＳＦ：Point Spread Function）を推定し、補正した上で画像を形成する。
　図３は、Ｘ線ＣＴ装置による画像形成処理の概要を示す説明図である。
　Ｘ線ＣＴ装置１００は、撮影（Ｆ１）によって得られた投影データ（Ｒａｗ　Ｄａｔａ）に対して、補正・再構成処理（Ｆ２）をおこない、再構成画像（Ｉｍａｇｅ）を生成する。なお、Ｆ２における補正・再構成処理とは、例えば、感度補正などの一般的な画像補正処理である。
　つづいて、Ｘ線ＣＴ装置１００の演算装置１０２は、補正・再構成処理（Ｆ２）により生成された再構成画像（Ｉｍａｇｅ）を基に、演算装置１０２（図１参照）上で仮想的に被写体３を模擬する（Ｆ３）。そして、Ｆ３で模擬した仮想被写体３Ａに対して、Ｘ線の詳細な物理相互作用を計算可能なモンテカルロシミュレーションを実施する（Ｆ４）。これにより、撮影ごとに変わる被写体３の構造を反映した点拡がり関数（ＰＳＦ）を精度良く推定することが可能となる。なお、散乱Ｘ線補正前の再構成画像（Ｉｍａｇｅ）には、散乱Ｘ線に起因する偽像が存在する可能性がある。しかし、偽像の影響が過度に深刻でなければ、再構成画像から被写体３の内部構造情報を読み取ることは可能である。
　演算装置１０２は、Ｆ４において推定した点拡がり関数（ＰＳＦ）を基に、撮影（Ｆ１）によって得られた投影データ（Ｒａｗ　Ｄａｔａ）を補正する（Ｆ５）。具体的には、投影データから散乱Ｘ線成分を除去する。そして、演算装置１０２は、散乱Ｘ線成分除去後の投影データを、あらためて補正・再構成処理する（Ｆ６）。これにより、演算装置１０２において、補正・再構成処理（Ｆ２）により生成された再構成画像（Ｉｍａｇｅ）と比較して、散乱Ｘ線の影響をより低減した良質な再構成画像（Ｎｅｗ　Ｉｍａｇｅ）が得られる。なお、Ｆ６における補正・再構成処理とは、Ｆ２と同様、感度補正などの一般的な画像補正処理である。
＜画像形成処理の詳細＞
　つづいて、図３に示した画像形成処理の各工程の詳細について説明する。
　なお、撮影（Ｆ１）、補正・再構成処理（Ｆ２）および補正・再構成処理（Ｆ６）は、一般的なＸ線ＣＴ装置１００の撮影および画像補正処理であり、説明を省略する。
（被写体３の模擬：Ｆ３）
　図３のＦ３において、Ｘ線ＣＴ装置１００の演算装置１０２は、再構成画像（Ｉｍａｇｅ）を基に、計算機上で仮想的に被写体３を模擬し仮想被写体３Ａを生成する。ここで、再構成画像（Ｉｍａｇｅ）は、被写体３内部のＸ線吸収係数の違いをＣＴ値として表現したものである。ＣＴ値は、水が０ＨＵ、空気が－１０００ＨＵになるように規格化されている。一方、モンテカルロシミュレーションで必要になるのは、物質の元素組成と密度であるが、ＣＴ値だけからは組成と密度は決定できない。
　そこで、Ｘ線ＣＴ装置１００では、組成、密度およびＣＴ値が既知であるいくつかの物質を被写体３を構成する物質としてあらかじめ定義しておき、それらの中間のＣＴ値をとる物質については、定義済み物質の混合物とみなすことによって、被写体３の内部構造（Ｘ線吸収係数分布）を模擬する。具体的には、例えば、再構成画像のピクセルごとに、定義済み物質またはそれらの混合物として構成物質を定義する。
　図４は、被写体３を模擬する物質の決定方法の一例を示す説明図である。図４（ａ）に示すように、例えば、空気、水、テフロン（登録商標、以下同じ）を、それぞれＣＴ値－１０００ＨＵ、０ＨＵ、１０００ＨＵをとる代表物質とみなし、任意のＣＴ値をとる物質を、これらの混合物として以下のように構成する。なお、図４において、横軸は画素のＣＴ値、縦軸は代表物質の混合割合である。
　空気：水：テフロン＝１００％：０％：０％（ＣＴ値≦－１０００ＨＵの場合）
　空気：水：テフロン＝Ａ％：（１００－Ａ）％：０％（－１０００ＨＵ＜ＣＴ値≦０ＨＵの場合）
　空気：水：テフロン＝０％：（１００－Ｂ）％：Ｂ％（０ＨＵ＜ＣＴ値≦１０００ＨＵの場合）
　空気：水：テフロン＝０％：０％：１００％（１０００ＨＵ＜ＣＴ値の場合）
　ここで、比率は体積比である。またＡおよびＢは、下記（式１），（式２）のように示される。
　Ａ　＝　－ＣＴ値［ＨＵ］／１０００×１００　・・・（式１）
　Ｂ　＝　ＣＴ値［ＨＵ］／１０００×１００　　・・・（式２）
　なお、被写体３の組成がある程度予測可能なら、被写体３に多量に存在する物質を代表物質として選択し、被写体３にほとんど存在しない物質は代表物質として選択しないことによって、精密化、高精度化することが可能である。ただし、任意のＣＴ値を補間可能とするため、代表物質のＣＴ値の範囲は広くすることが望ましい。
　また、反対に計算負荷低減のため、ＣＴ値が一定範囲にある領域を単一の物質で置き換えてもよい。例えば、図４（ｂ）に示すように、以下のような設定が考えられる。
　空気（ＣＴ値≦－５００ＨＵの場合）
　水（－５００ＨＵ＜ＣＴ値≦５００ＨＵの場合）
　テフロン（５００ＨＵ＜ＣＴ値の場合）
　つぎに、被写体３の形状の模擬について説明する。再構成画像（Ｉｍａｇｅ）は、例えば５１２×５１２画素のマトリックスとして表現される。これを反映して、演算装置１０２上でも１スライスあたり５１２×５１２のボクセルで被写体３を表わし、各ボクセルに上記方法で決定した構成物質情報を付加する。つまり、仮想被写体３Ａは、周囲の空気も含めて直方体状領域の集合として構成される。
　なお、被写体３の全体をスキャンする場合を除いて、Ｘ線ＣＴ装置１００で撮影する領域は、被写体３の一部分に限られる。その場合、再構成領域外の被写体３の情報は得られないが、散乱Ｘ線は再構成領域外でも相互作用をおこなって、再構成領域内に再帰する可能性もある。そのため、再構成領域外の被写体３の構造も模擬する必要がある。この場合、例えば、再構成領域の周辺画素が示す構造が、視野外にも一様に存在するとみなす方法がとれる。
（仮想被写体３Ａに対するモンテカルロシミュレーションの実行：Ｆ４）
　つづいて、演算装置１０２上に模擬した仮想被写体３Ａに対するモンテカルロシミュレーション（図３のＦ４）について説明する。Ｘ線ＣＴ装置１００の演算装置１０２は、モンテカルロシミュレーションによって、演算装置１０２上に再現した被写体３に対して実際の撮影と同様のＸ線投影をおこなった場合の投影データを求める。このとき、被写体３以外の条件、例えば、Ｘ線源（Ｘ線管球１、Ｘ線焦点９）やＸ線検出器４などの装置構造についても、実際の構造に準じて演算装置１０２上に模擬しておく。
　Ｘ線ＣＴ装置１００は、モンテカルロシミュレーションによって、ペンシル状Ｘ線ビームを模擬の被写体３に照射し、Ｘ線入射方向に対するＸ線検出器４面上でのＸ線強度分布（点拡がり関数ｐ）を求める。計測された投影データｇ（ｃｈ，ｓｌ）、散乱Ｘ線のない理想投影データｔ（ｃｈ，ｓｌ）、点拡がり関数ｐ（ｃｈ，ｓｌ）とし、これらのフーリエ変換をそれぞれＧ（ＣＨ，ＳＬ）、Ｔ（ＣＨ，ＳＬ）、Ｐ（ＣＨ，ＳＬ）とすると、下記（式３）および（式４）の関係が成り立つことが知られている。ここでＣＨ、ＳＬはそれぞれｃｈ方向、ｓｌ方向の周波数成分を表す。
　ｇ＝ｔ＊ｐ　・・・（式３）
　Ｇ＝Ｔ・Ｐ　・・・（式４）
ここで、＊は畳み込み積分、・は積を表す。なお変数表示は省略した。
　逆フーリエ変換をＦで表わせば、求める理想投影データｔは下記（式５）で与えられる。
　ｔ＝Ｆ［Ｇ／Ｐ］＝ｇ＊Ｆ［１／Ｐ］　・・・（式５）
　この方法を逆畳み込み積分法と呼ぶ。
　散乱Ｘ線の拡がり方は、被写体３の部位や撮影角度方向θ'に依存して変化する。ペンシルビームの入射方向にある検出素子６の位置を（ｃｈ'、ｓｌ'）とすると、Ｘ線入射方向（撮影角度方向）依存性をあらわにした点拡がり関数ｐは、下記（式６）であらわせる。
　ｐ（ｃｈ，ｓｌ，ｃｈ'，ｓｌ'，θ'）　・・・（式６）
　Ｘ線ＣＴ装置１００では、点拡がり関数は比較的緩やかに変化するものと仮定し、計算高速化のため、いくつかの代表投影角度・検出素子６位置（ｃｈ'，ｓｌ'，θ'）についてのみ点拡がり関数を求める。例えば、Ｘ線検出器４をチャネル方向に１０分割し、各領域の中央にある検出素子６位置（ｃｈ'，ｓｌ'）について点拡がり関数を計算する。求めた点拡がり関数と（式５）から、各領域における理想投影データｔが得られる（Ｆ５）。
　また、代表投影角度・位置の点拡がり関数ｐから、中間投影角度・位置の点拡がり関数を補間により推定することもできる。この場合、点拡がり関数ｐは、例えば、代表位置からの距離（｜ｃｈ｜）の関数として、Ｃ・ＥＸＰ（－Ｄ・｜ｃｈ｜）でフィッティングし、係数ＣとＤとして求めて、この係数に対してデータ補間を実施する。ここでＥＸＰは指数関数を表わす。
　なお、モンテカルロシミュレーションの統計量不足にともなう高周波成分を取り除くため、点拡がり関数の分布あるいはフィッティングパラメータ分布に、平滑化処理を実施してもよい。平滑化処理は、例えば、チャネル方向やスライス方向、投影角度方向へ、データの移動平均をとることによっておこなう。このような方法により、Ｘ線ＣＴ装置１００は、全投影角度方向について理想的投影ｔデータを求め、あらためて画像形成を実施する。これにより、散乱Ｘ線の影響のない良質な画像が得られる。
　以上説明したように、第１実施形態に係るＸ線ＣＴ装置１００では、被写体３の再構成画像を基に、被写体３に由来する散乱の点拡がり関数を推定するので、被写体３に応じた構造の変化を再現することができる。これにより、Ｘ線ＣＴ装置１００は、点拡がり関数を精度高く推定することができ、散乱Ｘ線に起因する画質低下を低減することができる。また、第１実施形態に係るＸ線ＣＴ装置１００によれば、物理相互作用を厳密に模擬可能なモンテカルロシミュレーションを用いるので、点拡がり関数をより精度高く推定することができる。
　また、第１実施形態に係るＸ線ＣＴ装置１００は、投影角度を代表角度方向にのみに限定し、Ｘ線検出素子６を代表する位置に限定してシミュレーションをおこなって点拡がり関数を求め、他の投影角度、他のＸ線検出素子６位置の点拡がり関数を補間によって求めることで、散乱Ｘ線補正に要する計算時間を大幅に短縮することができる。さらに、第１実施形態に係るＸ線ＣＴ装置１００は、計算によって得られた点拡がり関数の分布あるいはそのフィッティングパラメータ分布を平滑化処理するので、統計変動に起因する高周波ノイズを除去することができ、計算時間の短縮および過補正や誤補正の抑制を図ることができる。すなわち、第２実施形態に係るＸ線ＣＴ装置１００によれば、現実的な計算時間内で、散乱Ｘ線の影響のない良質なＣＴ画像を得ることができる。
≪第２実施形態≫
　第１実施形態では、代表点における点拡がり関数を推定することによって、散乱Ｘ線補正をおこなった。第２実施形態では、モンテカルロシミュレーションにより散乱Ｘ線の分布を推定して、散乱Ｘ線補正をおこなった。なお、第２実施形態において、Ｘ線ＣＴ装置１００の構成（図１、図２参照）については、第１実施形態と同様であるので詳細な説明を省略する。
＜Ｘ線ＣＴ装置１００による画像形成処理＞
　第２実施形態に係るＸ線ＣＴ装置１００による画像形成処理について説明する。Ｘ線ＣＴ装置１００によって撮影された投影データ（Ｘ線透過像データ）には、被写体３に起因する散乱Ｘ線が含まれている。このため、Ｘ線ＣＴ装置１００は、一般的な画像処理における画像補正処理の他、撮影ごとに変わる被写体３に起因する散乱Ｘ線の分布を推定し、補正した上で画像を形成する。
　図５は、第２実施形態に係るＸ線ＣＴ装置による画像形成処理の概要を示す説明図である。
　撮影（Ｆ１）、補正・再構成処理（Ｆ２）および仮想的な被写体３の模擬（Ｆ３）は、第１実施形態に係るＸ線ＣＴ装置１００による画像形成処理と同様であり、説明を省略する。
　Ｆ３で模擬した仮想被写体３Ａに対して、Ｘ線の詳細な物理相互作用を計算可能なモンテカルロシミュレーションを実施する（Ｆ７）。これにより、撮影ごとに変わる被写体３の構造を反映した散乱Ｘ線分布（Ｓｃａｔ．Ｄｉｓｔ．）を精度良く推定することが可能となる。なお、散乱Ｘ線補正前の再構成画像（Ｉｍａｇｅ）には、散乱Ｘ線に起因する偽像が存在する可能性がある。しかし、偽像の影響が過度に深刻でなければ、再構成画像から被写体３の内部構造情報を読み取ることは可能である。
　演算装置１０２は、Ｆ７において推定した散乱Ｘ線分布を基に、撮影（Ｆ１）によって得られた投影データ（Ｒａｗ　Ｄａｔａ）を補正する（Ｆ８）。具体的には、投影データから散乱Ｘ線成分を除去する。そして、演算装置１０２は、散乱Ｘ線成分除去後の投影データを、あらためて補正・再構成処理する（Ｆ９）。これにより、演算装置１０２において、補正・再構成処理（Ｆ２）により生成された再構成画像（Ｉｍａｇｅ）と比較して、散乱Ｘ線の影響をより低減した良質な再構成画像（Ｎｅｗ　Ｉｍａｇｅ）が得られる。なお、Ｆ９における補正・再構成処理とは、Ｆ２と同様、感度補正などの一般的な画像補正処理である。
＜画像形成処理の詳細＞
　つづいて、図５に示した画像形成処理の各工程の詳細について説明する。
　なお、撮影（Ｆ１）、補正・再構成処理（Ｆ２）および補正・再構成処理（Ｆ９）は、一般的なＸ線ＣＴ装置１００の撮影および画像補正処理であり、説明を省略する。また、被写体３の模擬（Ｆ３）は第１実施形態と同様であり、説明を省略する。
（仮想被写体３Ａに対するモンテカルロシミュレーションの実行：Ｆ７）
　つづいて、演算装置１０２上に模擬した仮想被写体３Ａに対するモンテカルロシミュレーション（図５のＦ７）について説明する。Ｘ線ＣＴ装置１００の演算装置１０２は、モンテカルロシミュレーションによって、演算装置１０２上に再現した被写体３に対して実際の撮影と同様のＸ線投影をおこなった場合の投影データを求める。このとき、被写体３以外の条件、例えば、Ｘ線源（Ｘ線管球１、Ｘ線焦点９）やＸ線検出器４などの装置構造についても、実際の構造に準じて演算装置１０２上に模擬しておく。
　モンテカルロシミュレーションでは、Ｘ線焦点９からＸ線を多数の光子として放射させる。Ｘ線のエネルギは、あらかじめ計測した各スキャン条件におけるエネルギスペクトルに基づいて各光子に分配される。各光子が受ける相互作用は、確率によって記述され、相互作用の有無は擬似乱数を用いて決定する。Ｘ線の関与する代表的な物理過程としては、例えば、コンプトン散乱、レーリー散乱、光電効果、特性Ｘ線放射過程などがある。
　モンテカルロシミュレーションにおいても、実際の撮影と同様に、仮想空間上でＸ線源（Ｘ線管球１、Ｘ線焦点９）とＸ線検出器４が被写体３の周りを回りながら複数の角度方向から撮影をおこない、計算によって投影データを求める。このとき、モンテカルロシミュレーションでは、Ｘ線と被写体３との相互作用の有無を知ることができるので、Ｘ線検出器４への直達Ｘ線の入射強度と散乱Ｘ線の入射強度とを、それぞれ別個に求めることができる。
　このように、仮想被写体３Ａに対してモンテカルロシミュレーションを実行することにより、被写体３を考慮した散乱Ｘ線分布（Ｓｃａｔ．Ｄｉｓｔ．）を精度良く推定することが可能となる。
　つぎに、モンテカルロシミュレーションにおける計算高速化方法について説明する。実際のＣＴ撮影では、スキャン条件にもよるが、１回転あたり１０００枚分程度の投影データを用いて画像を再構成している。すなわち、実際の撮影では、１回転あたり１０００回程度、異なる角度方向から撮影をおこなっている。
　一方、Ｘ線ＣＴ装置１００の演算装置１０２がおこなうシミュレーションでは、実際の撮影回数よりも少ない数の代表角度方向の投影データのみをシミュレーションで求めている。これは、散乱Ｘ線の分布は、直達Ｘ線の分布に比べて比較的変化が緩やかであるため、いくつかの代表角度方向の投影データに含まれる散乱Ｘ線分布から、他の角度方向の散乱Ｘ線分布を推定することが可能なためである。このように、代表角度方向のみに投影シミュレーションをおこなうことによって、シミュレーションに要する計算時間を短縮することができる。
　代表角度方向は、ランダムもしくは等角度間隔に選んでもよいが、散乱Ｘ線分布を限られた計算時間内で精度よく再現するために、例えば以下のように選択してもよい。
　図６は、モンテカルロシミュレーションの手順の一例を示すフローチャートである。
　ステップＳ１において、演算装置１０２は、あらかじめ決められた初期投影角度（例えば０度、９０度、１８０度、２７０度）について、上述したモンテカルロシミュレーションをおこない、散乱Ｘ線分布を計算する。
　ステップＳ２において、演算装置１０２は、散乱Ｘ線量の変化の大きさを、下記（式７）にしたがって計算する。
　下記（式７）において、ｄａｔａｓｃａｔ（ｃｈ，ｓｌ，θ）は、検出器１回転あたりに得られる散乱Ｘ線分布を表し、チャネル方向Ｘ線検出素子６の位置（ｃｈ）、スライス方向のＸ線検出素子６の位置（ｓｌ）、投影角度（θ）の関数として表わされる。ただし、ｃｈ、ｓｌ、θは離散的変数である。
Figure JPOXMLDOC01-appb-M000001
　ステップＳ３において、演算装置１０２は、散乱Ｘ線量の変化の大きさが最大となる区間を求める。
　ステップＳ４において、演算装置１０２は、ステップＳ３で求めた区間中に投影角度を設定し、上述したモンテカルロシミュレーションをおこない、散乱Ｘ線分布を計算する。
　例えば、投影角度θ１からθ２の間で上記（式７）の変化が最大値を示した場合（Ｓ３）、演算装置１０２は、つぎの投影角度として、例えば、（θ１＋θ２）／２をとり、再度モンテカルロシミュレーションをおこなって散乱Ｘ線分布を計算する（Ｓ４）。
　ステップＳ５において、演算装置１０２は、要求される統計精度を満たすまで、または、許容される計算時間を超えるまで（Ｓ５・Ｎｏ）、ステップＳ２に戻り、以降の処理をくり返す。そして、要求される統計精度を満たす、または、許容される計算時間を超えると（Ｓ５・Ｙｅｓ）、ステップＳ６に進む。
　ここで、モンテカルロシミュレーションに要する時間は、要求精度によって変わる。散乱Ｘ線補正を実施するに当たり、ユーザは、入力装置１０４を通して、モンテカルロシミュレーションに用いるＸ線統計量もしくは計算時間を任意に設定することができる。また、ユーザは、モンテカルロシミュレーションを任意のタイミングで中断することができる。ユーザによってシミュレーションが中断された場合、Ｘ線ＣＴ装置１００は、中断までに計算が完了しているデータを用いて散乱Ｘ線補正をおこなう。
　ステップＳ６において、演算装置１０２は、ステップＳ１およびステップＳ４で得られた散乱Ｘ線分布に対して補間処理および平滑化処理をおこない、本フローチャートによる処理を終了する。
　ステップＳ６の補間処理については、上記処理で得られた散乱Ｘ線分布は、投影角度方向に対して実測より少ないデータ数を持つため、中間角度データを補間するためにおこなう。データの補間方法としては、例えばスプライン補間法やラグランジュ補間法などが知られている。また、ステップＳ６の平滑化処理は、モンテカルロシミュレーションの統計量不足にともなう高周波成分を取り除くためにおこなう。平滑化処理は、例えば、チャネル方向やスライス方向、投影角度方向へ、データの移動平均をとることによっておこなう。
（散乱Ｘ線分布を用いた散乱Ｘ線補正：Ｆ８）
　つぎに、Ｆ７において推定した散乱Ｘ線分布を基に、撮影（Ｆ１）によって得られた投影データ（Ｒａｗ　Ｄａｔａ）を補正する散乱Ｘ線補正（図５のＦ８）の詳細について説明する。Ｘ線ＣＴ装置１００の演算装置１０２は、シミュレーションによって得られた散乱Ｘ線分布を用いて、撮影で得られた投影データから散乱Ｘ線の寄与を差し引く。
　具体的には、演算装置１０２は、下記（式８）に示す計算によって散乱Ｘ線補正をおこなう。下記（式８）において、ＤＡＴＡは撮影で得られた投影データ、ＤＡＴＡｒｅｆは撮影で得られた投影のリファレンスデータ（例えば、端部に位置するＸ線検出素子６の出力の平均データ）、ｄａｔａｓｃａｔは計算で得られた散乱Ｘ線分布データ、ｄａｔａｔｏｔａｌ ｒｅｆは計算で得られた全Ｘ線分布（＝直達Ｘ線分布＋散乱Ｘ線分布）のリファレンスデータであり、変数の表示は省略した。αは補正強度を表す定数である。なお計測データは、Ｘ線非照射時の出力（オフセット値）が０になるように補正されているものとする。
Figure JPOXMLDOC01-appb-M000002
　また、演算装置１０２は、下記（式９）に示す計算によってエアデータ（被写体３を空気とした場合の出力値）の散乱Ｘ線補正もおこなう。下記（式９）において、ＡＩＲは被写体３を空気として撮影した際の投影データである。また、各添え字の意味は上記（式８）と同様である。
Figure JPOXMLDOC01-appb-M000003
　上記（式８）および（式９）で算出した値を用いて、演算装置１０２は、Ｘ線検出素子６の感度補正（エア補正と呼ぶ）を下記（式１０）によっておこなう。
　ＮｅｗＤａｔａ／ＮｅｗＡｉｒ・・・（式１０）
　演算装置１０２は、上記（式１０）で示した散乱Ｘ線補正およびエア補正後のデータを基に、再び画像再構成をおこなう。これにより、Ｘ線ＣＴ装置１００は、散乱Ｘ線の影響のない良質な画像をユーザに提供することができる。
　なお、上記では散乱Ｘ線補正後にエア補正を実施した場合を示したが、エア補正後に散乱Ｘ線補正を実施することもできる。ただし、この場合には、エアデータに含まれる散乱Ｘ線は少ないものとして考慮しない。具体的には、下記（式１１）のような演算をおこなう。
Figure JPOXMLDOC01-appb-M000004
　以上説明したように、第２実施形態に係るＸ線ＣＴ装置１００では、被写体３の再構成画像（Ｉｍａｇｅ）を基に散乱Ｘ線分布（Ｓｃａｔ．Ｄｉｓｔ．）を推定するので、被写体３に応じた構造の変化を再現することができる。これにより、Ｘ線ＣＴ装置１００は、散乱Ｘ線分布を精度高く推定することができ、散乱Ｘ線に起因する画質低下を低減することができる。また、第２実施形態に係るＸ線ＣＴ装置１００によれば、物理相互作用を厳密に模擬可能なモンテカルロシミュレーションを用いるので、散乱Ｘ線分布をより精度高く推定することができる。
　また、第２実施形態に係るＸ線ＣＴ装置１００は、投影角度を代表角度方向にのみに限定してシミュレーションをおこなって散乱Ｘ線分布を求め、他の角度方向の散乱Ｘ線分布を補間によって求めるので、散乱Ｘ線分布の推定に要する計算時間を大幅に短縮することができる。さらに、第２実施形態に係るＸ線ＣＴ装置１００は、計算によって得られた散乱Ｘ線分布を平滑化処理するので、統計変動に起因する高周波ノイズを除去することができ、計算時間の短縮および過補正や誤補正の抑制を図ることができる。すなわち、第２実施形態に係るＸ線ＣＴ装置１００によれば、現実的な計算時間内で、散乱Ｘ線の影響のない良質なＣＴ画像を得ることができる。
≪第３実施形態≫
　第１実施形態および第２実施形態では、被写体３に対する計測（撮像：Ｆ１）をおこなう度に、仮想被写体３Ａ（Ｘ線吸収係数分布）に対するモンテカルロシミュレーションをおこない（Ｆ４）、散乱Ｘ線分布や点拡がり関数を求めた。第３実施形態では、事前のシミュレーションによって、被写体３を模擬したファントム（模擬被写体）を定義し、第１実施形態または第２実施形態で示した方法により、散乱Ｘ線分布または点拡がり関数の分布を詳細に求めておく。これにより、散乱Ｘ線補正に要する時間を大幅に短縮することができる。なお、第３実施形態において、Ｘ線ＣＴ装置１００の構成（図１、図２参照）および画像形成処理の流れ（図３，図５参照）については、第１実施形態または第２実施形態と同様であるので詳細な説明を省略する。
　第３実施形態において、Ｘ線ＣＴ装置１００は、例えば、ファントムとして楕円柱形の水を採用する。楕円柱であれば、変数は長径ａ、短径ｂのみである。以下、短径ｂの代わりに、扁平率ｆ＝１－ｂ／ａを変数として用いる。Ｘ線ＣＴ装置１００は、実際の被写体３の測定に先立って、さまざまな扁平率ｆを有する複数の楕円柱に対して、散乱Ｘ線分布または点拡がり関数の長径ａ依存性を求めておく。そして、Ｘ線ＣＴ装置１００は、任意の変数ａ，ｆに対する散乱Ｘ線分布や点拡がり関数の分布を、補間や関数フィッティングにより求め、データベース化しておく。
　実際に被写体３の撮影をおこなう場合、Ｘ線ＣＴ装置１００は、撮影によって得られた再構成画像を基に、被写体３の撮影部位形状を楕円状に近似し、その長径ａ、扁平率ｆを測定する。Ｘ線ＣＴ装置１００は、測定された長径ａ、扁平率ｆに近似する散乱Ｘ線分布または点拡がり関数の分布を、上述したデータベースから選択し、第１実施形態または２で示した方法によって、測定データに対して散乱Ｘ線補正を実施する。そして、Ｘ線ＣＴ装置１００は、散乱Ｘ線補正後の測定データを用いて画像再構成を実施する。これにより、散乱Ｘ線の影響のない良質な画像が得られる。
　以上説明したように、第３実施形態に係るＸ線ＣＴ装置１００によれば、第１実施形態または第２実施形態に係るＸ線ＣＴ装置１００と同様の効果に加え、つぎのような効果が得られる。すなわち、第３実施形態に係るＸ線ＣＴ装置１００によれば、計算に長時間を要するモンテカルロシミュレーションを被写体３の撮影前にあらかじめ実行しておくので、散乱Ｘ線補正に要する時間を大幅に短縮することができる。
　１　　　Ｘ線管球（Ｘ線源）
　２　　　開口部
　３　　　被写体
　３Ａ　　仮想被写体（模擬被写体、Ｘ線吸収係数分布）
　４　　　Ｘ線検出器
　５　　　散乱線防止グリッド
　６　　　セパレータ
　７　　　Ｘ線検出素子
　８　　　検出器モジュール
　９　　　Ｘ線焦点
　１００　Ｘ線ＣＴ装置
　１０１　記録装置
　１０２　演算装置
　１０３　制御装置
　１０４　入力装置
　１０５　出力装置
　Ｆ１　　撮影（撮影部)
　Ｆ２　　補正・再構成処理
　Ｆ３　　被写体の模擬（内部分布推定部）
　Ｆ４　　仮想被写体に対するモンテカルロシミュレーション（点拡がり関数推定部）
　Ｆ５　　投影データ補正（補正部）
　Ｆ６　　補正・再構成処理（画像化部）
　Ｆ７　　仮想被写体に対するモンテカルロシミュレーション（Ｘ線分布推定部）
　Ｆ８　　投影データ補正（補正部）
　Ｆ９　　補正・再構成処理（画像化部）

Claims (11)

  1. 　Ｘ線焦点からＸ線を発生するＸ線源と、前記Ｘ線を検出するためのＸ線検出素子が２次元配列されているＸ線検出器とが、被写体をはさんで対向配置された関係を保ちながら前記被写体の周りを回転し、複数の投影方向から前記被写体のＸ線透過像データを撮影する撮影部と、
    　前記撮影部によって撮影された前記Ｘ線透過像データに基づいて、前記被写体内部のＸ線吸収係数分布を推定する内部分布推定部と、
    　前記内部分布推定部によって推定された前記Ｘ線吸収係数分布を有する模擬被写体に対して、前記Ｘ線の物理相互作用を模擬するモンテカルロシミュレーションを実施して、前記被写体に由来する散乱の点拡がり関数を推定する点拡がり関数推定部と、
    　前記点拡がり関数推定部によって推定された前記点拡がり関数と、前記Ｘ線透過像データとを逆畳み込み積分して、前記Ｘ線透過像データを補正する補正部と、
    　前記補正部によって補正された前記Ｘ線透過像データを用いて、前記被写体のＸ線吸収係数分布画像を形成する画像化部と、を備える
    ことを特徴とするＸ線ＣＴ装置。
  2. 　前記点拡がり関数推定部は、
    　前記複数の投影方向よりも少ない数の投影方向である代表投影方向について前記モンテカルロシミュレーションを実施して、当該モンテカルロシミュレーションによって得られた前記代表投影方向間の前記点拡がり関数を補間することによって、前記複数の投影方向すべての前記点拡がり関数を推定する
    こと特徴とする請求項１に記載のＸ線ＣＴ装置。
  3. 　前記点拡がり関数推定部は、
    　前記２次元配列された前記Ｘ線検出素子のうち、所定数の素子の位置である代表素子位置において前記モンテカルロシミュレーションを実施して、当該モンテカルロシミュレーションによって得られた前記代表素子位置間の前記点拡がり関数を補間することによって、前記２次元配列された前記Ｘ線検出素子のすべての位置における前記点拡がり関数を推定する
    ことを特徴とする請求項１に記載のＸ線ＣＴ装置。
  4. 　Ｘ線焦点からＸ線を発生するＸ線源と、前記Ｘ線を検出するためのＸ線検出素子が２次元配列されているＸ線検出器とが、被写体をはさんで対向配置された関係を保ちながら前記被写体の周りを回転し、複数の投影方向から前記被写体のＸ線透過像データを撮影する撮影部と、
    　前記撮影部によって撮影された前記Ｘ線透過像データに基づいて、前記被写体内部のＸ線吸収係数分布を推定する内部分布推定部と、
    　前記内部分布推定部によって推定された前記Ｘ線吸収係数分布を有する模擬被写体に対して、前記Ｘ線の物理相互作用を模擬するモンテカルロシミュレーションを実施して、前記被写体に由来する散乱Ｘ線の分布を推定するＸ線分布推定部と、
    　前記Ｘ線分布推定部によって推定された前記分布に基づいて、前記Ｘ線透過像データから前記散乱Ｘ線の成分を除去する補正部と、
    　前記補正部によって前記散乱Ｘ線の成分が除去された前記Ｘ線透過像データを用いて、前記被写体のＸ線吸収係数分布画像を形成する画像化部と、を備える
    ことを特徴とするＸ線ＣＴ装置。
  5. 　前記Ｘ線分布推定部は、
    　前記複数の投影方向より少ない数の投影方向である代表投影方向について前記モンテカルロシミュレーションを実施して、当該モンテカルロシミュレーションによって得られた前記代表投影方向間の前記分布を補間することによって、前記複数の投影方向すべての前記分布を推定する
    こと特徴とする請求項４に記載のＸ線ＣＴ装置。
  6. 　前記Ｘ線分布推定部は、
    　前記代表投影方向の前記分布から前記散乱Ｘ線量の変化率が最大となる前記投影方向を推定し、前記変化率が最大となる前記投影方向の前記散乱Ｘ線の分布を前記モンテカルロシミュレーションによって求める手続きをくり返して、前記代表投影方向を逐次選択することを特徴とする請求項５に記載のＸ線ＣＴ装置。
  7. 　前記内部分布推定部は、
    　ＣＴ値が既知な複数の物質の混合物として前記被写体を模擬することによって前記Ｘ線吸収係数分布を推定する
    ことを特徴とする、請求項１に記載のＸ線ＣＴ装置。
  8. 　前記モンテカルロシミュレーションにおけるＸ線統計量の設定と、前記モンテカルロシミュレーションを実行する計算時間の設定と、前記モンテカルロシミュレーションの任意のタイミングでの中断と、のうち、少なくともいずれか１つを可能とするユーザインターフェースをさらに備える
    ことを特徴とする、請求項１に記載のＸ線ＣＴ装置。
  9. 　前記内部分布推定部は、
    　前記撮影部によって撮影された前記Ｘ線透過像データに基づいて、前記被写体内部のＸ線吸収係数分布を推定することに換えて、
    　前記被写体を模擬する模擬被写体内部のＸ線吸収係数分布を前記被写体内部のＸ線吸収係数分布として推定する
    ことを特徴とする、請求項１に記載のＸ線ＣＴ装置。
  10. 　Ｘ線焦点位置からＸ線を発生するＸ線源と、前記Ｘ線を検出するためのＸ線検出素子が２次元配列されているＸ線検出器とが、被写体をはさんで対向配置された関係を保ちながら前記被写体の周りを回転し、複数の投影方向から前記被写体のＸ線透過像データを撮影する撮影工程と、
    　前記撮影工程で撮影された前記Ｘ線透過像データに基づいて、前記被写体内部のＸ線吸収係数分布を推定する内部分布推定工程と、
    　前記内部分布推定工程で推定された前記Ｘ線吸収係数分布を有する模擬被写体に対して、前記Ｘ線の物理相互作用を模擬するモンテカルロシミュレーションを実施して、前記被写体に由来する散乱の点拡がり関数を推定する点拡がり関数推定工程と、
    　前記点拡がり関数推定工程で推定された前記点拡がり関数と、前記Ｘ線透過像データとを逆畳み込み積分して、前記Ｘ線透過像データを補正する補正工程と、
    　前記補正工程で補正された前記Ｘ線透過像データを用いて、前記被写体のＸ線吸収係数分布画像を形成する画像化工程と、
    　を含んだことを特徴とする散乱Ｘ線補正方法。
  11. 　Ｘ線焦点位置からＸ線を発生するＸ線源と、前記Ｘ線を検出するためのＸ線検出素子が２次元配列されているＸ線検出器とが、被写体をはさんで対向配置された関係を保ちながら前記被写体の周りを回転し、複数の投影方向から前記被写体のＸ線透過像データを撮影する撮影工程と、
    　前記撮影工程で撮影された前記Ｘ線透過像データに基づいて、前記被写体内部のＸ線吸収係数分布を推定する内部分布推定工程と、
    　前記内部分布推定工程で推定された前記Ｘ線吸収係数分布を有する模擬被写体に対して、前記Ｘ線の物理相互作用を模擬するモンテカルロシミュレーションを実施して、前記被写体に由来する散乱Ｘ線の分布を推定するＸ線分布推定工程と、
    　前記Ｘ線分布推定工程で推定された前記分布に基づいて、前記Ｘ線透過像データから前記散乱Ｘ線の成分を除去する補正工程と、
    　前記補正工程で前記散乱Ｘ線の成分が除去された前記Ｘ線透過像データを用いて、前記被写体のＸ線吸収係数分布画像を形成する画像化工程と、
    　を含んだことを特徴とする散乱Ｘ線補正方法。
PCT/JP2012/082251 2011-12-12 2012-12-12 Ｘ線ｃｔ装置および散乱ｘ線補正方法 WO2013089155A1 (ja)

Priority Applications (3)

Application Number Priority Date Filing Date Title
CN201280060843.XA CN103987320B (zh) 2011-12-12 2012-12-12 X线ct装置以及散射x线校正方法
US14/364,278 US9307949B2 (en) 2011-12-12 2012-12-12 X-ray CT device and method for correcting scattered X-rays
JP2013549295A JP5815048B2 (ja) 2011-12-12 2012-12-12 Ｘ線ｃｔ装置

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2011270845 2011-12-12
JP2011-270845 2011-12-12

Publications (2)

Publication Number Publication Date
WO2013089155A1 true WO2013089155A1 (ja) 2013-06-20
WO2013089155A8 WO2013089155A8 (ja) 2013-10-10

Family

ID=48612598

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2012/082251 WO2013089155A1 (ja) 2011-12-12 2012-12-12 Ｘ線ｃｔ装置および散乱ｘ線補正方法

Country Status (4)

Country Link
US (1) US9307949B2 (ja)
JP (1) JP5815048B2 (ja)
CN (1) CN103987320B (ja)
WO (1) WO2013089155A1 (ja)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2016028673A (ja) * 2014-07-14 2016-03-03 朝日レントゲン工業株式会社 画像処理装置、画像処理方法、及びｘ線撮影装置
JP2019111346A (ja) * 2017-12-22 2019-07-11 キヤノンメディカルシステムズ株式会社 医用処理装置及び放射線診断装置
WO2020184361A1 (ja) * 2019-03-13 2020-09-17 日本電産リード株式会社 検出値補正システム、係数算出方法、及び検出値補正方法
JP2020534929A (ja) * 2017-09-28 2020-12-03 コーニンクレッカ フィリップス エヌ ヴェＫｏｎｉｎｋｌｉｊｋｅ Ｐｈｉｌｉｐｓ Ｎ．Ｖ． 深層学習ベースの散乱補正
JP2021023761A (ja) * 2019-08-01 2021-02-22 恵一 中川 Ｘ線コーンビームｃｔ画像再構成方法

Families Citing this family (25)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103674979B (zh) * 2012-09-19 2016-12-21 同方威视技术股份有限公司 一种行李物品ct安检系统及其探测器装置
US9724056B2 (en) * 2013-11-28 2017-08-08 Toshiba Medical Systems Corporation Method and system for spectral computed tomography (CT) with inner ring geometry
JP6145889B2 (ja) * 2014-03-24 2017-06-14 富士フイルム株式会社 放射線画像処理装置および方法並びにプログラム
CN104783819B (zh) * 2014-08-27 2017-10-03 上海联影医疗科技有限公司 散射校正方法及装置
CN104840211B (zh) * 2015-05-18 2018-12-25 上海联影医疗科技有限公司 一种投影图像的散射矫正方法及装置
CN105574828B (zh) * 2015-12-22 2019-01-25 沈阳东软医疗系统有限公司 图像散射校正方法、装置及设备
CN105631910B (zh) * 2015-12-28 2019-02-05 沈阳东软医疗系统有限公司 一种pet图像的重建方法及装置
CN105894537B (zh) * 2016-04-01 2019-07-26 西安电子科技大学 一种发射式光学投影断层成像衰减与散射校正方法
CN107260194B (zh) * 2016-04-08 2020-08-28 山西锦地裕成医疗设备有限公司 用于正电子断层成像的方法和图像重建方法及系统
DE102016206559B3 (de) * 2016-04-19 2017-06-08 Siemens Healthcare Gmbh Verfahren zur Korrektur eines Röntgenbilds auf Effekte eines Streustrahlenrasters, Röntgeneinrichtung, Computerprogramm und elektronisch lesbarer Datenträger
CN108065950B (zh) * 2016-11-14 2021-05-11 通用电气公司 一种放射成像方法及其系统
CN106361367B (zh) * 2016-12-01 2019-10-08 上海联影医疗科技有限公司 一种检测器的校正方法和使用该校正方法的装置及设备
KR20180075220A (ko) 2016-12-26 2018-07-04 삼성전자주식회사 멀티미디어 신호의 프로세싱 방법, 장치 및 시스템
US11002866B2 (en) * 2017-07-25 2021-05-11 Shimadzu Corporation Scattering estimation method and image processor
CN107595314B (zh) * 2017-08-31 2020-12-25 上海联影医疗科技股份有限公司 一种校正杂散射线的方法
US10758201B2 (en) * 2017-12-13 2020-09-01 Carestream Health, Inc. Variance reduction for monte carlo-based scatter estimation
EP3731759A4 (en) 2017-12-29 2021-01-20 Shanghai United Imaging Healthcare Co., Ltd. SYSTEMS AND METHODS FOR SCRAP CORRECTION OF AN IMAGE
CN107928690B (zh) * 2017-12-29 2021-06-22 上海联影医疗科技股份有限公司 确定主射线比例方法及图像去散射方法、系统及存储介质
CN111096761B (zh) * 2018-10-29 2024-03-08 上海西门子医疗器械有限公司 修正楔形滤波器散射的方法、装置和相关设备
CN109655878A (zh) * 2019-01-31 2019-04-19 中国人民解放军火箭军工程大学 基于X射线源的X/γ剂量仪现场校准方法、辐射防护系统
CN111488701B (zh) * 2020-06-28 2020-10-16 中国人民解放军国防科技大学 一种基于点扩散函数的计算机辅助装调方法、装置及设备
CN111833410A (zh) * 2020-07-07 2020-10-27 南京革佳智能科技有限公司 基于深度学习的x射线散射抑制方法
CN113237903B (zh) * 2021-06-15 2022-06-10 清华大学 基于双层平板探测器锥束ct的散射校正方法及装置
CN114399564B (zh) * 2022-03-25 2022-08-12 康达洲际医疗器械有限公司 一种基于散射识别的锥束计算机断层扫描成像方法与系统
CN114942466B (zh) * 2022-07-22 2022-12-13 康达洲际医疗器械有限公司 一种基于气隙相关性的x射线散射估计方法与系统

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005058760A (ja) * 2003-07-31 2005-03-10 Toshiba Corp 画像データ処理装置及び画像データ処理方法

Family Cites Families (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP3069411B2 (ja) 1991-09-04 2000-07-24 科学技術振興事業団 導波路型非対称方向性光結合器
BE1007766A3 (nl) 1993-11-10 1995-10-17 Philips Electronics Nv Werkwijze en inrichting voor computer tomografie.
JP4218908B2 (ja) 1998-04-24 2009-02-04 株式会社東芝 Ｘ線ｃｔ装置
US7330594B2 (en) 2003-07-31 2008-02-12 Kabushiki Kaisha Toshiba Image enhancement or correction software, method, apparatus and system for substantially minimizing blur in the scanned image
JP5052281B2 (ja) 2007-10-02 2012-10-17 株式会社東芝 Ｘ線ｃｔにおける散乱線強度分布の推定方法およびｘ線ｃｔ装置
JP5093195B2 (ja) 2009-07-06 2012-12-05 カシオ計算機株式会社 液晶表示パネル
CN101987021B (zh) 2010-12-06 2012-10-17 中国科学院深圳先进技术研究院 Ct系统的散射校正方法及ct系统
US8433154B2 (en) * 2010-12-13 2013-04-30 Carestream Health, Inc. Enhanced contrast for scatter compensation in X-ray imaging

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005058760A (ja) * 2003-07-31 2005-03-10 Toshiba Corp 画像データ処理装置及び画像データ処理方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
GENEVIEVE JARRY: "Characterization of scattered radiation in kV CBCT images using Monte Carlo simulations", MEDICAL PHYSICS, vol. 33, no. ISSUE, 2006, pages 4320 - 4329, XP012091941, DOI: doi:10.1118/1.2358324 *
J. A. SEIBERT: "X-ray scatter removal by deconvolution", MEDICAL PHYSICS, vol. 15, no. ISSUE, 1988, pages 567 - 575, XP002340970, DOI: doi:10.1118/1.596208 *
TAKASHI SATO: "X-sen Kyushu Katei ni Okeru Ryoshi Yuragi no Denpa (Gazo Kogaku Noise Tokusei)", JAPANESE JOURNAL OF RADIOLOGICAL TECHNOLOGY, vol. 50, no. 8, 1 August 1994 (1994-08-01), pages 1022 *

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2016028673A (ja) * 2014-07-14 2016-03-03 朝日レントゲン工業株式会社 画像処理装置、画像処理方法、及びｘ線撮影装置
JP2020534929A (ja) * 2017-09-28 2020-12-03 コーニンクレッカ フィリップス エヌ ヴェＫｏｎｉｎｋｌｉｊｋｅ Ｐｈｉｌｉｐｓ Ｎ．Ｖ． 深層学習ベースの散乱補正
JP2019111346A (ja) * 2017-12-22 2019-07-11 キヤノンメディカルシステムズ株式会社 医用処理装置及び放射線診断装置
WO2020184361A1 (ja) * 2019-03-13 2020-09-17 日本電産リード株式会社 検出値補正システム、係数算出方法、及び検出値補正方法
JP2021023761A (ja) * 2019-08-01 2021-02-22 恵一 中川 Ｘ線コーンビームｃｔ画像再構成方法
JP7022268B2 (ja) 2019-08-01 2022-02-18 恵一 中川 Ｘ線コーンビームｃｔ画像再構成方法

Also Published As

Publication number Publication date
US20140328452A1 (en) 2014-11-06
WO2013089155A8 (ja) 2013-10-10
JP5815048B2 (ja) 2015-11-17
US9307949B2 (en) 2016-04-12
CN103987320A (zh) 2014-08-13
JPWO2013089155A1 (ja) 2015-04-27
CN103987320B (zh) 2016-04-13

Similar Documents

Publication Publication Date Title
JP5815048B2 (ja) Ｘ線ｃｔ装置
Sun et al. Improved scatter correction using adaptive scatter kernel superposition
US9265475B2 (en) Methods and apparatus for scatter correction for CBCT system and cone-beam image reconstruction
US8965078B2 (en) Projection-space denoising with bilateral filtering in computed tomography
Li et al. Scatter kernel estimation with an edge-spread function method for cone-beam computed tomography imaging
Desai et al. Practical evaluation of image quality in computed radiographic (CR) imaging systems
JP6567094B2 (ja) 放射線映像の処理方法及び放射線撮影システム
Mackenzie et al. Characterisation of noise and sharpness of images from four digital breast tomosynthesis systems for simulation of images for virtual clinical trials
US20200279410A1 (en) Scatter Correction Method and Apparatus for Dental Cone-Beam CT
Schörner et al. Comparison between beam-stop and beam-hole array scatter correction techniques for industrial X-ray cone-beam CT
US20140218362A1 (en) Monte carlo modeling of field angle-dependent spectra for radiographic imaging systems
JP2014236810A (ja) 画像処理装置、断層撮影装置、画像処理方法およびプログラム
JP2015118074A (ja) Ｘ線断層撮影装置およびｘ線断層撮影方法
Pauwels What is CBCT and how does it work?
US20060088140A1 (en) System and method for performing scatter measurement in volumetric CT
Marshall et al. Performance evaluation of digital breast tomosynthesis systems: physical methods and experimental data
Peterzol et al. A beam stop based correction procedure for high spatial frequency scatter in industrial cone-beam X-ray CT
JP6676337B2 (ja) 放射線撮像システム、放射線画像の情報処理装置、放射線画像の情報処理方法、及び、そのプログラム
JP4584550B2 (ja) Ｘ線計測装置
Zhao et al. Three-dimensional cascaded system analysis of a 50 µm pixel pitch wafer-scale CMOS active pixel sensor x-ray detector for digital breast tomosynthesis
Jang et al. Fourier analysis of noise characteristics in cone-beam microtomography laboratory scanners
JP4703221B2 (ja) Ｘ線ｃｔ装置
JP5753502B2 (ja) 画像処理装置および方法
WO2021157403A1 (ja) 画像処理装置および方法、プログラム
Kim et al. Spatial resolution and blurring artifacts in digital X-ray tomosynthesis

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: 12857330

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 2013549295

Country of ref document: JP

Kind code of ref document: A

WWE Wipo information: entry into national phase

Ref document number: 14364278

Country of ref document: US

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 12857330

Country of ref document: EP

Kind code of ref document: A1