Removal of Large-Scale Stripes Via Unidirectional Multiscale Decomposition

: Stripes are common in remote sensing imaging systems equipped with multichannel time delay integration charge-coupled devices (TDI CCDs) and have di ﬀ erent scale characteristics depending on their causes. Large-scale stripes appearing between channels are di ﬃ cult to process by most current methods. The framework of column-by-column nonuniformity correction (CCNUC) is introduced to eliminate large-scale stripes. However, the worst problem of CCNUC is the unavoidable cumulative error, which will cause an overall color cast. To eliminate large-scale stripes and suppress the cumulative error, we proposed a destriping method via unidirectional multiscale decomposition (DUMD). The striped image was decomposed by constructing a unidirectional pyramid and making di ﬀ erence maps layer by layer. The highest layer of the pyramid was processed by CCNUC to eliminate large-scale stripes, and multiple cumulative error suppression measures were performed to reduce overall color cast. The di ﬀ erence maps of the pyramid were processed by a designed ﬁlter to eliminate small-scale stripes. Experiments showed that DUMD had good destriping performance and was robust with respect to di ﬀ erent terrains.


Introduction
High-resolution optical remote sensing satellites are generally imaged by time delay integration charge-coupled devices (TDI CCD) sensors. A single CCD consists of multiple taps, and different taps have different levels of circuit noise. Moreover, a camera is often spliced by multiple CCDs to satisfy image width requirements. Therefore, a single linear camera has multiple channels. With a single channel, there is high-frequency response nonuniformity between different detectors, and this nonuniformity will produce small-scale stripes. With multiple channels, there is low-frequency response nonuniformity, which will produce large-scale stripes [1]. For example, the optical satellite Gaofen-1B (GF-1B) is equipped with two multispectral cameras, each camera has three CCDs, and each CCD has two taps. That is, a GF-1B multispectral image forms 12 channels. If the nonuniformity is not effectively corrected, the large-scale stripes will exist between these channels, and the small-scale stripes will exist in these channels.
Before and after the launch of a satellite, the staff will perform laboratory and on-orbit relative radiation calibration for the sensors to eliminate nonuniformity. However, errors can occur in the process of relative radiometric calibration, and the radiation response of the detectors will drift with

Conversion of the Destriping Problem
The no-stripe image is denoted L, the striped image is denoted L * , the destriped image is denoted L. Assuming that they are all of the size M × N, and that L * is affected by multiplicative and additive stripe noise simultaneously, where p * i,j is the pixel in L * , p i,j is the pixel in L, k j is the multiplicative stripe noise of the jth detector, b j is the additive stripe noise of the jth detector, and subscripts i, j are the row and column indexes, respectively.
The jth column of L is denoted P j = [p 1, j , p 2,j , · · · , p M,j ] T , the mean value of P j is denoted a j , and the standard deviation of P j is denoted s j . Similarly, the jth column of L * is denoted P * j = [p * 1, j , p * 2, j , · · · , p * M, j ] T , the mean value of P * j is denoted a * j , and the standard deviation of P * j is denoted s * j . These values have the following relationship: It can be seen from Equation (2) that standard deviation is only affected by the multiplicative noise, while the mean value is affected by both multiplicative and additive noise. For the jth column, the stripe noise can be eliminated if a j and s j can be determined. Moment matching assumes that a j and s j equal the mean and standard deviation of the overall image. The performance of moment matching depends on the image size and terrain uniformity of the imaged area. If the image size is too small or the terrain is not sufficiently uniform, moment matching will produce new stripes.
If the additive noise is ignored, then the destriped pixel can be expressed as: wherep i,j is a pixel inL andk j = a * j /a j is the approximate multiplicative noise. It can be seen from Equation (3) that if the additive noise is ignored, there will be an error term whose value is equal to Multiplicative noise is the main component of stripes; its main sources are the nonuniformity of responses, the gain amplification effect, and other factors. Additive noise is the secondary component of stripes; its main source is dark current whose intensity is relatively small. The error term caused by ignoring additive noise is relatively small. If the mean value of every column can be obtained, the no-stripe image can be reconstructed. The mean values of the columns can be written as a vector (the column mean vector, CMV). Therefore, the destriping problem can be transformed into estimating the optimal CMV of the no-stripe image.

Scale Characteristics of Stripes
Stripes have different scale characteristics depending on their sources. For a single channel (tap), the sources of the stripes include photo shot noise, reset noise, dark current, fixed pattern noise, and quantization noise [18,19]. The stripes in the channel appear as thin lines in the image and appear as small-scale impulse noise on the corresponding CMV. Figure 1a shows the Band4 of GF-1C as an example of small-scale stripes, and Figure 1c is its corresponding CMV.
The optical sensor in the satellite is generally formed by multiple CCDs, and each CCD is formed by multiple taps, causing a linear camera to have multiple channels. Different channels have different circuit noise, and the circuit noise is amplified by the gain, which causes the image to produce an obvious color cast. Stripes between channels are large-scale and of high-intensity and appear as large-scale steps in the corresponding CMV. Figure 1b shows the Band2 of GF-2 as an example of large-scale stripes, and Figure 1d is its corresponding CMV. This large-scale stripe appears at the junction of two CCDs. It is caused by the nonuniformity of the two CCDs, and this nonuniformity is not repaired by radiometric correction completely.
Small-scale stripes are easy to process, either by filtering the CMV or by using some variation-based methods. However, it is difficult to process large-scale stripes. To eliminate large-scale stripes, we introduced the framework of column-by-column nonuniformity correction (CCNUC), which is based on the spatial similarity of remote sensing images. Where (a) is small-scale stripe, (b) is large-scale stripe, (c) is the CMV of (a), (d) is the CMV of (b).

Overall Description of DUMD
When establishing the unidirectional pyramid of stripe image

Overall Description of DUMD
When establishing the unidirectional pyramid of stripe image L * ∼ L * 0 , L * 1 , · · · , L * n , the large-scale stripes are saved in the highest layer L * n , and the small-scale stripes are saved in difference maps. The difference map of the ith layer is written as D * i and is defined as: where L * i and L * i+1 both need to be resampled to have the same size as L * . The CMVs of L * 0 , L * 1 , · · · , L * n are denoted A * 0 , A * 1 , · · · , A * n , and the corresponding destriped CMVs are denoted {Â 0 ,Â 1 , · · · ,Â n }, respectively. The CMVs of D * 0 , D * 1 , · · · , D * n−1 are denoted E * 0 , E * 1 , · · · , E * n−1 , and the corresponding destriped CMVs are denoted {Ê 0 ,Ê 1 , · · · ,Ê n−1 }, respectively. Extending Equation (4), we obtained: Equation (5) illustrates that E * i can be obtained from A * i and A * i+1 and does not have to use D * i . Moreover, A * i can be obtained from the downsampling of A * (the CMV of L * ) and does not have to use L * i . Therefore, only A * needs to be calculated from two-dimensional data, while the other CMVs (A * i ,E * i ) are all calculated from one-dimensional data; this feature can greatly reduce the amount of calculation.
In Section 2.1, the destriping problem was converted into the CMV estimation of the no-stripe image. Therefore, the destriped imageL can be expressed as: In Equation (6), the destriping problem has been converted into the problem of seekingÂ n and E i . The small-scale stripes are saved in D * i , which is eliminated byÊ i . Small-scale stripes are easy to process, and we have designed a filter G(•) to remove them. Assuming 1 , e i,2 , · · · , e i,N ] (where h 1d is a one-dimensional smooth kernel), the filter can be expressed as:Ê where δ is a threshold used to judge small-scale stripe noise; δ was set at 1 in this paper. The large-scale stripes are saved in L * n , which is eliminated byÂ n . Large-scale stripes are difficult to process; this difficulty was the main issue addressed in this paper. The CCNUC framework was applied to obtainÂ; this process has been described more thoroughly in Section 3.2. However, the cumulative error produced by CCNUC cannot be ignored, and additional measures must be applied for suppression. These cumulative error suppression measures are described in Sections 3.3-3.5.
In summary, the overall process of the DUMD algorithm is shown in Figure 2, and mainly includes five steps: (1) A unidirectional pyramid of the stripe image is constructed, and the highest layer L * n is obtained. (2) CCNUC and cumulative error suppression measures are performed on L * n to obtain the destriped layerL n and the corresponding CMVÂ n .

The Framework of Column-by-Column Nonuniformity Correction (CCNUC)
Remote sensing images have spatial similarity, and any two random adjacent columns should have similar statistical characteristics. Because of this feature, the stripe image can be corrected by CCNUC with various methods. In Section 2.1, the destriping problem was converted into the problem of estimating the CMV of the no-stripe image, so the stripe noise can be removed by calculating the deviation between adjacent columns. The CCNUC framework can be expressed as: is a function that calculates the adjustment coefficient between two

The Framework of Column-by-Column Nonuniformity Correction (CCNUC)
Remote sensing images have spatial similarity, and any two random adjacent columns should have similar statistical characteristics. Because of this feature, the stripe image can be corrected by CCNUC with various methods. In Section 2.1, the destriping problem was converted into the problem of estimating the CMV of the no-stripe image, so the stripe noise can be removed by calculating the deviation between adjacent columns. The CCNUC framework can be expressed as: whereP j is the jth column of the destriped image, c k is the adjustment coefficient between P * k and P * k−1 , and F : (X, Y) → c is a function that calculates the adjustment coefficient between two columns. The adjustment coefficient c k can be divided into two parts: the accurate part (denoted γ k ) and the Remote Sens. 2019, 11, 2472 7 of 16 error term (denoted ε k ). There is a cumulative error ∆ between the estimated columnP j and the real column P j : In Equation (9), the role of P * j a * j is to normalize P * j . Therefore, the intensity of ∆ is positively correlated with the number of accumulations j and the magnitude of the single adjustment error ε k .
The cumulative error is an unavoidable problem of the CCNUC framework. In the image, the cumulative error appears as a color cast on the destriped image. The farther the error is from the initial column, the more serious the color cast. In the CMV, the cumulative error appears as an overall offset. The farther the error is from the initial point, the larger the offset. There are three steps that can be taken to compensate for the cumulative error produced by CCNUC: (1) Reduce the number of error accumulations. By Equation (9), the intensity of ∆ is positively correlated with j. Reducing the size of j is the most effective way to suppress the cumulative error. When the unidirectional pyramid of stripe image L * ∼ L * 0 , L * 1 , · · · , L * n is established, only the highest layer L * n is processed by CCNUC. In this way, the cumulative number can be reduced from N − 1 to N 2 n − 1. This step has been described in more detail in Section 3.3.
(2) Improve the accuracy of the single deviation estimation of adjacent columns (i.e., reduce the value of the single adjustment error ε k ). The cumulative error comes from the accumulation of multiple tiny errors, so improving the accuracy of every deviation estimation is also an effective suppression method. The size of ε k is related to the samples and the calculation method. The least-square method (LSM) is the most common method for calculating the linear relationship between two samples. However, its objective function is the root mean square error (RMSE), which is sensitive to extreme values (gross errors). If these extreme values cannot be excluded, the deviation estimation of the adjacent columns will be enlarged. Therefore, we proposed a deviation estimation method based on a normal distribution. This method can reduce the influence of extreme values and improve accuracy. This step has been described in more detail in Section 3.4.
(3) Design a compensation measure for cumulative error. The adjustment error ε k cannot be measured, but it has the property of accumulation. Using this property, we could detect it on a larger scale and compensate for it. Recording the initial corrected L * n asL (1) n , decomposingL (1) n to a larger scale, and making the second adjustment by the same method will compensate for cumulative error to some extent. This step has been described in more detail in Section 3.5.

Unidirectional Multiscale Decomposition
Because of their different causes, stripes have different scales. The stripes in the channel are small-scale and of low-intensity, and the stripes between channels are large-scale and of high-intensity. The unidirectional pyramid of stripe image L * ∼ L * 0 , L * 1 , · · · , L * n is established by mean filtering with a filter size 1 × 3. All layers of the pyramid are resampled to the same size, and the difference maps are obtained by Equation (4). Then, the stripe image can be decomposed into the multiscale sets L * → D * 0 , D * 1 , · · · , D * n−1 , L * n . Each layer of the unidirectional pyramid can only contain stripes of the same scale or larger scale. Therefore, each layer of the multiscale sets can only retain stripes of the same scale. In other words, stripes are completely decomposed into multiscale sets without overlap; L * n retains large-scale stripes, and D * 0 , D * 1 , · · · , D * n−1 retain small-scale stripes. Stripes appear as two-dimensional signals in the image, but they can be completely and losslessly compressed into one-dimensional signals in the corresponding CMV. Similarly, the CMVs of D * 0 , D * 1 , · · · , D * n−1 (Ê 0 ,Ê 1 , · · · ,Ê n−1 ) retain all the small-scale stripes, and the CMV of L * n (A * n ) retains all the large-scale stripes. Small-scale stripes can be eliminated by filteringÊ 0 ,Ê 1 , · · · ,Ê n−1 . Therefore, only L * n needs to be processed by CCNUC; this processing can reduce the error accumulation number from the original N − 1 to the current N 2 n − 1. In the CCNUC framework, the adjustment parameter c k is calculated by the statistical characteristics of the adjacent columns. The availability of the statistical characteristics depends on the number of samples. If the image is shrunk by the same factor in width and length, then the number of pixels in each column will be drastically reduced, and the reliability of the statistical characteristics cannot be guaranteed. Therefore, it is only necessary to shrink the image in a single direction (perpendicular to the stripe direction) in the process of establishing the pyramid.

Deviation Estimation of Adjacent Columns
The actual features of the remote sensing image change gradually, indicating that the remote sensing images have a spatial similarity. Therefore, the premise of the CCNUC framework is that two adjacent columns in the remote sensing image should have similar statistical characteristics. When the statistical characteristics of adjacent columns have large differences, it indicates that there are external factors causing systematic radiation deviation (e.g., existing stripe noise). The key issue of the CCNUC framework is how to measure the statistical characteristics of adjacent columns. If the mean and standard deviation are used as statistical indicators, CCNUC is equal to moment matching. If only the mean is used as a statistical indicator, CCNUC is equal to simplified moment matching, and its performance will be worse.
Stripe noise will cause a linear change in the columns. LSM is the most popular method for calculating the linear relationship between two sets of samples. Although the adjacent columns have similarities, they are not identical. There are a few differences between adjacent columns. These differences are caused by features and are especially noticeable in the high-gradient region. The object function of LSM is the RMSE of the two sets of samples, and this RMSE is susceptible to extreme values. Therefore, using LSM to estimate the deviation between adjacent columns will cause over-estimation (i.e., will increase the absolute value of ε k ) and enlarge the cumulative error. Based on this consideration, this paper proposed a deviation estimation method based on the normal distribution, which could reduce the influence of extreme values and improve the robustness of the measurement.
If the radiation responses of two adjacent detectors are consistent, then the difference between the adjacent columns of their images can be considered to be caused by features. In this situation, the difference between the two adjacent columns should satisfy the zero-mean normal distribution: where X k is the difference between adjacent columns, µ k is the mean of the normal distribution, which is equal to 0, and σ 2 is the variance of the normal distribution. If the radiation responses of two adjacent detectors are inconsistent, the fitting Gaussian curve will shift and µ k 0. Therefore, the adjustment parameter can be obtained by measuring the mean of the Gaussian fitting curve of X k : where F(•) is a function that calculates the adjustment parameter, P * k and P * k−1 are two adjacent columns, and η is the threshold introduced to reduce the algorithm error; it was set to the median of µ 2 , µ 3 · · · , µ N in this paper.

Cumulative Error Compensation
After L * n has been processed by CCNUC once, the large-scale stripes have been filtered out. At this time, the initial destriped imageL k−d ) contains two parts: the difference in features and the cumulative error. Therefore, F(P k−d ) has a certain cumulative error compensation effect. Since the two columns are not adjacent, the difference in features will be correspondingly larger, thereby reducing the compensation accuracy. Based on this consideration,L (1) n is downsampled again to obtainL (1) n+1 , which is then processed by CCNUC to obtainL (2) n+1 . Then, the final destriped layerL n can be expressed as: whereÂ n is the CMV ofL n ,Â (1) n is the CMV ofL (1) n , h 1d is a one-dimensional kernel, andÂ (2) n+1 is the CMV ofL (2) n+1 , which has been resampled to be the same size asÂ (1) n .

Experimental Data
This paper selected six multispectral images as experimental data. They were imaged by different satellites and included a variety of features, such as cities, mountains, waters, and clouds. The experimental data were sufficiently highly representative to evaluate the performances and characteristics of the proposed method and the comparison methods. The specific conditions of the experimental data are shown in Table 1. The original images had been radiometrically corrected during the data production. There were no stripes in the original images, so they could be treated as reference images to participate in subsequent evaluation. The striped images were obtained by adding random stripes with different scales to the original images and were used as the input data for testing. The striped images were processed by moment matching (MM) [11], improved moment matching (IMM) [12], spatially adaptive unidirectional total variation (SAUTV) [16], and the proposed method (DUMD) to obtain the final destriped images.
A combination of subjective and objective methods was applied to evaluate destriping performance. In the subjective evaluation, destriping performance was qualitatively evaluated by observing whether the destriped images have obvious stripes and color casts. In the objective evaluation, the root mean square error (RMSE) and signal-to-noise ratio (SNR) were applied to evaluate the experimental results; this evaluation could be done because of the existence of reference images.

Discussion
It can be seen from Figures 3-5 that after adding stripes with different scales, the striped images had a large number of high-frequency thin stripes on the microscopic surface, and multiple stripes with obvious color differences on the macroscopic surface. The small-scale stripes could be effectively removed by all four methods. However, for large-scale stripes, there were significant differences in the effects of these four methods. IMM and SAUTV were ineffective for large-scale stripes, and the color cast did not disappear. MM and DUMD could eliminate large-scale stripes, and the color cast was significantly improved. However, MM heavily relied on terrain uniformity. When there are terrains with large differences in grayscale, new stripes would appear in the image. There were obvious stripes in the MM-processed destriped images, such as the cloud area in Figure 3c and the water areas in Figures 3c and 4c. Compared with MM, DUMD was more robust with respect to terrain nonuniformity, and there were no obvious stripes in the corresponding areas.
In Figures 6-8, the CMVs of striped images had high-frequency impulse noise (small-scale stripe noise) and large-scale platform steps (large-scale stripe noise). After being processed by MM, IMM, SAUTV, and DUMD, the impulse noise was filtered out effectively. However, the overall differences among the four destriped CMVs were large, indicating that the processing results for large-scale stripe noise were quite different. MM assumed that the standard deviation and the mean of each column were the same for the whole image, so the column mean of the MM-processed destriped image appeared as a straight line parallel to the x-axis. IMM assumed that the standard deviation and the mean of each column were the same as for the block near the column, so the CMV of the IMM-processed destriped image appeared as a smooth curve. SAUTV eliminated stripes by constructing a model function that penalizes the partial derivative of the image in the x-direction. Therefore, the CMV of the SAUTV-processed destriped image was more consistent with the striped image. The DUMD eliminated large-scale stripes by CCNUC and multiple cumulative error suppression measures and eliminated small-scale stripes by filtering. However, a cumulative error was introduced by CCNUC. Therefore, the CMV of the DUMD-processed destriped image diverged from the intermediate point to both ends compared with the reference curve.
It can be seen from Table 2 that IMM and SAUTV had the best parameter performances, DUMD was second best, and MM was the worst. These findings occurred because both IMM and SAUTV dealt with the local part and did not have much influence on the grayscale distribution of the whole image, so the parameters were relatively good. However, since the two methods did not achieve the goal of removing large-scale stripes, the parameter performances did not have practical significance. MM adjusted the standard deviation and the mean of each column to be the same as the whole, so the absolute difference increased. Especially in the case of large variations in terrain, the performance of the parameters would be poorer. Compared with MM, DUMD retained the radiation description of the terrain to a certain extent. Although the grayscale distortion on both sides of the image was larger under the influence of cumulative error, the parameter performance of DUMD was still better than that of MM.

Conclusions
This paper analyzed the causes, characteristics, and manifestations of stripe noise in remote sensing images and illustrated the characteristics of stripe noise at different scales. In contrast with small-scale stripes, large-scale stripes are difficult to deal with using current methods. To address this problem, this paper proposed a destriping method called DUMD that is based on unidirectional multiscale decomposition. Large-scale stripes and small-scale stripes were separated by unidirectional multiscale decomposition. Large-scale stripes were eliminated by CCNUC and multiple cumulative error suppression measures, and small stripes were eliminated by a designed filter. Experiments with multiple sets of data showed that the proposed method had good destriping performance, could remove stripes with different scales, and was robust over varying terrains. There are still some shortcomings in this method that warrant future research and improvement: (1) Overexposure features, such as thick clouds, do not satisfy the linear assumption of stripe noise. How to improve robustness with respect to overexposure features is a worthy direction for future research. (2) It can be seen from the experiments that the cumulative error caused by DUMD cannot be eliminated and changes the radiation characteristics of the image itself. How to further reduce the cumulative error is also a question worth studying in the future.