A novel reduced-dimensional dictionary grid screening strategy for SR-STAP

The traditional space-time adaptive processing (STAP) methods require large training samples to estimate the clutter covariance matrix. Sparse recovery STAP (SR-STAP) can get a highly accurate estimation in the case of insufficient samples and it alleviates the above problem of conventional STAP methods, however, the price of SR-STAP is that it involves a huge computational load in optimization, especially when the input data size is very big. This paper proposed a novel SR-STAP dictionary construction method that can significantly improve computational efficiency. We introduce the noise grids as the supplement of the spectrum-aid methods and design a second grids screening method based on the prior spectrum information. After the process of the proposed method, the dimension of the space-time dictionary has a sharp decrease compared to the over-complete dictionary, moreover, the sparsity of variable can be maintained better than the common spectrum-aid methods. The numerical experiments verify that the proposed method not only can significantly reduce the optimization time to far smaller than the SR-STAP but also more robust than spectrum-aid methods even in the small input size.


Introduction
Space-time adaptive processing (STAP) is a mature technique which used in airborne radar for improving the signal-to-clutter-plus-noise ratio (SCNR) of slow-moving targets in severe clutter environments. It can adaptively construct two-dimensional space-time optimal filter to eliminate the ground clutter spreading [1]. However, it is unreal to achieve the optimal STAP in application due to the high computational load and unknown clutter-andnoise covariance matrix (CNCM) [2,3].
To overcome the aforementioned challenges of optimal STAP, sample covariance matrix (SCM) was developed to replace practical CNCM. Reed et al. proved that the requirement of independent and identically distributed (i.i.d.) training samples is at least twice degrees of freedom (DoFs) for the SCNR loss less than 3 dB [4]. Unfortunately, it is hard to obtain sufficient i.i.d. training samples during the short period of detection, especially in a heterogeneous environment [5]. Many algorithms have been designed to reduce sample limitation from different perspectives, which mainly consist of reduced-dimension (RD) and reduced-rank (RR). In RD-STAP, it reduces system DoFs by a linear transformation, which makes the number of training samples decrease to twice adaptive DoFs. Traditional RD-STAP algorithms include extended factored approach (EFA) [6], joint domain localized (JDL) [7], sum and difference beam-forming (SDBF) [8]. Among them, JDL-STAP and SDBF-STAP are sensitive to array error and they need to add additional channels to improve the robustness which leads to the increasing of requirement for training samples, EFA is more robust in the presence of spatial errors but it suffers from performance loss in low pulse repetition frequency (LPRF). As for RR-STAP, it projects the echo onto a lower-dimensional subspace by data-dependent transformation, thus reducing the sample requirement to twice the clutter rank. Typical RR-STAP methods include principal component (PC) [9], cross-spectral metric (CSM) [10], and multistage Wiener filter (MSWF) [11]. Despite their success to reduce the sample requirement, the performance of these methods is still constrained by secondary data in practical application [12].
STAP combined with the sparse recovery (SR) technique has attracted great interest in recent years. SR-STAP methods utilize the intrinsic sparsity of the clutter in angle-Doppler space to recover the clutter with advanced SR algorithms, and it can achieve a perfect estimation of CNCM under insufficient training samples condition, even in a single sample. Though the advantages of above SR-STAP methods have great potential, an unacceptable heavy computational burden is still required, which makes SR-STAP is unrealizable for application. Generally speaking, SR algorithms approximate the optimal solution by successive iterations, which causes a large time cost in repetitive operations, particularly in the region near to the ideal solution. Some studies have been devoted to reducing the computation complexity. Fast sparse Bayesian learning (SBL) turns sparse problem to maximum a posteriori estimation (MAP) and searches sparse coefficients according to Bayesian framework [13], it gets rid of the influence of regularization parameter. RD method transforms the global optimization problem to local optimization with RD techniques such as beam-space post-Doppler [14], which cuts down the dimension of variables and has a significant computational savings. Spectrum information can also be used to design a more suitable dictionary [15,16], this way is easy to achieve and reduces iterating time significantly, however, a relatively small number of atoms may destroy the stability of algorithms [17,18].
Consider the complexity limitation of hardware in practical application, the spectrum-aid method seems the most feasible way of the above fast methods, however, its performance and stability cannot be satisfied with the requirement of the application. This is because the sparsity of the residual grids after the spectrum-aid method is weaker than before. In this paper, the noise grids are introduced as additional information to improve the sparsity of the spectrum-aid methods and hence a novel RD dictionary is proposed, but different from the conventional overcomplete dictionary, only a few noise grids are used in the proposed method. It should be pointed out that our method is to focus on the improvement of the SR-STAP model but not on solving algorithms. Based on this core idea, we specially design a second grid screening structure to make the best use of spectrum information and get rid of the useless grids as much as possible. In the proposed algorithm, the region segmentation criterion is redesigned and the whole angle-Doppler space can be separated into three different areas according to the contribution toward the clutter. Then three different selection strategies are independently performed in their areas to select suitable space-time grids. Numerical experiment results verify the performance of the proposed method.
The remainder of this paper is arranged as follows: Sect. 2 establishes the basic signal mathematical model, which mainly includes the system working environment and the signal components of radar echo. The conventional STAP is also introduced in this section. Section 3 describes the method of SR-STAP and Sect. 4 makes a analysis of the effect of clutter grids and noise grids for performance. Section 5 gives a detailed description of the proposed RD dictionary construction algorithm. In Sect. 6, the numerical experiments are executed to verify the performance of the RD dictionary. Finally, Sect. 7 is the conclusion of the paper.

Basic signal model and conventional STAP principle
Consider a phased-array pulsed Doppler (PD) radar is placed in a moving platform which keeps a constant velocity v flying along with x-axis in the height of H as shown in Fig. 1. PD radar transmits M pulse with fixed PRF f r in each coherent processing interval (CPI), and the transmitted signal is a narrow-band modulating signal with the carrier Fig. 1 The geometry of airborne radar illumination frequency f c . Assume that the structure of the antenna array is a uniform linear array (ULA) which includes N omnidirectional ideal elements, and these antenna units are uniformly distributed in the array with the interval d (usually it is equal to half wavelength). Here, our analysis focuses on side-looking mode which means crab angle p = 0 • , and the array length L a is far smaller than the slantrange R which is the distance between the platform and scatter, hence all the elements can be viewed as having the same incidence direction with the reference element. , and are the cone angle, azimuth angle and elevation angle of the array relative to the ground unit, respectively.
Through the receiving process and matched filter, the multi-channel echo signal can be organized as a N × M × L three-dimensional data matrix, where N, M and L are the number of spatial elements, slow beat (pulse number) and snapshot (range samples), respectively. For a given range cell, the target detection problem can be demonstrated as a binary hypothesis testing problem as follow Hypothesis H 1 denotes the case of target presence and H 0 denotes the case of target absence. x c ∈ ℂ NM×1 is the clutter component of the echo signal, and x t ∈ ℂ NM×1 denotes the target signal component. n is NM × 1 white Gaussian noise vector. Note that the interference signal is not considered in Eq. (1) because it is not a key point in this paper. J.Ward clutter model [19] thought that the received signal could be expressed as the sum of reflection of small clutter patches when they are independent of each other. Based on this assumption, if no ambiguity appears in the current range cell, x c can be approximately written as where N c means the whole number of divided clutter patches. i is the complex scatter coefficient of the ith clut- antenna element and pulse. The formation of x t is similar to the x c and can be written as where i means the complex amplitude for ith moving target in the current range cell. The only difference between (4) and (2) is that N t and it means the moving targets count in CUT. Finally, the signal echo of H 1 is rewritten as and the NM × NM CNCM R of the received signal in the current snapshot is where | ⋅ | denotes the absolute value operation and (⋅) H is the conjugate transpose operation. 2 n represents the noise variance of the current test cell. I is NM × NM identity matrix in which diagonal elements are one. It should be noted that the targets are not included in CNCM. For the conventional STAP methods, the adaptive filter weight can be designed according to represents the normalization factor of weight. R −1 denotes the inverse of CNCM, which is mainly used to form the hollow of original clutter, s t ∈ ℂ NM×1 represents the expected target direction in the angle-Doppler domain, here the spatial angle is the azimuth angle of illumination and the Doppler frequency is the tested Doppler channel. In practice, R is unknown, and it is usually be replaced by its estimation R which is calculated from the secondary samples near the CUT by where E{⋅} is the expectation of range samples. Note that this estimation method supposes that training samples are satisfied with the condition of i.i.d, when in a nonideal environment, it would suffer from performance loss.

SR-STAP method
Different from the traditional STAP method, in SR-STAP operations, the continuous angle-Doppler plane is discretized into N d N s grid points ( N d points in Doppler axis, N s points in angle axis), usually N s = s N and N d = d M , where s and d are the coefficients of expansion in angle dimension and Doppler dimension. Note that grids number need to meet N d N s ≫ NM but it is not a strict condition. Then all the space-time steering vectors of grids can be used to construct the N d N s × L overcomplete dictionary set and it is where s n d ,n s is space-time steering vector corresponding to the grid with kth spatial angle and jth Doppler frequency. Then echo Eq. (2) can be rewritten according to (9) as follow with where is the complex amplitude of grids in the angle-Doppler domain. The SR-STAP problem can be modeled as where min (⋅) means to solve the minimum value of the object function. �� �� 1 = ∑ N d N s i=1 � i � denotes the l 1 -norm which is used to describe the sparsity of recovery, and || ⋅ || 2 is the l 2 -norm which measures the similarity between the recovery solution and original data. is noise threshold which is usually small value in optimization. One thing is worth pointing out that the original sparse recovery problem is l 0 -norm, here we take l 1 -norm relaxation [20]. Finally, the estimation R sr of SR-STAP can be easy calculated by where diag(⋅) is the diagonalization operation. The rest of the steps is the same as the conventional STAP filter design in Eq. (7).

The contribution of different type grids
The essence of (12) is to solve the underdetermined equation (N d N s ≫ NM) , and numerous solutions are existing in the SR-STAP problem. However, in the application, a recessive restraint (clutter ridge is only occupying a little area relative to the whole space-time plane [21]) from the real environment will be added in raw data, which means that is high sparsity and can be recovered by SR theory. Overcomplete dictionary plays a vital role in recovery, however observe Eq. (10), it should be known that the redundancy of is significant increasing with the rise of s and d , which leads to the large computation. In fact, based on the prior knowledge of clutter ridge, we know that the most atoms of has a little effect in performance, so whether there exists a way that reduce the number of atoms in while keeping the sparsity of signal. Generally speaking, it needs add more grids in clutter area and less grids in noise area as in [16]. However when the noise girds reduce to zero [15], it will lead to a poor recovery due to the decreasing of sparsity. This illustrates that clutter atoms number represents the accuracy of recovery and noise atoms are used to stay the sparsity of signal. Therefore, a more concrete selecting principle can be given, to increase the accuracy of recovery, it needs to select more clutter grids, but it also needs to add more noise grids to keep the sparsity of signal. Following RD dictionary design is mainly based on the above principle.

RD dictionary construction strategy
As shown in Fig. 2, the proposed structure of RD dictionary consists of four parts including spectrum extracting, region segmentation, grids screening and RD dictionary construction. The prior knowledge which hides in the echo signal is extracted in the step one and then is used to exclude useless grids in step three. Finally, in step four,the RD transformation matrix which is formed by the saving grids. Following is the detailed description for each part.

Extracting spectrum information
Spectrum analysis offers a set of spectrum transformation methods to convert the time signal to frequency domain signal, such as Fourier transform (FT), Capon transform, Multiple Signal Classification (MUSIC) and Estimation of Signal Parameters via Rotational Invariance Techniques (ESPRIT) [22]. In these algorithms, FT is the commonly used method in practical application due to the low-complexity. In the angle-Doppler domain, it is similar to the FT, the angle-Doppler transform can be achieved with the space-time steering vector basis as where R ∈ ℂ NM×NM is the estimation of the clutter covariance matrix with the current CUT and s k,j ∈ ℂ NM×1 denotes the space-time steering vector basis which is corresponding to the kth and jth angle-Doppler grid. Different from the conventional FT method, Eq. (14) can calculate any grids in the angle-Doppler plane, however, its precision is low. Capon transform is more accurate than FT, which can be expressed as follow The disadvantage of Capon method is that the big computation of the inverse operation, so it only suits to the situation which is small input data or weak real-time requirement. As for other high-accuracy algorithms, they all not suitable here due to the large running time.
Another alternative method is the least square (LS) algorithm which based on the optimal theory [23], and its spectrum estimate problem can be described as where P is the vectorization operation of power matrix P .
x represents the raw echo data. It can be solved directly according to general expression without searching. Let the minimum value of object function (16) is equal to zero, and then the simplified solution is where mat[⋅] is the transformation operation, which means the N d N s × 1 vector is converted to the N s × N d matrix.
After that, the threshold segmentation method is used to extract spectrum clutter distribution [24]. Consider that the main aim in this step is only to acquire the rough area of clutter rather than the exact distribution, so the mean of angle-Doppler is enough to use, and it is written as Some suggestions need to be mentioned. In a real application, it usually needs to process large data at the same time. For accelerating the transformation, (14) can be solved by fast Fourier transform (FFT), the price of this is that it only can calculate the uniform spacing grids. The most time cost step of (15) and (17) is the inverse operation which can be calculated according to LU factorization. On the other hand, a more accurate threshold method only has a little improvement in final result but it would spend much time on computation, so a complex segmentation method is not be recommended here.

Region segmentation
In this section, the slide-window process structure is chosen to achieve the region segmentation. As shown in Fig. 3, four vertices of the window are selected as test grids which are used to detect the type of area. The angle-Doppler plane can be separated into three types according to the test grids. First is the main clutter area, it represents that all the grids in this window are clutter grids of spectrum distribution, on the contrary, the noise area denotes that all the grids are noise components. As for the clutter edge area, it lies in the boundary between clutter and noise, so this region both includes clutter and noise components. In processing, the test window slides successively in the whole plane, and check the type of current area by test grids. The second division aims to shrink the clutter and noise area, and isolate weak effect grids from above region.

Grids screening
For the main clutter region, it includes the most effective clutter grids, in fact, the real ground clutter points mainly locates in the clutter ridge (slash in side-looking mode, curve in non-side-looking mode) but it would spread into neighbor girds due to the grid mismatch problem. The slight change of these grids will cause a severe influence on recovery. Therefore it best to save all the space-time grids in this region. The clutter edge belongs to the junction of noise area and main clutter area, and it occupies a relatively small area to the whole angle-Doppler space. This means that only the marginal clutter grids are included in this area, i.e., its clutter power is lower than the main clutter power. Therefore the clutter edge has a weak performance improvement in the final accuracy of recovery. Above idea is established in the assumption that the number of edge grids is far smaller than the number of main clutter grids. However when this assumption is not satisfied in the application, the importance of edge grids would enhance with the increase of grids number and it cannot be ignored in this case. Fortunately, the number of the discretized grid is usually large in the most radar system. Therefore, to reduced the computation complexity, the clutter edge grids are directly discarded in the proposed algorithm, even this operation would have a little performance but it further improves the efficiency of the proposed method.
Noise area is the most important region in the screening strategy because its grids number is far bigger than other regions, this also means that the process method of it finally decides the size of the RD dictionary. Noise girds are employed as auxiliary grids to improve the sparsity of the main area but the nature of them is useless grids, hence it is reasonable to only utilize a few grids of noise area. Assume that the sliding-window is small enough compared to the angle-Doppler domain, and slow change occurs in neighbor grids, then the performance of the window can be approximately represented by part of grids in the window, it also can be known the locations of grids have no effect in final result according to the assumption. Here, we define the deletion ratio of grids, and it can be expressed as where is deletion ratio. n o represents the deletion grids in the noise area. When window length a is constant, the lower boundary of is 1 a 2 which means that only one grid is used in the window. Combined with the previous analysis, the random method is designed to delete grids. First is to generate an index of n o deletion grids with random number generation (RNG) algorithm, and then exclude grids by index. In application, it could take a fixed index to replace the random index, the advantage of it is no need additional algorithm complexity but it will require more storage memory.

RD transform
Through the screening process, the residual grids are the atoms of RD dictionary. It can directly construct RD dictionary by transformation matrix as follow where T is the N d N s ×n logical reduced-dimensional transform matrix and n denotes the number of atoms of RD dictionary, for each column, the 1 only appears once in the index of saving grids (here each column only corresponds to one saving grids) and other locations are 0. However in practical application, the transformation matrix may be large due to the high-dimension of raw data. Here, we offer an easy way to solve this challenge, first is to prepare an additional matrix ̃ with the size N d N s ×n , then copy the column of the overcomplete dictionary to ̃ one by one according to the index of remaining grids. The advantage

Numerical experiments
In this section, the numerical experiments are executed to examine the performance of the proposed method. The simulation data is generated by adapting the J.Ward clutter model [19] and some key radar system parameters are summarized in Table 1. In order to simulate the real illumination, we only consider the azimuth angle from −90 • to90 • and set the number of clutter patch N c = 180 . In SR-STAP methods, the expansion factor is 6 and the noise level is 10 −6 . The number of i.i.d training samples is 2. For conciseness, the proposed method is distinguished from other SR-STAP methods by using the prefix RD, hence three methods are called RD-FT-SR-STAP, RD-LS-SR-STAP and RD-C-SR-STAP, respectively. Moreover, the deletion ration is fixed in 3 4 and L n is 2. CVX [25] is selected as the solving tool for SR-STAP in the following experiments. To evaluate the performance of SR-STAP, the signal-clutternoise ratio (SCNR) loss is used as the criterion and it can be calculated by 10 random trials for each algorithm are run to guarantee the accuracy of experiments.

The comparsion of SCNR loss
In the first experiment, we verify the SCNR loss of the proposed methods. As shown in Fig. 4, the notch of C-SR-STAP and FT-SR-STAP are far wider than other methods. This is because the residual grids of the above methods are insufficient to recover the CNCM. As for LS-SR-STAP, it gets better performance than C-SR-STAP and FT-SR-STAP for the reason that it includes more grids due to the low-resolution of the LS method. However compared with SR-STAP, it can be found that the above spectrum-aid SR-STAP methods have a severe degradation in performance. For the proposed methods, the SCNR loss of them only has a little bit down to the SR-STAP, but they still higher than the spectrum-aid methods. This result indicates the fact that the pure clutter grids cannot support the recovery of CNCM and utilizes the portion of noise grids as supplementary information can effectively improve the performance of spectrum-aid SR-STAP methods.

The effect of window size and detection ratio
In this experiment, we analyze the influence of window size and deletion ratio for the proposed method. As shown in Fig. 5, it is the SCNR loss of different window sizes. For three methods, their performance improve with the decreasing of the window size and they acquire the best result when L n is equal to 2. However, when L n increase to 5, they almost lost the ability due to the insufficient grids of the noise area. In these methods, the RD-LS-SR-STAP seems to have better performance than the other two methods but this is because the LS dictionary has already included more grids. As shown in Fig. 6, in order to evaluate the effect of the deletion ratio, we specially set L n = 3 and take the deletion ratio = [ 5 9 , 6 9 , 7 9 , 8 9 ] . There is almost no remarkable change for three RD dictionaries in small window size but it should be known that the effect of will be intensified with the increase of L n . Therefore, it can be found the window size is the main factor that decides the performance of the proposed methods. According to the above analysis, we recommend the combination of L n = 2 and = 3 4 and this can cut down the most noise grids in the application.  Fig. 4 The SCNR Loss of different dictionaries

The effect of spatial error
In this section, we focus on the stability of the proposed methods in existing spatial error. The Guerci error model [3] is used to generate a spatial error matrix and it can be described as following where f u ( a,n ) and f u ( p,n ) are the uniform distribution Probability Density Function (PDF) of amplitude error and phase error, respectively. In practical application, the gain error can be easily compensated by the preprocess, thereby this experiment only considers the case of phase error and the final spatial error can be expressed by As shown in Fig. 7, the phase error fluctuation p,n in here is set 5 • . It is not surprise that all the SR-STAP methods have a performance loss in the error case. In these methods, FT-SR-STAP and C-SR-STAP have a rapid decrease in performance, hence the spectrum-aid methods are instability to the error condition. In contrary to the above two methods, the proposed methods only have a slight decrease compared to the no error case. This demonstrates that the proposed methods are more robust than spectrum-aid methods, moreover, it also suggests that the noise grids have a great influence in improving stability.

Real measurement data experiment
In order to verify the effectiveness of the RD dictionaries in real environment, we execute SR-STAP in the MCARM data [26]. In the 1990s, the Rome Lab and Northrop Grumman Company developed the MCARM plan, an L-band radar was installed in the BAC-111, and it recorded different ground types echo, these radar data can reveal the property of the real environment. The antenna array has 16 columns and each column has 8 elements, every four elements are organized into one channel. The index of using data is RL050575. The crab angle of the flight is   raw data as shown in Fig. 8, and the red circles represent the location of simulated moving targets.
The suppression of different SR-STAP methods in MCARM data is shown in Fig. 9. It should be pointed out that the spectrum-aid methods fail to recover the CCM in = 2 . The range-Doppler spectrum of SR-STAP is shown in Fig. 9a and the other three RD-SR-STAP methods are shown in Fig. 9b-d, respectively. It seems that there is no visible difference between the SR-STAP and the proposed methods, which illustrates that only adds a few noise grids to spectrum-aid methods can realize the performance of an overcomplete dictionary. Furthermore, we display the output power of the different range bin in normalized Doppler frequency 0.223 as shown in Fig. 10. It can be easy to find that the moving target emerges in the range cell 313 and the surrounded clutter is eliminated after STAP. For three methods, the RD-FT-SR-STAP and RD-LS-SR-STAP are better than RD-C-SR-STAP. This means that fewer grids cannot effectively recover the CCM even it is combined with noise grids. Therefore, the capon method is not suitable in application because of its computational load and poor suppression.
The computational complexity of the proposed methods is evaluated in the last experiment and the expansion factor is used as the variable to change the size of input data. As shown in Fig. 11, the red bar areas represent the addition grid number of SR-STAP compared to the other methods. This means that the SR-STAP needs large grids to recover CCM perfectly. As for the proposed methods, the grid number only has a slight increase compared to spectrum-aid methods but far smaller than the SR-STAP. Next is the running time as shown in Fig. 12, and the symbol 'X' represents that a failed recovery in current . The spectrum-aid methods often lost efficacy in the small input size, especially the C-SR-STAP cannot normally work in such conditions. Therefore, pure spectrum-aid methods lack practical significance. The proposed methods can still stay the effective recovery in all the input size, which means they are more robust than spectrum-aid methods in real data. The only price is a slight time extension, but the proposed methods can be easier to accept in using than SR-STAP and spectrumaid methods.

Conclusion
For conventional SR-STAP methods, the complete dictionary leads to a huge computational load, and this limits the application of SR-STAP methods. The spectrum-aid methods have a great effect in reducing the computation of SR-STAP, but they are sensitive to the input data size. Except that, the spectrum-aid methods have poor performance because they discard lots of grids in recovery. In order to solve the above problems, the proposed methods first introduce the noise grids as the auxiliary information to enhance the sparsity of the recovery variable. Furthermore, a second grids screening method is developed to reject useless grids as much as possible. The exhaustive experiments demonstrate that it is feasible to get a similar recovery as the overcomplete dictionary by adding a few noise grids. Although we give three methods to extract spectrum information, the experiments show that the RD-FT-SR-STAP is more suitable than other two methods because of its balance between performance and computational load. In a word, the proposed method has a significant improvement and its computation is only a slight increase compared to the spectrum-aid methods, hence the proposed methods are more realistic significance. However, the proposed method only adopts the simple random sampling method to eliminate the useless noise grid points, hence a more appropriate screening scheme of noise grids is worth exploring. Furthermore, extending this method to other existing solving algorithms is also one of the main works in the future.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creat iveco mmons .org/licen ses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.

Data availability statement
The MCARM data used to support the findings of this study are included within the supplementary information file(s).

Compliance with ethical standards
Conflict of interest On behalf of all authors, the corresponding author states that there is no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.