Sequence set design for waveform-agile coherent radar systems

With increased degrees of freedom of the transmitter, coherent waveform-agile radar system can change its transmission on-the-fly in response to target detection’s requirement. This approach can provide better performance than a single waveform. In this paper, we consider unimodular sequence set design (USSD) problem based on summed range-Doppler ambiguity function (SRDAF) for coherent waveform-agile radar system. We present two algorithms for constructing unimodular sequence sets under the constant module constraint with desired minimized sidelobes on a predefined area of range-Doppler plane. The proposed algorithms are constructed based on singular value decomposition (SVD) and cyclic algorithm (CA), respectively. Numerical examples show the effectiveness of the proposed algorithms.


Introduction
In the last decade, radar systems have been able to take full advantage of degrees of freedom of the receiver to improve performance. The emergence of prototype radar systems equipped with highly agile, software-driven waveform generators has provided the ability to change the transmit waveform at each time step to match environments and sensing objectives, such as increased signal-to-noise ratio (SNR), reduced estimation errors, or increased collection of information. Since there is an infinity of possible waveforms, it becomes critical to select or optimize the transit waveform at each time step. In coherent radar detection applications, the waveform-agile transmitter optimizes the waveform on a pulse-to-pulse basis and the radar system can obtain the information by optimally improving matched filter's response of the target.
In the literature, extensive research has been investigated single waveform design from the viewpoint of autocorrelation function and ambiguity function . Some studies have also been focused on the choice of agile waveform in sensing. In the works [26][27][28][29][30][31][32][33], the dynamic selection of waveforms for target tracking was considered and the optimal waveform parameters were derived for tracking target motion using a linear/nonlinear observations model with detection in a clutter/clutter-free environment. These approaches can provide performance improvements over waveform optimization (2020) 2020: 31 Page 2 of 17 or improved tracking algorithm. Agile waveform can also be designed from the viewpoint of ambiguity function. The authors of [34] demonstrated that suitably transmitted and processed radar waveforms based on Golay sequences provide new primitives for adaptive waveform transmission. The adaptive transmission enables improved detection and finer resolution, while managing computational complexity at the receiver. Actually, Golay complementary sequences [35] can be seen as the simplest sequences in waveform-agile coherent radar system, and there is extensive literature because of the importance of such sequences in communications, coding theory, cryptology, and radar [36][37][38][39][40][41].
It should be noted that unimodular (i.e., constant modulus) sequences with good autocorrelation properties are useful in several areas, including communications and radar. The integrated sidelobe level (ISL) of the correlation function is often used to express the goodness of the correlation properties of a given sequence. In the references [3,9,25], several cyclic algorithms for the local minimization of ISL-related metrics are presented. The power-like method are used for synthesizing sequence with good correlation and can also be extended for synthesizing sequences with good autocorrelation and cross-correlation functions in multiple-input multiple-output (MIMO) radar. However, these methods can not be directly utilized for optimizing sequences in waveform-agile coherent radar because they can only optimize range sidelobes level for correlation function.
In the aforementioned works, we have only considered ambiguity function synthesis for a single waveform, and little attention was paid to waveform design of coherent waveformagile radar based on summed range-Doppler ambiguity function (SRDAF). In addition, sidelobe shape control has also not been considered. In this paper, we investigate sequence set design problem with constant module, i.e., unimodular sequence set design (USSD) problem by SRDAF mathematical tool. Furthermore, we also consider that minimizing sidelobes on range-Doppler plane around the origin can improve the detection performance for closely spaced targets. We present two algorithms for constructing unimodular sequence sets under the constant module constraint with desired minimized sidelobes on a predefined area of range-Doppler plane. The proposed algorithms are constructed based on singular value decomposition (SVD) and cyclic algorithm (CA), respectively, and can be convergent by several hundreds of iterations.
In [42][43][44], coordinate descent (CD) algorithm is an important search algorithm and can solve optimization problems by successively performing approximate minimization along coordinate directions or coordinate hyperplanes. It has been used in applications for many years, and their popularity continues to grow because of their usefulness in data analysis, machine learning, and other areas of current interest. CD approach has been successfully applied for discrete-phase sequence design and minimum ISL code sequence design. In this paper, the performance of CD approach is compared with the proposed algorithms.
The rest of this work is organized as follows. Section 2 discusses the mathematical model and ambiguity function of coherent waveform-agile radar system. Section 3 proposes two optimization algorithms for optimizing the shape of SRDAF based on SVD and CA, i.e., SVD-USSD and CA-USSD algorithms. Several numerical examples are presented in Section 4. Finally, concluding remarks and directions for future research are presented in Section 5.

Mathematical model and ambiguity function
A radar system transmits a coherent pulse train consisting of N pulses, which can be expressed as the following general form: where s n (t) is the complex envelope of nth transmitted pulse and T r is the pulse repetition interval (PRI). We assume that s n (t) = 0, ∀t [ 0, T], where T is the pulse duration. s n (t) is modulated by the following way where M is the subpluse length, x(n, m), n = 0, 1, ..., N − 1, m = 0, 1, ..., M − 1 is the modulating code sequence that is to be designed, and is an ideal rectangular shaping pulse of time length t p ( and thus T = Mt p ). It is usually desired to transmit unimodular sequences where {φ n,m } are the phases. In common pulsed Doppler radar, the pulses with the same modulation, i.e., s n (t) = s 0 (t), ∀n, are transmitted consecutively in coherent processing interval (CPI). In this paper, we investigate the ambiguity function properties of waveform-agile pulsed Doppler radar. The first processing step of pulsed Doppler is to implement matched filtering of the received echo. For a point target for nth pulse, the matched filter output can be given by where τ and ν are time delay and Doppler frequency shift of the target, respectively. We show that if, for all relevant Doppler frequency shift −1/2T r ≤ ν ≤ 1/2T r , the pulse duration T is small enough such that |2πνT| < π/5, which suggests T ≤ T r /5, then the phase change in the received echo that is caused by the Doppler frequency of target motion can be neglected. Additionally, consider the time grid τ = pt p , p = −(M − 1), ..., 0, ..., (M − 1) whose points are integer multiples of the subpulse length t p . It is not difficult to calculate χ n (τ , ν) at τ = pt p , which can be expressed as Note that T r t p , and the Doppler unambiguity constraint −1/2T r ≤ ν ≤ 1/2T r , we obtain The second step for pulsed Doppler radar is to process the coherent pulse train, which is based on discrete Fourier transformation (DFT). To perform coherent processing of a particular delay, we must collect N samples at a certain delay from last N pulses. Therefore, we construct a two dimensional (2M − 1) × N array of matched filter output: first dimension is the delay, second dimension is the pulse number. The IDFT operation is then performed on N samples of a row of the matched filter output, corresponding to a certain delay. The pth IDFT processor output can be expressed as Consider the frequency grid v = q/NT r , where q = 0, 1, ..., N − 1. For convenience of expression, we ignore the constant in (8) and rewrite (8) as It is no surprising that an ideal ambiguity function should have a high narrow peak in the origin and zero sidelobes everywhere else. Note that the target Doppler frequency ν only affect the high peak position in Doppler frequency axis, and the sidelobe energy is constant with all relevant ν; thus, we can set ν = 0 and define a range-Doppler ambiguity function of pulse train, i.e. SRDAF, as where ρ n (p) = M−1 m=0 x * (n, m)x(n, m − p) denotes the auto-correlation coefficient of the modulation code of the nth pulse.
It is easy to verify ρ n (0) = M, ∀n. Thus, we have ξ(0, 0) = NM, and ξ(0, q) = 0, q = 0. Note that the shape of ξ(p, q) is uncorrelated with the modulation codes {x(n, m)} and ξ(0, 0), which indicates the peak value of SRDAF, is determined by the modulation code energy. Additionally, for p = 0, the shape of ξ(p, q) is determined by the modulation codes {x(n, m)} of N pulses. Therefore, we can control the sidelobe of a a predefined area on range-Doppler plane by optimizing the modulation codes {x(n, m)}. For instance, the range and Doppler frequency of the observed target are a prior known by prescan or a prior information, and we can minimize the sidelobe of ξ(p, q) in that area. Assume I is a subset including the positions of (p, q) on range-Doppler plane for minimizing sidelobe and excluding positions on zero delay (i.e., p = 0), the partial sidelobe minimization of ξ(p, q) can be described by Considering the objective function in (11) is a quartic form, which is difficult to tackle, we try to transform the above minimization problem into a quadratic form.

Optimization method
The optimization problem requires N transmitted waveform {s n (t)} to be optimized. It is hard for us to optimize them simultaneously. By utilizing the alternating direction method, we can optimize the nth waveform s n (t) or its modulation code {x(n, m)} M−1 m=0 while the others are fixed. Hence, the optimization problem in (11) can be transformed to N optimization problems, and the nth problem can be given by where x n =[ x(n, 0) x(n, 1) ... x(n, M−1)] T , and (·) T denotes the transpose operation of a vector/matrix. In the following subsection, we will introduce an unimodular sequence set design method based on singular eigenvalue decomposition, namely SVD-USSD.

SVD-USSD method
By substituting (10) into (12), we optimize (12) by its autocorrelation coefficient {ρ n (p)} instead of x n . The corresponding optimization problem can be written as Suppose that the sidelobe of AF shape on a region of range-Doppler plane is to be minimized, and q(p) ⊂[ q p,min , q p,max ], where q p,min and q p,max denote the bounds of the selected Doppler bins of the pth range bin.
By further utilizing alternating direction method, the optimization problem P n can be split into multiple minimization problems, and the pth sub optimization problem can be described by where a p = e j2πnq p,min /N ... e j2πnq p,max /N T , and b p = −[ N−1 k=0,k =n e j2πkq p,min /N ρ k (p) ...  where (·) H denotes the conjugate transpose operation, N p is the number of the selected Doppler bins of the pth range bin. Let With the obtained ρ n (p), we can construct a matrix, which is given by Because X H n X n = G n indicates that the modulation code {x(n, m)} M−1 m=0 is the optimal solution for (14), we can think of designing {x(n, m)} M−1 m=0 by minimizing the following criterion over the set of unimodular sequences. This optimization problem can be approximated by a simpler criterion where U is a (2M − 1) × M unitary matrix, i.e., U H U = I. The design problem associated with (19) can be stated as follows Regarding the minimization problem in (20), we note the following facts. For given X n , let Then, the solution of U for fixed X n is given by Let V = UG Considering the unimodular constraint, the minimizer of (20) can be given by The SVD-USSD method is summarized in Table 1.

Remark 1
In SVD-USSD optimization procedure, iterations are required to be performed and convergence can be obtained after a number of iterations. Note that eigenvalue decomposition (EVD) and SVD operations are performed once in every iteration step and hence the computational complexity of this algorithm is a little higher than expected. In practical occasions, where low complexity algorithms are more popular, algorithms without EVD and SVD operations can be preferable by radar system designer. In the next subsection, we will introduce a low complexity algorithm, i.e., CA-USSD algorithm.

CA-USSD method
The objective function in (12) is a quartic form, which is relatively difficult to obtain the global optimum by the analytical expression or the optimization method. In this section, we expect to find the local optimum for the problem in (28) and propose a computationally efficient approach.
Let U (n) p,q = e j2πnq/N T p , where T p is a M × M shift matrix and defined as Table 1 The SVD-USSD algorithm Step 1: set k = 0, initialize modulation codes {x(n, m)}, n = 0, 1, ..., N − 1, m = 0, 1, ..., M − 1 using randomly generate codes; Step 2: for n=0:N-1 Step 2.1: solve P1 Step 2.2: solve P2 n ; for fixed X n , perform SVD operation on X n G  The objective function in (11) can be expressed as where and The optimization problem in (12) can be transformed to where the M × M Hermitian matrix Q n is positive definite.
Considering that the diagonal loading does not change the optimum solution of (29), this minimization problem can also be transformed to where Q n = Q n − λ c I M is a negative definite matrix if λ c > λ max (Q n ) + 1. Here, λ max (·) denotes the largest eigenvalue, and we have Q n < −I.
In the following, we try to optimize (29) in a cyclic way. Let x (t) n and Q (t) n be the obtained sequence and matrix at the tth iteration.

Lemma 1 The following two inequalities
can be obtained if the following condition is satisfied The proof of Lemma can be seen in [13].
Proof Because of the negative definiteness of the matrices Q (t) n and Q (t+1) n , we can have and With the above two inequalities, we can obtain and If Re x n are satisfied at the same time, we can have and Note that n , therefore we can conclude

Remark 2 Lemma 1 indicates the requirement for the two inequalities, and it also states that the objective function x
n is convergent and approaches the minimum value if (32) is satisfied. In the next step, the relationship between the above two inequalities in (31) is shown in Lemma 2.

Lemma 2 The two inequalities x
n establish simultaneously if the sequence x (t+1) n can be obtained by the following minimization problem min Proof By expanding the objective function in (40), we obtain is the minimizer of Re n . Hence, the second requirement can be satisfied by solving the minimization problem in (41). The solution for this problem can be easily given by .
Considering the solution x Considering the constant modulus constraint, the following operation can be added after (42) x where arg(·) denotes the argument of a complex number. The obtained unimodular sequence consists of different phase values, arbitrarily generated by [ 0, 2π), and it is a suboptimal solution to (40).
The CA-USSD method is summarized in Table 2.

Convergence analysis:
The CA-USSD algorithm given in Table 2 is based on the general cyclic scheme; thus, according to Section 3.2, we know that the sequence of objective values evaluated at generated by the algorithm is nonincreasing. And it is easy to see that the objective value is bounded below by 0; thus, the sequence of objective values is guaranteed to converge to a finite value.
In this part, we will further analyze the convergence property of the sequence generated by the CA-USSD algorithm and show the convergence.
According to Lemma 1, x n is satisfied.

end
Step 3: set k = k + 1, repeat the Step (2)   We can further note that and we have Q (t+1) n < Q (t) n < −I. We can also note that and verify the negative definiteness of the matrix Q n . The results of (45) and (46) highlight that if the equality x (t+1) n are correspondingly obtained. Therefore, the objective function n in iteration procedure. Due to iterative calculations, the obtained sequence x (t) n usually does not satisfy the constant modulus constraint.

Computational complexity
In the SVD-USSD procedure, the computational complexity of G where K is the total number of iterations. In the CA-USSD procedure, the overall computational complexity is O(KL(4M 2 + 2M)), where L is the number of the selected range and Doppler bins. It is obvious that the computational complexity of CA-USSD algorithm is obviously lower than that of SVD-USSD algorithm.

Simulation results
In this section, we provide several simulation examples to demonstrate the performance of the proposed methods, SVD-USSD and CA-USSD. In the following examples, it is assumed that the radar pulse width T is 10 us, and the PRI equals to 1 ms to ensure the Doppler frequency insensitivity of the transmit waveform. The coherent radar system transmits N = 64 pulses in a CPI and M = 100 subpulses in a pulse. To evaluate the shape of SRDAF of the waveform-agile radar system, a point-like target is considered in the simulations. In SRDAF image, the time delay axis is normalized by subpulse duration t p , and the Doppler frequency shift axis is normalized by T r . The convergence of the proposed algorithms will be tested by using randomly generated sequences in the initialization.
In this example, the shape of SRDAF is desired to have sidelobes in an interested area as low as possible. Therefore, we suppose that = {(p, q)||q| ≤ N/4, (p, q) = (0, 0)} is the (2020) 2020:31 Page 12 of 17 area on range-Doppler plane, which is close to the origin but excludes the origin. With randomly generated sequences in the initialization, the proposed algorithms in Section 3 is utilized to minimized the ISL of the SRDAF of the synthesized sequences. The shape of SRDAF of sequence set optimized by SVD-USSD and CA-USSD algorithms are shown in Fig. 1a and b. The sidelobes of the SRDAFs of the sequence sets synthesized by SVD-USSD and CA-USSD algorithms are suppressed to about − 50 dB in the interested area on range-Doppler plane . In Fig. 1c, the curves of the objective functions vs. iteration number in SVD-USSD and CA-USSD algorithms are given. The convergence performance of these two algorithms are almost the same and can converge to a stable point after 300 iterations. However, SVD-USSD algorithm shows faster convergence speed.
To evaluate the sidelobe synthesis performances of the proposed SVD-USSD and CA-USSD algorithms , we define the average sidelobe level and peak sidelobe level by , and PSL = 20log 10 ξ(0, 0) max (p,q)⊂I |ξ(p, q)| (dB) Figure 2 plots the curves of average sidelobe level vs. sequence length and sequence number of the sequence set optimized by SVD-USSD and CA-USSD algorithms after 200 iterations, i.e., average sidelobe level vs. sequence length M, average sidelobe level vs. sequence number N, peak sidelobe level vs. sequence length M, and peak sidelobe level vs. sequence number N, respectively. As sequence length M increases from 20 to 100 and sequence number N equals to 64, ASL is decreased from − 40 dB to − 48 dB for CA-USSD algorithm and decreased from − 44 dB to − 51 dB for SVD-USSD algorithm. Meanwhile, PSL is decreased from − 57 dB to − 68 dB for CA-USSD algorithm and decreased from − 62 dB to − 75 dB for SVD-USSD algorithm.
As sequence number M increases from 10 to 130 and sequence length M equals to 100, ASL is decreased from − 36 dB to − 48 dB for CA-USSD algorithm, and ASL is decreased from − 40 dB to − 55 dB for SVD-USSD algorithm. Meanwhile, PSL is decreased from − 52 dB to − 72 dB for CA-USSD algorithm and decreased from − 59 dB to − 79 dB for SVD-USSD algorithm.
In the simulation result of Fig. 2, CD approach is utilized for synthesizing unimodular sequence set with 4PSK, 16PSK, and 64PSK modulation when the difference of objective function between two adjacent iterations is small enough. In the iteration process, the parameter is set to be 10 −4 . Due to the limited number of discrete phases, the ISL and PSL performance of CD-4PSK, CD-16PSK, and CD-64PSK is higher than that of SVD-USSD and CA-USSD algorithms. However, with increased number of discrete phases, the ISL and PSL performance is improved. The performance of CD-64PSK is close to that of SVD-USSD and CA-USSD algorithms. In Fig. 2b and d, we can also find that the performance of CD-64PSK is even better with small sequence number N.
In Tables 3 and 4, runtime of SVD-USSD and CA-USSD is given after 200 iterations. Runtime of the CD-USSD algorithm is also given with 16PSK modulation. CA-USSD algorithm has the shortest running time, CD-USSD has the longest running time and SVD-USSD algorithm has the middle running time. We also note that SVD-USSD algorithm can give better optimization results than CA-USSD algorithm because of faster convergence speed. Although the proposed two algorithms can achieve almost the same optimization performance after 300 iterations, SVD-USSD algorithm can obtain better performance under fewer iteration condition. Figure 3 shows the range-Doppler processing results by using Golomb sequence, random sequences and sequence sets optimized by SVD-USSD and CA-USSD algorithms. It is assumed that there exists five scattering points, which are located at (0, 0), (20, 7.5), (20, −7.5), (−20, 7.5),and(−20, −7.5) on range-Doppler plane. As we can see, two scattering points are emerged by the sidelobes of the range-Doppler processing result of Golomb sequence. The sidelobes of the range-Doppler processing result of random sequences is lower than that of Golomb sequence. With sequence sets optimized  by SVD-USSD and CA-USSD algorithms, sidelobes is obviously suppressed and the best detection performance can be obtained.

Discussions
It should be noted that CA-USSD algorithm avoids EVD and SVD operations in iteration procedure, and the computational complexity of this algorithm is obviously lower than that of SVD-USSD algorithm. However, it does not mean that CA-USSD algorithm outperforms SVD-USSD algorithm in every aspect. In the above section, we compare the two algorithms thoroughly. SVD-USSD algorithm has higher computational complexity and better optimization result. The computational complexity of CA-USSD algorithm is lower, but the optimization result of CA-USSD algorithm is not as good as that of SVD-USSD algorithm.

Conclusion
In this paper, two algorithms, i.e., SVD-USSD and CA-USSD, for constructing unimodular sequence sets under constant module constraint with desired minimized sidelobes on a predefined area of range-Doppler plane are presented. SVD-USSD and CA-USSD algorithms show their advantages in different aspect. SVD-USSD algorithm has higher computational complexity and better optimization result. The computational complexity of CA-USSD algorithm is lower, but the optimization result of CA-USSD is not as good as that of SVD-USSD. We further note that the computational efficiency of SVD-USSD limited by SVD operation. This algorithm is better for the sequences of the length no longer than 10 4 . Although the convergence speed of CA-USSD is slower than that of SVD-USSD, computational efficiency of CA-USSD is better and more suitable for the sequences with longer length.