JP2011203160A - Method and program for reconstruction of x-ray ct image - Google Patents

Abstract

PROBLEM TO BE SOLVED: To reduce an artifact of a radial pattern occurring in image reconstruction by X-ray absorption difference of a sample.SOLUTION: An X-ray continuous spectrum is stored, X-ray projection data g of a single spectrum is generated, and the correction amount is initialized. Correction X-ray projection data g' is acquired that corresponds to the case where the X-ray projection data is corrected using the correction amounts h and t, and X-rays of the single spectrum is radiated to a photographing object in consideration of the energy dependence of the X-ray absorption coefficients of two substances included in the photographing object. A CT image is reconstructed by reverse projection operation using an image reconstruction algorithm for single wavelength. The region of each substance is extracted from the CT image using each of thresholds TH1, TH2, or the like for distinguishing each substance, and the length of the X-ray path of this region is calculated. The X-ray projection data f of the continuous spectrum is calculated using the energy of an occurring X-ray continuous spectrum, the length of the region of each substance, and the X-ray absorption coefficient of each substance. The correction amounts h and t are calculated based on a difference between the X-rays projection data g of the single spectrum and the X-ray projection data f of the continuous spectrum.

Description

  The present invention relates to an X-ray CT (Computed Tomography) image reconstruction method and an X-ray CT image reconstruction program.

  In general, in the X-ray CT image reconstruction method, X-rays having a normal wavelength distribution are projected from a number of directions around an object to be photographed, and projection data as X-ray transmission energy is obtained. Is backprojected using a single wavelength image reconstruction algorithm such as a filter-corrected backprojection method or a successive approximation method to reconstruct an image (CT image).

  In the X-ray CT image reconstruction method described above, if a substance region having an extremely large X-ray absorption coefficient, such as a metal, is present on the object to be imaged, the radiation quality of the X-ray CT image having the wavelength distribution is transmitted through the substance region. Due to the beam hardening phenomenon that changes, a radial noise pattern called an image artifact is generated in the CT image. This image artifact hinders discrimination of the reconstructed image at a medical site, a nondestructive inspection site, or the like.

  Next, the beam hardening phenomenon will be described with reference to FIGS. 40, 41, 42, and 43.

FIG. 40 is a diagram showing a general X-ray CT image reconstruction apparatus.

  In FIG. 40, 1 is an X-ray generator, 2 is an X-ray detector, 3 is an object to be imaged, 4 is a control circuit for controlling the X-ray generator 1, the X-ray detector 2 and the object to be imaged 3. . The control circuit 4 is configured by a computer and includes a CPU, a ROM, a RAM, and the like.

  FIG. 41 shows the X-ray absorption coefficient μ of the substance of the imaging object 3 in FIG.

  As shown in FIG. 41, not only the distribution of the X-ray absorption coefficient μ differs for each substance, but the X-ray absorption coefficient μ varies depending on the X-ray wavelength even for the same substance, that is, it has energy dependency.

  42 shows an X-ray energy spectrum of the CT image reconstruction apparatus of FIG. 40, (A) shows a continuous energy spectrum (hereinafter simply referred to as a continuous spectrum) of X-rays generated by the X-ray generator 1, and (B). Indicates a continuous spectrum of X-rays detected by the X-ray detector 2 after passing through the imaging object 3 having a size L and an X-ray absorption coefficient μ (E).

  As shown in FIG. 42, the continuous spectrum of X-rays differs before and after transmission through the object 3 to be imaged.

  Incidentally, in the CT image reconstruction calculation in the X-ray CT image reconstruction apparatus of FIG. 40, the energy dependence of the X-ray absorption amount is not taken into consideration. That is, it is premised on X-ray energy irradiation of a single energy spectrum with constant μ (hereinafter simply referred to as single spectrum), that is, single wavelength irradiation. This is because the X-ray continuous spectrum depends on the apparatus, the type of X-rays to be irradiated, imaging conditions, etc., and it is difficult to cope with this. However, actual CT imaging has a wavelength distribution for normal X-rays. As a result, a contradiction occurs due to the calculation at the time of reconstruction of the single wavelength CT image, and an image artifact as shown in FIG. 43 occurs.

  A first conventional X-ray CT image reproduction method for reducing image artifacts is to extract a material region having a high X-ray absorption coefficient, for example, a metal region from a sinogram of projection data, and modify the projection data of the material region by virtual projection. Corrected projection data is calculated, and the corrected projection data is backprojected to reconstruct a reconstructed image (CT image) (see Patent Document 1).

In the first conventional X-ray CT image construction method described above, severe image artifacts can be reduced, but it is difficult to extract a soft tissue portion and correction of projection data is insufficient.

A second conventional X-ray CT image construction method for reducing image artifacts is to change the projection conditions and perform multiple CT scans to obtain a plurality of projection data to limit the types of X-ray continuous spectrum, The absorption amount is obtained (see: Patent Document 2, Non-Patent Document 1).

JP 2006-167161 A JP 2008-241376 A

Asahi Lee, Michihiko Koseki, Hitoshi Kimura, Norio Inoue: Research on image quality improvement of X-ray CT images (proposal of image reconstruction method considering X-ray spectrum), Society of Instrument and Control Engineers 24th Sensing Forum Data, pp.107- 111,2007

  However, in the second conventional X-ray CT image reconstruction method described above, since a plurality of X-ray continuous spectra are required, the number of CT imaging increases, and different X-ray continuation by changing the imaging conditions. There is a problem that it is difficult to obtain a spectrum.

  In order to solve the above-described problems, the present invention is an X-ray CT image reconstruction method, wherein a generated X-ray continuous energy spectrum is stored as known information, and X-ray projection data of a single energy spectrum of an imaging object is generated. The correction amount is initialized, the X-ray projection data of the single energy spectrum is corrected with the correction amount, and the energy dependency of the X-ray absorption coefficients of at least two substances included in the imaging target is taken into consideration. Corrected X-ray projection data corresponding to irradiation with X-rays with a single energy spectrum is obtained, and the CT image is reconstructed by backprojecting the corrected X-ray projection data using a single wavelength image reconstruction algorithm. The region of each substance is extracted from the CT image using at least one threshold value for distinguishing each substance, the length of the X-ray path of these areas is calculated, and the generated X-ray continuous The X-ray projection data of the continuous energy spectrum is calculated using the energy of the energy spectrum, the length of each substance region, and the X-ray absorption coefficient of each substance, and the X-ray projection data of the single energy spectrum and the X of the continuous energy spectrum are calculated. The correction amount is calculated based on the difference from the line projection data. That is, when a CT image (original CT image) obtained using a conventional image reconstruction algorithm for a single wavelength is irradiated with X-rays having a single energy spectrum on an estimated imaging target, X Estimate line projection data. For this estimation, at least two substance regions are extracted from the original CT image using at least one threshold, and the length of each X-ray path is calculated. Next, the X-ray projection data of the continuous energy spectrum actually irradiated is corrected using the energy of the X-ray continuous energy spectrum of the known information, the length of the region of the substance, and the X-ray absorption coefficient, and the single energy spectrum A value close to the X-ray projection data is obtained. In this way, the degree of beam hardening exerted by at least two substances classified by at least one threshold is individually converted and expressed as a correction amount.

  In addition, the correction amount calculation is repeated from the above correction so that the correction amount converges within a predetermined range. In this case, if the threshold value is unknown, the threshold value is estimated from the histogram of the CT image to increase the number of threshold values. In this way, the correction amount is updated and the image artifacts are gradually reduced so that the length of the X-ray path of the material region becomes accurate.

  According to the present invention, since the X-ray projection data of the single energy spectrum is not changed, CT imaging is only required once, and image artifacts can be reduced by repeating the correction.

It is a flowchart for demonstrating 1st Embodiment of the X-ray CT image reconstruction method which concerns on this invention. It is a graph which shows the example of the X-ray continuous spectrum of step 101 of FIG. It is a figure explaining the threshold value of step 106 of FIG. 1, (A) is a CT image, (B) is a pixel value histogram. It is a graph for demonstrating the correction amount of step 108 of FIG. It is a table | surface which shows the continuous spectrum of the X-ray | X_line generate | occur | produced of an X-ray generator. It is a table | surface which shows the X-ray absorption coefficient of a living bone. It is a figure which shows the bone | column cylinder model as a imaging | photography target object used for the 1st simulation example of FIG. The bone cylinder model before correction | amendment of FIG. 7 is shown, (A) is a CT image, (B) is a graph which shows (A) CT value on a center line. FIG. 9 is a CT image showing extraction of a region of a substance from the CT image in FIG. The bone cylinder model of the 1st correction | amendment of FIG. 7 is shown, (A) is a CT image, (B) is a graph which shows CT value on the centerline of (A). It is a table | surface which shows the X ray absorption coefficient of iron. It is a table | surface which shows the X-ray absorption coefficient of an acrylic resin. FIG. 4 is a diagram illustrating a plurality of same-type pin models as photographing objects used in the second simulation example of FIG. 1, where (A) is a perspective image and (B) is a top image. FIG. 14 shows a plurality of same-type pin models before correction in FIG. 13, (A) is a CT image, and (B) is a graph showing (A) CT values on the center line. 14A shows a region extraction of the substance from the CT image of FIG. 14A, FIG. 14A shows a CT image showing the iron region extraction, and FIG. 14B shows a CT image showing the acrylic resin region extraction. FIG. 14 shows a plurality of same-type pin models for the first correction in FIG. 13, (A) is a CT image, and (B) is a graph showing CT values on the center line of (A). It is a table | surface which shows the X-ray absorption coefficient of aluminum. FIG. 6 is a top image showing a plurality of different pin models as photographing objects used in the third simulation example of FIG. 1. FIG. 18 shows a plurality of different pin models before correction in FIG. 18, (A) is a CT image, and (B) is a graph showing (A) CT values on the center line. It is a graph which shows the pixel value histogram of FIG. FIG. 19A shows a region extraction of a substance from the CT image of FIG. 19, (A) is a CT image showing the extraction of the iron region, (B) is a CT image showing the region extraction of aluminum, (C) is an acrylic resin It is CT image which shows area | region extraction. FIG. 19 shows a plurality of different pin models for the first correction in FIG. 18, (A) is a CT image, and (B) is a graph showing a pixel value histogram of (A). FIG. 19 shows a plurality of different pin models for the second correction in FIG. 18, (A) is a CT image, and (B) is a graph showing a pixel value histogram of (A). FIG. 19 is a graph showing a plurality of different pin models for the third correction in FIG. 18, (A) is a CT image, and (B) is a pixel value histogram of (A). It is a graph which shows convergence of the corrected amount of the 3rd simulation of FIG. FIG. 1 shows a first experimental example in FIG. 1, where (A) is a CT image before correction, and (B) is a CT image after correction. The 2nd experiment example of FIG. 1 is shown, (A) is a CT image before correction | amendment, (B) is a CT image after correction | amendment. FIG. 3 shows a third experimental example of FIG. 1, (A) is a perspective view of an iron nail-rice grain model, and (B) is a graph showing an X-ray spectrum at the time of photographing. The 3rd example of an experiment of Drawing 1 is shown, (A) is a CT image before amendment, and (B) is a CT image after amendment. It is a CT image after correction by extracting only iron in the third experimental example. It is a flowchart for demonstrating 2nd Embodiment of the X-ray CT image reconstruction method which concerns on this invention. It is a graph which shows the example of the X-ray absorption coefficient of substance iron in the steps 202 and 203 of FIG. 31, aluminum, and an acrylic resin. FIG. 32 is a CT image after the first correction in the first simulation example of FIG. 31. FIG. FIG. 31 shows a second simulation example of FIG. 31, where (A) is a CT image after the first correction, and (B) is a CT image after correction. It is a flowchart for demonstrating 3rd Embodiment of the X-ray CT image reconstruction method which concerns on this invention. FIG. 35 shows a first simulation example, and (A), (B), (C), (D), (E), and (F) are 3 before correction, after the first correction, after the second correction, CT images after the fourth correction, after the fourth correction, and after the fifth correction. FIG. 35 shows a second simulation example, and (A), (B), (C), (D), and (E) are before correction, after the first correction, after the second correction, and after the third correction, It is a CT image after the fourth correction. It is a figure which shows the biological tissue model as a imaging | photography object used for the 1st simulation example of FIG. The CT image of FIG. 38 is shown, (A) is a CT image before correction, and (B) is a CT image after correction. It is a figure which shows a general X-ray CT image reconstruction apparatus. It is a graph which shows the X-ray absorption coefficient of the substance of the imaging target object of FIG. It is a graph which shows the X-ray continuous spectrum of FIG. It is CT image which shows the example of an image artifact.

  FIG. 1 is a flowchart for explaining a first embodiment of an X-ray CT image reconstruction method according to the present invention. The flowchart in FIG. 1 is executed by the X-ray CT image reconstruction apparatus in FIG.

First, in step 101, the X-ray continuous spectrum of the X-ray generator 1 is stored in advance. The X-ray continuous spectrum determined by the type of X-ray tube and the imaging conditions of the X-ray generator 1 can be measured in advance with a dedicated device. The X-ray continuous spectrum is actually stored as discrete values as E _i and its photon number n (E _i ).

Next, in step 102, a correction amount _ht and a correction count t described later are initialized. In other words,
h ₀ ← 0
t ← 0
Here, since h ₀ exists as many as the number of projection data, h ₀ represents the generic name.

Next, in step 103, single spectrum projection data (X-ray absorption amount) g in a single spectrum, that is, wavelength i is generated from the detection signal of the X-ray detector 2. The wavelength i employs a representative value calculated from the X-ray tube voltage.

When the X-ray absorption amount g in step 103 is theoretically considered, the X-ray absorption amount g is expressed by the equation (1).
However, I _in : Energy at the time of incidence of substance k
I _out : Energy after permeation of substance k
L _k : Length of the X-ray path of the substance k μ _ik : X-ray absorption coefficient depending on the wavelength i and the substance k If the substance through which the X-ray passes is composed of n kinds of substances, the formula 1 is Equation 2
When Expression 2 is expressed logarithmically, Expression 3 is obtained.

Next, in step 104, it obtains the correction projection data g 'by subtracting the correction amount h _t in a single spectral projection data g.
g '= g-h _t
In this case, since g ′ exists as many as the number of projection data, it represents the generic name.

  Next, in step 105, the corrected projection data g 'is backprojected by a filter-corrected backprojection method or the like to reconstruct a CT image. The backprojection method in this case uses a single wavelength image reconstruction algorithm for reconstructing conventional single spectrum X-ray projection data.

Next, at step 106, the length Lk of the region of the X-ray path of the substance _k is calculated from the CT image using the known threshold values TH1, TH2,. For example, if the CT image shown in FIG. 3A includes three kinds of substances, that is, iron, aluminum, and acrylic resin as shown in FIG. 3B, in step 106, TH1 is set to a pixel value of 14000. , TH2 is set in advance to the pixel value 600. If the target object 3 contains a plurality of unknown substances, the pixel value histogram of FIG. 3B is calculated in step 106, and the pixel value concentration point, for example, (( B) The pixel values −900, 1800, 8000, and 20000 are discriminated, and the threshold value is estimated and set between them. This threshold value estimation algorithm can be easily obtained.

  Next, in step 107, continuous spectrum projection data f corresponding to each single spectrum projection data g is calculated.

The continuous spectrum projection data f includes the energy E _i of the X-ray continuous spectrum of the X-ray generator 1 stored in Step 101, the X-ray absorption coefficient μ _ik of the substance k, and the X-ray of the substance k calculated in Step 106. Calculation is performed based on the path length L _k . That is, the X-ray energy I _in at the time of incidence of the substance k is expressed by the equation (4).
Further, the X-ray energy I _out after the transmission of the substance k is expressed by the equation (5).
Expressions 6 and 6 are obtained from Expressions 4 and 5.
When Expression 6 is expressed logarithmically, continuous spectrum projection data (X-ray absorption amount) f is expressed by Expression 7.

Next, in step 108,
t ← t + 1
And

Next, in step 109, the t-th correction amount _ht is calculated. That is,
h _t = g-f

That is, as shown in FIG. 4, even with the same length L _k of the X-ray path of the material k, single spectral projection data (X-ray absorption) g and continuous spectrum projection data (X-ray absorption) the value is in the f are different from each other by the correction amount h _t. Since the correction amount h _t is different, the correction amount h _t appears as an image artifact. There is also a slight difference in the length L _k of the X-ray path of the substance k from the projection data g ′ where the image artifact in step 106 exists. Accordingly, in step 104 can be made close to a single spectral projection data g in the continuous spectrum projection data f by correcting by subtracting this correction amount h _t, also stably converge length L _k, along with this, An image artifact reduction effect occurs.

Step 110 repeats the flow of steps 104 to 109 until the correction amount h _t converges to a predetermined range and | h _t −h _t−1 | becomes less than ε. When | h _t −h _t−1 | becomes less than ε, the routine ends in step 111. In this case as well, | h _t -h _t-1 | indicates the total sum of projection data.

The important flow chart of FIG. 1, a single spectral projection data g obtained by one CT imaging unchanged, by updating only the correction amount h _t, gradually reduce image artifacts by this update The calculation of the length L _k of the X-ray path of the substance k is also becoming accurate.

  A simulation example according to the flowchart of FIG. 1 will be described. In other words, in this case, a virtual projection method in which CT imaging is performed on the control circuit (computer) 4 is used, and the X-ray continuous spectrum of the X-ray generator 1 is approximated by nine discrete values shown in FIG.

The first simulation example is applied to a bone cylinder model shown in FIG. 7 using a substance having an X-ray absorption coefficient μ _ik shown in FIG. 6, that is, a living bone. That is, the substance k and its X-ray absorption coefficient μ _ik are assumed to be known. In this case, since the bone cylinder model is uniform and simple in shape, the influence of beam hardening is small, and therefore, it is considered that there are few image artifacts.

In the first simulation example, in step 101, the X-ray continuous spectrum of FIG. 5 and the X-ray absorption coefficient μ _ik of FIG. 6 are stored. In step 103, single spectrum projection data g is generated based on these, and in step 105, a CT image is reconstructed from the single spectrum projection data g. In this case, the CT image shown in FIG. 8A is obtained. At first glance, no image artifact appears, but as shown in FIG. 8B, the CT values on the center line of the CT image are plotted. Then, it can be seen that the luminance increases at the pixel at the center of the cylinder and decreases at the pixels at the periphery of the cylinder. That is, the CT value of the pixels in the periphery of the cylinder is larger than the true value 1240. Thus, even with a simple shape, an appropriate amount of X-ray absorption cannot be obtained due to the effect of beam hardening, and the quantitativeness of the CT value is lost.

Next, the first correction is performed on the CT image of FIG. That is, in step 106, the region of the material living bone shown in FIG. 9 is extracted using the known threshold value TH1, for example, 500 for the pixel of the CT image of FIG. As shown in FIG. 9, the region of the material living bone is accurately extracted. The length _Lk is calculated for each X-ray path of this extraction region. Next, in step 107, the continuous spectrum projection data f is calculated using the length L _k , the energy E _i of the X-ray continuous spectrum in FIG. 5 and the X-ray absorption coefficient μ _ik in FIG. The first correction amount h ₁ is calculated, the single spectrum projection data g is corrected with the correction amount h ₁ , and the CT image shown in FIG. 10A is reconstructed from the corrected projection data.

  When the pixels on the center line are plotted with respect to the CT image of the first correction in FIG. 10A, as shown in FIG. 10B, the luminance of the pixels around the cylinder rises and becomes uniform inside the cylinder. It can be seen that the CT true value 1240 of 50 keV is approaching.

The second simulation example is applied to a plurality of similar pin models shown in FIG. 13 using a substance having an X-ray absorption coefficient μ _ik shown in FIGS. 11 and 12, that is, iron and acrylic resin. As a similar model, when non-destructive inspection of electronic parts and semiconductors, if metal parts are embedded in the resin, image artifacts occur around the metal parts and the shape of the resin part around the metal parts, The composition may not be discriminated. Further, when the implant metal is embedded in the human body, an image artifact may occur around the implant metal, and the shape and composition of the muscle portion around the implant metal may not be discriminated. In this case as well, the substance k and its X-ray absorption coefficient μ _ik are assumed to be known. This multiple-type pin model has a slightly complicated shape, so that the effect of beam hardening is large, and therefore image artifacts are considered to be large to some extent.

In the second simulation example, in step 101, the energy E _i of the X-ray continuous spectrum in FIG. 5 and the X-ray absorption coefficient μ _ik in FIGS. 11 and 12 are stored. In step 103, single spectrum projection data g is generated based on these, and in step 105, a CT image is reconstructed from the single spectrum projection data g. In this case, the CT image shown in FIG. 14A is obtained, and a strong image artifact appears. As shown in FIG. 14B, when the CT values on the center line of the CT image are plotted, the brightness rapidly increases at the pixel at the center of the pin, and the brightness rapidly decreases at the pixel at the periphery of the pin. I understand. That is, the CT value of the pixel around the pin is larger than the true value 20000. Thus, in a complicated shape, an appropriate amount of X-ray absorption cannot be obtained due to the influence of strong beam hardening, and the quantitative property of the CT value is lost.

Next, the first correction is performed on the CT image of FIG. That is, in step 106, each region of the substance iron and the acrylic resin shown in FIG. 15 using a known threshold value TH1 such as 1300 and threshold value TH2 such as 900 for the pixel of the CT image of FIG. To extract. In this case, as shown in FIG. 15A, the region of the substance iron is accurately extracted, but as shown in FIG. 15B, the region of the substance acrylic resin is affected by artifacts caused by the substance iron. Is inaccurate. The length L _k is calculated for each X-ray path of these extraction regions. Next, in step 107, the continuous spectrum projection data f is calculated using these lengths L _k , the energy E _i of the X-ray continuous spectrum of FIG. 5, and the X-ray absorption coefficient μ _ik of FIGS. In Step 109, the first correction amount h ₁ is calculated, the single spectrum projection data g is corrected with the correction amount h ₁ , and the CT image shown in FIG. 16A is obtained from the corrected projection data g ′. Reconfigure.

  When the pixels on the center line are plotted with respect to the first-corrected CT image of FIG. 16A, the brightness of the pixels around the pin increases and becomes uniform inside the pin as shown in FIG. It can be seen that the CT true value is approaching 20000. At the same time, it can be seen that the luminance of the pixels in the material acrylic resin region becomes uniform and approaches the CT true value of -680 of 50 keV.

The third simulation example is applied to a plurality of different pin models shown in FIG. 18 using a substance having an X-ray absorption coefficient μ _ik shown in FIGS. 11, 12, and 17, that is, iron, acrylic resin, and aluminum. In this case as well, the substance k and its X-ray absorption coefficient μ _ik are assumed to be known. This multi-species pin model includes a material having a plurality of high X-ray absorption coefficients, so that the effect of beam hardening is very large, and therefore, the image artifact is considered to be large.

In the third simulation example, in step 101, the energy E _i of the X-ray continuous spectrum in FIG. 5 and the X-ray absorption coefficient μ _ik in FIGS. 11, 12, and 17 are stored. In step 103, single spectrum projection data g is generated based on these, and in step 105, a CT image is reconstructed from the single spectrum projection data g. In this case, the CT image shown in FIG. 19A is obtained, and strong image artifacts due to the material iron appear. As shown in FIG. 19B, when the CT values on the center line of the CT image are plotted, the brightness rapidly increases at the pixel at the center of the pin, and the brightness rapidly decreases at the pixel at the periphery of the pin. I understand. That is, the CT value of the pixel around the pin is larger than the true value 20000 or 1800. When the pixel histogram at this time is shown in FIG. 20, the distribution of the substance iron is at a place of about 23000 away from the true value 20000. Accordingly, the threshold value TH1 in step 106 in FIG. 1 is estimated and set to 8000, and the region of the material iron can be extracted. In this case, even if the material iron is unknown, one threshold value TH1 is estimated and set. As a result, as shown in FIG. 21A, the region of the substance iron can be accurately extracted. However, as shown in FIGS. 21B and 21C, the regions of the substance aluminum and the acrylic resin are changed to the substance iron. It cannot be accurately extracted due to the effect of the artifact.

Next, the first correction is performed on the CT image of FIG. That is, the length L _k is calculated for each X-ray path of the extraction region. Next, the continuous spectrum projection data f is calculated using these lengths L _k , the energy E _i of the X-ray continuous spectrum of FIG. 5, and the X-ray absorption coefficient μ _ik of FIGS. 11, 12, and 17, The first correction amount h ₁ is calculated, the single spectrum projection data g is corrected with the correction amount h ₁ , and the CT image shown in FIG. 22A is reconstructed from the corrected projection data g ′. In the CT image of FIG. 22A, image artifacts remain and correction is insufficient. The reason for this is considered that the pixel value of the image artifact caused by the material iron is almost the same as the pixel value of the material aluminum, so that the shape of the image artifact is reconstructed assuming that the material aluminum exists. When the pixel histogram at this time is shown in FIG. 22B, the peak values of the materials aluminum and acrylic resin appear at their true values −680 and 1800. Therefore, if the threshold value TH2 in step 106 in FIG. 1 is set between -680 and 1800, for example, 1200, the respective regions of the substance aluminum and the acrylic resin can be accurately extracted. In this case also, the threshold value TH2 is estimated and set even if the materials aluminum and acrylic resin are unknown. Therefore, when n substances contained in the object to be imaged are unknown, n-1 threshold values are estimated and set by repeating steps 104 to 110 in FIG. As a result, an unknown substance region can be accurately extracted.

Next, the second correction is performed on the CT image of FIG. That is, in step 106, each region of the substance iron, aluminum, and acrylic resin is extracted using the estimated threshold value TH1 such as 8000 and threshold value TH2 such as 1200 for the CT image pixel of FIG. To do. Next, in step 107, the continuous spectrum projection data f is calculated using these lengths L _k , the energy E _i of the X-ray continuous spectrum of FIG. 5, and the X-ray absorption coefficient μ _ik of FIGS. and, then, at step 109, a second correction amount h ₂ is calculated, the correction amount h ₂ by correcting the single spectral projection data g, indicating a corrected projection data g 'in (a) of FIG. 23 Reconstruct the CT image.

Next, the third correction is performed on the CT image in FIG. That is, in step 106, each region of the substance iron, aluminum, and acrylic resin is extracted using the estimated threshold value TH1 such as 8000 and threshold value TH2 such as 1200 for the CT image pixel of FIG. To do. Next, in step 107, the continuous spectrum projection data f is calculated using these lengths L _k , the energy E _i of the X-ray continuous spectrum of FIG. 5, and the X-ray absorption coefficient μ _ik of FIGS. and, then, at step 109, it calculates the correction amount h ₃ of the third correction amount h ₃ by correcting the single spectral projection data g, indicating a corrected projection data g 'in (a) of FIG. 23 Reconstruct the CT image.

  When the pixels on the center line are plotted with respect to the CT image of the third correction in FIG. 24A, as shown in FIG. 24B, the luminance of the pixels around the pin is increased and becomes uniform inside the pin, 50 keV. It can be seen that the CT true value is approaching 20000 or 1800. At the same time, it can be seen that the luminance of the pixels in the material acrylic resin region becomes uniform and approaches the CT true value of -680 of 50 keV.

In the third simulation example, even if the CT image is subjected to the fourth correction, the fifth correction,..., The CT image hardly changes. That is, as the number of corrections increases, image artifacts are reduced, and the regions of the substance iron, aluminum, and acrylic resin can be extracted almost accurately. Therefore, as shown in FIG. 25, when the number of corrections exceeds 3 and the total | h _t -h _t-1 | becomes ε or less, each substance approaches its true value, and the error of the CT image is reduced. It can be confirmed that the effect of the correction is improved. At this time, the difference | h _t −h _t−1 | between the correction amount and the previous correction amount also approaches 0, indicating that the correction amount converges to a constant value.

  Next, an experimental example according to the flowchart of FIG. 1 will be described.

FIG. 26 is a first experimental example performed corresponding to the second simulation example described above. Since the continuous spectrum of the generated X-rays at the time of imaging was unknown, it was approximated by the discrete values of the continuous spectrum of FIG. Further, the values shown in FIGS. 11 and 12 were used as known values for the X-ray absorption coefficient μ _ik of the substances iron and acrylic resin. Furthermore, an appropriate value was used as a known threshold value for extracting each region of the substance iron and acrylic resin. As a result, the CT image before correction was as shown in FIG. 26A, but as a result of performing correction based on the flowchart of FIG. 1, the CT image shown in FIG. 26B was obtained. In other words, image artifacts were greatly reduced, and a correction effect was recognized.

FIG. 27 shows a second experimental example performed on the iron cross shape model. Also in this case, since the continuous spectrum of generated X-rays at the time of imaging was unknown, it was approximated by discrete values of the continuous spectrum of FIG. Further, the values shown in FIGS. 11 and 12 were used as known values for the X-ray absorption coefficient μ _ik of the substances iron and acrylic resin. Furthermore, an appropriate value was used as a known threshold value for extracting each region of the substance iron and acrylic resin. As a result, the CT image before correction was as shown in FIG. 27A, but as a result of performing correction based on the flowchart of FIG. 1, the CT image shown in FIG. 27B was obtained. Again, image artifacts were greatly reduced, and a correction effect was observed.

This is performed on a model in which nine iron nails shown in FIG. 28A are embedded in a pile of rice grains. The continuous spectrum of generated X-rays with a tube voltage of 130 kV at the time of imaging was approximated with discrete values of the continuous spectrum shown in FIG. Further, the values shown in FIGS. 11 and 12 were used as known values for the X-ray absorption coefficient μ _ik of the substance iron and rice grains. That is, the rice grain was regarded as an acrylic resin. Furthermore, 13500 and 750 were used as known thresholds for the extraction of each region of the substance iron and rice grains. As a result, the CT image before correction was as shown in FIG. That is, strong image artifacts are generated in the line connecting the iron nails. On the other hand, as a result of performing correction based on the flowchart of FIG. 1, a CT image shown in FIG. 29B was obtained. Again, image artifacts were greatly reduced, and a correction effect was observed.

  In the third experimental example, when only the threshold value 13500 was set and an area of only iron nails (iron) was extracted, rice grains could not be extracted accurately as shown in FIG. This is because it is considered that excessive absorption of the X-ray absorption amount occurs not only in iron nails (iron) but also in rice grains due to the continuous spectrum.

  Thus, according to the experimental example, it can be seen that the present invention is effective in individual modeling using a CT image and reverse engineering nondestructive inspection.

FIG. 31 is a flowchart for explaining a second embodiment of the X-ray CT image reconstruction method according to the present invention, in which the number of substances k constituting the object to be imaged 3 is small, and the X of these substances k. This is applied when the linear absorption coefficient μ _ik is unknown.

In FIG. 31, step 101 in FIG. 1 is replaced with step 201, and steps 202 and 203 are added. That is, in step 201, only the energy E _i of the X-ray continuous spectrum of FIG. 5 is stored, and in steps 202 and 203, based on the energy E _i of the X-ray continuous spectrum and the length L _k of the X-ray path of the substance k. To calculate the X-ray absorption coefficient μ _ik of the substance k.

The basis of steps 202 and 203 will be described. The projection data D at a certain position of the X-ray detector 2 is expressed by the equation (8).
Here, since the energy E _i of the X-ray continuous spectrum is known, the projection data of the formula 8 is a function of the X-ray absorption coefficient μ _ik and the length L _k determined in step 106. On the other hand, the X-ray absorption coefficient μ _ik is expressed by the equation (9).
As shown in FIG. 32, the validity of the formula 9 is the X-ray absorption obtained by the formula 9 when a _k = 134000, 7000, 200 for the materials iron, aluminum, and acrylic resin. The coefficient μ _ik can be confirmed from the fact that it substantially matches the actual X-ray absorption coefficient μ _ik .
When the formula of the formula 8 is expressed by the formula of the formula 9, the formula of the formula 10 is obtained.
In the equation (10), a _k can be solved as an n-order equation with an unknown number of types of substance k. Since the data actually detected can be obtained for the number of columns of the X-ray detector 2, for example, 512 × number of rotation directions of imaging, for example, 800, the number of types n of the substance k can be solved at most. That is, since projection data D ′ is obtained for each X-ray beam irradiated from a plurality of directions, the number of Expressions 10 is established. In this case, a _k can be _obtained if the number of Expression 10 is n or more. However, since the area length L _k has an error due to image artifacts, a _k is _obtained by the least square method. At this time, the substance k can be estimated based on the obtained a _k .

Accordingly, in step 202, a _k in the equation (10) is calculated by the least square method. Next, in step 203, the X-ray absorption coefficient μ _ik of the substance k is calculated by the equation (9).

  A simulation example according to the flowchart of FIG. 31 will be described. In other words, in this case, a virtual projection method in which CT imaging is performed on the control circuit (computer) 4 is used, and the X-ray continuous spectrum of the X-ray generator 1 is approximated by nine discrete values shown in FIG.

The first simulation example is applied to a plurality of same pin models shown in FIG. 13 using two substances, namely iron and acrylic resin. However, the X-ray absorption coefficient μ _ik of the materials iron and acrylic resin is unknown.

In the first simulation example, in step 201, the energy E _i of the X-ray continuous spectrum of FIG. 5 is stored. In step 103, single spectrum projection data g is generated based on these, and in step 105, a CT image is reconstructed from the single spectrum projection data g.

Next, the first correction is performed on the CT image. That is, in step 106, the regions of the substance iron and the acrylic resin are extracted using the known threshold value TH1 such as 1300 and threshold value TH2 such as 900 for the pixels of the CT image. The length L _k is calculated for each X-ray path of these extraction regions. Next, in step 202, a _k in the equation (10) is calculated by the method of least squares. as a result,
a ₁ = 124651
a ₂ = 41
was gotten. In addition, since the true value of a (iron) and a (acrylic resin) is 130,000 and 200, it turns out that the above-mentioned result value is close to a true value. Therefore, it is estimated that a ₁ = 124651 corresponds to iron and a ₂ = 41 corresponds to acrylic resin. Next, in step 203, the X-ray absorption coefficient μ _ik corresponding to the substance iron and the acrylic resin is calculated using the equation (9). Next, in step 107, the continuous spectrum projection data f is calculated using the length L _k described above, the energy E _i of the X-ray continuous spectrum in FIG. 5 and the X-ray absorption coefficient μ _ik obtained in step 203, then, in step 109, it calculates the first correction amount h _1, the correction amount h ₁ in correcting a single-spectrum projection data g, to reconstruct the CT image shown in FIG. 33 from the corrected projection data g ' . Thus, even when the X-ray absorption coefficient μ _ik is unknown, the image artifact is greatly reduced.

The second simulation example is applied to the multiple different pin model shown in FIG. 18 using three substances, that is, iron, aluminum and acrylic resin. However, the X-ray absorption coefficient μ _ik of the substances iron, aluminum and acrylic resin is unknown.

In the second simulation example, in step 201, the energy E _i of the X-ray continuous spectrum of FIG. 5 is stored. In step 103, single spectrum projection data g is generated based on these, and in step 105, a CT image is reconstructed from the single spectrum projection data g.

Next, the first correction is performed on the CT image. That is, in step 106, the regions of the substance iron, aluminum, and acrylic resin are extracted using the known threshold value TH1 such as 1300 and threshold value TH2 such as 900 for the pixels of the CT image. The length L _k is calculated for each X-ray path of these extraction regions. Next, in step 202, a _k in the equation (10) is calculated by the method of least squares. as a result,
a ₁ = 127682
a ₂ = -1702
a ₃ = 416
was gotten. In this case, the region equivalent to aluminum is inappropriate due to the influence of image artifacts, and thus the true value 7000 of a (aluminum) has been separated. Also, since the X-ray absorption coefficient is not negative,
a ₁ = 127682
a ₂ = 0
a ₃ = 416
And Next, in step 203, the X-ray absorption coefficient μ _ik corresponding to the substance iron, aluminum, and acrylic resin is calculated using the equation (9). Next, in step 107, the continuous spectrum projection data f is calculated using the length L _k described above, the energy E _i of the X-ray continuous spectrum in FIG. 5 and the X-ray absorption coefficient μ _ik obtained in step 203, then, in step 109, calculates the first correction amount h _1, the correction amount to correct a single-spectrum projection data g by h _1, the corrected CT image from the projection data g 'shown in (a) of FIG. 34 Reconfigure. As described above, the first correction reduces the image artifacts although it is insufficient.

Next, the second correction is performed on the CT image of FIG. That is, in step 106, the regions of the substance iron, aluminum, and acrylic resin are extracted using the known threshold value TH1 such as 1300 and threshold value TH2 such as 900 for the pixels of the CT image. The length L _k is calculated for each X-ray path of these extraction regions. Next, in step 202, a _k in the equation (10) is calculated by the method of least squares. as a result,
a ₁ = 123476
a ₂ = 6962
a ₃ = 21
was gotten. In addition, since the true value of a (iron), a (aluminum), and a (acrylic resin) is 130,000, 7000, and 200, it turns out that the above-mentioned result value is close to a true value. Therefore, it is estimated that a ₁ = 123476 corresponds to iron, a ₂ = 6962 corresponds to aluminum, and a ₃ = 21 corresponds to acrylic resin. Next, in step 203, the X-ray absorption coefficient μ _ik corresponding to the substance iron, aluminum, and acrylic resin is calculated using the equation (9). Next, in step 107, the continuous spectrum projection data f is calculated using the length L _k described above, the energy E _i of the X-ray continuous spectrum in FIG. 5 and the X-ray absorption coefficient μ _ik obtained in step 203, then, in step 109, calculates the first correction amount h _1, the correction amount h ₁ in correcting a single-spectrum projection data g, corrected CT image from the projection data g 'shown in (B) of FIG. 34 Reconfigure. Thus, image artifacts are greatly reduced by the second correction.

FIG. 35 is a flowchart for explaining a third embodiment of the X-ray CT image reconstruction method according to the present invention, in which the number of substances k constituting the object to be imaged 3 is large and the X of these substances k. This is applied when the linear absorption coefficient μ _ik is unknown.

In FIG. 35, step 202 is replaced with step 202A in FIG. That is, in step 201A, a _k in the flowchart of FIG.
a _k = 9.873 ・ TH _n + 9873
= 9.873 ・ (TH _n + 1000)
Where TH _n (n = 1,2, ...) is a provisional threshold,
Calculate by Here, values 9.873 and 9873 are linearly defined so that a (air) = 0 when the CT value of air is −1000, and a (iron) = 134000 when the CT value of iron is 12572. It is. Here, the iron CT value 12572 is the maximum CT value of an actual CT image. However, 9.873 and 9873 can be set to other values α _k and β _k so that two given points are linear. Also in this case, the substance k can be estimated based on the obtained a _k .

  A simulation example according to the flowchart of FIG. 35 will be described. In other words, in this case, a virtual projection method in which CT imaging is performed on the control circuit (computer) 4 is used, and the X-ray continuous spectrum of the X-ray generator 1 is approximated by nine discrete values shown in FIG.

The first simulation example is applied to a plurality of same pin models shown in FIG. 13 using two substances, namely iron and acrylic resin. However, the X-ray absorption coefficient μ _ik of the materials iron and acrylic resin is unknown.

In the first simulation example, in step 201, the energy E _i of the X-ray continuous spectrum of FIG. 5 is stored. Then, in step 103, single spectrum projection data g is generated based on these, and in step 105, a CT image shown in FIG. 36A is reconstructed from the single spectrum projection data g.

Next, the first correction is performed on the CT image. That is, in step 106, the regions of the substance iron and the acrylic resin are extracted using the known threshold value TH1 such as 1300 and threshold value TH2 such as 900 for the pixels of the CT image. The length L _k is calculated for each X-ray path of these extraction regions. Next, in step 202A, a _k is calculated. as a result,
a (equivalent to iron) = 128349
a (acrylic equivalent) = 8885
was gotten. Next, in step 203, the X-ray absorption coefficient μ _ik corresponding to the substance iron and the acrylic resin is calculated using the equation (9). Next, in step 107, the continuous spectrum projection data f is calculated using the length L _k described above, the energy E _i of the X-ray continuous spectrum in FIG. 5 and the X-ray absorption coefficient μ _ik obtained in step 203, then, in step 109, calculates the first correction amount h _1, the correction amount to correct a single-spectrum projection data g by h _1, the corrected CT image from the projection data g 'shown in (B) of FIG. 36 Reconfigure. Similarly, when the second correction, the third correction, the fourth correction, and the fifth correction are performed, CT images shown in (C), (D), (E), and (F) of FIG. 36 are obtained. Although the number of calculations increases, image artifacts gradually decrease.

The second simulation example is applied to a plurality of different pin models shown in FIG. 18 using three substances, that is, iron, aluminum, and acrylic resin. However, the X-ray absorption coefficient μ _ik of the substances iron, aluminum, and acrylic resin is unknown.

In the second simulation example, in step 201, the energy E _i of the X-ray continuous spectrum of FIG. 5 is stored. Then, in step 103, single spectrum projection data g is generated based on these, and in step 105, a CT image shown in FIG. 37A is reconstructed from the single spectrum projection data g.

Next, the first correction is performed on the CT image. That is, in step 106, the regions of the substance iron and the acrylic resin are extracted using the known threshold value TH1 such as 1300 and threshold value TH2 such as 900 for the pixels of the CT image. The length L _k is calculated for each X-ray path of these extraction regions. Next, in step 202A, a _k is calculated. as a result,
a (equivalent to iron) = 128349
a (equivalent to aluminum) = 8885
a (acrylic equivalent) = 8885
was gotten. Next, in step 203, the X-ray absorption coefficient μ _ik corresponding to the substance iron and the acrylic resin is calculated using the equation (9). Next, in step 107, the continuous spectrum projection data f is calculated using the length L _k described above, the energy E _i of the X-ray continuous spectrum in FIG. 5 and the X-ray absorption coefficient μ _ik obtained in step 203, then, in step 109, calculates the first correction amount h _1, the correction amount h ₁ in correcting a single-spectrum projection data g, corrected CT image from the projection data g 'shown in (B) of FIG. 37 Reconfigure. Similarly, when the second correction, the third correction, the fourth correction, and the fifth correction are performed, CT images shown in (C), (D), and (E) of FIG. 37 are obtained. The number of calculations increases, but again the image artifacts gradually decrease.

The third simulation example is applied to the biological tissue model shown in FIG. However, the X-ray absorption coefficient μ _ik of the living tissue is unknown.

In the third simulation example, in step 201, the energy E _i of the X-ray continuous spectrum of FIG. 5 is stored. Then, in step 103, single spectrum projection data g is generated based on these, and in step 105, a CT image shown in FIG. 39A is reconstructed from the single spectrum projection data g.

Next, the first correction is performed on the CT image shown in FIG. That is, in step 106, the regions of the substance iron and the acrylic resin are extracted using the known threshold value TH1 such as 1300 and threshold value TH2 such as 900 for the pixels of the CT image. The length L _k is calculated for each X-ray path of these extraction regions. Next, in step 202A, a _k is calculated. as a result,
a (equivalent to iron) = 128349
a (acrylic equivalent) = 8885
was gotten. Next, in step 203, the X-ray absorption coefficient μ _ik corresponding to the substance iron and the acrylic resin is calculated using the equation (9). Next, in step 107, the continuous spectrum projection data f is calculated using the length L _k described above, the energy E _i of the X-ray continuous spectrum in FIG. 5 and the X-ray absorption coefficient μ _ik obtained in step 203, then, in step 109, it calculates the first correction amount h _1, the correction amount h ₁ in correcting a single-spectrum projection data g, to reconstruct the CT image from the corrected projection data g '. Similarly, when correction is performed a plurality of times, a CT image shown in FIG. 39B is obtained. The number of calculations increases, but again the image artifacts gradually decrease.

1: X-ray generator 2: X-ray detector 3: Target object 4: Control circuit
I _in : X-ray energy at incidence
I _out : X-ray energy after transmission μ _ik : X-ray absorption coefficient
g: Single spectrum projection data
f: Continuous spectrum projection data
h _t : Correction amount

Claims (16)

Storing a generated X-ray continuous energy spectrum as known information;
An X-ray projection data generation stage for generating X-ray projection data of a single energy spectrum of an imaging object;
An initialization stage for initializing the correction amount;
The X-ray projection data of the single energy spectrum is corrected with the correction amount, and the single energy is applied to the imaging object in consideration of the energy dependency of the X-ray absorption coefficients of at least two substances contained in the imaging object. A correction stage for obtaining corrected X-ray projection data corresponding to the case of irradiating spectral X-rays;
A CT image reconstruction step for reprojecting the corrected X-ray projection data using a single wavelength image reconstruction algorithm to reconstruct a CT image;
A length calculating step of extracting an area of each substance from the CT image using at least one threshold value for distinguishing each substance and calculating a length of an X-ray path of the area;
An X-ray projection data calculation step of calculating X-ray projection data of the continuous energy spectrum using the energy of the generated X-ray continuous energy spectrum, the length of the region of each substance and the X-ray absorption coefficient of each substance;
An X-ray CT image reconstruction method comprising: a correction amount calculation step of calculating the correction amount based on a difference between the X-ray projection data of the single energy spectrum and the X-ray projection data of the continuous energy spectrum.   The substance is known, and the correction amount converges to a predetermined range by repeating the correction step, the CT image reconstruction step, the length calculation step, the X-ray projection data calculation step, and the correction amount calculation step. The X-ray CT image reconstruction method according to claim 1, further comprising the step of:   The X-ray CT image reconstruction method according to claim 1, further comprising a threshold estimation step of estimating the threshold from a histogram of the CT image, wherein the substance is unknown.   Further, the correction amount converges to a predetermined range by repeating the correction step, the CT image reconstruction step, the threshold value estimation step, the length calculation step, the X-ray projection data calculation step, and the correction amount calculation step. The X-ray CT image reconstruction method according to claim 3, further comprising the step of:   The X-ray CT image reconstruction method according to claim 1, wherein an X-ray absorption coefficient of each substance is set in advance. further,
The coefficient a _k of each substance k
A step of calculating using
Using the calculated coefficient a _k and the energy E _i of the generated X-ray continuous energy spectrum, the X-ray absorption coefficient μ _ik of the substance is obtained.
The X-ray CT image reconstruction method according to claim 1, further comprising the step of: further,
The coefficient a _k of each substance k
a _k = α _k・ TH _n (n = 1,2,…) + β
Where TH _n is the threshold value, α _k and β _k are calculated by constants;
Using the calculated coefficient a _k and the energy E _i of the generated continuous X-ray energy spectrum, the X-ray absorption coefficient μ _ik of each substance k is calculated.
The X-ray CT image reconstruction method according to claim 1, further comprising the step of: The X-ray CT image reconstruction method according to claim 6 or 7, further comprising a step of estimating the substance _k by the coefficient a _k . A storage procedure for storing the generated X-ray continuous energy spectrum as known information;
X-ray projection data generation procedure for generating X-ray projection data of a single energy spectrum of an imaging object;
An initialization procedure for initializing the correction amount;
The X-ray projection data of the single energy spectrum is corrected with the correction amount, and the single energy is applied to the imaging object in consideration of the energy dependency of the X-ray absorption coefficients of at least two substances contained in the imaging object. A correction procedure for obtaining corrected X-ray projection data corresponding to the case of irradiating spectral X-rays;
A CT image reconstruction procedure for reprojecting the corrected X-ray projection data using a single wavelength image reconstruction algorithm to reconstruct a CT image;
A length calculation procedure for extracting an area of each substance from the CT image using at least one threshold value for distinguishing each substance and calculating a length of an X-ray path of the area;
X-ray projection data calculation procedure for calculating the X-ray projection data of the continuous energy spectrum using the energy of the generated X-ray continuous energy spectrum, the length of each substance region, and the X-ray absorption coefficient of each substance;
An X-ray CT image reconstruction program comprising: a correction amount calculation procedure for calculating the correction amount based on a difference between the X-ray projection data of the single energy spectrum and the X-ray projection data of the continuous energy spectrum.   The substance is known, and the correction amount converges to a predetermined range by repeating the correction procedure, the CT image reconstruction procedure, the length calculation procedure, the X-ray projection data calculation procedure, and the correction amount calculation procedure. The X-ray CT image reconstruction program according to claim 9, further comprising a procedure for performing the operation.   The X-ray CT image reconstruction program according to claim 9, further comprising a threshold estimation procedure for estimating the threshold from a histogram of the CT image when the substance is unknown.   Further, the correction amount converges to a predetermined range by repeating the correction procedure, the CT image reconstruction procedure, the threshold value estimation procedure, the length calculation procedure, the X-ray projection data calculation procedure, and the correction amount calculation procedure. The X-ray CT image reconstruction program according to claim 11, further comprising: a procedure for performing the operation.   The X-ray CT image reconstruction program according to claim 9, wherein an X-ray absorption coefficient of each substance is set in advance. further,
The coefficient a _{k of the} substance k
A procedure for computing using
Using the calculated coefficient a _k and the energy E _i of the generated X-ray continuous energy spectrum, the X-ray absorption coefficient μ _ik of the substance k is obtained.
The X-ray CT image reconstruction program according to claim 9, further comprising: further,
The coefficient a _{k of the} substance k
a _k = α _k・ TH _n (n = 1,2,…) + β _k
However, TH _n is the threshold value, α _k , β _k are a procedure to calculate by a constant,
Using the calculated coefficient a _k and the energy E _i of the generated continuous X-ray energy spectrum, the X-ray absorption coefficient μ _ik of each substance k is calculated.
The X-ray CT image reconstruction program according to claim 7, further comprising: The X-ray CT image reconstruction program according to claim 14, further comprising a step of estimating the substance _k by the coefficient a _k .