CN111190181A - Real-time imaging processing method for unmanned aerial vehicle SAR (synthetic aperture radar) of bumpy platform - Google Patents

Abstract

The invention belongs to the field of unmanned aerial vehicle remote sensing mapping, and discloses a bumpy platform unmanned airborne SAR real-time imaging processing method, wherein the unmanned airborne SAR real-time imaging processing based on a heterogeneous multi-core chip utilizes the advantage of high computation efficiency of an FPGA to execute digital down conversion, distance pulse pressure and azimuth pre-filtering processing which mainly comprise multiplication and accumulation operation, and the FPGA performs data distribution and task allocation according to DSP resources; by utilizing the advantage of high programming flexibility of the DSP, the method executes the imaging processing of the bumping platform SAR based on multi-level motion compensation, the DSP end performs motion compensation processing based on inertial navigation information, the DSP end performs motion compensation processing based on data, the DSP end performs image domain self-focusing processing and other multiple motion error compensation technologies, realizes high-resolution imaging of the bumping platform represented by an unmanned aerial vehicle, greatly improves the utilization rate of hardware resources, and has good focusing effect and strong robustness on the real-time imaging processing of the bumping platform SAR.

Description

Real-time imaging processing method for unmanned aerial vehicle SAR (synthetic aperture radar) of bumpy platform

Technical Field

The invention relates to the field of unmanned aerial vehicle remote sensing mapping, in particular to a bumpy platform unmanned aerial vehicle-mounted SAR real-time imaging processing method.

Background

The Synthetic Aperture Radar (SAR) can acquire a two-dimensional high-resolution image of a ground scene, and has a wide application prospect in the fields of remote sensing detection and topographic mapping. Meanwhile, the unmanned aerial vehicle platform has the unique advantages of low cost, small size and the like, and the SAR carried on the unmanned aerial vehicle platform can be organically combined with the advantages of the SAR and the unmanned aerial vehicle platform, so that the unmanned aerial vehicle SAR becomes a research focus in recent years.

The unmanned aerial vehicle platform is light in weight, and its flight path is easily influenced by weather factors such as crosswind, air current to introduce great motion error, belong to typical jolting platform. Most of the existing solutions are combined with an Inertial Navigation System (INS) to perform motion compensation processing, most of motion errors along the direction of a beam line of sight (LOS) are eliminated, and imaging requirements can be met when residual motion errors are smaller than half of a resolution unit. However, for high-resolution SAR imaging, the current situation that the INS accuracy is insufficient is limited, the residual motion error still has a certain influence on focusing, and the image has phenomena such as defocusing and tailing to a certain extent. At the same time, the bumpiness of the stage can result in non-uniform position sampling, which can lead to phase errors in the direction perpendicular to the LOS, again with a non-negligible effect on focusing.

With the rapid development of semiconductor technology, the operational efficiency of hardware platforms such as Field Programmable Gate Arrays (FPGAs), Digital Signal Processors (DSPs) and the like is further improved, and the method is more applied to the field of SAR real-time imaging. In engineering practice, cost and power consumption are two major factors that trade off hardware platforms. Most of the existing SAR real-time processors adopt a chip stacking mode, imaging real-time performance is met, meanwhile, a large amount of performance redundancy is caused, and control cost and power consumption reduction are very unfavorable. Therefore, the engineering scheme focuses more on how to fully utilize hardware resources, and exerts the advantages of each platform while performing real-time processing, so as to achieve the balance between algorithm performance and cost and power consumption.

Disclosure of Invention

Aiming at the problems in the prior art, the invention aims to provide a bumpy platform unmanned aerial vehicle-mounted SAR real-time imaging processing method, which is based on the unmanned aerial vehicle-mounted SAR real-time imaging processing of a heterogeneous multi-core chip, realizes high-resolution imaging of a bumpy platform represented by an unmanned aerial vehicle, greatly improves the utilization rate of hardware resources, and has good focusing effect on the bumpy platform SAR real-time imaging processing and strong robustness.

The technical idea of the invention is as follows:

(1) improve traditional formation of image flow, eliminate the influence of jolting platform motion error: (1a) correcting the deviation of the platform from an ideal track according to INS parameters, and mainly eliminating phase errors caused by motion errors along the LOS direction of the beam; (1b) aiming at the problem that the Doppler modulation frequency is greatly influenced by the motion parameters of the bumping platform, estimating and compensating the residual Doppler modulation frequency by using an image bias (MD) algorithm, thereby greatly improving the focusing effect of the image; (1c) and performing phase gradient self-focusing (PGA) processing to further compensate residual phase errors caused by acceleration in an image domain.

(2) The real-time processing software architecture based on the heterogeneous multi-core chip is designed by combining the characteristics of an FPGA chip XC7A100T and a DSP chip TMS320C 6678: (2a) the advantage of high operation efficiency of the FPGA is fully exerted, and Digital Down Conversion (DDC), distance pulse voltage (PC) and pre-filtering (PF) processing which mainly comprise multiply-accumulate operation are executed at the FPGA end; (2b) partitioning the data matrix along the distance direction, and distributing the data to a plurality of DSP chips through the FPGA; (2c) and (3) the advantage of high flexibility of DSP programming is fully played, and the improved pitch platform SAR imaging algorithm in the step (1) is executed at the DSP end.

In order to achieve the purpose, the invention is realized by adopting the following technical scheme.

A bumpy platform unmanned aerial vehicle carries SAR real-time imaging processing method, including the following steps:

step 1, receiving and transmitting an LFM pulse signal by using an unmanned airborne SAR system to obtain an intermediate frequency echo real signal s_r(t_r) (ii) a Using analog-to-digital converter pairs s_r(t_r) After sampling and quantization, a discrete intermediate frequency echo real signal s is obtained_r(k) FPGA end pair s_r(k) Performing digital down-conversion processing to generate zero intermediate frequency orthogonal dual-path signal s_I(k) And s_Q(k) And integrated as time domain echoes ss (t)_r,t_a；R_S)＝s_I(k)+js_Q(k) (ii) a Wherein, t_rAs a fast time variable of distance, t_aFor azimuthal slow time variation, R_SThe center slant distance of the scene is shown, and j is an imaginary number unit;

step 2, the FPGA end judges the bandwidth B of the transmitting signal_rAnd ADC sampling frequency f_sThe relation between the time domain echo ss (t) and the time domain echo t is selected to be a direct sampling mode or a de-skew mode_r,t_a；R_S) Distance direction pulse compression processing is carried out to obtain a distance direction pulse compression signal ss_PC,R(t_r,t_a；R_S)；

Step 3, FPGA end-to-end distance direction pulse compression signal ss_PC,R(t_r,t_a；R_S) Performing azimuth pre-filtering to obtain an azimuth filtered signal ss_PF,R(t_r,t_a；R_S)；

Step 4, the FPGA end enables the azimuth-filtered signal ss_PF,R(t_r,t_a；R_S) Partitioning along the distance direction to obtain a plurality of data blocks, wherein the number of overlapping points of adjacent data blocks is G, and each data block is sent to a DSP end; wherein each data block is denoted as ss_n,R(t_r,t_a；R_S) N is the number of data blocks;

step 5, the DSP end obtains the phase error

Using said phase error

For ss_n,R(t_r,t_a；R_S) Performing INS-based motion compensation to obtain an INS-based motion compensation signal ss_INSn,R(t_r,t_a；R_S)；

Step 6, the DSP end pair is used for motion compensation signal ss based on INS_INSn,R(t_r,t_a；R_S) Matrix transposition is carried out to obtain a transposed signal ss_INSn,A(t_r,t_a；R_S)；

Step 7, DSP terminalSignals ss after inversion_INSn,A(t_r,t_a；R_S) Sequentially carrying out distance Fourier transform and azimuth Fourier transform to obtain two-dimensional frequency domain echo signal SS_INSn,A(f_r,f_a；R_S) (ii) a Constructing a range migration correction function H_RCMC(f_r,f_a；R_S) Two-dimensional frequency domain echo signal SS_INSn,A(f_r,f_a；R_S) And H_RCMC(f_r,f_a；R_S) Multiplying and performing two-dimensional inverse Fourier transform to obtain range migration correction signal ss_RCMCn,A(t_r,t_a；R_S)；

Step 8, the DSP end constructs a phase compensation function H_MDCorrection of range migration signal ss_RCMCn,A(t_r,t_a；R_S) And H_MDMultiplying to obtain a data-based motion compensated signal ss_MDn,A(t_r,t_a；R_S)；

Step 9, the DSP port pair motion compensation signal ss based on data_MDn,A(t_r,t_a；R_S) Performing direction Fourier transform to obtain a signal sS in a range-Doppler domain_MDn,A(t_r,f_a；R_S) (ii) a Constructing an orientation matching filter function, and adopting the orientation matching filter function to sS_MDn,A(t_r,f_a；R_S) Performing azimuth compression to obtain an azimuth compression signal ss_ACn,A(t_r,t_a)；

Step 10, the DSP end obtains an estimation value of the azimuth phase error by utilizing a phase gradient self-focusing algorithm, and an azimuth compression signal ss is obtained according to the estimation value of the azimuth phase error_ACn,A(t_r,t_a) Corrected to obtain accurate range block focusing signal ss_PGAn,A(t_r,t_a)；

Step 11, obtaining accurate range block focusing signal ss from multiple DSP terminals_PGAn,A(t_r,t_a) And carrying out image splicing according to the number G of the overlapped points of the adjacent data blocks, and assembling into a complete SAR image.

Compared with the prior art, the invention has the beneficial effects that:

(1) the traditional unmanned aerial vehicle SAR imaging algorithm is low in motion error containment degree, and when the platform jolts greatly, the imaging result is prone to defocusing, tailing and the like. The invention fully considers the flight characteristics of a jolting platform and the influence of the motion error on the imaging performance, and makes a plurality of effective improvements aiming at the traditional strip mode imaging algorithm, such as a plurality of motion error compensation technologies of performing motion compensation processing based on inertial navigation information at a DSP end, performing motion compensation processing based on data at the DSP end, performing image domain self-focusing processing at the DSP end, and the like. Compared with the existing bumpy platform unmanned aerial vehicle SAR imaging algorithm, the method has the advantages of better focusing effect and stronger robustness.

(2) Most of the existing real-time processing technologies are designed in a chip stacking mode, so that part of computing capacity and storage resources are idle and wasted. The advantages of heterogeneous multi-core chips are fully exerted, the advantage of high computational efficiency of the FPGA is utilized, digital down conversion, distance pulse pressure and azimuth pre-filtering processing which take multiply-accumulate operation as a main body are executed, and the FPGA performs data distribution and task allocation according to DSP resources; and (3) executing bump platform SAR imaging processing based on multilevel motion compensation by utilizing the advantage of high flexibility of DSP programming, such as phase estimation, phase compensation, matrix transposition and the like. Compared with the existing scheme, the method effectively avoids the idleness and redundancy of hardware resources, and has higher resource utilization rate and lower power consumption and cost.

(3) The DSP terminal processes the INS-based motion compensation signal ss between the distance direction processing and the azimuth direction processing_INSn,R(t_r,t_a；R_S) The data is subjected to matrix transposition, and the matching of a data storage mode and a module processing mode is ensured, so that the bus occupancy rate of reading and writing data is greatly reduced, the data reading and writing efficiency is improved, and the real-time performance of system imaging is effectively improved.

Drawings

The invention is described in further detail below with reference to the figures and specific embodiments.

FIG. 1 is a flow chart of the real-time imaging processing method of the unmanned airborne SAR of the bumpy platform of the present invention;

FIG. 2 is an unmanned airborne SAR front side view stripe imaging geometry model;

FIG. 3 is a flow chart of an FPGA digital down conversion process;

FIG. 4 is a flow chart of the FPGA performing range pulse pressure processing; wherein, fig. 4(a) is a flow chart of direct sampling pulse pressure; FIG. 4(b) is a flow chart of declivity pulse pressure;

FIG. 5 is a schematic diagram of data distribution and task allocation by the FPGA to the DSP1 and the DSP 2;

FIG. 6 is a diagram of DSP matrix transposition;

FIG. 7 is a flow chart of a data motion compensation method combining DSP with MD idea;

FIG. 8 is a schematic diagram of dividing azimuth sub-apertures in the MD algorithm;

fig. 9 is a flowchart of the DSP executing PGA;

fig. 10 is an imaging result diagram of the tossing platform unmanned aerial vehicle-mounted SAR real-time imaging processing method in the invention.

Detailed Description

Embodiments of the present invention will be described in detail below with reference to examples, but it will be understood by those skilled in the art that the following examples are only illustrative of the present invention and should not be construed as limiting the scope of the present invention.

The invention provides a real-time imaging processing method of an unmanned airborne SAR (synthetic aperture radar) of a bumpy platform, and a specific processing flow chart is shown in figure 1. The FPGA end executes digital down conversion, distance pulse pressure, azimuth pre-filtering and data distribution, then the DSP1 and the DSP2 execute INS-based motion compensation, matrix transposition, distance migration correction, data-based motion compensation, azimuth compression and PGA, and finally image splicing is completed at the DSP1 end and is sent to a system board through a PCIe interface. The present invention will be described in further detail with reference to the accompanying drawings.

Step 1, receiving and transmitting an LFM pulse signal by using an unmanned airborne SAR system to obtain an intermediate frequency echo real signal s_r(t_r) (ii) a Using analog-to-digital converter (ADC) pairs s_r(t_r) After the sampling and the quantization are carried out,obtaining discrete intermediate frequency echo real signal s_r(k) FPGA end pair s_r(k) Performing digital down-conversion (DDC) processing to generate zero-intermediate-frequency orthogonal dual-path signal s_I(k) And s_Q(k) And integrated as time domain echoes ss (t)_r,t_a；R_S)＝s_I(k)+js_Q(k)。

Specifically, step 1 comprises the following substeps:

substep 1.1, referring to fig. 2, the unmanned airborne SAR operates in a front side view stripe mode, ideally with the platform flying at a horizontal velocity v along the Y-axis direction and an azimuth beam width θ_BWThe flying height is H. Point B is the scene center point, R_SIs the center slant distance of the scene. The point P is any point distributed along the azimuth direction in the scene, and the distance between the point P and the point B is X_nI.e. X_nIs the distance between the center point of the scene and any point in the scene that is distributed along the azimuth direction. According to the geometric relationship, the instantaneous slope distance between the platform and the point target is obtained as follows:

wherein, t_aIs an azimuth slow time variable. Based on the instantaneous slope distance between the platform and the point target, an expression of the intermediate frequency echo real signal can be obtained:

wherein, t_rAs a function of distance and fast time, c denotes the speed of light, f_IRepresenting the intermediate frequency, k, before ADC sampling_rIndicating the tuning frequency of the LFM signal.

Substep 1.2, discrete intermediate frequency echo real signal s_r(k) Comprises the following steps:

wherein k represents a discrete time variable, N_aIndicating the number of azimuth sample points.

Substep 1.3, referring to fig. 3, the FPGA receives the discrete intermediate frequency real echo signal s_r(k) Then storing the data into a Read Only Memory (ROM), and utilizing the floating point intellectual property core (IP core) to respectively correspond to cos (2 pi f)_Ik △ t) and sin (2 π f)_Ik △ t), wherein △ t represents the time interval between adjacent sampling points, then a low-pass filter f (k) is generated by means of an fdatool tool in MATLAB software and is preset in a ROM of the FPGA, the two signals are respectively shifted, multiplied and accumulated by using the f (k), the low-pass filtering is completed, and a zero intermediate frequency orthogonal two-way signal s is generated_I(k) And s_Q(k) In that respect The above process can be expressed as:

wherein the content of the first and second substances,

which represents the calculation of the convolution,

indicating that a summation operation is performed from M-0 to M-1, where M indicates the number of total sample points and an increase in the summation variable M indicates a shift in the vector.

Substep 1.4, quadrature two-way signal s of zero intermediate frequency_I(k) And s_Q(k) Respectively as the real part and the imaginary part of the complex signal to form a time domain echo ss (t)_r,t_a；R_S) Comprises the following steps:

wherein, w_r(t_r) Representing a distance-time-domain window function, w_a(t_a) Representing the azimuth time domain window function, j is an imaginary unit, and λ is a wavelength.

Step 2, judging the bandwidth of the transmitting signal and the ADC sampling frequency f by the FPGA terminal_sThe relation between the time domain echo ss (t) and the time domain echo t (t) is selected to be in a direct sampling mode or a de-skew mode_r,t_a；R_S) Performing distance direction pulse compression to obtain distanceOff-directional pulse compression signal ss_PC,R(t_r,t_a；R_S) (ii) a Where the subscript PC indicates the distance-wise pulse compression and the subscript R indicates that the data is continuous along the distance-wise (Range).

Specifically, when the system bandwidth is low, the sampling frequency of the ADC meets the Nyquist sampling theorem, and the FPGA end can obtain the range-wise pulse compression result ss by adopting a direct pulse pressure sampling mode_PC,R(t_r,t_a；R_S) Wherein the subscript R indicates that the data is continuous along the distance direction (Range); when the system bandwidth is high, the maximum sampling frequency of the ADC cannot meet the Nyquist sampling theorem, and the FPGA end needs to obtain ss in a way of removing the oblique pulse pressure_PC,R(t_r,t_a；R_S) This approach equivalently reduces the system sampling frequency. Step 2 comprises the following substeps:

substep 2.1, the FPGA end adopts a frequency domain matching filtering mode to carry out distance pulse pressure, firstly, an IP core of Fourier transform (FFT) is called to lead time domain echo ss (t)_r,t_a；R_S) And (3) transforming to a distance frequency domain to obtain a distance echo:

Ss(f_r,t_a；R_S)＝FFT_R[ss(t_r,t_a；R_S)]

wherein, FFT_RRepresenting a distance fourier transform.

Substep 2.2, based on the transmission signal bandwidth B_rAnd ADC sampling frequency f_sThe relationship between the two is that the distance direction pulse compression form is selected, and the specific form is as follows:

(a) when B is present_r<f_sAnd the FPGA executes direct sampling pulse pressure processing. According to the frequency modulation rate k of LFM signal_rConstructing a frequency domain matched filtering reference function H_PCAnd presetting the FPGA in a ROM of the FPGA:

wherein f is_rRepresenting the distance versus frequency variation. Referring to FIG. 4(a), the multiplier is invoked to convert Ss (f)_r,t_a；R_S) And H_PCMultiplying; and (3) calling an IP core of inverse Fourier transform (IFFT) to transform the multiplication result to a two-dimensional time domain to obtain a range direction pulse compression signal:

ss_PC,R(t_r,t_a；R_S)＝IFFT_R[Ss(f_r,t_a；R_S)·H_PC]

wherein the IFFT_RRepresenting the inverse fourier transform of the distance.

(b) When B is present_r≥f_sAnd the FPGA executes the deskew pulse pressure processing. Constructing a deskew function S that eliminates Residual Video Phase (RVP) terms_c(f_i) And is preset in ROM of FPGA

Wherein f is_iIs a baseband range frequency variable. Referring to FIG. 4(b), the multiplier is invoked to convert Ss (f)_r,t_a；R_S) And S_c(f_i) Multiplying, namely calling the IP core of IFFT to check the product result and perform inverse Fourier transform of distance direction, then calling a multiplier and the IP core of FFT to check the time domain signal and perform windowed Fourier transform to obtain a distance direction pulse compression signal:

ss_PC,R(t_r,t_a；R_S)＝FFT_R[w_r(t_r)·IFFT_R[Ss(f_r,t_a；R_S)·S_c(f_i)]]

wherein, w_r(t_r) Representing a distance-time-domain window function, FFT_RRepresenting a distance fourier transform.

Step 3, FPGA end-to-end distance direction pulse compression signal ss_PC,R(t_r,t_a；R_S) Performing azimuth pre-filtering to obtain an azimuth filtered signal ss_PF,R(t_r,t_a；R_S)。

Specifically, the step 3 is as follows: first, the Doppler bandwidth B is calculated according to the radar wavelength, the beam width and the flight speed_a：

Wherein D is_aIs the antenna aperture, theta_BWIs the beam width, λ is the wavelength; subsequently, the low-pass filter coefficients are generated by means of the fdatool tool in MATLAB software and preset in the ROM of the FPGA, wherein the pass-band width of the low-pass filter is the Doppler bandwidth B_aThe sampling frequency is the Pulse Repetition Frequency (PRF) of the system. In azimuth direction ss_PC,R(t_r,t_a；R_S) Performing a filtering operation, and storing the output result as a matrix, i.e. the azimuth-filtered signal ss_PF,R(t_r,t_a；R_S)ss_PF,R(t_r,t_a；R_S). The influence of most of moving targets is filtered through the azimuth pre-filtering treatment, so that the imaging effect is improved.

Step 4, the FPGA end enables the azimuth-filtered signal ss_PF,R(t_r,t_a；R_S)ss_PF,R(t_r,t_a；R_S) Partitioning along the distance direction to obtain a plurality of data blocks, wherein the number of overlapping points of adjacent data blocks is G, and each data block is sent to a DSP end; wherein each data block is denoted as ss_n,R(t_r,t_a；R_S) And n is the number of data blocks.

Specifically, referring to fig. 5, the FPGA performs task control and data distribution on a plurality of DSP chips. Will ss_PF,R(t_r,t_a；R_S) Divided into two data blocks along the distance direction, denoted as ss respectively_1,R(t_r,t_a；R_S) And ss_2,R(t_r,t_a；R_S) N is 1 and 2, and the number of overlapped points G of the two data blocks in the distance direction is 512; and distributing the data blocks to the DSP1 and the DSP2 through a serial Rapid I/O (SRIO) interface, and issuing an imaging task instruction.

Step 5, the DSP end obtains the phase error

Using said phase error

For ss_n,R(t_r,t_a；R_S) Performing INS-based motion compensation to obtain an INS-based motion compensation signal ss_INSn,R(t_r,t_a；R_S)。

Specifically, step 5 comprises the following substeps:

substep 5.1, the DSP terminal analyzes the SRIO frame according to the interrupt protocol and analyzes the eastern speed v in the geographic coordinate system transmitted by the INS_eastNorth velocity v_northVelocity v_skyInformation, combined with track angle theta_yawCalculating the along-track velocity v under the imaging coordinate system_alongVertical track velocity v_crossAnd a lifting speed v_eleva：

Step 5.2, integrating inertial navigation speed information by the DSP end, and combining an ideal track to obtain a motion error of the platform along the LOS direction; conversion of motion error into phase error, ss_n,R(t_r,t_a；R_S) Performing INS-based motion compensation to obtain ss_INSn,R(t_r,t_a；R_S) Wherein n is 1, 2.

Specifically, for the three-dimensional velocity component (along-track velocity v) in the imaging coordinate system_alongVertical track velocity v_crossAnd a lifting speed v_eleva) Performing linear integration to obtain three-dimensional motion error e under the coordinate system_along、 e_cross、e_elevaThen the motion error e in the LOS direction is calculated from the pitch angle β of the beam_LOSAnd corresponding phase error

Substep 5.3. mixing

And ss_n,R(t_r,t_a；R_S) Multiplying to obtain an INS-based motion compensation signal:

step 6, the DSP end pair is used for motion compensation signal ss based on INS_INSn,R(t_r,t_a；R_S) Matrix transposition is carried out to obtain a transposed signal ss_INSn,A(t_r,t_a；R_S) (ii) a Where subscript a indicates that the data is continuous in the Azimuth direction (Azimuth).

Specifically, the steps 1 to 5 are all distance direction processing, and require data to be continuous along the distance direction. The subsequent processing is carried out estimation and compensation aiming at the azimuth signal, and in order to improve the data read-write speed and reduce the hardware resource loss, the ss is subjected to_INSn,R(t_r,t_a；R_S) A matrix transposition operation is performed. Referring to FIG. 6, the DSP performs matrix transposition to obtain ss that is continuous along the azimuth direction_INSn,A(t_r,t_a；R_S)。

Step 7, the DSP end contra-rotating signal ss_INSn,A(t_r,t_a；R_S) Sequentially carrying out distance Fourier transform and azimuth Fourier transform to obtain two-dimensional frequency domain echo signal SS_INSn,A(f_r,f_a；R_S) (ii) a Constructing a range migration correction function H_RCMC(f_r,f_a；R_S) Two-dimensional frequency domain echo signal SS_INSn,A(f_r,f_a；R_S) And H_RCMC(f_r,f_a；R_S) Multiplying and performing two-dimensional inverse Fourier transform to obtain range migration correction signal ss_RCMCn,A(t_r,t_a；R_S)。

In particular, the meterCalculating ss_INSn,A(t_r,t_a；R_S) The middle-order range bending amount and the high-order range migration amount are used as a phase generation range migration correction function H_RCMC(f_r,f_a；R_S) (ii) a Then to ss_INSn,A(t_r,t_a；R_S) Compensating to obtain range migration correction signal ss_RCMCn,A(t_r,t_a；R_S) Wherein n is 1, 2.

Step 7 comprises the following substeps:

substep 7.1, DSP port pair ss_INSn,A(t_r,t_a；R_S) Distance direction FFT is carried out to obtain a distance frequency domain expression

Wherein, W_r(f_r) As a function of a distance frequency domain window, f_cCarrier frequency of radar system, f_rRepresenting the distance versus frequency variation.

Substep 7.2, DSP port pair Ss_INSn,A(f_r,t_a；R_S) Performing azimuth FFT to obtain a two-dimensional frequency domain echo signal:

wherein, W_a(f_a) Representing the azimuth frequency domain window function, f_aRepresenting the azimuth frequency variation.

Substep 7.3, relating the second exponential term in the above equation to the variable f_rThe low-order craulin is unfolded:

from which f is found_rAnd f_aAnd constructing a range migration correction function:

substep 7.4, the DSP side will SS_INSn,A(f_r,f_a；R_S) And H_RCMC(f_r,f_a；R_S) Multiplying and performing two-dimensional IFFT to obtain a range migration correction signal:

ss_RCMCn,A(t_r,t_a；R_S)＝IFFT₂[SS_INSn,A(f_r,f_a；R_S)·H_RCMC(f_r,f_a；R_S)]

wherein the IFFT₂[·]Representing a two-dimensional inverse fourier transform.

Step 8, the DSP end constructs a phase compensation function H_MDCorrection of range migration signal ss_RCMCn,A(t_r,t_a；R_S) And H_MDMultiplying to obtain a data-based motion compensated signal ss_MDn,A(t_r,t_a；R_S)。

Specifically, step 8 comprises the following substeps:

substep 8.1, referring to fig. 7, the DSP performs data-based motion compensation using the MD algorithm; referring to FIG. 8, the DSP maps the azimuth signal vector s (t)_a) Divided into two sub-apertures, respectively denoted as s₁(t_a) And s₂(t_a) I.e. s₁(t_a) And s₂(t_a) Front and back subaperture signals, respectively:

wherein T represents the synthetic aperture time, k_aRepresenting the azimuth modulation frequency, a (-) represents the signal envelope. DSP will s₁(t_a) And s₂(t_a) Transforming to azimuth frequency domain to obtain S₁(f_a) And S₂(f_a) I.e. S₁(f_a) And S₂(f_a) Respectively signals transformed into the azimuth frequency domain：

Wherein

Substeps 8.2, for S₁(f_a) And S₂(f_a) Correlation (Correlation) is performed to obtain the relative movement amount △ n between the two view images, and the estimated value of the Doppler modulation frequency is calculated

Where PRF denotes the pulse repetition frequency of the system.

Substep 8.2, estimation value of Doppler frequency modulation by DSP end

The phase error of the Doppler FM estimation is obtained by double integration and is recorded as

Wherein ^ integral ^ (·) represents double integral, dt_aRepresenting time t along azimuth_aAnd (6) performing integral operation.

Substep 8.3 phase error estimation from doppler frequency modulation

Construction of the phase Compensation function H_MD：

Substep 8.4, the DSP end will ss_RCMCn,A(t_r,t_a；R_S) And H_MDMultiplying to obtain a data-based motion compensation signal:

ss_MDn,A(t_r,t_a；R_S)＝ss_RCMCn,A(t_r,t_a；R_S)·H_MD

step 9, the DSP port pair motion compensation signal ss based on data_MDn,A(t_r,t_a；R_S) Performing direction Fourier transform to obtain a signal sS in a range-Doppler domain_MDn,A(t_r,f_a；R_S) (ii) a Constructing an orientation matching filter function, and adopting the orientation matching filter function to sS_MDn,A(t_r,f_a；R_S) Performing azimuth compression to obtain an azimuth compression signal ss_ACn,A(t_r,t_a)。

Specifically, the DSP end calculates the azimuth frequency modulation rate to generate an azimuth compression reference function, and the ss is subjected to the comparison_MDn,A(t_r,t_a；R_S) Performing primary focusing to obtain ss_ACn,A(t_r,t_a) Wherein n is 1, 2. Step 9 comprises the following sub-steps:

substep 9.1, DSP port pair ss_MDn,A(t_r,t_a；R_S) And performing azimuth FFT to obtain an expression of the signals in the range-Doppler domain after estimation and compensation of the residual Doppler frequency modulation rate:

wherein sinc (·) represents sinc function, the first exponential term reflects the azimuth position of the target, and the second exponential term is azimuth frequency modulation term.

Substep 9.2, constructing an orientation-matched filter function based on the second exponential term in the above equation:

substep 9.3, DSP compares sS_MDn,A(t_r,f_a；R_S) And H_AC(f_a) Multiplying and performing an azimuth IFFT to obtain an azimuth compression signal ss_ACn,A(t_r,t_a)：

Wherein the IFFT_A[·]Representing an azimuthal inverse fourier transform.

Step 10, the DSP end obtains an estimation value of the azimuth phase error by utilizing a phase gradient self-focusing algorithm, and an azimuth compression signal ss is obtained according to the estimation value of the azimuth phase error_ACn,A(t_r,t_a) Corrected to obtain accurate range block focusing signal ss_PGAn,A(t_r,t_a)。

Specifically, PGA processing is performed for a high-order phase error in the image domain, and ss is processed_ACn,A(t_r,t_a) Further compensating the residual phase error to obtain a more accurate range-blocked focus signal ss_PGAn,A(t_r,t_a) Wherein n is 1, 2.

Step 10 comprises the following substeps:

substep 10.1, referring to fig. 9, first the DSP takes out the azimuth compressed signal ss_ACn,A(t_r,t_a) The maximum amplitude of each range gate in the azimuth direction and its index value (azimuth sampling position) are determined for a number of range cells with strong energy in the matrix.

Substep 10.2, the saliency points in the image are all located at the centre of the image by a circumferential shift.

Substep 10.3, a segment of complex data centered at the maximum is truncated by windowing and the first derivative of the segment of complex data is calculated.

Substep 10.4, estimating a first derivative of the azimuth phase error, and integrating the first derivative to obtain an estimated value of the azimuth phase error;

and a substep 10.5 of correcting the azimuth signal using the estimate of the azimuth phase error. And repeating the above processes, wherein the window length is reduced by half every cycle until the judgment condition is met, namely the window length is reduced to 1-3 azimuth units, and then all phase errors are corrected. Since PGA processing takes the idea of transform domain, multiple azimuthal FFTs and IFFTs appear in fig. 9.

Step 11, obtaining accurate range block focusing signal ss from multiple DSP terminals_PGAn,A(t_r,t_a) And carrying out image splicing according to the number G of the overlapped points of the adjacent data blocks, and assembling into a complete SAR image.

Specifically, the DSP2 will ss_PGA2,A(t_r,t_a) Sending to DSP1, DSP1 dividing ss according to the overlapping point number 512 of distance blocks_PGA1,A(t_r,t_a) And ss_PGA2,A(t_r,t_a) Splicing the images, and gathering the images into a complete SAR image; and then the DSP1 sends the complete SAR image to a system board through a PCIe interface, and the unmanned airborne SAR real-time imaging processing design is completed.

The improvement effect of the present invention can be further explained by the following actual measurement test.

And carrying out an SAR system on the unmanned aerial vehicle platform, and carrying out imaging test on the ground cooperation area. Since the unmanned platform is light in weight, the air flow in the air will cause the platform to jolt violently. The unmanned aerial vehicle-mounted SAR real-time imaging processing software framework provided by the invention is loaded on a bumpy unmanned aerial vehicle platform to carry out SAR imaging processing, and the obtained imaging result is shown in figure 10.

By comparing the two images, the software architecture provided by the invention can be summarized to bring the following improvements to the SAR imaging of the bumpy platform.

(1) The invention eliminates defocusing and tailing phenomena caused by the motion error of the jolting flying platform;

(2) according to the invention, through good focusing on the special display points, the number of pixel value saturation points in the image is effectively reduced, so that the quantization coefficient is controlled in a normal range, and the problem of integral darkening of the image is solved;

(3) the method carries out estimation and compensation aiming at the influence of the motion error on the beam pointing, and effectively inhibits the influence of the out-of-band clutter on the imaging algorithm.

The motion compensation method described by the invention can well focus the irradiation area, and no tailing phenomenon occurs at the periphery of the strong scatterer. The real-time processing software framework designed by the invention completes the imaging task quickly and efficiently, meets the requirement of the real-time performance of the system, and can continuously detect and image a large scene.

Although the present invention has been described in detail in this specification with reference to specific embodiments and illustrative embodiments, it will be apparent to those skilled in the art that modifications and improvements can be made thereto based on the present invention. Accordingly, such modifications and improvements are intended to be within the scope of the invention as claimed.

Claims (7)

1. A bumpy platform unmanned aerial vehicle carries SAR real-time imaging processing method is characterized by comprising the following steps:

step 1, receiving and transmitting an LFM pulse signal by using an unmanned airborne SAR system to obtain an intermediate frequency echo real signal s_r(t_r) (ii) a Using analog-to-digital converter pairs s_r(t_r) After sampling and quantization, a discrete intermediate frequency echo real signal s is obtained_r(k) FPGA end pair s_r(k) Performing digital down-conversion processing to generate zero intermediate frequency orthogonal dual-path signal s_I(k) And s_Q(k) And integrated as time domain echoes ss (t)_r,t_a；R_S)＝s_I(k)+js_Q(k) (ii) a Wherein, t_rAs a fast time variable of distance, t_aFor azimuthal slow time variation, R_SThe center slant distance of the scene is shown, and j is an imaginary number unit;

step 2, the FPGA end judges the bandwidth B of the transmitting signal_rAnd ADC sampling frequency f_sThe relation between the time domain echo ss (t) and the time domain echo t (t) is selected to be in a direct sampling mode or a de-skew mode_r,t_a；R_S) Performing a radial pulseCompressing to obtain a range-wise pulse compression signal ss_PC,R(t_r,t_a；R_S)；

Step 3, FPGA end-to-end distance direction pulse compression signal ss_PC,R(t_r,t_a；R_S) Performing azimuth pre-filtering to obtain an azimuth filtered signal ss_PF,R(t_r,t_a；R_S)；

Step 4, the FPGA end enables the azimuth-filtered signal ss_PF,R(t_r,t_a；R_S) Partitioning along the distance direction to obtain a plurality of data blocks, wherein the number of overlapping points of adjacent data blocks is G, and each data block is sent to a DSP end; wherein each data block is denoted as ss_n,R(t_r,t_a；R_S) N is the number of data blocks;

step 5, the DSP end obtains the phase error

Using said phase error

For ss_n,R(t_r,t_a；R_S) Performing INS-based motion compensation to obtain an INS-based motion compensation signal ss_INSn,R(t_r,t_a；R_S)；

Step 6, the DSP end pair is used for motion compensation signal ss based on INS_INSn,R(t_r,t_a；R_S) Matrix transposition is carried out to obtain a transposed signal ss_INSn,A(t_r,t_a；R_S)；

Step 7, the DSP end contra-rotating signal ss_INSn,A(t_r,t_a；R_S) Sequentially carrying out distance Fourier transform and azimuth Fourier transform to obtain two-dimensional frequency domain echo signal SS_INSn,A(f_r,f_a；R_S) (ii) a Constructing a range migration correction function H_RCMC(f_r,f_a；R_S) Two-dimensional frequency domain echo signal SS_INSn,A(f_r,f_a；R_S) And H_RCMC(f_r,f_a；R_S) Multiplying and performing two-dimensional inverse Fourier transform to obtain range migration correction signal ss_RCMCn,A(t_r,t_a；R_S)；

Step 8, the DSP end constructs a phase compensation function H_MDCorrection of range migration signal ss_RCMCn,A(t_r,t_a；R_S) And H_MDMultiplying to obtain a data-based motion compensated signal ss_MDn,A(t_r,t_a；R_S)；

Step 9, the DSP port pair motion compensation signal ss based on data_MDn,A(t_r,t_a；R_S) Performing direction Fourier transform to obtain a signal sS in a range-Doppler domain_MDn,A(t_r,f_a；R_S) (ii) a Constructing an orientation matching filter function, and adopting the orientation matching filter function to sS_MDn,A(t_r,f_a；R_S) Performing azimuth compression to obtain an azimuth compression signal ss_ACn,A(t_r,t_a)；

Step 10, the DSP end obtains an estimation value of the azimuth phase error by utilizing a phase gradient self-focusing algorithm, and an azimuth compression signal ss is obtained according to the estimation value of the azimuth phase error_ACn,A(t_r,t_a) Corrected to obtain accurate range block focusing signal ss_PGAn,A(t_r,t_a)；

Step 11, obtaining accurate range block focusing signal ss from multiple DSP terminals_PGAn,A(t_r,t_a) And carrying out image splicing according to the number G of the overlapped points of the adjacent data blocks, and assembling into a complete SAR image.

2. The toss platform unmanned aerial vehicle SAR real-time imaging processing method according to claim 1, characterized in that step 1 comprises the following substeps:

substep 1.1, intermediate frequency echo real signal s_r(t_r) Comprises the following steps:

wherein, t_rAs a function of distance and fast time, c denotes the speed of light, f_IRepresenting the intermediate frequency, k, before ADC sampling_rIndicating the modulation frequency of the LFM signal;

v is horizontal velocity, X_nThe distance between the center point of the scene and any point distributed along the azimuth direction in the scene;

substep 1.2, discrete intermediate frequency echo real signal s_r(k) Comprises the following steps:

wherein k represents a discrete time variable, N_aRepresenting the number of azimuth sampling points;

substep 1.3, discrete intermediate frequency echo real signal s_r(k) Respectively with cos (2 π f)_Ik △ t) and sin (2 π f)_Ik △ t), then the two signals are respectively shifted, multiplied and accumulated by using low-pass filters f (k), and a zero intermediate frequency orthogonal two-path signal s is generated_I(k) And s_Q(k)：

Where △ t represents the time interval between adjacent sample points,

which represents the calculation of the convolution,

means that a summation operation is performed from M to 0 to M-1, M representing the number of total sampling points, and an increase in the summation variable M representing a shift of the vector;

substep 1.4, quadrature two-way signal s of zero intermediate frequency_I(k) And s_Q(k) Respectively doForm time domain echo ss (t) for real part and imaginary part of complex signal_r,t_a；R_S) Comprises the following steps:

wherein, w_r(t_r) Representing a distance-time-domain window function, w_a(t_a) Representing an azimuthal time domain window function, and λ is the wavelength.

3. The toss platform unmanned aerial vehicle SAR real-time imaging processing method according to claim 1, wherein the step 2 specifically comprises the following substeps:

substep 2.1, FPGA end-to-time domain echo ss (t)_r,t_a；R_S) Fourier transform is carried out to a distance frequency domain to obtain a distance echo:

Ss(f_r,t_a；R_S)＝FFT_R[ss(t_r,t_a；R_S)]

wherein, FFT_RRepresents a distance-to-fourier transform;

substep 2.2, based on the transmission signal bandwidth B_rAnd ADC sampling frequency f_sThe relationship between the distance and the pulse compression is selected as follows:

when B is present_r<f_sThen, the FPGA executes direct sampling pulse pressure processing; the method specifically comprises the following steps: according to the frequency modulation rate k of LFM signal_rConstructing a frequency domain matched filtering reference function H_PC：

Wherein f is_rRepresenting a distance direction frequency variable;

will Ss (f)_r,t_a；R_S) And H_PCMultiplying, and performing inverse Fourier transform on the multiplied result to obtain a distance direction pulse compression signal:

ss_PC,R(t_r,t_a；R_S)＝IFFT_R[Ss(f_r,t_a；R_S)·H_PC]

wherein the IFFT_RRepresenting an inverse fourier transform of the distance;

when B is present_r≥f_sWhen the pulse voltage is not detected, the FPGA executes oblique pulse voltage removing processing; the method specifically comprises the following steps: constructing a declivity function S_c(f_i)：

Wherein f is_iIs a baseband distance frequency variable;

will Ss (f)_r,t_a；R_S) And S_c(f_i) Multiplying, namely sequentially performing distance inverse Fourier transform and windowed Fourier transform on the multiplied result to obtain a distance pulse compression signal:

ss_PC,R(t_r,t_a；R_S)＝FFT_R[w_r(t_r)·IFFT_R[Ss(f_r,t_a；R_S)·S_c(f_i)]]

wherein, w_r(t_r) Representing a distance-time-domain window function, FFT_RRepresenting a distance fourier transform.

4. The toss platform unmanned aerial vehicle SAR real-time imaging processing method according to claim 1, characterized in that step 5 comprises the following substeps:

substep 5.1, the DSP terminal analyzes the east velocity v under the geographic coordinate system transmitted by the INS_eastNorth velocity v_northVelocity v_skyInformation, combined with track angle theta_yawCalculating the along-track velocity v under the imaging coordinate system_alongVertical track velocity v_crossAnd a lifting speed v_eleva：

Substep 5.2, DSP end pair imaging stationTrack-following velocity v under a coordinate system_alongVertical track velocity v_crossAnd a lifting speed v_elevaPerforming linear integration to obtain three-dimensional motion error e under the coordinate system_along、e_cross、e_elevaCalculating the motion error e in the LOS direction according to the pitch angle β of the beam_LOSAnd corresponding phase error

Wherein λ is the wavelength;

substep 5.3. mixing

And ss_n,R(t_r,t_a；R_S) Multiplying to obtain an INS-based motion compensation signal:

5. the toss platform unmanned aerial vehicle SAR real-time imaging processing method according to claim 2, characterized in that step 7 comprises the following substeps:

substep 7.1, DSP port pair ss_INSn,A(t_r,t_a；R_S) Performing distance-to-Fourier transform to obtain a distance frequency domain expression:

wherein, W_r(f_r) As a function of the distance frequency domain window, w_a(t_a) Representing the azimuthal time-domain window function, f_cCarrier frequency of radar system, f_rRepresenting a distance direction frequency variable;

substep 7.2, DSP port pair Ss_INSn,A(f_r,t_a；R_S) Carrying out azimuth Fourier change to obtain a two-dimensional frequency domain echo signal:

wherein, W_a(f_a) Representing the azimuth frequency domain window function, f_aRepresenting an azimuth frequency variable;

substep 7.3, relating the second exponential term in the above equation to the variable f_rLow-order mullin expansion:

from which f is found_rAnd f_aAnd constructing a range migration correction function:

substep 7.4, the DSP side will SS_INSn,A(f_r,f_a；R_S) And H_RCMC(f_r,f_a；R_S) Multiplying and performing two-dimensional inverse Fourier transform to obtain a range migration correction signal:

ss_RCMCn,A(t_r,t_a；R_S)＝IFFT₂[SS_INSn,A(f_r,f_a；R_S)·H_RCMC(f_r,f_a；R_S)]

wherein the IFFT₂[·]Representing a two-dimensional inverse fourier transform.

6. The toss platform unmanned aerial vehicle SAR real-time imaging processing method according to claim 1, characterized in that step 8 comprises the following substeps:

substep 8.1, DSP end will orient signal vector s (t)_a) Is divided into a front part and a rear partSub-apertures, respectively denoted as s₁(t_a) And s₂(t_a)：

Wherein T represents the synthetic aperture time, k_aRepresenting the azimuth modulation frequency, a (-) represents the signal envelope; DSP will s₁(t_a) And s₂(t_a) Transforming to azimuth frequency domain to obtain S₁(f_a) And S₂(f_a)：

Wherein

Wherein f is_aRepresenting an azimuth frequency variable;

substeps 8.2, for S₁(f_a) And S₂(f_a) Performing correlation operation to obtain relative shift amount △ n between two view images, and calculating Doppler frequency modulation estimated value according to relative shift amount △ n

Wherein N is_aThe number of sampling points in the azimuth direction is represented, and the PRF represents the pulse repetition frequency of the system;

substep 8.2 estimation of the Doppler modulation frequency

The phase error of the Doppler FM estimation is obtained by double integration and is recorded as

Wherein ^ integral ^ (·) represents double integral, dt_aRepresenting time t along azimuth_aPerforming integral operation;

substep 8.3 phase error estimation from doppler frequency modulation

Construction of the phase Compensation function H_MD：

Substep 8.4, let ss_RCMCn,A(t_r,t_a；R_S) And H_MDMultiplying to obtain a data-based motion compensated signal:

ss_MDn,A(t_r,t_a；R_S)＝ss_RCMCn,A(t_r,t_a；R_S)·H_MD。

7. the toss platform unmanned aerial vehicle SAR real-time imaging processing method according to claim 2, characterized in that step 9 comprises the following substeps:

substep 9.1, DSP port pair ss_MDn,A(t_r,t_a；R_S) And (3) carrying out azimuth Fourier change to obtain an expression of signals in a range Doppler domain after estimation and compensation of residual Doppler frequency modulation rate:

wherein sinc (·) represents sinc function, the first exponential term reflects the azimuth position of the target, and the second exponential term reflects the azimuth position of the targetThe plurality of terms are azimuth frequency modulation terms; w_a(f_a) Representing the azimuth frequency domain window function, f_aRepresenting the azimuthal frequency variation, f_cIs the carrier frequency of the radar system;

substep 9.2, constructing an orientation-matched filter function based on the second exponential term in the above equation:

substep 9.3. converting sS_MDn,A(t_r,f_a；R_S) And H_AC(f_a) Multiplying and performing azimuth inverse Fourier transform to obtain an azimuth compression signal ss_ACn,A(t_r,t_a)：

Wherein the IFFT_A[·]Representing an azimuthal inverse fourier transform.

Images