A method for oil spill detection in SAR images based on wavelet threshold transform

Abstract

The present invention discloses a SAR image oil spill detection method based on wavelet threshold transformation. First, the SAR image is transformed by logarithmic operator, then sharpened, and then wavelet decomposition with a scale of 3 is performed to obtain wavelet coefficients of different frequency sub-bands, the threshold corresponding to each layer of sub-bands is calculated, and the wavelet coefficients of each layer of sub-bands are threshold processed by a threshold function, and then the SAR image is reconstructed in sequence by two-dimensional discrete wavelet inverse transformation, and finally the reconstructed SAR image is segmented, and the pixel-by-pixel accurate judgment is performed with the visual interpretation result to calculate the oil spill detection accuracy. The SAR of the present invention is an effective means for monitoring marine oil spills due to its all-day and all-weather characteristics. The present invention improves the threshold estimation method and threshold function in the wavelet threshold transformation, improves the segmentation accuracy of the SAR image, and has a good application prospect in SAR image oil spill detection.

Description

Technical Field

[0001] The invention belongs to the technical field of synthetic aperture radar oil spill image detection, and in particular relates to a SAR image oil spill detection method based on wavelet threshold transformation. Background Art

[0002] Oil spill accidents often occur under adverse weather conditions. Affected by natural wind, waves, currents and other conditions, oil spills will quickly spread, penetrate, evaporate and other phenomena, and the electromagnetic signals detected by sensors will also be affected. Therefore, remote sensing monitoring methods such as visible light and infrared sensors and ultraviolet sensors cannot complete comprehensive monitoring. Synthetic aperture radar (SAR) technology has been proven to be an effective means of marine oil spill monitoring due to its all-day and all-weather characteristics.

[0003] Using effective extraction algorithms to extract oil spill information quickly and accurately from oil spill SAR images is of great significance for oil pollution detection and other processing. At present, there are many methods for denoising, such as wavelet transform-based methods, non-local means algorithms, sparse coding shrinkage, Bayesian methods, principal component analysis, compressed sensing theory, bilateral filtering, etc. These technologies can smooth and denoise images and can maintain useful details and edge characteristics to a certain extent. However, these methods have certain limitations in terms of image quality or computational efficiency. For example, the non-local means method requires a relatively high computational cost when performing a global search for similar grayscales; the bilateral filtering method is easy to implement and has good edge preservation characteristics, but it is not suitable for denoising images containing complex information in high frequencies.

[0004] Due to the sparsity and scale decomposition characteristics of wavelet transform coefficients, it has attracted many researchers to study the wavelet domain. These characteristics of the wavelet domain provide flexibility for signal representation. The main energy of the representative signal is represented by several large coefficients, and the remaining energy is represented by many small coefficients. Because most of the noise power appears in small coefficients, these coefficients can be modified according to certain rules to achieve image denoising. Wavelet threshold denoising includes two key issues: threshold function and threshold selection. The modification rule for wavelet coefficients is the threshold function, and the threshold is the threshold value of the threshold function in modifying the wavelet coefficients. It defines which wavelet coefficients need to be processed in what way. In order to obtain a better wavelet threshold denoising effect, it is necessary to design the optimal threshold function and threshold.

[0005] SAR image denoising has experienced a development from spatial domain filtering to frequency domain filtering, during which different types of algorithms have emerged, such as anisotropic denoising algorithms, block matching denoising algorithms, and hybrid denoising algorithms. After studying and analyzing typical denoising algorithms, it is found that these algorithms have the following main disadvantages:

[0006] (1) Spatial domain filtering algorithms are usually prone to damaging image details, resulting in poor denoising results. Frequency domain filtering has made great progress compared to spatial domain filtering, but the denoised image is generally blurry and cannot achieve a good balance between denoising and edge preservation.

[0007] (2) The denoising results obtained by using frequency domain denoising algorithms or some hybrid denoising algorithms will produce pseudo-Gibbs phenomenon, resulting in artificial texture or local distortion in the denoised image. Since the traditional wavelet threshold transform does not perform well in SAR image denoising, the final segmentation accuracy can only reach about 85% at most. Summary of the invention

[0008] The purpose of the present invention is to solve the above technical problems and to provide a SAR image oil spill detection method based on wavelet threshold transformation.

[0009] The method of the present invention adopts a new threshold estimation method and threshold function to replace the traditional general threshold and hard and soft threshold functions, which greatly eliminates multiplicative noise and retains the edge details of the original picture to the greatest extent possible, thereby improving the detection accuracy.

[0010] In order to achieve the above object, the present invention provides the following technical solutions:

[0011] The SAR image oil spill detection method based on wavelet threshold transform includes the following steps:

[0012] Step 1: The speckle noise with multiplicative characteristics on the SAR image is first transformed into additive noise through logarithmic operator transformation, and then the transformed SAR image is sharpened;

[0013] Step 2: Perform a two-dimensional discrete wavelet transform on the SAR image processed in step 1, perform a wavelet decomposition of the SAR image with a scale of 3, decompose the SAR image into sub-bands of different frequencies, and obtain the discretized wavelet coefficients of each sub-band;

[0014] Step 3: Set a threshold λ, process each layer of sub-bands separately, and calculate the threshold λ corresponding to each layer of sub-bands j ;

[0015] Step 4: Use the threshold function to perform threshold processing on the discretized wavelet coefficients of each layer sub-band, retain or shrink the wavelet coefficients greater than the threshold set in step 3, and do not blindly eliminate the wavelet coefficients less than this threshold, but also retain and shrink them to a certain extent;

[0016] Step 5: Restore the final reconstructed SAR image in the order of decomposition by inverse two-dimensional discrete wavelet transform;

[0017] Step 6: Segment the reconstructed SAR image and make precise pixel-by-pixel judgments based on the expert visual interpretation results to calculate the segmentation accuracy.

[0018] The present invention is further improved, step 1 also includes using the USM sharpening enhancement algorithm for the edge blur phenomenon after the SAR image is logarithmically transformed. In order to control the sharpening parameters, the sharpening amount is 200%, the radius is 1.2, and the threshold is 3.

[0019] The present invention is further improved. When the wavelet basis is used to perform wavelet decomposition on the SAR image in step 2, the selected wavelet basis is db5, taking into account that the selected wavelet basis must satisfy the principles of orthogonality, symmetry, regularity, compact support, vanishing moment, etc.

[0020] The present invention is further improved, the threshold λ set in step 3 is a hierarchical adaptive threshold, and the estimation formula of the threshold is as follows:

[0021]

[0022] Where: σ represents the standard deviation of the noise image, N is the length of the signal after wavelet transformation, j represents the decomposition scale, and when j is equal to 1, it is the universal threshold. During wavelet decomposition, as the scale j increases, the image signal amplitude and the noise amplitude change in opposite directions, and the threshold λ j It will become smaller accordingly.

[0023] The present invention is further improved, the threshold function used in step 4 is added with an adjustment factor, has bidirectional contraction, and is continuously differentiable at the threshold. The threshold function is as follows:

[0024]

[0025] Where m and k are adjustment factors, m∈[0,1], k is a positive integer. The threshold function proposed in the present invention is continuous and controllable at point λ. The parameters m and k in the threshold function proposed in the present invention are adjustment factors; when k≥|λ| and m→0, the function has a hard threshold, and when 0<k<|λ| and m→1, the function is infinitely close to the soft threshold; in the function, m adjusts the change of the wavelet coefficient, and k adjusts the threshold function at |w j，k ≤λ| for smoothness.

[0026] The present invention is further improved. The specific method for segmenting the reconstructed SAR image described in step 6 is: adding a weight coefficient to the inter-class variance of the background and the target, and introducing the proportion of the target pixel points in the whole image and the proportion of the background pixel points in the whole image as weights into the target equation.

[0027] The present invention is further improved, and the formula of the weighted maximum inter-class variance method described in step 6 is as follows:

[0028] u(t) 2 =ωω0(t)u0(t) 2 +ω2ω1(t)u1(t) 2 (3)

[0029] T=Argmax 0＜L≤L-1 (u(t) 2 ) (4)

[0030] Where: ω and ω2 are the ratios of target pixels to background pixels, and ω=ω0(t), ω2=ω1(t); u0(t) and u1(t) are the target mean and background mean, respectively.

[0031] Beneficial effects of the present invention:

[0032] After the threshold estimation method and threshold function proposed in the present invention are applied to the wavelet threshold transformation, the image details are well preserved on the basis of image smoothing. The weighted maximum inter-class variance method proposed in the present invention is used to perform threshold segmentation on the pre-processed image, which greatly improves the accuracy of SAR oil spill information extraction by wavelet threshold transformation. The performance in MSE and PSNR evaluation indicators is better than that of traditional wavelet threshold transformation, and the final segmentation accuracy is as high as about 95%. BRIEF DESCRIPTION OF THE DRAWINGS

[0033] Figure 1 is a flow chart of the method of the present invention;

[0034] Figure 2 A schematic diagram showing a comparison between the threshold function of the present invention and a traditional threshold function;

[0035] Figure 3 is the original SAR image;

[0036] Figure 4 The results of the present invention on the oil spill detection of SAR images;

[0037] Figure 5 The results of oil spill detection on SAR images using traditional wavelet threshold transform. DETAILED DESCRIPTION

[0038] In order to deepen the understanding of the present invention, the present invention will be further described in detail below in conjunction with embodiments. The embodiments are only used to explain the present invention and do not limit the protection scope of the present invention.

[0039] Example: Figure 1 As shown, a SAR image oil spill detection method based on wavelet threshold transformation of the present invention comprises the following steps:

[0040] Step 1: The speckle noise with multiplicative characteristics on the SAR image is first transformed into additive noise through logarithmic operator transformation, and then the transformed SAR image is sharpened. The specific implementation method is as follows.

[0041] The SAR image is logarithmically transformed to convert the multiplicative noise into additive noise, which is convenient for subsequent removal. The logarithmic operator transformation formula is:

[0042] u0(x,y)=u(x,y)q(x,y) (1)

[0043] logu0(x,y)=logu(x,y)+logq(x,y) (2)

[0044] Where: u0(x, y) is the observed image, u(x, y) is the original image, q(x, y) is the noise, and the distribution of the noise q(x, y) is unknown.

[0045] Aiming at the edge blur phenomenon after logarithmic transformation of SAR images, the USM sharpening enhancement algorithm is adopted. In order to control the sharpening parameters, the USM sharpening amount is 200%, the radius is 1.2, and the threshold is 3.

[0046] Step 2: Perform a two-dimensional discrete wavelet transform on the SAR image processed in step 1, perform a wavelet decomposition of scale 3 on the SAR image, decompose the SAR image into sub-bands of different frequencies, and obtain the discretized wavelet coefficients of each sub-band. The specific implementation method is as follows.

[0047] (1) Considering that the selected wavelet basis must satisfy the principles of orthogonality, symmetry, regularity, compact support, and vanishing moments, db5 is selected as the wavelet basis.

[0048] (2) Through the two-dimensional discrete wavelet transform (DWT), after the first scale decomposition, the SAR image is decomposed into four sub-bands with different frequencies, represented by LL1, LH1, HL1, and HH1. Among them, LL1 is a low-frequency sub-band, LH1, HL1, and HH1 are high-frequency sub-bands, and the HH1 sub-band has the highest frequency and contains the most noise. The LL1 sub-band corresponds to the details of the image, and the LL1 sub-band can be further recursively decomposed.

[0049] (3) After the second scale decomposition, only the low-frequency sub-band LL1 is decomposed, and it can be decomposed into four sub-bands, represented by LL2, LH2, HL2, and HH2 respectively.

[0050] (4) Similarly, after the third scale decomposition, four sub-bands with different frequencies can be obtained, which are represented by LL3, LH3, HL3, and HH3 respectively.

[0051] The formula for discrete two-dimensional wavelet transform is as follows:

[0052]

[0053] Where a0 is the discrete scale, τ0 and τ1 are the discrete displacement steps, and j is the decomposition scale.

[0054] Step 3: Set a threshold λ, process each layer of sub-bands separately, and calculate the threshold λ corresponding to each layer of sub-bands j The specific implementation method is as follows.

[0055] The threshold λ of each layer is calculated by the threshold estimation formula j , the specific threshold estimation formula is as follows:

[0056]

[0057] Set the threshold λ j is the hierarchical adaptive threshold, where σ represents the standard deviation of the noise image, N is the length of the signal after wavelet transform, j represents the decomposition scale, and when j is equal to 1, it is the universal threshold. During wavelet decomposition, as the scale j increases, the image signal amplitude and the noise amplitude change in opposite directions, and the threshold λ j It will become smaller accordingly.

[0058] Step 4: Use the threshold function to perform threshold processing on the discretized wavelet coefficients of each layer subband, retain or shrink the wavelet coefficients greater than the threshold set in step 3, and do not blindly eliminate the wavelet coefficients less than this threshold, but also retain and shrink them to a certain extent. The specific implementation method is as follows.

[0059] Through the processing of the threshold function, the wavelet coefficients of each layer after wavelet decomposition are retained and shrunk. The specific threshold function formula is as follows:

[0060]

[0061] Where: m and k are adjustment factors, m∈[0,1], k is a positive integer. The threshold function proposed in the present invention is continuous and controllable at point λ. The parameters m and k in the threshold function proposed in the present invention are adjustment factors; when k≥|λ| and m→0, the function has a hard threshold, and when 0<k<|λ| and m→1, the function is infinitely close to the soft threshold; in the function, m adjusts the change of the wavelet coefficient, and k adjusts the threshold function at |w j，k ≤λ| for smoothness.

[0062] Step 5: Through the two-dimensional discrete wavelet inverse transform, the final reconstructed SAR image is restored in sequence according to the order of decomposition.

[0063] In this embodiment, the specific implementation method of step 5 is as follows.

[0064] (1) First, LL3, LH3, HL3, and HH3 are reconstructed to obtain a new low-frequency subband LL2.

[0065] (2) LL2, LH2, HL2, and HH2 are reconstructed to obtain a new low-frequency subband LL1.

[0066] (3) Finally, LL1, LH1, HL1, and HH1 are reconstructed to obtain the final reconstructed image.

[0067] Step 6: Segment the reconstructed SAR image, and make accurate pixel-by-pixel judgments based on the expert visual interpretation results to calculate the segmentation accuracy. The specific implementation method is as follows.

[0068] The ratio of target pixels to the entire image and the ratio of background pixels to the entire image are introduced as weights into the target equation of the maximum inter-class variance method, and the weighted maximum inter-class variance method is used to perform threshold segmentation on the denoised SAR image. The formula of the weighted maximum inter-class variance method is as follows:

[0069] u(t) 2 =ωω0(t)u0(t) 2 +ω2ω1(t)u1(t) 2 (6)

[0070] T=Argmax 0＜L≤L-1 (u(t) 2 ) (7)

[0072] Where: ω and ω2 are the ratios of target pixels to background pixels, and ω=ω0(t), ω2=ω1(t); u0(t) and u1(t) are the target mean and background mean, respectively.

[0073] The final reconstructed image is segmented, and the segmentation result is as follows Figure 4 and Figure 5 The accuracy of segmentation is calculated by comparing the results of the expert visual interpretation and making accurate pixel-by-pixel judgments.

[0074] The effect of the present invention is further described below in conjunction with experiments:

[0075] The experiment used data from an offshore oil spill experiment conducted by the Norwegian Ocean Cleaning Association (NOFO) in the North Sea of ​​Europe from June 6 to 9, 2011 ( Figure 3 ), three substances were splashed in the experiment, namely crude oil, emulsified oil and bio-oil, with an interval of about 12 hours. At the end of the experiment, in order to protect the marine ecological environment, chemical decomposition was used to decompose and recycle the substances splashed in the experiment.

[0076] Table 1 Image denoising performance evaluation

[0077]

[0078] It can be clearly seen from Table 1 that the threshold function and threshold calculation method proposed in this paper have obvious denoising effects on SAR images.

[0079] Table 2 Detection accuracy of oil spill images

[0080]

[0081] It can be clearly seen from Table 2 that the threshold function and threshold calculation method proposed in the present invention have significantly improved the accuracy of SAR image oil spill detection, and have good application prospects in oil spill detection.

[0082] The above is only an embodiment of the present invention, and the common sense such as the known specific structure and characteristics in the scheme is not described in detail here. For those skilled in the art, it is obvious that the present invention is not limited to the details of the above exemplary embodiments, and the present invention can be implemented in other specific forms without departing from the spirit or basic characteristics of the present invention. Therefore, no matter from which point of view, the embodiments should be regarded as exemplary and non-restrictive. The scope of the present invention is defined by the attached claims rather than the above description, and it is intended to include all changes within the meaning and scope of the equivalent elements of the claims. Any mark in the claims should not be regarded as limiting the claims involved.

Claims

1. A SAR image oil spill detection method based on wavelet threshold transform, characterized in that: The steps include: Step 1: The speckle noise with multiplicative characteristics on the SAR image is first transformed into additive noise through logarithmic operator transformation, and then the transformed SAR image is sharpened; Step 2: Perform a two-dimensional discrete wavelet transform on the SAR image processed in step 1, decompose the SAR image into sub-bands of different frequencies by wavelet decomposition with a scale of 3, and obtain the discretized wavelet coefficients of each sub-band; Step 3: Set a threshold λ, process each layer of sub-bands separately, and calculate the threshold λ corresponding to each layer of sub-bands j , threshold λ j For the hierarchical adaptive threshold, the estimation formula is as follows: Where: σ represents the standard deviation of the noise image, N is the length of the signal after wavelet transformation, j represents the decomposition scale, and when j is equal to 1, it is the universal threshold. During wavelet decomposition, as the scale j increases, the image signal amplitude and the noise amplitude change in opposite directions, and the threshold λ j will become smaller accordingly; Step 4: Use the threshold function to perform threshold processing on the discretized wavelet coefficients of each layer subband, and set the coefficients greater than the threshold λ set in step 3. j The wavelet coefficients are retained or shrunk, and those smaller than this threshold λ j The wavelet coefficients of are also retained and shrunk; the threshold function adds an adjustment factor, has bidirectional shrinkage, and is continuously differentiable at the threshold. The threshold function is as follows: Where m and k are adjustment factors, m∈[0,1], k is a positive integer; Step 5: Through the two-dimensional discrete wavelet inverse transform, the final reconstructed SAR image is restored in sequence according to the order of decomposition; Step 6: Segment the reconstructed SAR image and make precise pixel-by-pixel judgments based on the expert visual interpretation results to calculate the segmentation accuracy.

2. A SAR image oil spill detection method based on wavelet threshold transform according to claim 1, characterized in that: Step 1 also includes using the USM sharpening enhancement algorithm to deal with the edge blur after the SAR image is logarithmically transformed, and the sharpening parameters are set to 200%, 1.2, and 3.

3. The SAR image oil spill detection method based on wavelet threshold transform according to claim 1 is characterized in that: In step 2, when the processed SAR image is decomposed by wavelet basis, the wavelet basis selected is db5.

4. The SAR image oil spill detection method based on wavelet threshold transform according to claim 1 is characterized in that: The specific method of step 6 for segmenting the reconstructed SAR image is: adding weight coefficients to the inter-class variances of the background and the target, and introducing the proportion of the target pixels in the entire image and the proportion of the background pixels in the entire image as weights into the objective equation of the maximum inter-class variance method.

5. The SAR image oil spill detection method based on wavelet threshold transform according to claim 4, characterized in that: The objective equation of the maximum inter-class variance method is as follows: u(t) 2 =ωω0(t) u 0(t) 2 +ω2ω1(t)u1(t) 2 (3) T=Argmax 0＜L≤L-1 (u(t) 2 ) (4) Where: ω and ω2 are the ratios of target pixels to background pixels, and ω=ω0(t), ω2=ω1(t); u0(t) and u1(t) are the target mean and background mean, respectively.

Images