CN109961411B - Non-downsampling shear wave transformation medical CT image denoising method - Google Patents
Non-downsampling shear wave transformation medical CT image denoising method Download PDFInfo
- Publication number
- CN109961411B CN109961411B CN201910186481.2A CN201910186481A CN109961411B CN 109961411 B CN109961411 B CN 109961411B CN 201910186481 A CN201910186481 A CN 201910186481A CN 109961411 B CN109961411 B CN 109961411B
- Authority
- CN
- China
- Prior art keywords
- image
- decomposition
- noise
- geometric
- low
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 71
- 230000009466 transformation Effects 0.000 title claims abstract description 27
- 238000000354 decomposition reaction Methods 0.000 claims abstract description 73
- 238000002591 computed tomography Methods 0.000 claims abstract description 66
- 238000012545 processing Methods 0.000 claims abstract description 26
- 230000009467 reduction Effects 0.000 claims abstract description 26
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 21
- 238000001914 filtration Methods 0.000 claims abstract description 17
- 230000000694 effects Effects 0.000 claims description 21
- 230000008569 process Effects 0.000 claims description 18
- 230000003044 adaptive effect Effects 0.000 claims description 15
- 239000000654 additive Substances 0.000 claims description 15
- 230000000996 additive effect Effects 0.000 claims description 15
- 230000002146 bilateral effect Effects 0.000 claims description 9
- 238000004364 calculation method Methods 0.000 claims description 9
- 230000008901 benefit Effects 0.000 claims description 7
- 230000010355 oscillation Effects 0.000 claims description 6
- 230000003534 oscillatory effect Effects 0.000 claims description 6
- 238000000926 separation method Methods 0.000 claims description 6
- 230000008602 contraction Effects 0.000 claims description 5
- 230000007547 defect Effects 0.000 claims description 5
- 238000009826 distribution Methods 0.000 claims description 5
- 230000004927 fusion Effects 0.000 claims description 3
- 230000000717 retained effect Effects 0.000 claims description 3
- 238000005303 weighing Methods 0.000 claims description 3
- 238000010008 shearing Methods 0.000 claims 2
- 239000000284 extract Substances 0.000 claims 1
- 238000003745 diagnosis Methods 0.000 abstract description 5
- 238000004458 analytical method Methods 0.000 abstract description 4
- 230000009286 beneficial effect Effects 0.000 abstract 1
- 238000010586 diagram Methods 0.000 description 7
- 238000013170 computed tomography imaging Methods 0.000 description 6
- 238000005516 engineering process Methods 0.000 description 5
- 238000003384 imaging method Methods 0.000 description 4
- 230000003902 lesion Effects 0.000 description 4
- 238000011961 computed axial tomography Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 238000000605 extraction Methods 0.000 description 2
- 230000008859 change Effects 0.000 description 1
- 238000003759 clinical diagnosis Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000002059 diagnostic imaging Methods 0.000 description 1
- 201000010099 disease Diseases 0.000 description 1
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 238000002595 magnetic resonance imaging Methods 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- 210000000056 organ Anatomy 0.000 description 1
- 230000005855 radiation Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000012285 ultrasound imaging Methods 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
- G16H50/00—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics
- G16H50/20—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics for computer-aided diagnosis, e.g. based on medical expert systems
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10081—Computed x-ray tomography [CT]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20024—Filtering details
Landscapes
- Engineering & Computer Science (AREA)
- Health & Medical Sciences (AREA)
- Medical Informatics (AREA)
- Biomedical Technology (AREA)
- Public Health (AREA)
- Epidemiology (AREA)
- Pathology (AREA)
- Databases & Information Systems (AREA)
- Data Mining & Analysis (AREA)
- General Health & Medical Sciences (AREA)
- Primary Health Care (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Image Processing (AREA)
- Apparatus For Radiation Diagnosis (AREA)
Abstract
The non-subsampled shear wave transformation medical CT image denoising method comprises the following steps: step 1) establishing a medical CT scanning image model; step 2) carrying out NSST multi-scale decomposition and multi-directional decomposition on the CT image to obtain a low-frequency sub-band and a plurality of high-frequency sub-bands; step 3) decomposing the original image into a smooth geometric part and a texture part containing noise by using rapid geometric texture decomposition, and fusing components of the geometric part with low-frequency sub-bands; step 4), carrying out noise reduction treatment on the fused new low-frequency component by adopting a trilateral filtering method to obtain a new low-frequency sub-band; step 5) processing the high-frequency sub-band coefficient after shear wave transformation by adopting a self-adaptive threshold shrinkage method; and 6) combining the processed low-frequency component and the processed high-frequency component to carry out NSST reconstruction to obtain the denoised medical CT image. Compared with the traditional denoising field algorithm, the invention is effectively applied to the medical CT denoising field and can be better beneficial to the analysis and diagnosis of doctors.
Description
Technical Field
The invention relates to a non-downsampling shear wave medical CT image denoising method.
Background
With the development of science and technology, in the field of medical imaging, imaging technologies such as ultrasound imaging, ct (computed tomography), mri (magnetic Resonance imaging), and the like are generally applied to medical clinical diagnosis. Computed Tomography (CT) is an examination technique for imaging diagnostics. This technique has been referred to as Computed Axial Tomography (Computed Axial Tomography). Computed tomography, using a computer to process a number of scanned objects of a particular area across a cross-section resulting from combined X-ray measurements from different directional angles, allows a physician to see the interior of the diagnosed object. The CT imaging technology is used for detecting cross section imaging, tissues or organs in any directions can be displayed through image reconstruction, lesion parts can be displayed more comprehensively, and omission can be effectively prevented. CT imaging techniques have high density resolution and can show changes in the density of subtle lesions, thereby defining the nature of the lesions.
Due to the mechanism of CT imaging, even the CT values of each pixel in an image of an object with uniform internal distribution are uneven, and the image has graininess and affects the density resolution, that is, CT noise. There are various causes of CT noise, among which are detector causes, such as: the sensitivity of the detector, the pixel size, the scanning layer thickness, the amount of X-rays, etc. There are mechanical and electronic reasons, and CT noise is also caused by the method of image reconstruction and scattered radiation. Since the quality of CT imaging is better when the CT noise is less, and the quality of CT imaging is worse when the CT noise is more, the observation is affected, it is necessary to know the cause of the generation of the CT noise and suppress the generation as much as possible. The present invention uses medical CT images as the object of study. CT imaging is susceptible to a number of objective factors that produce speckle noise, which can degrade the imaging quality of CT images, resulting in poor quality medical diagnostic images. From the viewpoint of observation, CT noise is small spots with different correlation shapes in a spatial domain, and is called speckle noise, and the existence of the speckle noise affects the medical image features with small gray scale differences. The presence of such speckle noise can therefore greatly affect the accurate diagnosis of lesions by physicians, particularly young physicians who are not yet highly experienced, in the case of clinical diagnosticians. Therefore, considering the requirement of practical clinical application, there is a need to develop a method for reducing speckle noise of medical CT images, which helps the diagnostician from the technical level to ensure more accurate diagnosis of the disease. Therefore, it can be very meaningful to study the noise reduction method of the medical CT image.
Disclosure of Invention
The invention aims to overcome the defects in the prior art and provides a non-downsampling shear wave medical CT image denoising method based on rapid geometric texture decomposition and trilateral filtering
The invention discloses a medical CT image denoising algorithm based on image non-subsampled shear wave transformation, which is used for decomposing an image and processing coefficients on a decomposition coefficient layer so as to realize a method for reducing CT noise. In the image decomposition transformation technology adopted in the prior art, the wavelet transformation can well process one-dimensional singular points, but cannot effectively deal with the problem of abrupt change of straight lines and curves. The ridgelet transform overcomes this effectively in capturing the singularity of the line relative to the wavelet, but is only very effective in capturing the singularity of the curve, but not very effective in capturing the line. The non-subsampled shear wave applied by the invention is an improvement on the basis of shear wave denoising in recent years, and the medical image denoising technology has great span from the proposal of shear wave to the application of a discrete shear wave algorithm on a medical image. The shear wave transformation has many advantages, the multi-resolution, the directivity, the locality and the anisotropy are the most sparse representation modes of the image, and the characteristics of the shear wave transformation have very effective effects when being applied to the field of medical denoising. In addition, the invention also fuses the geometric part obtained by decomposition and the low-frequency distribution image obtained by decomposing the CT image by non-down-sampling shear wave transformation in the medical CT image decomposition by using the rapid geometric texture image decomposition method, generates a new low-frequency part and further carries out the next processing. The invention uses a novel trilateral filter with very strong denoising capability to process denoising in the processing of a new low-frequency part, invents a better denoising method for the non-subsampled shear wave medical CT image, and finally verifies the feasibility and the optimization effect of the method through simulation.
Compared with the prior art, the invention has the creativity and advantages that: the invention applies the rapid geometric texture image decomposition method to the decomposition of the de-noised image. And then fusing the geometric part obtained by decomposition and the low-frequency component obtained by non-subsampled shear wave decomposition to form a new low-frequency component. And reducing the noise of a low-frequency part contained in the CT image by using a trilateral filter denoising method. The method overcomes the defect that wavelet transformation can only process one-dimensional singularity and the defect that ridge transformation can only effectively solve straight line but not curve singularity by utilizing the characteristic of non-subsampled shear wave, and makes up the bad defect of gradient inversion of shear wave. The application of the invention can bring more convenient and more accurate diagnosis to the diagnosis of the clinician.
In order to make the objects, technical solutions and advantages of the present invention clearer, the following further describes the technical solutions of the present invention, and the method for denoising a non-downsampling shear wave transform medical CT image of the present invention includes the following steps:
step 1) establishing a medical CT scanning image model;
the medical CT image is computed tomography, X-rays scan the fixed part of a human body from a plurality of different directions and angles and send the fixed part to different cross sections obtained by computer processing scanning to establish images, so that doctors, patients and the like can see the scanning objects in a specific examination area, and further medical judgment is carried out. However, the emission current with too low intensity may generate a large amount of gaussian noise, so that the image quality is reduced, and the observation and determination result may be affected.
The model of the CT image can be built in two parts, which are respectively the human tissue reflection signal required for medical observation and the noise signal obstructing the medical observation, wherein the noise signal can also be divided into multiplicative noise and additive noise, and the additive noise has very little influence on the CT image compared with the multiplicative noise from the viewpoint of influencing the observation, so the additive noise is generally ignored in the process. The general model of the CT electrical signal is thus expressed as:
o(x，y)＝p(x，y)q(x，y) (1)
in the formula, x represents the abscissa of the CT image, y represents the ordinate of the image, p (x, y) represents a noise-free signal, and q (x, y) represents multiplicative noise.
Since the additive noise model is easier to separate than the model multiplied by the noise, the model of equation (1) above is logarithmically transformed into an additive model, which is expressed as:
log(o(x，y))＝log(p(x，y))+log(q(x，y)) (2)
step 2) performing NSST multi-scale decomposition and multi-direction decomposition on the CT image;
firstly, the logarithmically transformed CT image which is convenient for noise separation is subjected to multi-scale decomposition, and one low-frequency CT image component with the same size as the original image and a plurality of high-frequency CT image components with the same size as the original image are obtained after the multi-scale decomposition. The low frequency components are not processed in this step, which is done in step 3), and the high frequency component subbands resulting from the scaling are processed using a shear filter bank, i.e. each subband is directionally decomposed using a shear filter bank.
Step 3) decomposing the original image into a smooth geometric part and a texture part containing noise by using rapid geometric texture decomposition;
fast geometric texture decomposition decomposes the original image into smooth geometric parts and texture parts containing noise. Fast geometric texture decomposition allows efficient texture extraction by using a high-pass filter while preserving the main features of the image by creating a non-linear filter containing local indicators to determine whether the image is local to the texture or the geometric part. The main differences between the texture part and the geometric part are: the main feature of the texture region is the high total variation due to its oscillation, whereas the main feature of the geometric local region is that its total variation is not affected by the low-pass filtering. Its local total variation can be expressed as the following equation:
LTVσ(f)(x):＝Lσ*|Df|(x) (3)
wherein σ represents texture scale, | Df | represents variation, L σ | Df | is variation of local range, and its relative local reduction rate in the process of noise reduction using local total variation can be represented as:
the meaning reflected on the image can be understood as the local oscillation behavior of the function, when the reduction ratio is close to 0:
that is, the low-pass filter has very little effect on the local total variation reduction when the reduction ratio is close to 0, and if the reduction ratio is close to 1, the currently considered point should be attributed to the texture part, so that a fast low-pass filter and high-pass filter pair can be calculated from the relative reduction of the local total variation using a weighted average of f and L σ f. The calculation method can be expressed as:
u(x)＝w(λσ(x))(Lσ*f)(x)+(1-w(λσ(x)))f(x) (6)
v(x)＝f(x)-u(x) (7)
wherein the threshold function w (x) is an increasing function whose value is constantly equal to 0 around 0 and constantly equal to 1 around 1, the calculation of the geometric decomposition part is u (x), and the texture decomposition part is v (x), where the definition of the threshold function can be represented by the following piecewise function:
wherein the parameter a1And a2The size of (a) varies from the separation result to 0.25 and 0.5, respectively, and if λ σ (x) is very small, the function f is non-oscillatory around x, so that f (x) can be attributed to u (x), i.e. to the geometric part, around this range of x values. Conversely, if λ σ (x) is very large, then the function f is locally oscillatory near x, so that f (x) can be attributed to v (x), i.e., texture portion, near this range of x values.
In this step, the texture part resulting from the decomposition is discarded, and only the geometric part resulting from the decomposition is taken as the geometric component involved in the calculation. Respectively calculating the entropy of the geometric component and the two-dimensional image of the low-frequency component generated by the image in the step 2) through non-subsampled shear wave transformation, performing fusion operation according to the entropy, and combining the processed image with the high-frequency part processed in the same way after the next step of processing.
Step 4), carrying out noise reduction treatment on the fused new low-frequency component by adopting a trilateral filtering method;
the new low-frequency component obtained by fusing the low-frequency components subjected to the geometric component and the non-subsampled shear wave decomposition still contains certain noise, so that a further denoising process is necessary before the denoising image is reconstructed with the high-frequency coefficient subjected to the adaptive threshold shrinkage process. The bilateral filtering method that has been commonly used can effectively remove noise in the spatial domain and the value domain and can retain the edge information of the image to a certain extent, but has the disadvantage that the gradient distortion phenomenon often occurs. Therefore, the trilateral filtering method is selected to be used for denoising. The trilateral filter is characterized in that a new weight is added on the basis of an original bilateral filter, the original bilateral filter only performs weighing processing on distance weight on an airspace level and gray scale weight on a value domain level, and the weight added by the trilateral filter is quality weight which represents image edge information to a certain extent, so that the phenomenon of gradient distortion can be effectively reduced by adding the quality weight, and the integral denoising effect is improved. The trilateral filter can be expressed as follows:
wherein m isS(x, xi) denotes a spatial filter, mR(x, xi) denotes a value domain filter, mIAnd (ξ) represents the quality weight component. The spatial filter part considers the influence of pixel points near the neighborhood of the current pixel point on the gray value of the current pixel point, the closer the weight is to the Gaussian distance, the larger the weight occupied by the weighted average is, the value domain filter considers the weighted average of the points close to the gray value of the current pixel point, the occupied weight is related to the gray difference of the pixel point, the larger the gray difference is, the smaller the weight occupied by the weighted average is, otherwise, the larger the weight occupation ratio is.
Step 5) processing the high-frequency sub-band coefficient after shear wave transformation by adopting a self-adaptive threshold shrinkage method;
in the traditional filtering process, soft threshold method and hard threshold method are generally adopted for decomposing coefficients, but when the decomposing method is shear wave transformation with multi-scale and multi-directivity, other methods are required for processing the coefficients, and an adaptive threshold shrinking algorithm is adopted. The advantage of choosing to use adaptive threshold puncturing is to avoid that too much noise is retained while ensuring that no useful information is filtered because the threshold is chosen too small, and to avoid that some useful detail information is filtered because the threshold is too large. Therefore, the invention selects the following self-adaptive threshold function to perform contraction processing on the high-frequency sub-band coefficient obtained by the shear wave transformation decomposition:
wherein σnIs the standard deviation of the noise, M corresponds to the total number of transform coefficients in the transform domain, tjThe adaptive parameters representing the j layers have different distributions of coefficients of different decomposition layers after the shear wave decomposition, so that tjAnd selecting based on different j layers and specific experimental targets. In the present invention, the CT graph of 512X 512 is decomposed into a low frequency part and 4 high frequency parts, and the threshold value T is approximately equal to 5 sigmanFor more detailed decomposition subbands T ≈ 4 σ is usednOr T ≈ 3 σn。
Step 6) combining the processed low-frequency component and the high-frequency component to carry out NSST reconstruction;
the low-frequency component of the denoised shear wave coefficient is processed by a trilateral filter, the high-frequency component is obtained by adaptive threshold shrinkage, and finally, the denoised CT image which can enable a doctor to diagnose a focus more conveniently and more accurately is obtained by combining the high-frequency part and the low-frequency part and performing NSST reconstruction.
The invention has the following advantages:
1. the method decomposes the fast geometric texture decomposition CT denoised image, fuses the geometric components obtained by decomposition and the low-frequency components obtained by non-downsampling shear wave decomposition to obtain new low-frequency components.
2. The invention reduces the noise of the low-frequency part contained in the CT image by using a trilateral filter denoising method.
3. The invention adopts a self-adaptive threshold contraction algorithm to effectively process the coefficients of the molecular band of the high-frequency part of the shear wave.
4. The method has definite steps and simple structure and has perfect theoretical support.
Drawings
Fig. 1 is a fused image of a geometric portion of a head CT image obtained by fast geometric texture decomposition and a low-frequency portion obtained by non-downsampling shear wave decomposition.
FIG. 2 is a flow chart of the overall steps of the present invention.
Fig. 3 is an overall flow of case analysis.
Fig. 4a to 4e are experimental effect diagrams of various algorithms applied to a lens diagram (σ ═ 40), and an original and noise diagram, fig. 4a is the original, fig. 4b is the noise diagram, fig. 4c is the FDCT algorithm effect diagram, fig. 4d is the FFST algorithm effect diagram, and fig. 4e is the algorithm effect diagram of the present invention.
Fig. 5a to 5d are an experimental effect graph and a noise graph of various algorithms applied to a head CT graph (σ ═ 40), fig. 5a is a noise graph, fig. 5b is an FDCT algorithm effect graph, fig. 5c is an FFST algorithm effect graph, and fig. 5d is an algorithm effect graph of the present invention.
Detailed Description
The invention is further described below with reference to the accompanying drawings.
The invention discloses a denoising method for a non-subsampled shear wave transformation medical CT image, which comprises the following steps of:
step 1) establishing a medical CT scanning image model;
the medical CT image is computed tomography, X-rays scan the fixed part of a human body from a plurality of different directions and angles and send the fixed part to different cross sections obtained by computer processing scanning to establish images, so that doctors, patients and the like can see the scanning objects in a specific examination area, and further medical judgment is carried out. However, the emission current with too low intensity may generate a large amount of gaussian noise, so that the image quality is reduced, and the observation and determination result may be affected.
The model of the CT image can be built in two parts, which are respectively the human tissue reflection signal required for medical observation and the noise signal obstructing the medical observation, wherein the noise signal can also be divided into multiplicative noise and additive noise, and the additive noise has very little influence on the CT image compared with the multiplicative noise from the viewpoint of influencing the observation, so the additive noise is generally ignored in the process. The general model of the CT electrical signal is thus expressed as:
o(x，y)＝p(x，y)q(x，y) (1)
in the formula, x represents the abscissa of the CT image, y represents the ordinate of the image, p (x, y) represents a noise-free signal, and q (x, y) represents multiplicative noise. Since the additive noise model is easier to separate than the model multiplied by the noise, the model of equation (1) above is logarithmically transformed into an additive model, which is expressed as:
log(o(x，y))＝log(p(x，y))+log(q(x，y)) (2)
step 2) performing NSST multi-scale decomposition and multi-direction decomposition on the CT image;
firstly, the logarithmically transformed CT image which is convenient for noise separation is subjected to multi-scale decomposition, and one low-frequency CT image component with the same size as the original image and a plurality of high-frequency CT image components with the same size as the original image are obtained after the multi-scale decomposition. The low frequency components are not processed in this step, which is done in step 3), and the high frequency component subbands resulting from the scaling are processed using a shear filter bank, i.e. each subband is directionally decomposed using a shear filter bank.
Step 3) decomposing the original image into a smooth geometric part and a texture part containing noise by using rapid geometric texture decomposition;
fast geometric texture decomposition decomposes the original image into smooth geometric parts and texture parts containing noise. Fast geometric texture decomposition allows efficient texture extraction by using a high-pass filter while preserving the main features of the image by creating a non-linear filter containing local indicators to determine whether the image is local to the texture or the geometric part. The main differences between the texture part and the geometric part are: the main feature of the texture region is the high total variation due to its oscillation, whereas the main feature of the geometric local region is that its total variation is not affected by the low-pass filtering. Its local total variation can be expressed as the following equation:
LTVσ(f)(x):＝Lσ*|Df|(x) (3)
wherein σ represents texture scale, | Df | represents variation, L σ | Df | is variation of local range, and its relative local reduction rate in the process of noise reduction using local total variation can be represented as:
the meaning reflected on the image can be understood as the local oscillation behavior of the function, when the reduction ratio is close to 0:
that is, the low-pass filter has very little effect on the local total variation reduction when the reduction ratio is close to 0, and if the reduction ratio is close to 1, the currently considered point should be attributed to the texture part, so that a fast low-pass filter and high-pass filter pair can be calculated from the relative reduction of the local total variation using a weighted average of f and L σ f. The calculation method can be expressed as:
u(x)＝w(λσ(x))(Lσ*f)(x)+(1-w(λσ(x)))f(x) (6)
v(x)＝f(x)-u(x) (7)
wherein the threshold function w (x) is an increasing function whose value is constantly equal to 0 around 0 and constantly equal to 1 around 1, the calculation of the geometric decomposition part is u (x), and the texture decomposition part is v (x), where the definition of the threshold function can be represented by the following piecewise function:
wherein the parameter a1And a2The size of (d) is chosen to vary the separation result, here 0.25 and 0.5, respectively, if λ σ (x) is very small, this functionThe number f is non-oscillatory around x, so that f (x) can be classified as u (x), i.e. as a geometric part, around the value range of x. Conversely, if λ σ (x) is very large, then the function f is locally oscillatory near x, so that f (x) can be attributed to v (x), i.e., texture portion, near this range of x values.
In this step, the texture part resulting from the decomposition is discarded, and only the geometric part resulting from the decomposition is taken as the geometric component involved in the calculation. Respectively calculating the entropy of the geometric component and the two-dimensional image of the low-frequency component generated by the image in the step 2) through non-subsampled shear wave transformation, performing fusion operation according to the entropy, and combining the processed image with the high-frequency part processed in the same way after the next step of processing. Fig. 1 is a fused image of a geometric portion of a head CT image subjected to a fast geometric texture decomposition and a low frequency portion subjected to a non-downsampling shear wave decomposition process.
Step 4), carrying out noise reduction treatment on the fused new low-frequency component by adopting a trilateral filtering method;
the new low-frequency component obtained by fusing the low-frequency components subjected to the geometric component and the non-subsampled shear wave decomposition still contains certain noise, so that a further denoising process is necessary before the denoising image is reconstructed with the high-frequency coefficient subjected to the adaptive threshold shrinkage process. The bilateral filtering method that has been commonly used can effectively remove noise in the spatial domain and the value domain and can retain the edge information of the image to a certain extent, but has the disadvantage that the gradient distortion phenomenon often occurs. Therefore, the trilateral filtering method is selected to be used for denoising. The trilateral filter is characterized in that a new weight is added on the basis of an original bilateral filter, the original bilateral filter only performs weighing processing on distance weight on an airspace level and gray scale weight on a value domain level, and the weight added by the trilateral filter is quality weight which represents image edge information to a certain extent, so that the phenomenon of gradient distortion can be effectively reduced by adding the quality weight, and the integral denoising effect is improved. The trilateral filter can be expressed as follows:
wherein m isS(x, xi) denotes a spatial filter, mR(x, xi) denotes a value domain filter, mIAnd (ξ) represents the quality weight component. The spatial filter part considers the influence of pixel points near the neighborhood of the current pixel point on the gray value of the current pixel point, the closer the weight is to the Gaussian distance, the larger the weight occupied by the weighted average is, the value domain filter considers the weighted average of the points close to the gray value of the current pixel point, the occupied weight is related to the gray difference of the pixel point, the larger the gray difference is, the smaller the weight occupied by the weighted average is, otherwise, the larger the weight occupation ratio is.
Step 5) processing the high-frequency sub-band coefficient after shear wave transformation by adopting a self-adaptive threshold shrinkage method;
in the traditional filtering process, soft threshold method and hard threshold method are generally adopted for decomposing coefficients, but when the decomposing method is shear wave transformation with multi-scale and multi-directivity, other methods are required for processing the coefficients, and an adaptive threshold shrinking algorithm is adopted. The advantage of choosing to use adaptive threshold puncturing is to avoid that too much noise is retained while ensuring that no useful information is filtered because the threshold is chosen too small, and to avoid that some useful detail information is filtered because the threshold is too large. Therefore, the invention selects the following self-adaptive threshold function to perform contraction processing on the high-frequency sub-band coefficient obtained by the shear wave transformation decomposition:
wherein σnIs the standard deviation of the noise, M corresponds to the total number of transform coefficients in the transform domain, tjThe adaptive parameters representing the j layers have different distributions of coefficients of different decomposition layers after the shear wave decomposition, so that tjAnd selecting based on different j layers and specific experimental targets. In the inventionIn 512 x 512, the CT map is decomposed into a low frequency part and 4 high frequency parts, and the threshold T ≈ 5 σnFor more detailed decomposition subbands T ≈ 4 σ is usednOr T ≈ 3 σn。
Step 6) combining the processed low-frequency component and the high-frequency component to carry out NSST reconstruction;
the low-frequency component of the denoised shear wave coefficient is processed by a trilateral filter, the high-frequency component is obtained by adaptive threshold shrinkage, and finally, the denoised CT image which can enable a doctor to diagnose a focus more conveniently and more accurately is obtained by combining the high-frequency part and the low-frequency part and performing NSST reconstruction.
The overall steps of the present invention are shown in flow chart form in fig. 2.
Case analysis
The invention takes a specific head medical CT image as an object, uses a non-downsampling shear wave medical CT image denoising method based on rapid geometric texture decomposition and trilateral filtering to perform a test, processes a high-frequency subband coefficient by using a self-adaptive threshold algorithm, takes a PSNR (peak signal-to-noise ratio) value as an evaluation index to perform verification (the larger the PSNR value is, the higher the denoising effect is), and simultaneously embodies the advancement of the invention by comparing with the prior art. The experiment is carried out by taking a classical Lena image (the image size is 512 multiplied by 512) of a medical CT noise image as input data and bringing the input data into different denoising methods, and the overall flow chart of case analysis is shown in FIG. 3. Experiments were performed by comparing FDCT (fast discrete curvelet transform), FFST (fast finite shear transform). Various algorithms are applied to the experimental effect map of the lens map and the original and noise maps as shown in fig. 4, and various algorithms are applied to the experimental effect map of the head CT map and the original and noise maps as shown in fig. 5.
As can be seen from tables 1 and 2, the higher the noise variance is, the higher the requirement on the denoising algorithm is, as can be seen from the experimental data of the classical image Lena and the medical CT image. The invention can keep relatively stable denoising effect under the condition of increasing noise variance. The effect is better than FDCT and FFST on the same noise variance.
Table 1: PSNR/dB values of different noise of different denoising algorithms of Lena image
Table 2: PSNR/dB value of different noise reduction algorithms of medical CT image
While the invention has been described in connection with specific embodiments thereof, it will be understood that these should not be construed as limiting the scope of the invention, which is defined in the following claims, and any variations which fall within the scope of the claims are intended to be embraced thereby.
Claims (1)
1. The non-subsampled shear wave transformation medical CT image denoising method comprises the following steps:
step 1) establishing a medical CT scanning image model;
the X-ray scans the fixed part of the human body from a plurality of different directions and angles and sends the fixed part to different cross sections obtained by computer processing scanning to establish images, so that doctors and patients can see the scanning object in a specific examination area to further carry out medical judgment; however, the emission current with too low intensity can generate a large amount of Gaussian noise, so that the image quality is reduced, and the observation and judgment result can be influenced;
the method comprises the steps of establishing a model of the CT image in two parts, wherein the two parts are human tissue reflection signals required by medical observation and noise signals obstructing the medical observation respectively, the noise signals are also divided into multiplicative noise and additive noise, and the additive noise has very little influence on the CT image compared with the multiplicative noise from the viewpoint of influencing the observation, so the additive noise is ignored in the processing; the general model of the CT electrical signal is thus expressed as:
o(x，y)＝p(x，y)q(x，y) (1)
wherein x represents the abscissa of the CT image, y represents the ordinate of the image, p (x, y) represents a noise-free signal, and q (x, y) represents multiplicative noise;
since the additive noise model is easier to separate than the model multiplied by the noise, the model of equation (1) above is logarithmically transformed into an additive model, which is expressed as:
log(o(x，y))＝log(p(x，y))+log(q(x，y)) (2)
step 2) performing NSST multi-scale decomposition and multi-direction decomposition on the CT image;
firstly, carrying out multi-scale decomposition on a CT image which is subjected to logarithmic transformation and is convenient for noise separation, and obtaining a low-frequency CT image component with the same size as an original image and a plurality of high-frequency CT image components with the same size as the original image after the multi-scale decomposition; the low-frequency component is not processed in the step, the processing is carried out in the step 3), and each high-frequency component subband obtained by the scale decomposition is processed by using a shearing filter bank, namely, each subband is subjected to direction decomposition by using a shearing filter bank;
step 3) decomposing the original image into a smooth geometric part and a texture part containing noise by using rapid geometric texture decomposition;
the fast geometric texture decomposition decomposes the original image into a smooth geometric part and a texture part containing noise; the fast geometric texture decomposition effectively extracts textures on the basis of retaining main characteristics of an image by using a high-pass filter, and determines whether a part of the image belongs to a texture part or a geometric part by establishing a nonlinear filter containing a local indicator; the main differences between the texture part and the geometric part are: the texture region is mainly characterized by high total variation caused by oscillation, and the geometric local region is mainly characterized in that the total variation is not influenced by low-pass filtering; the local total variation is expressed as the following formula:
LTVσ(f)(x)：＝Lσ*|Df|(x) (3)
wherein σ represents texture scale, | Df | represents variation difference, L σ | Df | is variation difference of local scope, and relative local reduction rate is represented as:
local oscillation behavior as a function of the sense reflected on the image, when the reduction ratio is close to 0:
that is, the low-pass filter has very little effect on the local total variation reduction when the reduction ratio is close to 0, and if the reduction ratio is close to 1, the currently considered point is attributed to the texture part, whereby a fast low-pass filter and high-pass filter pair is calculated from the relative reduction of the local total variation using a weighted average of f and L σ f; the calculation method is represented as follows:
u(x)＝w(λσ(x))(Lσ*f)(x)+(1-w(λσ(x)))f(x) (6)
v(x)＝f(x)-u(x) (7)
wherein the threshold function w (x) is an increasing function whose value is constantly equal to 0 around 0 and constantly equal to 1 around 1, the calculation of the geometric decomposition part is u (x), the texture decomposition part is v (x), and the definition of the threshold function is represented by the following piecewise function:
wherein the parameter a1And a2The size of (a) varies from the separation result to 0.25 and 0.5, respectively, and if λ σ (x) is very small, the function f is non-oscillatory around x, so that f (x) is reduced to u (x) around the value range of x, and is classified as a geometric part; conversely, if λ σ (x) is very large, then the function f is locally oscillatory near x, so that f (x) is reduced to v (x) near the range of x values, i.e., to texture;
in the step, texture parts obtained by decomposition are abandoned, and only the geometric parts generated by decomposition are taken as the geometric components participating in calculation; respectively calculating the entropy of the geometric component and the entropy of the two-dimensional image of the low-frequency component generated by the image in the step 2) through non-subsampled shear wave transformation, performing fusion operation according to the entropy, and combining the processed image with the high-frequency part processed in the same way after the next step of processing;
step 4), carrying out noise reduction treatment on the fused new low-frequency component by adopting a trilateral filtering method;
the new low-frequency component obtained by fusing the low-frequency components subjected to geometric component and non-subsampled shear wave decomposition still contains certain noise, so that the new low-frequency component is subjected to further denoising processing before reconstructing a denoised image with the high-frequency coefficient subjected to adaptive threshold shrinkage processing; the bilateral filtering method which is often used at one time can effectively remove noise in a space domain and a value domain angle and can keep the edge information of an image to a certain extent, but has the defect that the gradient distortion phenomenon often occurs; therefore, the trilateral filtering method is selected to be used for denoising; the trilateral filter is characterized in that a new weight is added on the basis of the original bilateral filter, the original bilateral filter only performs weighing processing on the distance weight on an airspace level and the gray scale weight on a value domain level, and the weight added by the trilateral filter is the quality weight which represents image edge information to a certain extent, so that the phenomenon of gradient distortion can be effectively reduced by adding the quality weight, and the integral denoising effect is improved; the trilateral filter is represented by the following equation:
wherein m isS(x, xi) denotes a spatial filter, mR(x, xi) denotes a value domain filter, m1(ξ) represents the quality weight component; the spatial filter part considers the influence of the pixel points near the neighborhood of the current pixel point on the gray value of the current pixel point, the closer the weight is to the Gaussian distance, the larger the weight occupied by the weighted average is, and the value domain filter considers the point with the gray value close to the gray value of the current pixel pointThe weight occupied by the weighted average is related to the gray level difference of the pixel points, the larger the gray level difference is, the smaller the weight occupied by the weighted average is, and otherwise, the larger the weight proportion is;
step 5) processing the high-frequency sub-band coefficient after shear wave transformation by adopting a self-adaptive threshold shrinkage method;
the traditional filter processing decomposition coefficient adopts a soft threshold method and a hard threshold method, but when the decomposition method is shear wave transformation with multi-scale and multi-directivity, the other method is required to process the coefficient, and an adaptive threshold contraction algorithm is adopted; the advantage of choosing to use adaptive threshold shrinking is to avoid that too much noise is retained while ensuring that no useful information is filtered because the threshold is chosen too small, and to avoid that some useful detail information is filtered because the threshold is too large; therefore, the following adaptive threshold function is selected to perform contraction processing on the high-frequency sub-band coefficient obtained by shear wave transformation decomposition:
wherein σnIs the standard deviation of the noise, M corresponds to the total number of transform coefficients in the transform domain, tjThe adaptive parameters representing the j layers have different distributions of coefficients of different decomposition layers after the shear wave decomposition, so that tjSelecting based on different j layers and specific experimental targets; the 512 x 512 CT image is decomposed into a low frequency part and 4 high frequency parts, and the threshold value T ≈ 5 σnFor more detailed decomposition subbands T ≈ 4 σ is usednOr T ≈ 3 σn；
Step 6) combining the processed low-frequency component and the high-frequency component to carry out NSST reconstruction;
and finally, combining the high-frequency part and the low-frequency part and carrying out NSST reconstruction to obtain a denoised CT image which can enable a doctor to diagnose a focus more conveniently and more accurately.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910186481.2A CN109961411B (en) | 2019-03-12 | 2019-03-12 | Non-downsampling shear wave transformation medical CT image denoising method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910186481.2A CN109961411B (en) | 2019-03-12 | 2019-03-12 | Non-downsampling shear wave transformation medical CT image denoising method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109961411A CN109961411A (en) | 2019-07-02 |
| CN109961411B true CN109961411B (en) | 2021-02-02 |
Family
ID=67024295
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910186481.2A Active CN109961411B (en) | 2019-03-12 | 2019-03-12 | Non-downsampling shear wave transformation medical CT image denoising method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109961411B (en) |
Families Citing this family (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110517198B (en) * | 2019-08-22 | 2022-12-27 | 太原科技大学 | High-frequency sensitive GAN network for denoising LDCT image |
| CN110517199B (en) * | 2019-08-26 | 2022-03-08 | 电子科技大学 | Image rain removing method convenient for intelligent vehicle driving |
| CN110634112B (en) * | 2019-10-15 | 2022-02-22 | 中国矿业大学(北京) | Method for enhancing noise-containing image under mine by double-domain decomposition |
| CN111242853B (en) * | 2019-12-31 | 2023-08-08 | 浙江工业大学之江学院 | Medical CT image denoising method based on optical flow processing |
| CN111340726B (en) * | 2020-02-26 | 2022-08-02 | 青海民族大学 | Image auxiliary denoising method based on supervised machine learning |
| CN112241759B (en) * | 2020-07-22 | 2023-06-09 | 西安交通大学 | Pyramid decomposition method and system for multi-scale geometric analysis |
| CN112819728A (en) * | 2021-02-19 | 2021-05-18 | 燕山大学 | Method for removing radiation image Poisson noise based on shear wave |
| CN115063302B (en) * | 2022-05-10 | 2024-03-29 | 华南理工大学 | Effective removing method for spiced salt noise of fingerprint image |
| CN114782414B (en) * | 2022-06-15 | 2022-12-13 | 深圳市迈捷生命科学有限公司 | Artificial bone data analysis method based on image data processing |
Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| EP2618302A2 (en) * | 2012-01-20 | 2013-07-24 | Kabushiki Kaisha Toshiba | Method and system for image denoising using discrete total variation (TV) minimization with one-direction condition |
| CN104268835A (en) * | 2014-09-23 | 2015-01-07 | 西安电子科技大学 | Weak and small infrared target strengthening method based on image geometric separation |
| CN105631820A (en) * | 2015-12-25 | 2016-06-01 | 浙江工业大学 | Medical ultrasonic image denoising method based on wavelet transform and trilateral filter |
| CN106127711A (en) * | 2016-06-23 | 2016-11-16 | 浙江工业大学之江学院 | Shearlet conversion and quick two-sided filter image de-noising method |
| CN107845079A (en) * | 2017-11-15 | 2018-03-27 | 浙江工业大学之江学院 | 3D shearlet medicine CT video denoising methods based on compact schemes |
| CN108846813A (en) * | 2018-06-04 | 2018-11-20 | 浙江工业大学之江学院 | The medicine CT image denoising method of frame and NSST is decomposed based on MFDF |
| CN109003232A (en) * | 2018-06-15 | 2018-12-14 | 浙江医院 | Medical MRI image de-noising method based on the smooth Shearlet of frequency domain scale |
-
2019
- 2019-03-12 CN CN201910186481.2A patent/CN109961411B/en active Active
Patent Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| EP2618302A2 (en) * | 2012-01-20 | 2013-07-24 | Kabushiki Kaisha Toshiba | Method and system for image denoising using discrete total variation (TV) minimization with one-direction condition |
| CN104268835A (en) * | 2014-09-23 | 2015-01-07 | 西安电子科技大学 | Weak and small infrared target strengthening method based on image geometric separation |
| CN105631820A (en) * | 2015-12-25 | 2016-06-01 | 浙江工业大学 | Medical ultrasonic image denoising method based on wavelet transform and trilateral filter |
| CN106127711A (en) * | 2016-06-23 | 2016-11-16 | 浙江工业大学之江学院 | Shearlet conversion and quick two-sided filter image de-noising method |
| CN107845079A (en) * | 2017-11-15 | 2018-03-27 | 浙江工业大学之江学院 | 3D shearlet medicine CT video denoising methods based on compact schemes |
| CN108846813A (en) * | 2018-06-04 | 2018-11-20 | 浙江工业大学之江学院 | The medicine CT image denoising method of frame and NSST is decomposed based on MFDF |
| CN109003232A (en) * | 2018-06-15 | 2018-12-14 | 浙江医院 | Medical MRI image de-noising method based on the smooth Shearlet of frequency domain scale |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109961411A (en) | 2019-07-02 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109961411B (en) | Non-downsampling shear wave transformation medical CT image denoising method | |
| JP5453401B2 (en) | System and method for discriminating tissue properties by homomorphic deconvolution analysis of backscattered ultrasound | |
| CN109598680B (en) | Shear wave transformation medical CT image denoising method based on rapid non-local mean value and TV-L1 model | |
| Bhuiyan et al. | Spatially adaptive thresholding in wavelet domain for despeckling of ultrasound images | |
| CN106023200A (en) | Poisson model-based X-ray chest image rib inhibition method | |
| Dolui et al. | A new similarity measure for non-local means filtering of MRI images | |
| CN106725612B (en) | Four-dimensional ultrasonic image optimization method and system | |
| CN109003232B (en) | Medical MRI image denoising method based on frequency domain scale smoothing Shearlet | |
| Bhateja et al. | An improved medical image fusion approach using PCA and complex wavelets | |
| JP5818592B2 (en) | Ultrasonic diagnostic apparatus, medical image processing apparatus, and medical image processing method | |
| CN109035156B (en) | DNST (deep depth transform) -based medical CT (computed tomography) image denoising method | |
| Raj et al. | Ultrasound medical image denoising using hybrid bilateral filtering | |
| CN111242853B (en) | Medical CT image denoising method based on optical flow processing | |
| CN108846813A (en) | The medicine CT image denoising method of frame and NSST is decomposed based on MFDF | |
| Zhang et al. | Despeckling Methods for Medical Ultrasound Images | |
| Tsantis et al. | Inter-scale wavelet analysis for speckle reduction in thyroid ultrasound images | |
| Chinnaswamy et al. | Performance evaluation of filters for de-noising the intravascular ultrasound (IVUS) images | |
| CN109584322B (en) | Shearlet medical PET image denoising method based on frequency domain direction smoothing | |
| CN111815527A (en) | Mixed high-order variational ultrasonic image denoising method based on Weibull distribution | |
| Paul et al. | Preprocessing techniques with medical ultrasound common carotid artery images | |
| Nicolae et al. | Image analysis of kidney using wavelet transform | |
| Degirmenci et al. | High dynamic range ultrasound imaging | |
| Radhi et al. | Anisotropic Diffusion Method for Speckle Noise Reduction in Breast Ultrasound Images. | |
| AU2020103375A4 (en) | Speckle Denoising System for Ultrasound Images with Framelet Transform and Gaussian Filter | |
| Kumar et al. | Denoising of Iris image using stationary wavelet transform |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |