Flying target inter-pulse and inter-frame combined imaging method based on satellite-borne platform

Abstract

The invention discloses a flight target inter-pulse and inter-frame combined imaging method based on a satellite-borne platform, and relates to the field of ISAR imaging, in particular to a flight target inter-pulse and inter-frame combined imaging method based on the satellite-borne platform. The objective of the invention is to solve the problem that ISAR imaging quality is reduced in a flying target imaging process of an existing spaceborne platform. The method comprises the following steps: 1, dividing observation time; 2, constructing a signal model of imaging of a flying target by the satellite-borne platform after translation compensation in an imaging time period; 3, obtaining an optimal solution estimation value of the equivalent rotation speed and an optimal solution estimation value of the equivalent rotation acceleration; 4, a phase compensation function is constructed, a focused ISAR image is obtained, and calibration is completed; 5, repeating the steps 3 and 4 to obtain focused and calibrated ISAR images in different imaging time periods; and carrying out image registration on the focused and calibrated ISAR image, and carrying out wavelet transform on the registered image to realize image fusion.

Description

Technical Field

[0001] The invention relates to the field of ISAR imaging, and in particular to a method for inter-pulse and inter-frame joint imaging of a flying target based on a satellite-borne platform. Background Art

[0002] Imaging of flying targets has always been a hot topic in the field of ISAR imaging. Traditional ground-based imaging platforms are limited by the platform location and geographical terrain, and airborne and shipborne imaging platforms are limited by the range, making it difficult to achieve continuous and stable detection. Therefore, consider imaging flying targets with satellite platforms equipped with ISAR. However, during long observation times, the satellite-borne platform and the aerial target move along their respective orbits. The complex motion between the platform and the target will cause the image projection surface to be non-stationary, resulting in rotational range migration and two-dimensional spatial phase error. In the ISAR imaging model, phase error cannot be ignored. Among the current ISAR compensation methods for complex motion situations, the non-parametric method has a degraded algorithm performance when the target has high maneuverability; the parametric method regards the signal in each range unit as an independent part, ignoring the integrity of the target's rotational motion, thereby reducing the accuracy of parameter estimation. At the same time, when detecting and imaging the target, some scattering points on the target are blocked, resulting in the loss of some scattering points during ISAR imaging. Summary of the invention

[0003] The purpose of the present invention is to solve the problem of reduced ISAR imaging quality caused by rotational range migration, two-dimensional spatial phase error, and missing scattering points during the imaging of flying targets by existing satellite platforms, and to propose an inter-pulse and inter-frame joint ISAR imaging method based on a satellite platform.

[0004] A method for inter-pulse and inter-frame joint imaging of a flying target based on a satellite-borne platform has the following specific processes:

[0005] Step 1: Estimate the three-dimensional rotation angle of the flying target using narrowband tracking data, and calculate the equivalent rotation speed ω of the target based on the three-dimensional rotation angle of the flying target e and the equivalent rotational acceleration α e , based on the equivalent rotation speed ω of the target e and the equivalent rotational acceleration α e Calculate the distance traveled R caused by rotation Rot,p , the observation time is divided into time periods based on the condition that the distance traveled does not exceed one distance unit;

[0006] Step 2: Using the result of dividing the observation time obtained in step 1, a signal model of the imaging of the flying target by the satellite platform after translation compensation within the imaging time period is constructed;

[0007] Step 3: Obtain image entropy based on the signal model constructed in step 2, and obtain the optimal solution estimate of the parameter vector M by minimizing the image entropy, that is, the optimal solution estimate of the equivalent rotation speed. and the optimal solution estimate of the equivalent rotational acceleration

[0008] Step 4: Use the optimal solution estimate of the equivalent rotation speed obtained in step 3 and the optimal solution estimate of the equivalent rotational acceleration Construct phase compensation function, obtain focused ISAR image, and complete ISAR image calibration;

[0009] Step 5: Repeat steps 3 and 4 to obtain focused and calibrated ISAR images at different imaging time periods; perform image registration on the focused and calibrated ISAR images, and perform wavelet transform on the registered images to achieve image fusion.

[0010] The beneficial effects of the present invention are:

[0011] The present invention proposes an imaging method based on ISAR of a satellite platform. The present invention relates to the field of ISAR imaging, and in particular to an inter-pulse and inter-frame joint ISAR imaging method based on a satellite platform. The present invention mainly aims at the problem that when a satellite platform images a flying target, there is a non-stationary imaging plane caused by complex relative motion, and the target itself occludes and causes the missing of scattering points. The main contents of the present invention include: first, the long observation time is divided into imaging time periods using narrow-band tracking data, and a satellite-borne ISAR signal model after the imaging time period is divided is constructed, and a phase error model is further derived; then, a parameter optimization algorithm based on minimum image entropy search is proposed; further, a phase compensation function is constructed using the estimated motion parameters to obtain a focused and calibrated ISAR image in the current imaging time period; finally, the ISAR images in different imaging time periods are fused to obtain a well-focused and complete ISAR image.

[0012] In order to solve the above-mentioned problems such as rotational range migration, two-dimensional spatial phase error, missing scattering points, etc., the present invention proposes a pulse-to-pulse and frame-to-frame joint ISAR imaging method for a satellite platform through long observation time segmentation, parameter optimization based on minimum image entropy search, and image fusion.

[0013] The present invention proposes an inter-pulse and inter-frame joint imaging algorithm for a flying target based on a satellite-borne platform. Firstly, the long observation time is divided into imaging time periods using narrow-band tracking data, and a satellite-borne ISAR signal model after the imaging time period division is constructed, and a phase error model is further derived; then, a parameter optimization algorithm based on minimum image entropy search is proposed, and an ISAR image that is focused and calibrated in the current imaging time period is obtained; finally, ISAR images in different imaging time periods are fused to obtain an ISAR image with good focus and complete targets. BRIEF DESCRIPTION OF THE DRAWINGS

[0014] Figure 1 It is a flow chart of the present invention;

[0015] Figure 2 The imaging results of different algorithms are shown in Figure 2. (a) RD, (b) STFT, (c) DSA, (d) PGA, (e) IGCPF, and (f) the algorithm proposed by the present invention.

[0016] Figure 3 This is the ISAR image calibration result diagram;

[0017] Figure 4 The calibration results of ISAR images in different imaging time periods, (a) reference image, (b) image to be fused;

[0018] Figure 5 The image to be fused is the registration result map;

[0019] Figure 6 To fuse the images. DETAILED DESCRIPTION

[0020] Specific implementation method 1: This implementation method is based on a satellite-based platform to conduct inter-pulse and inter-frame joint imaging of a flying target. The specific process is as follows:

[0021] Step 1: Estimate the three-dimensional rotation angle of the flying target using narrowband tracking data, and calculate the equivalent rotation speed ω of the target based on the three-dimensional rotation angle of the flying target e and the equivalent rotational acceleration α e , based on the equivalent rotation speed ω of the target e and the equivalent rotational acceleration α e Calculate the distance traveled R caused by rotation Rot,p , the long observation time is divided by taking the distance traveled not exceeding one distance unit as the time period division condition;

[0022] Step 2: Using the result of dividing the long observation time obtained in step 1, a signal model of the imaging of the flying target by the spaceborne platform after translation compensation within the imaging time period is constructed;

[0023] Step 3: Obtain image entropy based on the signal model constructed in step 2, and obtain the optimal solution estimate of the parameter vector M by minimizing the image entropy, that is, the optimal solution estimate of the equivalent rotation speed. and the optimal solution estimate of the equivalent rotational acceleration

[0024] Step 4: Use the optimal solution estimate of the equivalent rotation speed obtained in step 3 and the optimal solution estimate of the equivalent rotational acceleration Construct phase compensation function, obtain focused ISAR image, and complete ISAR image calibration;

[0025] Step 5: Repeat steps 3 and 4 to obtain well-focused and calibrated ISAR images at different imaging time periods; perform image registration on the well-focused and calibrated ISAR images, and perform wavelet transform on the registered images to achieve image fusion.

[0026] Specific implementation method 2: The difference between this implementation method and specific implementation method 1 is that in step 1, the three-dimensional rotation angle of the flying target is estimated using narrowband tracking data, and the equivalent rotation speed ω of the target is calculated based on the three-dimensional rotation angle of the flying target. e and the equivalent rotational acceleration α e , based on the equivalent rotation speed ω of the target e and the equivalent rotational acceleration α e Calculate the distance traveled R caused by rotation Rot,p , the long observation time is divided into time periods with the distance traveled not exceeding one distance unit as the time period division condition; the specific process is:

[0027] Step 1: In order to effectively track the target, the radar system also has a narrowband ranging function, which can measure the position of the target at each moment in the radar coordinate system. r (t m ), t m is the azimuth time domain, i.e. the slow time dimension; t m =mT r , T r represents the pulse repetition period, m is the mth accumulated pulse, m=0,1,…,M-1, M is the number of accumulated pulses;

[0028] Radar line of sight unit vector i los (t m ) is expressed as

[0029]

[0030] Among them, | |2 represents the second-order norm;

[0031] The velocity vector of the target is expressed as

[0032]

[0033] Among them, v n (t m ) represents the speed of the target along the n-axis of the radar coordinate system, v u (t m ) represents the velocity of the target along the u-axis of the radar coordinate system, v w (t m ) represents the speed of the target along the w-axis of the radar coordinate system;

[0034] The origin of the radar coordinate system is located at the center of the radar, the u-axis of the radar coordinate system is parallel to the movement direction of the radar satellite platform, the w-axis of the radar coordinate system points away from the center of the earth, and the n-axis of the radar coordinate system is obtained by the right-hand rule;

[0035] The acceleration vector of the target is expressed as

[0036]

[0037] Among them, a n (t m ) represents the acceleration of the target along the n-axis of the radar coordinate system, a u (t m ) represents the acceleration of the target along the u-axis of the radar coordinate system, a w (t m ) represents the acceleration of the target along the w-axis of the radar coordinate system;

[0038] The three-dimensional rotation angle is obtained based on the target's velocity vector and the target's acceleration vector. The three-dimensional rotation angle is expressed as

[0039]

[0040]

[0041]

[0042] Among them, θ w (t m ) represents the rotation angle of the target around the w-axis of the radar coordinate system, θ u (t m ) represents the rotation angle of the target around the u-axis of the radar coordinate system, θ n (t m ) represents the rotation angle of the target around the n-axis of the radar coordinate system, g is the gravity of the earth’s surface, R(t m ) is the curvature radius of the target turn;

[0043]

[0044] Step 1 and 2: Use the extended Kalman filter to calculate the three-dimensional rotation angle θ w (t m ), θ u (t m ), θ n (t m ) to obtain accurate three-dimensional rotation angle

[0045] Based on accurate 3D rotation angles Get the rotation vector ω caused by the target rotation s (t m ),ω s (t m ) is expressed as

[0046]

[0047] Among them, T means to find the transpose;

[0048] The rotation vector ω caused by the change in the radar line of sight direction l (t m ) is expressed as

[0049]

[0050] Among them, × represents vector cross product;

[0051] The rotation vector ω caused by the target rotation s (t m ) and the rotation vector ω caused by the change in the radar line of sight direction l (t m ), get the actual rotation vector ω i (t m ),ω i (t m ) is expressed as

[0052] ω i (t m )=ω l (t m )+ω s (t m )

[0053] The equivalent rotation speed vector ω(t m ) is the actual rotation vector ω i (t m ) is perpendicular to the radar line of sight unit vector i los (t m ) is calculated by the following formula

[0054]

[0055] ω e (t m )=||ω(t m )||2 represents t m The equivalent rotation speed of the target on the radar projection plane at the moment, || ||2 represents the second-order norm of the vector;

[0056] The equivalent rotational acceleration of the target on the projection plane is expressed as

[0057]

[0058] The equivalent rotation angle of the target within the imaging time is expressed as

[0059]

[0060] Where, t is time;

[0061] Then, during the imaging time, the distance change between a scattering point P on the target and the radar caused by rotation is:

[0062] R Rot,p (t m )≈y P cos(θ R (t m ))+x P sin(θ R (t m ))

[0063] Among them, (x P ,y P ) is the coordinate of the scattering point P on the radar projection plane;

[0064] The condition for dividing the time period is that the distance migration of the scattering point within the time period is less than one distance resolution unit, that is,

[0065] R Rot,p (t m )≤Δr

[0066] Where Δr represents the distance resolution unit, c represents the speed of light, and B represents the signal bandwidth;

[0067] If the above conditions are met, it can be considered that the scattering point does not migrate within the imaging time period. This condition can be used to segment the long observation time. The total number of segments divided into time periods is Q, and the length of each time period after division is T. q ,q=1,2,...,Q.

[0068] The other steps and parameters are the same as those in the first embodiment.

[0069] Specific implementation method three: This implementation method is different from specific implementation method one or two in that: in step two, the result of dividing the long observation time obtained in step one is used to construct a signal model S (f) of the satellite-borne platform imaging the flying target after translation compensation within the imaging time period. r ,t m ); the specific process is:

[0070] Step 21: Compensate the target echo for translational motion within the divided time period. At this time, the equivalent rotation speed within the target imaging time period is assumed to be ω e (Unknown quantity, to be estimated later), the equivalent rotation acceleration during the target imaging period is assumed to be α e (Unknown quantity, to be estimated later);

[0071] Then the rotation angle change of the scattering point P is expressed as

[0072]

[0073] Where θ0 is the initial rotation angle of the scattering point P;

[0074] Step 2. The distance traveled caused by the rotational motion is expressed as

[0075]

[0076] The trigonometric function of the above equation is approximated by using the Maclaurin formula, which is expressed as

[0077]

[0078]

[0079] The distance traveled due to the rotational motion is simplified to

[0080]

[0081] Step 2: The received signal spectrum S (f r ,t m )for

[0082]

[0083] in,

[0084] ·Represents a multiplication sign;

[0085] σ p is the scattering intensity at scattering point p;

[0086] f r represents the fast time dimension frequency;

[0087] λ is the wavelength, λ=c / f c , f c is the carrier frequency;

[0088] j is the imaginary unit, j 2 = -1;

[0089] After time division, the rotation distance migration can be ignored, so the received signal spectrum S(f r ,t m ) can be ignored, and S(f r ,t m ) is simplified to S'(f r ,t m );

[0090] S'(f r ,t m ) is divided into two parts, one part is the signal without phase error The remainder is the phase error caused by the rotational motion

[0091]

[0092]

[0093] Among them, S'(f r ,t m ) represents the spectrum of the received signal after ignoring the rotation distance migration, A signal indicating the absence of phase error;

[0094] Step 24: S'(f r ,t m )and Performing a discrete inverse Fourier transform in the distance dimension, we get

[0095]

[0096]

[0097] Among them, s(k,n) represents the pair S'(f r ,t m ) is the result of discrete inverse Fourier transform in the distance dimension;

[0098] Express The result after performing a discrete inverse Fourier transform in the distance dimension;

[0099] n represents t mDiscrete form, n = 0, 1, ..., N-1, N represents the number of pulses in the time period after the time period is divided (m is the index before the time period is divided, and the current n is the index after the time period is divided);

[0100] k represents the distance dimension unit index, k = 0, 1, ..., K-1, K represents the number of distance units;

[0101] Step 25: s(k,n) and Performing a discrete inverse Fourier transform in the slow time dimension, we obtain

[0102]

[0103]

[0104] Among them, g(k,h) represents the result of discrete inverse Fourier transform of s(k,n) in the slow time dimension;

[0105] Express The result after performing a discrete inverse Fourier transform in the slow time dimension;

[0106] h represents the discrete sequence index after the slow-time discrete inverse Fourier transform, h = 0, 1, ..., N-1;

[0107] (k,h) represents the pixel position in the image domain;

[0108] (y P ,ω1x p ) represents the position of the scattering point P in the image;

[0109] The distance resolution and azimuth resolution are obtained by P ,ω1x p ), expressed as

[0110]

[0111]

[0112] Substituting the above formula into g(k,h) we get

[0113]

[0114] The other steps and parameters are the same as those in the first or second embodiment.

[0115] Specific implementation method 4: This implementation method is different from any one of specific implementation methods 1 to 3 in that: in step 3, the image entropy is obtained based on the signal model constructed in step 2, and the optimal solution estimate of the parameter vector M is obtained by minimizing the image entropy, that is, the optimal solution estimate of the equivalent rotation speed and the optimal solution estimate of the equivalent rotational acceleration

[0116] The specific process is:

[0117] Step 3: In order to obtain a focused ISAR image, a phase compensation term needs to be constructed to compensate the echo signal s(k,n). The parameter vector required for compensation is defined as M = (ω e ,α e );

[0118] The phase compensation term is expressed as

[0119]

[0120] The signal after phase compensation is expressed as

[0121]

[0122] Image entropy E g (M) is defined as a function of M, expressed as

[0123]

[0124]

[0125] Among them, S g is the image density; || 2 represents the square of the absolute value;

[0126] Step 32: Obtain the optimal solution estimate of the parameter vector M by minimizing the image entropy, that is, the optimal solution estimate of the equivalent rotation speed and the optimal solution estimate of the equivalent rotational acceleration

[0127] The other steps and parameters are the same as those in Specific Embodiments 1 to 3.

[0128] Specific implementation method 5: This implementation method is different from specific implementation methods 1 to 4 in that: in step 32, the optimal solution estimate of the parameter vector M is obtained by minimizing the image entropy, that is, the optimal solution estimate of the equivalent rotation speed and the optimal solution estimate of the equivalent rotational acceleration The specific process is:

[0129] The optimization function is constructed as

[0130]

[0131] in, represents the optimal solution estimate of the equivalent rotation speed, represents the estimated value of the equivalent rotational acceleration, represents the optimal solution estimate of the parameter vector M;

[0132] The parameter vector M = (ω e ,α e ) is used as the initial input value of the optimization function, and the quasi-Newton method BFGS (Broyden-Fletcher-Goldfarb-Shanno) is used to find the optimal solution estimate of the equivalent rotation speed and the optimal solution estimate of the equivalent rotational acceleration

[0133] After the above analysis, the idea of ​​phase error estimation based on minimum entropy optimization is established. By minimizing the image entropy, the optimal parameters can be obtained. The estimated parameters are used to achieve phase compensation, which ultimately significantly improves the clarity of the image. In order to improve computational efficiency, the BFGS algorithm is used to find the optimal result. The initial value of the optimization algorithm can be obtained from step one.

[0134] The other steps and parameters are the same as those in Specific Embodiments 1 to 4.

[0135] Specific implementation method 6: This implementation method is different from any one of the specific implementation methods 1 to 5 in that: the optimal solution estimate of the equivalent rotation speed obtained in step 3 is used in step 4 and the optimal solution estimate of the equivalent rotational acceleration Construct phase compensation function, obtain focused ISAR image, and complete ISAR image calibration;

[0136] The specific process is:

[0137] Based on the optimal solution estimate of the equivalent rotation speed obtained in step 3 and the optimal solution estimate of the equivalent rotational acceleration Construct a phase compensation function, which is expressed as

[0138]

[0139] The focused ISAR image is expressed as

[0140]

[0141] Among them, the lateral calibration factor δ of the ISAR image is a And the distance dimension calibration factor δ r Respectively expressed as

[0142]

[0143]

[0144] At this point, a well-focused and calibrated ISAR image is obtained.

[0145] The other steps and parameters are the same as those in Specific Implementation Methods 1 to 5-1.

[0146] Specific implementation method 7: This implementation method is different from any one of specific implementation methods 1 to 6 in that: in step 5, step 3 and step 4 are repeated to obtain well-focused and calibrated ISAR images at different imaging time periods; image registration is performed on the well-focused and calibrated ISAR images, and wavelet transform is performed on the registered images to achieve image fusion;

[0147] The specific process is:

[0148] Step 5.1: The focused and calibrated ISAR images at different imaging time periods are expressed as RD q ,q=1,2,...,Q;

[0149] The optimal solution estimate of the equivalent rotation speed corresponding to each imaging time period in step 52 and step 3 is expressed as The optimal solution estimate of the equivalent rotational acceleration is expressed as

[0150] Then the rotation angle of the target in each imaging time period is expressed as

[0151]

[0152] The lateral calibration factor of the ISAR image in each imaging time period is expressed as

[0153]

[0154] The range dimension calibration factor of the ISAR image in each imaging time period is expressed as

[0155]

[0156] Step 5.3: Next, image registration is performed. The focused and calibrated ISAR image obtained in the first imaging time period is considered as the reference image. The rotation angle of the focused and calibrated ISAR image obtained in the qth imaging time period relative to the reference image is expressed as

[0157]

[0158] Step 54: Given the rotation angle between the reference image and the image to be registered, and the two-dimensional resolution δ of the image a,q ,q=1,2,...,Q and δ r , based on the reference image RD1 to be registered image RD q,q=2,3,...,Q performs rotation interpolation processing to obtain the registered image RD' q ,q=2,3,...,Q;

[0159] The reference image is the focused and calibrated ISAR image RD obtained in the first imaging time period. q ,q=1;

[0160] The image to be registered is the focused and calibrated ISAR image RD obtained in the 2nd, 3rd, ..., Qth imaging time period q ,q=2,3,...,Q;

[0161] Step 55: Compare the reference image RD1 and the registered image RD' q ,q=2,3,...,Q is used to perform wavelet transform to achieve image fusion.

[0162] The other steps and parameters are the same as those in Specific Implementations 1 to 6-1.

[0163] Specific implementation eight: This implementation differs from any one of specific implementations one to seven in that: in step five, the reference image RD1 and the registered image RD' q ,q=2,3,...,Q to perform wavelet transform to realize image fusion; the specific process is:

[0164] Step 551: Compare the reference image RD1 and the registered image RD' q ,q=2,3,...,Q performs wavelet transform and decomposes it into transform coefficients in different frequency domains;

[0165] Step 552: Select the transformation coefficient according to the principle of taking the maximum absolute value;

[0166] Step 553: Perform inverse wavelet transform on the selected transformation coefficients to obtain a fused image.

[0167] The other steps and parameters are the same as those in Specific Implementation 1 to 7-1.

[0168] Simulation analysis:

[0169] The inter-pulse and inter-frame joint ISAR imaging method based on a space-borne platform designed in the present invention is verified by simulation, and the results are explained.

[0170] Simulation 1 is used to verify the effectiveness of the present invention for ISAR imaging focusing. The relative motion parameters and simulation signal parameters are shown in Table 1. The imaging results corresponding to different algorithms are shown in Table 1. Figure 2 As shown, Figure 2 Medium (a) to Figure 2(f) shows the imaging results of RD algorithm, short-time Fourier transform (STFT), distinctive point algorithm (DSA), phase gradient autofocus (PGA), iterative generalized cubic phase function (IGCPF) and the method proposed in this invention, and the corresponding image entropy is shown in Table 2. Figure 2 It can be seen from Table 2 that after the image is focused, the algorithm proposed in the present invention can obtain an image with the highest focusing degree and the smallest image entropy.

[0171] Table 1 Simulation parameters

[0172]

[0173]

[0174] Table 2 Image entropy

[0175]

[0176] Simulation 2 is used to verify the effectiveness of the present invention for ISAR image calibration. The simulation parameters are shown in Table 1. Figure 3 This is the imaging result after calibration using the method proposed in the present invention, where the red dot is the actual projection of the target scattering point on the image projection plane, which is highly consistent with the image after calibration by the present invention, verifying the effectiveness of the calibration method proposed in the present invention.

[0177] Simulation 3 is used to verify the effectiveness of the present invention for ISAR image fusion. Figure 4 Two target images to be fused are shown, Figure 4 (a) is the reference image. Figure 4 (b) is the image to be fused. Figure 5 Give the result after the image to be fused is registered. Figure 6 It is the image after the two images are fused. It can be seen that the two images before fusion are missing the left wing and the right wing respectively. After the image fusion, the missing scattering points are completed, which verifies the effectiveness of the image fusion method proposed in the present invention.

[0178] The present invention may also have many other embodiments. Without departing from the spirit and essence of the present invention, those skilled in the art may make various corresponding changes and modifications based on the present invention, but these corresponding changes and modifications should all fall within the scope of protection of the claims attached to the present invention.

Claims

1. A method for inter-pulse and inter-frame joint imaging of a flying target based on a satellite-borne platform, characterized in that: The specific process of the method is: Step 1: Estimate the three-dimensional rotation angle of the flying target using narrowband tracking data, and calculate the equivalent rotation speed ω of the target based on the three-dimensional rotation angle of the flying target e and the equivalent rotational acceleration α e , based on the equivalent rotation speed ω of the target e and the equivalent rotational acceleration α e Calculate the distance traveled R caused by rotation Rot,p , the observation time is divided into time periods based on the condition that the distance traveled does not exceed one distance unit; Step 2: Using the result of dividing the observation time obtained in step 1, a signal model of the imaging of the flying target by the satellite platform after translation compensation within the imaging time period is constructed; Step 3: Obtain image entropy based on the signal model constructed in step 2, and obtain the optimal solution estimate of the parameter vector M by minimizing the image entropy, that is, the optimal solution estimate of the equivalent rotation speed. and the optimal solution estimate of the equivalent rotational acceleration Step 4: Use the optimal solution estimate of the equivalent rotation speed obtained in step 3 and the optimal solution estimate of the equivalent rotational acceleration Construct phase compensation function, obtain focused ISAR image, and complete ISAR image calibration; Step 5: Repeat steps 3 and 4 to obtain focused and calibrated ISAR images at different imaging time periods; perform image registration on the focused and calibrated ISAR images, and perform wavelet transform on the registered images to achieve image fusion.

2. The method for combined pulse-to-pulse and frame-to-frame imaging of a flying target based on a satellite-borne platform according to claim 1, characterized in that: In the step 1, the three-dimensional rotation angle of the flying target is estimated by using the narrowband tracking data, and the equivalent rotation speed ω of the flying target is calculated based on the three-dimensional rotation angle of the flying target. e and the equivalent rotational acceleration α e , based on the equivalent rotation speed ω of the target e and the equivalent rotational acceleration α e Calculate the distance traveled R caused by rotation Rot,p , the observation time is divided into time periods based on the condition that the distance traveled does not exceed one distance unit; The specific process is: Step 1: The radar system measures the position r of the target at each moment in the radar coordinate system r (t m ), t m is the azimuth time domain, i.e. the slow time dimension; t m =mT r , T r represents the pulse repetition period, m is the mth accumulated pulse, m=0,1,…,M-1, M is the number of accumulated pulses; Radar line of sight unit vector i los (t m ) is expressed as Among them, 2 represents the second-order norm; The velocity vector of the target is expressed as Among them, v n (t m ) represents the speed of the target along the n-axis of the radar coordinate system, v u (t m ) represents the velocity of the target along the u-axis of the radar coordinate system, v w (t m ) represents the speed of the target along the w-axis of the radar coordinate system; The origin of the radar coordinate system is located at the center of the radar, the u-axis of the radar coordinate system is parallel to the movement direction of the radar satellite platform, the w-axis of the radar coordinate system points away from the center of the earth, and the n-axis of the radar coordinate system is obtained by the right-hand rule; The acceleration vector of the target is expressed as Among them, a n (t m ) represents the acceleration of the target along the n-axis of the radar coordinate system, a u (t m ) represents the acceleration of the target along the u-axis of the radar coordinate system, a w (t m ) represents the acceleration of the target along the w-axis of the radar coordinate system; The three-dimensional rotation angle is obtained based on the target's velocity vector and the target's acceleration vector. The three-dimensional rotation angle is expressed as Among them, θ w (t m ) represents the rotation angle of the target around the w-axis of the radar coordinate system, θ u (t m ) represents the rotation angle of the target around the u-axis of the radar coordinate system, θ n (t m ) represents the rotation angle of the target around the n-axis of the radar coordinate system, g is the gravity of the earth’s surface, R(t m ) is the curvature radius of the target turn; Step 1 and 2: Use the extended Kalman filter to calculate the three-dimensional rotation angle θ w (t m ), θ u (t m ), θ n (t m ) to obtain the three-dimensional rotation angle Based on 3D rotation angle Get the rotation vector ω caused by the target rotation s (t m ),ω s (t m ) is expressed as Among them, T means to find the transpose; The rotation vector ω caused by the change in the radar line of sight direction l (t m ) is expressed as Among them, × represents vector cross product; The rotation vector ω caused by the target rotation s (t m ) and the rotation vector ω caused by the change in the radar line of sight direction l (t m ), get the actual rotation vector ω i (t m ),ω i (t m ) is expressed as oh i (t m )=ω l (t m )+ω s (t m ) The equivalent rotation speed vector ω(t m ) is the actual rotation vector ω i (t m ) is perpendicular to the radar line of sight unit vector i los (t m ) is calculated by the following formula ω e (t m )=ω(t m )2 represents t m The equivalent rotation speed of the target on the radar projection plane at the moment, 2 represents the second-order norm of the vector; The equivalent rotational acceleration of the target on the projection plane is expressed as The equivalent rotation angle of the target within the imaging time is expressed as Where, t is time; Then, during the imaging time, the distance change between a scattering point P on the target and the radar caused by rotation is R Rot,p (t m )≈y P cos(θ R (t m ))+x P sin(θ R (t m )) Among them, (x P ,y P ) is the coordinate of the scattering point P on the radar projection plane; The condition for dividing the time period is that the distance migration of the scattering point within the time period is less than one distance resolution unit, that is, R Rot,p (t m )≤Δr Where Δr represents the distance resolution unit, c represents the speed of light, and B represents the signal bandwidth; The total number of segments divided into time periods is Q, and the length of each time period after division is T q ,q=1,2,...,Q.

3. The method for combined pulse-to-pulse and frame-to-frame imaging of a flying target based on a satellite-borne platform according to claim 2, characterized in that: In the step 2, the result of dividing the observation time obtained in the step 1 is used to construct a signal model of the imaging of the flying target by the satellite-borne platform after translation compensation within the imaging time period; The specific process is: Step 21: Compensate the target echo for translational motion within the divided time period. At this time, the equivalent rotation speed within the target imaging time period is assumed to be ω e , the equivalent rotation acceleration during the target imaging period is assumed to be α e ; Then the rotation angle change of the scattering point P is expressed as Where θ0 is the initial rotation angle of the scattering point P; Step 2. The distance traveled caused by the rotational motion is expressed as The trigonometric function of the above equation is approximated by using the Maclaurin formula, which is expressed as The distance traveled due to the rotational motion is simplified to Step 2: The received signal spectrum S (f r ,t m )for in, ·Represents the multiplication sign; σ p is the scattering intensity of the scattering point p; f r represents the fast time dimension frequency; λ is the wavelength, λ=c / f c , f c is the carrier frequency; j is the imaginary unit, j 2 = -1; Received signal spectrum S(f r ,t m ) can be ignored, and S(f r ,t m ) is simplified to S'(f r ,t m ); S'(f r ,t m ) is divided into two parts, one part is the signal without phase error The remainder is the phase error caused by the rotational motion Among them, S'(f r ,t m ) represents the spectrum of the received signal after ignoring the rotation distance migration, A signal indicating the absence of phase error; Step 24: S'(f r ,t m )and Performing a discrete inverse Fourier transform in the distance dimension, we get Among them, s(k,n) represents the pair S'(f r ,t m ) is the result of discrete inverse Fourier transform in the distance dimension; Express The result after performing a discrete inverse Fourier transform in the distance dimension; n means t m Discrete form, n = 0, 1, ..., N-1, N represents the number of pulses in the time period after the time period is divided; k represents the distance dimension unit index, k = 0, 1, ..., K-1, K represents the number of distance units; Step 25: s(k,n) and Performing a discrete inverse Fourier transform in the slow time dimension, we obtain Among them, g(k,h) represents the result of discrete inverse Fourier transform of s(k,n) in the slow time dimension; Express The result after discrete inverse Fourier transform in the slow time dimension; h represents the discrete sequence index after the slow-time discrete inverse Fourier transform, h = 0, 1, ..., N-1; (k,h) represents the pixel position in the image domain; (y P ,ω1x p ) represents the position of the scattering point P in the image; The distance resolution and azimuth resolution are obtained by P ,ω1x p ), expressed as Substituting the above formula into g(k,h) we get 4. The method for combined pulse-to-pulse and frame-to-frame imaging of a flying target based on a satellite-borne platform according to claim 3, characterized in that: In step 3, the image entropy is obtained based on the signal model constructed in step 2, and the optimal solution estimate of the parameter vector M is obtained by minimizing the image entropy, that is, the optimal solution estimate of the equivalent rotation speed. and the optimal solution estimate of the equivalent rotational acceleration The specific process is: Step 3: Construct a phase compensation term to compensate the echo signal s(k,n). The parameter vector required for compensation is defined as M = (ω e ,α e ); The phase compensation term is expressed as The signal after phase compensation is expressed as Image entropy E g (M) is defined as a function of M, expressed as Among them, S g is the image density; 2 represents the square of the absolute value; Step 32: Obtain the optimal solution estimate of the parameter vector M by minimizing the image entropy, that is, the optimal solution estimate of the equivalent rotation speed and the optimal solution estimate of the equivalent rotational acceleration 5. The method for combined inter-pulse and inter-frame imaging of a flying target based on a satellite-borne platform according to claim 4, characterized in that: In step 32, the optimal solution estimate of the parameter vector M is obtained by minimizing the image entropy, that is, the optimal solution estimate of the equivalent rotation speed and the optimal solution estimate of the equivalent rotational acceleration The specific process is: The optimization function is constructed as in, represents the optimal solution estimate of the equivalent rotation speed, represents the estimated value of the equivalent rotational acceleration, represents the optimal solution estimate of the parameter vector M; The parameter vector M = (ω e ,α e ) is used as the initial input value of the optimization function, and the quasi-Newton method BFGS is used to find the optimal solution estimate of the equivalent rotation speed and the optimal solution estimate of the equivalent rotational acceleration 6. The method for combined inter-pulse and inter-frame imaging of a flying target based on a satellite-borne platform according to claim 5, characterized in that: The optimal solution estimate of the equivalent rotation speed obtained in step 3 is used in step 4. and the optimal solution estimate of the equivalent rotational acceleration Construct phase compensation function, obtain focused ISAR image, and complete ISAR image calibration; The specific process is: Based on the optimal solution estimate of the equivalent rotation speed obtained in step 3 and the optimal solution estimate of the equivalent rotational acceleration Construct a phase compensation function, which is expressed as The focused ISAR image is expressed as Among them, the lateral calibration factor δ of the ISAR image is a And the distance dimension calibration factor δ r Respectively expressed as At this point, a focused and calibrated ISAR image is obtained.

7. The method for combined pulse-to-pulse and frame-to-frame imaging of a flying target based on a satellite-borne platform according to claim 6, characterized in that: In the step 5, step 3 and step 4 are repeated to obtain focused and calibrated ISAR images at different imaging time periods; image registration is performed on the focused and calibrated ISAR images, and wavelet transform is performed on the registered images to realize image fusion; The specific process is: Step 5.1: The focused and calibrated ISAR images at different imaging time periods are expressed as RD q ,q=1,2,...,Q; The optimal solution estimate of the equivalent rotation speed corresponding to each imaging time period in step 52 and step 3 is expressed as The optimal solution estimate of the equivalent rotational acceleration is expressed as Then the rotation angle of the target in each imaging time period is expressed as The lateral calibration factor of the ISAR image in each imaging time period is expressed as The range dimension calibration factor of the ISAR image in each imaging time period is expressed as Step 5.3: The focused and calibrated ISAR image obtained in the first imaging time period is used as the reference image. The rotation angle of the focused and calibrated ISAR image obtained in the qth imaging time period relative to the reference image is expressed as: Step 54: The rotation angle between the reference image and the image to be registered, and δ a,q ,q=1,2,...,Q and δ r , based on the reference image RD1 to be registered image RD q ,q=2,3,...,Q performs rotation interpolation processing to obtain the registered image RD' q ,q=2,3,...,Q; Step 55: Compare the reference image RD1 and the registered image RD' q ,q=2,3,...,Q is used to perform wavelet transform to achieve image fusion.

8. The method for combined pulse-to-pulse and frame-to-frame imaging of a flying target based on a satellite-borne platform according to claim 7, characterized in that: In step 55, the reference image RD1 and the registered image RD′ are compared. q ,q=2,3,...,Q performs wavelet transform to realize image fusion; The specific process is: Step 551: Compare the reference image RD1 and the registered image RD' q ,q=2,3,...,Q performs wavelet transform and decomposes it into transform coefficients in different frequency domains; Step 552: Select the transformation coefficient according to the principle of taking the maximum absolute value; Step 553: Perform inverse wavelet transform on the selected transformation coefficients to obtain a fused image.

Images