Research Article  Open Access
Focusing of Spotlight TandemConfiguration Bistatic Data with Frequency Scaling Algorithm
Abstract
A frequency scaling (FS) imaging algorithm is proposed for spotlight bistatic SAR data processing. Range cell migration correction (RCMC) is realized through phase multiplication. The proposed algorithm is insensitive to the length of the baseline due to the high precision of the point target (PT) spectrum that we are based on. It is capable of handling bistatic SAR data with a large baseline to range ratio. The algorithms suitable for small and high squint angles are both discussed according to whether the range dependence of the second range compression (SRC) can be neglected or not. Simulated experiments validate the effectiveness of the proposed algorithm.
1. Introduction
Bistatic synthetic aperture radar (SAR) imaging has been widely discussed in recent years [1–4]. Besides the synchronization problems [5, 6], fast imaging algorithms are still in search. However, we cannot acquire accurate bistatic point target (PT) spectrum through the principle of stationary phase (POSP) because of the doublesquareroot term existing in the range history of bistatic SAR. Although the time domain method is the best in theory, it accomplishes optimum reconstruction with heavier computational burden than the frequency domain imaging algorithms [7–9]. The numerical methods [10–12] can handle the bistatic SAR focusing well; however, they also suffer the computational burden [13–15]. The accurate analytical PT spectrum is the precondition of designing most fast bistatic imaging algorithms in the frequency domain, and some effective spectra in the frequency domain have been presented recently [13–17], and many imaging algorithms have been proposed based on them. Tandem bistatic SAR has a simple formation and is relatively easy to accomplish in engineering, which has a bright future in multiple applications [18, 19]. Wu et al. [20, 21] deduced an exact analytical PT spectrum of bistatic SAR in stripmap mode based on the concept of instantaneous Doppler wavenumber (IDW) [22] in tandem configuration.
We focus on the wellknown tandem configuration here. A frequency scaling (FS) imaging algorithm is proposed in this paper for bistatic SAR imaging in spotlight mode based on the spectrum presented in [20]. The process is implemented in the frequency domain, and no interpolation is needed. The azimuth spectral folding effect is solved by the combination of the subaperture approach and the deramping process, when the azimuth bandwidth is larger than the pulse repetition frequency (PRF). Due to the high precision of the spectrum, the proposed algorithm is capable of handling bistatic SAR data with a large baseline.
The paper is organized as follows. In Section 2, the geometry and the signal model of spotlight bistatic SAR are introduced. The proposed FS algorithm is discussed in detail in Section 3 for both small and high squint angles. And the performance of the proposed algorithm is given in Section 4 with simulations. Conclusions are drawn in Section 5.
2. Geometry Relationship and Signal Model of Spotlight Bistatic SAR
Figure 1 shows the geometry of spotlight bistatic SAR, the antenna is constantly steering to the imaging spot during the whole synthetic aperture time to realize a finer azimuth resolution. and represent the transmitter and the receiver, respectively. They travel along the same track with the same velocity , and and are the slant distances from the target to the transmitter and the receiver, respectively. and are the corresponding squint angles. is the closest distance from the target to the flight track, is the half bistatic angle, and is half the length of the baseline.
In spotlight mode, suppose that radar transmits the linear frequencymodulated (LFM) pulses, and the echo signal after dechirp on receiver can be written aswhere stands for the bistatic backscattering coefficient,represents the rectangular window, is the shape of the antenna illumination, is the displacement of the platform in azimuth, and and are the fast and dwell time, respectively. is the pulse duration, is the chirp rate, is the reference range for dechirp, is the closest distance from the scene center to the flight track, is the speed of light, is the carrier frequency, and .
Transforming the signal into the range wavenumber domain, in case of large time bandwidth product [23], the signal can be expressed aswhere is the range wavenumber, is the center of , is the variation of the range wavenumber, , and “” means the convolution operation.
3. Bistatic FS Algorithm
3.1. FS Imaging Algorithm for Tandem Bistatic SAR with Small Squint Angles
Transforming the signal into the twodimensional (2D) wavenumber domain based on the IDW concept [22], we havewhere represents the Doppler wavenumber and indicates the azimuth position where the imaged target is located.
The proposed FS imaging algorithm is based on this equation. Wu et al. deduced an exact analytical expression of the half bistatic angle for tandem bistatic SAR [21]. After substituting it into (4) and making Taylor series expansion up to the cubic term with respect to , the signal comes to
The detailed expressions of , , , , and are given in Appendix. The expressions may seem rather complicated due to the complex formation of the spectrum. However, they will turn to be the familiar azimuth modulation term , the range cell migration (RCM) factor , and the second range compression term in monostatic case when , .
In the following, we fit a straight line to the RCM factorwhere , , and . To see the physical implications of the parameters clearly, again, we can get and , which is just the familiar scale factor in monostatic case when , . The profile error caused by (6) is much smaller than the range resolution and can be neglected.
Substituting (6) into (5), the signal arrives at
To equalize the RCM of all the ranges to the one of the scene center, resembling the monostatic case [23], the frequency scaling function is given aswhere and is the center of . and are the squint angles from the scene center to the transmitter and the receiver at zero time. The reason for introducing into the scale factor is to reduce the amount of shift in range caused by the scaling operation when the squint angles are high, ensuring that the signal after scaling is still in the bandwidth that we are dealing with. When the system works with small squint angles, can be chosen accordingly. Details about this issue are discussed in [24] in monostatic case.
The multiplication between (7) and (8) results in
From the comparison between (7) and (9), it can be seen that the spacevariant term in the RCM term has been eliminated. The linear phase term that indicates the range position of the target and the bulk range shift term have emerged. After an inverse fast Fourier transform (IFFT) in range, is multiplied by the signal to correct the range time shift caused by the term called the residual video phase (RVP), which results from the dechirp operation on receiverwhere indicates the range summation of the transmitter and the receiver of bistatic SAR and and are Fourier transform pairs.
Again, transforming the signal into the 2D wavenumber domain, the signal moves to
Then, the inverse frequency scaling function is introduced to correct the secondorder range phase error caused by the frequency scaling operation
The following operation is the bulk range shift function for the RCMC:
The range dependence of SRC can be neglected when the bistatic SAR system works with small squint angles, and it can be compensated with the parameters of the scene center. The SRC function is given aswhere and . The residual phase error caused by the approximation will be smaller than , which will not give a crucial influence to the final imaging quality. However, when the squint angles are high enough, neglecting the influence of the approximation will no longer be available. In this situation, we can extend the nonlinear chirp scaling algorithm (NCSA) [24] of the monostatic case into the bistatic one.
The operations in range are finished after transforming the signal into the domain.
If the azimuth bandwidth is smaller than PRF, the azimuth compression filter can be directly given as
Note that, due to the change of the image in range as is shown in (11), we have to accommodate the range parameter into , so are the rangedependent variables and . At the end of the whole imaging process, an IFFT is performed in azimuth to transform the focused data into the complex image domain.
However, the PRF is usually smaller than the whole azimuth bandwidth in practice in spotlight mode, especially in spaceborne case; thus the imaging algorithms cannot be directly applied into the whole aperture because of the spectral folding effect. One way to solve the problem is the subaperture method and combining it with the deramping process. The data is divided into several subapertures according to the azimuth time, ensuring that the bandwidth within each subaperture is smaller than PRF. The FS algorithm is then implemented in each subaperture.
In the following, we come to get the highorder phase compensation function to transform the bistatic azimuth phase history into a purely quadratic one as the monostatic case [23]. Making Taylor series expansion of the phase term in (15) with respect to
contains the highorder terms of . After the compensation of , the purely quadratic term emerges. On analog of the monostatic case [23], the deramping function is given after an IFFT with respect to where and are the start and the end time of the th subaperture and is the central time of the th subaperture. In the end, the interested image is obtained with a full aperture FFT in azimuth after the subapertures are recombined.
To summarize and make the proposed FS algorithm more clearly, the workflows of the proposed FS algorithm are given as follows: (1) divide the bistatic SAR data into subapertures; (2) perform azimuth FFT and implement the frequency scaling by (8); (3) perform range IFFT and implement residual video phase correction by (10); (4) perform range FFT and implement the inverse frequency scaling by (12), the bulk range shift for RCMC by (13), and the second range compression by (14); (5) perform range IFFT and implement highorder terms compensation with respect to the Doppler wavenumber; (6) perform azimuth IFFT and implement deramping by (17); and (7) recombine the subapertures and perform full aperture azimuth FFT to obtain the final bistatic SAR image. The block diagram of the proposed FS imaging algorithm is shown in Figure 2.
3.2. FS Imaging Algorithm for Tandem Bistatic SAR with High Squint Angles
In some cases, the bistatic SAR system works with high squint angles, such as observing the interested object which is squintdepended or inspecting the front situations. The range dependence of SRC has to be taken into account in such cases. Taking (14) to compensate the SRC would be inappropriate, the residual phase error will be larger than , and ideal focusing quality cannot be guaranteed. Facing the problem, one instinctive way is to extend the NCSA in monostatic case into the bistatic one. Although the problem can also be solved by processing small range blocks, it requires a large overlap of blocks, which will make the process inefficient [24].
To eliminate the range dependence of SRC, firstly, a thirdorder filter term is introduced
The coefficient will be given subsequently. Similar to the monostatic case, we believe in the approximation that the thirdorder term of SRC can be compensated by the parameters of the scene center. In the same time a new thirdorder phase term is introduced, and thus the secondorder term in the SRC is preserved for the scaling operation. The SRC term shown by (14) is replaced by (18) accordingly. And the signal arrives atwhere . Using the POSP and neglecting the influence of the thirdorder phase to the stationary point, the signal is subsequently transformed into the domain
Multiply the signal by the nonlinear scaling function which is constructed as
A range FFT is performed to the signal after the multiplication between (20) and (21). We take a close look at the phase term
Due to the range dependence of , we fit it by a straight linewhere and , provided that the following approximations are met:
Substituting (23) and (24) into (22), the phase term comes towhere
Making and to eliminate the terms that are spacevariant, we get , , and . The scaling factor , and is chosen to be here, and the value can be adjusted for satisfying imaging results [24]. Substituting the results into (25), the signal arrives at
As can be seen, the range dependence of SRC has been eliminated. The function for range compression is given as
And, then, the phase correction factor is given as
To construct accurate azimuth compression filter, we have to accommodate into according to (27), so are the rangedependent parameters and :
An azimuth IFFT is performed to transform the focused data into the complex image domain in the end.
Similarly, we summarize the proposed nonlinear FS algorithm that is suitable for high squint angle bistatic SAR. The corresponding workflows are given as follows: (1) perform azimuth FFT and implement frequency scaling by (8); (2) perform range IFFT and implement residual video phase correction by (10); (3) perform range FFT and implement the inverse frequency scaling by (12), the bulk range shift for RCMC by (13), and the thirdorder filter by (18); (4) perform range IFFT and implement nonlinear frequency scaling by (21); (5) perform range FFT and implement the secondorder range compression by (28); (6) perform range IFFT and implement phase correction by (29) and azimuth compression by (30); and (7) perform azimuth IFFT to obtain the final bistatic SAR image. The block diagram is shown in Figure 3.
From the above discussion, we can tell that the range history of bistatic SAR is the sum of two square roots, which results in different bistatic spectrum from the monostatic case. The different bistatic spectrum leads to different range cell migration. Thus the conventional FS algorithm cannot be directly utilized into the bistatic case. In this paper, we design a new frequency scaling algorithm which is suitable for bistatic case. The proposed bistatic frequency scaling function is established based on an exact analytical bistatic spectrum to realize rangedependent RCM correction. And the corresponding functions for bistatic FS algorithm are also established, such as the inverse frequency scaling function, the residual video phase correction function, the range migration correction function, and the second range compression function. We construct all these new functions based on the bistatic spectrum. Besides, the rangedependent parameters need to be updated for bistatic case. In addition, a nonlinear frequency scaling algorithm is proposed to deal with the high squint angle bistatic SAR based on the analytical bistatic spectrum. We also establish the functions of the nonlinear FS algorithm for bistatic SAR accordingly. Although these functions for the bistatic FS algorithm are much more complex than the monostatic one, fortunately, we can evaluate the correctness of the algorithm by degenerating the conditions into the monostatic case. In practice, we can choose the suitable algorithm according to the squint angles.
4. Simulation Results
4.1. Simulations with Small Squint Angles
Experiments are carried out to validate the effectiveness of the proposed algorithm. Table 1 illustrates the main parameters of the experiment in three cases in tandem configuration. We put nine point targets in the imaging scene with the center one located at . The distance between any two neighboring targets is 1500 m in both range and azimuth. Taking Case I for illustration, the profile error caused by the linear approximation of the RCM (6) is shown in Figure 4, which is far less than the range resolution. The influence can be neglected. The residual phase error of SRC for the edge target compensated by the local (14) is shown in Figure 5. As can be seen, the error is tolerable. The final imaging results of the center and the edge target by using the proposed FS algorithm are shown in Figure 6. We can see that they are quite well focused; the detailed values of the impulseresponse width (IRW), the peak sidelobe ratio (PSLR), and the integrated sidelobe ratio (ISLR) are shown in Table 2 with theoretic values (IRW, PSLR, and ISLR) 0.75 m in range, 1 m in azimuth, −13.3 dB, and −10 dB, respectively.


(a)
(b)
To show the advantage of the proposed algorithm, we compare it with the MSR and DMO spectra based algorithms. Satisfying focusing quality can be obtained by using all the algorithms for Case I except the DMO based one. For Case II, the baseline increases to 420 km, and the imaging results obtained by using different algorithms are shown in Figures 7–9, respectively. As can be seen, the center target is still well focused by using both the proposed and the MSR based algorithms, while the imaging result obtained by using the DMO spectrum based algorithm is much worse than the others. As for the edge target, the focusing qualities differ much under different algorithms. As can be seen, the one obtained by using the proposed algorithm is the best.
(a)
(b)
(a)
(b)
(a)
(b)
In the following, we come to Case III with the baseline 600 km; corresponding imaging results by using the algorithms are shown in Figures 10–12, respectively. Again, ideal focusing results can be obtained for the center target with both the proposed and the MSR spectrum based algorithms but not the DMO spectrum based one. The focusing quality degrades dramatically with the increasing length of the baseline. As for the edge target, the focusing qualities are intolerable by using all the algorithms except the proposed one. As can be seen, the focusing quality for the edge target decreases by using the DMO and the MSR based algorithms in this case. The reason for the phenomenon lies in the fact that, with the increasing length of the baseline, the precision of the DMO and MSR spectrum degrades, and the corresponding imaging results based on them degrade accordingly. The detailed values of the IRW, PSLR, and ISLR of the center target under different algorithms for Case II and Case III are shown in Table 3, and the corresponding descriptions of the edge target are given in Table 4.


(a)
(b)
(a)
(b)
(a)
(b)
We can tell that the proposed algorithm is capable of dealing with spotlight bistatic SAR data with a large baseline. Note that, for the MSR spectrum based algorithm, we only expand the Taylor series expansion up to the third term for comparison here. For better focusing quality of the MSR based algorithm, higher orders can be expanded.
4.2. Simulations with High Squint Angles
In the following, we take a look at the case with high squint angels. The main parameters are shown in Table 5. Nine point targets are put in the imaging scene with the center one located at . The distance between any two neighboring point targets is 400 m in both range and azimuth. The distances from the center target to the transmitter and the receiver are 17.47 km and 10.99 km, respectively. Corresponding squint angles are 62.75° and 43.28°, respectively. The PRF is larger than the azimuth bandwidth in this case, so the proposed FS imaging algorithm can be applied onto the whole aperture.

If the parameters of the scene center are still chosen to compensate the rangedependent SRC, the residual phase error will be larger than , as shown in Figure 13. Imaging results for the center and the edge targets are shown in Figure 14 if the algorithm suitable for small squint angles is adopted. It can be seen that the center target can be well focused, but the focusing quality of the edge target degrades dramatically. The results obtained by using the proposed algorithm suitable for high squint angles are shown in Figure 15. We can tell that the focusing quality for the edge target has improved explicitly from the comparison between Figures 14(b) and 15(b). It implies that the range dependence of SRC must be taken into consideration in this case. To show the advantage of the proposed algorithm, the MSR based algorithm is carried out for comparison. Corresponding imaging results are shown in Figure 16. Detailed descriptions of the focusing quality of the edge target under different algorithms are shown in Table 6.

(a)
(b)
(a)
(b)
(a)
(b)
In practice, when the data is obtained, we can determine the SRC to test the amplitude of the residual phase error and choose the proper imaging algorithm subsequently.
5. Conclusion
An FS algorithm suitable for tandem bistatic SAR in spotlight mode is proposed, which is insensitive to the baseline to range ratio. The subaperture approach and the deramping process are combined to handle the problem of the azimuth spectral folding effect like the monostatic case. Ideal focusing results are obtained in the frequency domain without interpolation. A nonlinear FS algorithm is also discussed to deal with the tandem bistatic data with high squint angles, in which situation the range dependence of SRC must be taken into consideration. Satisfactory results are obtained with simulated experiments.
However, we have to note that the imaging algorithm presented is based on ideal conditions without considering motion errors, which cannot be avoided in practice [25]. Motion error analysis and compensation algorithms deserve further studying.
Appendix
Here, we give the detailed formations of some parameters. Wu et al. deduced an exact expression for tandem bistatic SAR [20], which is given aswhere
Substituting it into (4) and making Taylor series expansion up to the cubic term with respect to , the signal will come to (5). The detailed expressions of , , , and are
In (A.3)
As can be seen, the expressions seem rather complicated due to the complex formation of the half bistatic angle; however, they will turn into the familiar ones under the monostatic condition when and .
Competing Interests
The authors declare that they have no competing interests.
Acknowledgments
This work was supported by the Natural Science Foundation of China (Grant 61222108) and the Fundamental Research Funds for the Central Universities (Grant GK201603089).
References
 M. Antoniou, M. Cherniakov, and H. Ma, “Spacesurface bistatic synthetic aperture radar with navigation satellite transmissions: a review,” Science China Information Sciences, vol. 58, no. 6, pp. 1–20, 2015. View at: Publisher Site  Google Scholar
 S. Chen, M. Xing, S. Zhou, L. Zhang, and Z. Bao, “Focusing of tandem bistatic SAR data using the chirpscaling algorithm Emerging radar techniques,” EURASIP Journal on Advances in Signal Processing, vol. 2013, no. 1, article 38, 2013. View at: Publisher Site  Google Scholar
 S. Cherouat, F. Soltani, F. Schmitt, and F. Daout, “Using fractal dimension to target detection in bistatic SAR data,” Signal, Image and Video Processing, vol. 9, no. 2, pp. 365–371, 2013. View at: Publisher Site  Google Scholar
 Q. Zhang, W. Chang, and X. Li, “An integrative synchronization and imaging approach for bistatic spaceborne/stratospheric SAR with a fixed receiver,” Eurasip Journal on Advances in Signal Processing, vol. 2013, no. 1, article 165, 2013. View at: Publisher Site  Google Scholar
 M. Weiß, “Time and frequency synchronization aspects for bistatic SAR systems,” in Proceedings of the European Conference on Synthetic Aperture Radar (EUSAR '04), pp. 395–399, Ulm, Germany, 2004. View at: Google Scholar
 G. Krieger and A. Moreira, “Spaceborne bi and multistatic SAR: potentials and challenges,” IET Radar, Sonar and Navigation, vol. 153, no. 3, pp. 184–198, 2006. View at: Publisher Site  Google Scholar
 B. Barber, “Theory of digital imaging from orbit synthetic aperture radar,” International Journal of Remote Sensing, vol. 6, no. 6, pp. 1009–1057, 1985. View at: Publisher Site  Google Scholar
 M. RodriguezCassola, P. Prats, G. Krieger, and A. Moreira, “Efficient timedomain image formation with precise topography accommodation for general bistatic SAR configurations,” IEEE Transactions on Aerospace and Electronic Systems, vol. 47, no. 4, pp. 2949–2966, 2011. View at: Publisher Site  Google Scholar
 Y. F. Shao, R. Wang, Y. K. Deng et al., “Fast backprojection algorithm for bistatic SAR imaging,” IEEE Geoscience and Remote Sensing Letters, vol. 10, no. 5, pp. 1080–1084, 2013. View at: Publisher Site  Google Scholar
 V. Giroux, H. Cantalloube, and F. Daout, “An OmegaK algorithm for SAR bistatic systems,” in Proceedings of the IEEE International Geoscience and Remote Sensing Symposium (IGARSS '05), pp. 1060–1063, Seoul, South Korea, July 2005. View at: Publisher Site  Google Scholar
 I. Walterscheid, J. H. G. Ender, A. R. Brenner, and O. Loffeld, “Bistatic SAR processing and experiments,” IEEE Transactions on Geoscience and Remote Sensing, vol. 44, no. 10, pp. 2710–2717, 2006. View at: Publisher Site  Google Scholar
 R. Bamler, F. Meyer, and W. Liebhart, “Processing of bistatic SAR data from quasistationary configurations,” IEEE Transactions on Geoscience and Remote Sensing, vol. 45, no. 11, pp. 3350–3358, 2007. View at: Publisher Site  Google Scholar
 Y. L. Neo, F. Wong, and I. G. Cumming, “A twodimensional spectrum for bistatic SAR processing using series reversion,” IEEE Geoscience and Remote Sensing Letters, vol. 4, no. 1, pp. 93–96, 2007. View at: Publisher Site  Google Scholar
 Y. L. Neo, F. Wong, and I. G. Cumming, “A comparison of point target spectra derived for bistatic SAR processing,” IEEE Transaction on Geoscience and Remote Sensing, vol. 46, no. 9, pp. 93–96, 2008. View at: Google Scholar
 R. Wang, Y. K. Deng, O. Loffeld et al., “Processing the azimuthvariant bistatic SAR data by using monostatic imaging algorithms based on twodimensional principle of stationary phase,” IEEE Transactions on Geoscience and Remote Sensing, vol. 49, no. 10, pp. 3504–3520, 2011. View at: Publisher Site  Google Scholar
 D. D' Aria, A. M. Guarnieri, and F. Rocca, “Focusing bistatic synthetic aperture radar using dip move out,” IEEE Transactions on Geoscience and Remote Sensing, vol. 42, no. 7, pp. 1362–1376, 2004. View at: Publisher Site  Google Scholar
 O. Loffeld, H. Nies, V. Peters, and S. Knedlik, “Models and useful relations for bistatic SAR processing,” IEEE Transactions on Geoscience and Remote Sensing, vol. 42, no. 10, pp. 2031–2038, 2004. View at: Publisher Site  Google Scholar
 M. RodriguezCassola, P. Prats, D. Schulze et al., “First bistatic spaceborne SAR experiments with TanDEMX,” IEEE Geoscience and Remote Sensing Letters, vol. 9, no. 1, pp. 33–37, 2012. View at: Publisher Site  Google Scholar
 A. Jäggi, O. Montenbruck, Y. Moon et al., “Interagency comparison of TanDEMX baseline solutions,” Advances in Space Research, vol. 50, no. 2, pp. 260–271, 2012. View at: Publisher Site  Google Scholar
 Q. Wu, M. Xing, H. Shi, X. Hu, and Z. Bao, “Exact analytical twodimensional spectrum for bistatic synthetic aperture radar in tandem configuration,” IET Radar, Sonar and Navigation, vol. 5, no. 3, pp. 349–360, 2011. View at: Publisher Site  Google Scholar
 Q. Wu, Y. Liang, M. Xing, C. Qiu, Z. Bao, and T.S. Yeo, “Focusing of tandem bistaticconfiguration data with range migration algorithm,” IEEE Geoscience and Remote Sensing Letters, vol. 8, no. 1, pp. 88–92, 2011. View at: Publisher Site  Google Scholar
 Z. Zhenhua, X. Mengdao, D. Jinshan, and B. Zheng, “Focusing parallel bistatic SAR data using the analytic transfer function in the wavenumber domain,” IEEE Transactions on Geoscience and Remote Sensing, vol. 45, no. 11, pp. 3633–3645, 2007. View at: Publisher Site  Google Scholar
 J. Mittermayer, A. Moreira, and O. Loffeld, “Spotlight SAR data processing using the frequency scaling algorithm,” IEEE Transactions on Geoscience and Remote Sensing, vol. 37, no. 5, pp. 2198–2214, 1999. View at: Publisher Site  Google Scholar
 W. Davidson, I. G. Cumming, and M. R. Ito, “A chirp scaling approach for processing squint model SAR data,” IEEE Transaction on Aerospace and Electronic Systems, vol. 32, no. 1, pp. 121–133, 1996. View at: Publisher Site  Google Scholar
 B. D. Rigling and R. L. Moses, “Motion measurement errors and autofocus in bistatic SAR,” IEEE Transactions on Image Processing, vol. 15, no. 4, pp. 1008–1016, 2006. View at: Publisher Site  Google Scholar
Copyright
Copyright © 2016 Shichao Chen et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.