WO2023015623A1 - 一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法 - Google Patents

一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法 Download PDF

Info

Publication number
WO2023015623A1
WO2023015623A1 PCT/CN2021/115882 CN2021115882W WO2023015623A1 WO 2023015623 A1 WO2023015623 A1 WO 2023015623A1 CN 2021115882 W CN2021115882 W CN 2021115882W WO 2023015623 A1 WO2023015623 A1 WO 2023015623A1
Authority
WO
WIPO (PCT)
Prior art keywords
segment
platform
distance
imaging
estimated
Prior art date
Application number
PCT/CN2021/115882
Other languages
English (en)
French (fr)
Inventor
雒梅逸香
徐丰
Original Assignee
复旦大学
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 复旦大学 filed Critical 复旦大学
Priority to US18/576,667 priority Critical patent/US20240319364A1/en
Publication of WO2023015623A1 publication Critical patent/WO2023015623A1/zh

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/904—SAR modes
    • G01S13/9041—Squint mode
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9004—SAR image acquisition techniques
    • G01S13/9011—SAR image acquisition techniques with frequency domain processing of the SAR signals in azimuth
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9021—SAR image post-processing techniques
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/904—SAR modes
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/904—SAR modes
    • G01S13/9064—Inverse SAR [ISAR]

Definitions

  • the invention belongs to the technical field of signal processing, and in particular relates to a segmented aperture imaging and positioning method of a multi-rotor UAV-borne synthetic aperture radar.
  • the airborne imaging algorithms adopted by research institutions at home and abroad mainly rely on the data obtained by the high-precision inertial navigation system for motion compensation and imaging, and because the servo system carried by the large-scale platform can compensate the impact of the platform attitude angle change in real time, There is no need to consider the platform angle change, so it cannot solve three problems in the multi-rotor UAV SAR system: 1) the trajectory of the multi-rotor platform is usually very unstable because it is easily affected by environmental factors such as wind, It may cause severe shaking and sharp turns; 2) The unique flight principle of multi-rotor will bring high-frequency error and incident angle error caused by rotor vibration and fuselage tilting operation; 3) High maneuvering trajectory and low-cost inertial navigation system combined, resulting in rather poor motion state and attitude angle data for motion compensation.
  • the positioning accuracy of the drone's navigation system will also decrease. If the aircraft's flight trajectory is reversed based on the imaging results during the flight, it can be The navigation system provides the basis for more intelligent track planning. Therefore, it is of great significance to study the imaging and positioning algorithm of multi-rotor UAV-borne synthetic aperture radar that does not rely on inertial navigation equipment.
  • the present invention aims at the shortcomings and deficiencies of the above-mentioned existing airborne imaging algorithms, and aims to provide a multi-rotor UAV-borne synthetic aperture radar segmented aperture imaging (SAI) and positioning method, which is suitable for non-inertial navigation equipment Or the synthetic aperture radar system of low-precision inertial navigation equipment to improve the focusing effect and imaging efficiency of radar imaging.
  • SAI segmented aperture imaging
  • the present invention provides a multi-rotor UAV-borne synthetic aperture radar segmented aperture imaging method, which is used to perform segmented aperture imaging according to the original echo signal of the synthetic aperture radar. It has such characteristics and includes the following steps: Step 1, Carry out range pulse compression on the original echo signal s(t, ⁇ ) to obtain the range pulse pressure signal s RC (t, ⁇ ), where t is the fast time for the distance, and ⁇ is the slow time for the azimuth.
  • the range pulse pressure signal s RC (t , ⁇ ) the phase history of the strong scattering point Computing the estimated velocity of the maneuvering platform and the estimated beam center oblique angle Step 2, according to the estimated speed
  • Step 3 according to the estimated velocity and the estimated beam center oblique angle Calculate the phase compensation amount of the segment pulse pressure signal s RC corresponding to each segment, i (t, ⁇ )
  • the segment pulse pressure signal s RC, i (t, ⁇ ) is multiplied by the motion error compensation filter
  • step 1 the phase of the strong scattering point in the distance pulse pressure signal s RC (t, ⁇ ) course Perform second-order fitting to obtain the phase history of the strong scattering point
  • step 2 the phase of the strong scattering point in the distance pulse pressure signal s RC (t, ⁇ ) course
  • step 2 Perform second-order fitting to obtain the phase history of the strong scattering point
  • t is the distance fast time
  • ⁇ is the azimuth slow time
  • o( ⁇ ) is the high-order phase error
  • ⁇ is the wavelength of the system transmitting signal
  • the segmented aperture imaging method of multi-rotor unmanned aerial vehicle-borne synthetic aperture radar provided by the present invention, it can also have such a feature: wherein, in step 3, the phase compensation amount
  • R 0 is the The mean value of , ⁇ 0 is the oblique angle of N estimated beam centers mean value.
  • the azimuth compression filter In the formula, f is the frequency corresponding to the fast distance time t, f d is the Doppler frequency corresponding to the azimuth slow time ⁇ , f c is the carrier frequency of the system transmitting signal, and c is the speed of light.
  • the segmented aperture imaging method of multi-rotor UAV-borne synthetic aperture radar provided by the present invention, it can also have such a feature: wherein, in step 6, the corresponding imaging result s IMG of the adjacent described segment, i ( t, ⁇ ) for geometric correction to obtain corrected imaging results The corrected imaging results are then to rotate degrees, to obtain corrected imaging results with the slant distance perpendicular to the trajectory of the maneuvering platform Then the corrected imaging results corresponding to the adjacent segments are sequentially Alignment and coherent accumulation of distance-dimension envelopes where strong focal points in overlapping areas are located, and splicing non-overlapping areas to obtain the final imaging result S all .
  • step 6 carries out geometric correction to the imaging result s IMG, and i is carried out in the following sub-steps: step 6 -1, carry out azimuth Fourier transform to imaging result s IMG, i (t, ⁇ ), obtain distance Doppler domain map s IMG, i (t, f d ); Step 6-2, according to Fourier Transform features and object space geometry to construct an expression for a skew correction filter that corrects image skew: In the formula, l is the range scale of a typical building; step 6-3, multiply the range-Doppler domain map s IMG,i (t, f d ) by the tilt correction filter H GC-1 to obtain the tilt-corrected frequency domain plot Step 6-4, for the tilt-corrected frequency domain map Perform azimuth inverse Fourier transform to obtain the tilt-corrected time-domain map Step 6
  • step 1 the estimated speed and estimated beam center oblique angle Calculate according to the following steps: step 1-1, according to the definition of Doppler FM slope K a and Doppler center f dc , get the calculation formula of Ka and f dc :
  • step 1-1 the space formed by the original echo signal is called the signal space, and the phase history of the strong scattering point in the distance pulse pressure signal s RC (t, ⁇ )
  • step ⁇ is the second-order item coefficient
  • ⁇ is the first-order item coefficient
  • o( ⁇ ) is the high-order phase error
  • step 1-3 the phase history Substituting the expressions of K a
  • step 4 the two-dimensional spectrum s MC is adopted by the series inversion method, i (f , f d ) decompose to obtain the azimuth compression filter expression H AC, the process of i is carried out according to the following sub-steps: Step 4-2-1, obtain the two-dimensional spectrum s MC corresponding to the segment according to the stationary phase point method, i (f , f d ) expression:
  • fc represents the carrier frequency of the signal transmitted by the system, and c represents the speed of light
  • step 4-2-2 the two-dimensional spectrum s MC, i (f, f The expression of d ) is decomposed, and the two-dimensional spectrum expression of eliminating the coupling term of f and f d is obtained: Step 4-2-3, constructing an ideal phase filter expression according to the two-dimensional spectrum expression of the coupling term of
  • the present invention also provides a multi-rotor unmanned aerial vehicle-borne synthetic aperture radar segmentation aperture imaging method, which is used to calculate the flight trajectory of the unmanned aerial vehicle according to the original echo signal of the unmanned aerial vehicle-borne synthetic aperture radar. It is characterized in that it has such feature, comprising the following steps: Step S1, performing distance pulse compression on the original echo signal s (t, ⁇ ) to obtain a distance pulse pressure signal s RC (t, ⁇ ), according to the distance pulse pressure signal s RC (t, Phase history of strong scattering points in ⁇ ) Computing the estimated velocity of the maneuvering platform and the estimated beam center oblique angle Step S2, according to the estimated speed
  • the platform motion error is accurately compensated by estimating the platform motion parameters based on the phase history of the echo signal and segmenting the echo signal, It achieves good imaging focusing effect and high imaging success rate, improves the efficiency of multi-rotor UAV platform data acquisition, and can effectively perform high-resolution imaging on the multi-rotor UAV platform synthetic aperture radar system, and can Effectively estimate the platform trajectory, obtain the coordinate information of the platform relative to the measurement area, and achieve the effect of platform positioning.
  • the present invention has the following advantages:
  • the present invention takes into account the relationship between the platform motion and the signal phase, and is applicable to the synthetic aperture radar imaging system without inertial navigation equipment or low-precision inertial navigation equipment for imaging;
  • the present invention takes into account the influence of squint caused by the change of attitude angle of the UAV platform, compensates the coupling phase of distance and azimuth, improves the image focusing effect, and is applicable to antennaless Servo synthetic aperture radar imaging system;
  • the present invention adopts the method of phase filter multiplication instead of the interpolation method, and improves the imaging speed of single imaging in each segment;
  • the invention adopts the method of parallel imaging of each segment and then splicing into a complete image to obtain images, thereby improving the imaging speed of the complete image.
  • the invention can estimate the position of the platform while imaging, achieve the effect of positioning the platform, and facilitate the establishment of an automatic, intelligent and integrated UAV detection system.
  • FIG. 2 schematic diagram of target space and signal space of the airborne synthetic aperture radar system in an embodiment of the present invention
  • Fig. 3 is a flow chart of geometry correction in an embodiment of the present invention.
  • Fig. 4 is a schematic diagram of trajectory estimation in an embodiment of the present invention.
  • Fig. 6 is a schematic diagram of the relationship between the platform trajectory and the spatial position of the lattice target in the embodiment of the present invention.
  • Fig. 7 is a comparison diagram of the imaging results of each segment, the full aperture imaging results and the imaging results of the traditional imaging method in the embodiment of the present invention.
  • Figure 8 is the comparison result of the point spread function between the present invention and the traditional imaging method in the embodiment of the present invention.
  • Fig. 9 is a comparison diagram of the traditional imaging method imaging results, the imaging results of this method, and the optical model in the embodiment of the present invention.
  • Fig. 10 is the comparison and error curve of the estimated trajectory curve of the platform and the inertial navigation measurement trajectory curve in the embodiment of the present invention.
  • the synthetic aperture radar carried by the unmanned aerial vehicle transmits a chirp signal with a carrier frequency fc through the transmitting antenna, and the transmitted signal is received by the radar through the receiving antenna after being scattered by the target, and the received signal is
  • Synthetic aperture length L s R ⁇ BW , where R is the reference distance, and ⁇ BW is the azimuth beam width.
  • the segmented aperture imaging method of the multi-rotor unmanned aerial vehicle-borne synthetic aperture radar provided by the present invention is used to perform segmented aperture imaging according to the original echo signal of the unmanned aerial vehicle-borne synthetic aperture radar, comprising the following steps:
  • Step 1 the distance pulse compression is carried out to the original echo signal s (t, ⁇ ) to obtain the range pulse pressure signal s RC (t, ⁇ ), where t is the fast time of the distance, and ⁇ is the slow time of the azimuth, according to the distance pulse pressure signal s Phase history of strong scattering points in RC (t, ⁇ ) Computing the estimated velocity of the maneuvering platform and the estimated beam center oblique angle
  • Step 2 according to the estimated speed
  • Step 3 according to the estimated speed and the estimated beam center oblique angle Calculate the phase compensation amount of the segment pulse pressure signal s RC corresponding to each segment, i (t, ⁇ )
  • the segment pulse pressure signal s RC, i (t, ⁇ ) is multiplied by the motion error compensation filter Obtain N compensated signals as s MC, i (t, ⁇ ), imaginary number
  • Step 4 carry out two-dimensional Fourier transform to the compensated signal s MC, i (t, ⁇ ) to obtain N two-dimensional spectrum s MC, i (f, f d ), and use the series inversion method to analyze the two-dimensional spectrum s MC,i (f, f d ) is decomposed to construct the azimuth compression filter H AC,i , where f is the frequency corresponding to the fast time t of the range, and f d is the Doppler frequency corresponding to the slow time ⁇ of the azimuth ;
  • Step 5 multiply the two-dimensional frequency spectrum s MC, i (f, f d ) by the azimuth compression filter H AC, i , and then perform two-dimensional inverse Fourier transform to obtain N imaging results s IMG, i (t, ⁇ ) ;
  • Step 6 sequentially align and coherently accumulate the distance-dimension envelopes of the strong focal points in the overlapping regions of the imaging results s IMG, i (t, ⁇ ) corresponding to the adjacent segments, and splicing the non-overlapping regions to obtain the final imaging result S all .
  • step 1 the phase history of the strong scattering point in the distance pulse pressure signal s RC (t, ⁇ ) Perform second-order fitting to obtain the phase history of the strong scattering point Calculate the estimated velocity of the motorized platform from the second-order term coefficient ⁇ and the first-order term coefficient ⁇ and estimated beam center oblique angle Among them, t is the distance fast time, ⁇ is the azimuth slow time, o( ⁇ ) is the high-order phase error, is the constant phase term, ⁇ is the wavelength of the system transmitting signal, is the estimated value of the system reference distance.
  • step 3 the amount of phase compensation
  • R 0 is the The mean value of
  • ⁇ 0 is the oblique angle of N estimated beam centers mean value.
  • step 4 the azimuth compression filter
  • f is the frequency corresponding to the fast distance time t
  • f d is the Doppler frequency corresponding to the azimuth slow time ⁇
  • f c is the carrier frequency of the system transmitting signal
  • c is the speed of light.
  • step 6 geometrically correct the imaging results sIMG, i (t, ⁇ ) corresponding to the adjacent segments to obtain the corrected imaging results
  • the corrected imaging results are then to rotate degrees, to obtain corrected imaging results with the slant distance perpendicular to the trajectory of the maneuvering platform
  • the corrected imaging results corresponding to the adjacent segments are sequentially Alignment and coherent accumulation of distance-dimension envelopes where strong focal points in overlapping areas are located, and splicing non-overlapping areas to obtain the final imaging result S all .
  • the present invention also provides a multi-rotor unmanned aerial vehicle-borne synthetic aperture radar segmentation aperture imaging method, which is used to calculate the flight trajectory of the unmanned aerial vehicle according to the original echo signal of the unmanned aerial vehicle-borne synthetic aperture radar. It is characterized in that it has such feature, comprising the following steps: Step S1, performing distance pulse compression on the original echo signal s (t, ⁇ ) to obtain a distance pulse pressure signal s RC (t, ⁇ ), according to the distance pulse pressure signal s RC (t, Phase history of strong scattering points in ⁇ ) Computing the estimated velocity of the maneuvering platform and the estimated beam center oblique angle
  • Step S2 according to the estimated speed
  • ⁇ in represents the angle between the UAV-borne SAR beam and the ground normal
  • Step S3 in the adjacent areas in the i-th and i-1th segments, extract three strong scattering points, and the coordinates of the strong scattering points in the i-th segment are respectively The coordinates of the strong scattering points in the i-1th segment are respectively
  • Step S4 according to the coordinates of strong scattering points and Compute the rotation matrix ⁇ for the i-th and i-1-th segments,
  • Step S5 the platform trajectory coordinates of the i-th segment Rotate to the coordinate system of the platform trajectory of the i-1th segment, so that the platform trajectory coordinates of adjacent segments are aligned, that is:
  • Step S6 carry out coherent accumulation of platform track coordinates in the overlapping area of the i-th and i-1th segments, and splicing the platform track coordinates of the non-overlapping areas to obtain the spliced i-th and i-1th segments trajectory coordinates [P x , P y , P z ],
  • Step S7 repeat steps S3-S6 until the spliced trajectory coordinates [P x , P y , P z ] of all adjacent segments are obtained, so as to obtain the final trajectory coordinates of the platform
  • the mobile platform refers to a multi-rotor UAV
  • the airborne synthetic aperture radar is a Ku-band chirp continuous wave radar
  • the target refers to the ground area where synthetic aperture imaging needs to be performed.
  • the multi-rotor UAV-borne synthetic aperture radar segmented aperture imaging method includes the following steps:
  • Step 1 After the mobile platform takes off, the airborne synthetic aperture radar transmits a chirp signal with a carrier frequency fc through the transmitting antenna.
  • the transmitted signal is scattered by the target and then received by the radar through the receiving antenna.
  • the received signal is the original echo signal s (t, ⁇ ), where t is the fast time for distance, and ⁇ is the slow time for azimuth.
  • Perform range pulse compression on the original echo signal to obtain the pulse compression signal s RC (t, ⁇ ), according to the phase history of the strong scattering point in the pulse compression signal s RC (t, ⁇ )
  • Computing the estimated velocity of the maneuvering platform and the estimated beam center oblique angle Specifically follow the following sub-steps:
  • Step 1-1 according to the definition of the Doppler FM slope K a and the Doppler center f dc , the calculation formulas of K a and f dc are obtained:
  • Step 1-2 with reference to Fig. 2, the space formed by the original echo signal is called the signal space, and the original echo signal in the signal space is subjected to distance pulse compression to obtain the pulse compression signal s RC (t, ⁇ ), Among them, t is the fast distance time, the phase history of the strong scattering point in the pulse compression signal s RC (t, ⁇ ) Perform second-order fitting to obtain the expression of the strong scattering point in the signal space:
  • ⁇ is the second-order term coefficient of azimuth slow time
  • ⁇ is the first-order term coefficient of azimuth slow time
  • o( ⁇ ) is the high-order phase error, is a constant phase term
  • Steps 1-4 the space formed by the target and the actual position of the maneuvering platform is called the target space, and according to the geometric relationship in the target space and the cosine law, the expression of the distance R between the target and the platform is obtained:
  • v represents the movement speed of the maneuvering platform
  • ⁇ represents the oblique angle of the beam center caused by the movement of the maneuvering platform
  • R0 represents the shortest distance between the target and the maneuvering platform
  • phase history in the target space is obtained expression for:
  • ⁇ is the wavelength of the signal transmitted by the system
  • Steps 1-5 substituting the above formula ⁇ 8> into formula ⁇ 1> and formula ⁇ 2> respectively, to obtain the expressions of K a and f dc in the target space:
  • Steps 1-6 respectively compare the above formula ⁇ 4> with formula ⁇ 9>, formula ⁇ 5> and the right side of the equation of formula ⁇ 10>, and get the expressions of estimated velocity, estimated oblique angle of beam center, and estimated distance, respectively for:
  • Step 2 according to the estimated speed
  • Step 3 according to the estimated speed and the estimated beam center oblique angle Calculate the phase compensation amount of the segment compression signal s RC corresponding to each segment, i (t, ⁇ )
  • the segment compressed signal sRC, i (t, ⁇ ) is multiplied by the motion error compensation filter
  • N compensated signals as s MC, i (t, ⁇ ), imaginary number Realize motion compensation for each segment compressed signal s RC,i (t, ⁇ ).
  • Step 3-1 calculate the influence of the slope distance change ⁇ R on the phase caused by the change of the platform motion state in each section, the expression is:
  • v 0 is The mean value of , ⁇ 0 is the mean value of
  • Step 3-2 construct the corresponding segment compression signal s RC of each segment, the expression of the phase compensation amount of i (t, ⁇ ):
  • Step 3-3 construct motion error compensation filter expression:
  • j represents an imaginary number
  • Step 3-4 the segment compression signal s RC, i (t, ⁇ ) is multiplied by the motion error compensation filter shown in the above-mentioned formula ⁇ 14>, obtains the compensated signal s MC, i (t, ⁇ ), completes the Motion compensation of the compressed signal per segment.
  • Step 4 carry out two-dimensional Fourier transform to the compensated signal s MC, i (t, ⁇ ) to obtain N two-dimensional spectrum s MC, i (f, f d ), and use the series inversion method to analyze the two-dimensional spectrum s MC,i (f, f d ) is decomposed to construct the azimuth compression filter H AC,i , where f is the frequency corresponding to the fast time t of the range, and f d is the Doppler frequency corresponding to the slow time ⁇ of the azimuth .
  • H AC azimuth compression filter
  • Step 4-1 carry out two-dimensional Fourier transform to the compensated signal s MC, i (t, ⁇ ) to obtain the corresponding two-dimensional frequency spectrum s MC of each segment, i (f, f d ), wherein, f is the distance The frequency corresponding to the fast time t, f d is the Doppler frequency corresponding to the azimuth slow time ⁇ ;
  • Step 4-2-1 obtain the expression of each two-dimensional spectrum s MC, i (f, f d ) according to the stationary phase point method:
  • f c represents the carrier frequency of the signal transmitted by the system, and c represents the speed of light;
  • Step 4-2-2 using the series inversion method to decompose the above formula ⁇ 17>, and obtain the two-dimensional spectrum expression that eliminates the coupling term of f and f d :
  • Step 4-2-3 construct the ideal phase filter expression according to the above formula ⁇ 15>:
  • Step 4-2-4 will estimate the velocity Estimated Beam Center Slant Angle distance estimate Substituting into the above formula ⁇ 19>, the phase filter expression is obtained:
  • Step 5 multiply the two-dimensional frequency spectrum s MC, i (f, f d ) by the azimuth compression filter H AC, i , and then perform two-dimensional inverse Fourier transform to obtain N imaging results s IMG, i (t, ⁇ ) .
  • Step 6 sequentially align and coherently accumulate the distance-dimension envelope of the strong focal point in the overlapping area of the imaging results s IMG corresponding to the adjacent segments, i (t, ⁇ ), splicing the non-overlapping area,
  • the final imaging result S all is obtained. Specifically follow the following sub-steps:
  • Step 6-1 performing azimuth Fourier transform on the imaging result s IMG, i (t, ⁇ ), to obtain the range-Doppler domain map s IMG, i (t, f d );
  • Step 6-2 according to the Fourier transform characteristics and the geometric structure of the target space, construct the expression of the tilt correction filter for correcting the tilt of the image:
  • l is the distance scale of typical buildings
  • Step 6-3 multiply the range Doppler domain map s IMG,i (t, f d ) by the tilt correction filter H GC-1 to obtain the tilt corrected frequency domain map
  • Step 6-4 for the tilt-corrected frequency domain map Perform azimuth inverse Fourier transform to obtain the tilt-corrected time-domain map
  • Step 6-5 according to the geometric structure of the target space, get the stretch/compression factor expression:
  • Step 6-6 substituting the above formula ⁇ 22> into the tilt-corrected time-domain diagram Among them, the deformation-corrected time-domain map is obtained
  • Steps 6-7 correcting the deformed time domain map Perform distance-to-Fourier transform to obtain the deformed corrected frequency domain map
  • Steps 6-8 according to the Fourier transform characteristics and the target spatial geometry, construct the expression of the two-position correction filter for correcting image translation:
  • Steps 6-9 the deformed corrected frequency domain map Multiply by the position correction filter H GC-3 in the above formula ⁇ 23> to get the geometrically corrected frequency domain map
  • Steps 6-10 the geometrically corrected frequency domain map Carry out the distance inverse Fourier transform to obtain the geometrically corrected time domain diagram of each segment that has completed the geometric correction operation
  • Steps 6-11 the geometrically corrected time domain map Rotate counterclockwise degrees, to obtain the time-domain map to be spliced with the slant distance perpendicular to the trajectory of the maneuvering platform
  • Steps 6-12 sequentially perform adjacent time-domain images to be spliced Alignment of distance-dimension envelopes where strongly focused points in overlapping regions are located;
  • steps 6-13 perform coherent accumulation on overlapping regions, sequentially connect non-overlapping regions, complete sub-aperture splicing, and obtain a full-aperture imaging result S all .
  • the positioning part of the rotary-wing UAV-borne synthetic aperture radar segmented aperture imaging and positioning method includes the following steps:
  • Step S1 performing range pulse compression on the original echo signal s(t, ⁇ ) to obtain the pulse compressed signal s RC (t, ⁇ ), according to the phase history of the strong scattering point in the pulse compressed signal s RC (t, ⁇ ) Computing the estimated velocity of the maneuvering platform and the estimated beam center oblique angle
  • This step is exactly the same as step 1 of the imaging method part, and will not be repeated here.
  • Step S2 according to the estimated speed
  • M is the azimuth length of each section.
  • ⁇ in represents the angle between the UAV-borne SAR beam and the ground normal.
  • Step S3 in the adjacent areas in the i-th and i-1th segments, extract three strong scattering points, and the coordinates of the strong scattering points in the i-th segment are respectively The coordinates of the strong scattering points in the i-1th segment are respectively
  • Step S4 according to the coordinates of strong scattering points and Calculate the rotation matrix for the i-th and i-1-th segments
  • Step S5 the platform trajectory coordinates of the i-th segment Rotate to the coordinate system of the platform trajectory of the i-1th segment, so that the platform trajectory coordinates of adjacent segments are aligned
  • Step S6 coherently accumulate the platform trajectory coordinates in the overlapping area of the i-th and i-1th segments, splice the platform trajectory coordinates in the non-overlapping area, and obtain the spliced trajectory coordinates of the i-th and i-1th segments [P x , P y , P z ].
  • Step S7 repeat steps S3 to S6 N-1 times until the spliced trajectory coordinates [P x , P y , P z ] of all adjacent segments are obtained, so as to obtain the final trajectory coordinates of the platform
  • Verification 1 Under the conditions in Table 1, the original echo signal s r (t, ⁇ ) of the dot matrix target obtained by simulating the measured trajectory, the measured trajectory comes from the measured data collected by a certain test inertial navigation, and its trajectory is shown in Figure 5 (a ), (b) shown.
  • Figure 5(a) shows the change curve of the attitude angle of the maneuvering platform during flight and the movement of the maneuvering platform recorded by the inertial navigation system.
  • the three figures from top to bottom are the roll angle, pitch angle, and yaw angle relative to the movement of the platform.
  • the change curve of , the abscissa in the figure represents the flight path of the maneuvering platform from -200m to 200m.
  • Figure 5(b) shows the curves of the motion trajectory of the maneuvering platform recorded by inertial navigation in three dimensions, and the three figures from top to bottom are the change curves of the azimuth dimension, distance dimension and height dimension respectively.
  • the lattice targets are 7 ⁇ 7 point targets with distance and azimuth intervals of 25m respectively.
  • the spatial position relationship between the lattice targets and the platform trajectory is shown in Figure 6.
  • the coordinate systems X, Y, and Z shown in Figure 6 correspond to the dimensions of azimuth, distance, and height respectively.
  • R ref represents the reference distance of the system at the current incident angle ⁇ in
  • H represents the average flight height of the maneuvering platform
  • ⁇ a represents The beam width of the azimuth dimension of the system
  • ⁇ bw represents the beam width of the system distance dimension
  • S represents the width of the area to be measured
  • v represents the movement speed of the maneuvering platform
  • L s represents the length of a synthetic aperture.
  • the dot matrix on the right side of the figure represents the positional relationship of the simulation target, and the center of the dot matrix is located at the intersection of R ref and the Y axis on the left side of the picture.
  • the length and width of the lattice are uniformly arranged with 7 points, and the point spacing is 25m. Refer to Table 1 for the values of other parameters.
  • the s RD (t, ⁇ ) is obtained after s r (t, ⁇ ) distance pulse compression, and the SAI algorithm and the traditional imaging method are respectively used for imaging comparison of s RD (t, ⁇ ), and the results are shown in Fig. 7 and Fig. 8 .
  • the path is divided into three sections (section 1, section 2, and section 3) by the SAI algorithm proposed in this paper, and (a), (b), and (c) in Fig. 7 are section 1 and section 2 in this embodiment, respectively. 1.
  • the imaging result after geometric correction of segment 3 since each segment can only image part of the area, the imaging result of each segment only includes a part of the target area.
  • FIG. 7 is the full-aperture imaging result after splicing each segment in this embodiment
  • FIG. 7 is the imaging result of the traditional imaging method. Comparing Figure 7(d) and (e), it can be seen that the focus points have different outlines and different degrees of lightness and darkness, among which the contrast of several points framed by boxes is the most obvious.
  • FIG. 8 is a comparison result of the point spread function between this embodiment and the traditional imaging method. The upper figure in Fig.
  • the final full-aperture image of SAI greatly reduces the peak amplitude loss due to the complementary relationship of each segment; after measuring the point spread function curve index of a certain line of imaging results, the amount of information contained in the complete image is judged by image entropy.
  • the SAI method proposed in this paper has higher resolution, larger signal-to-noise ratio, and peak The amplitude loss is equivalent to the original method, and it contains more information. Therefore, the SAI proposed in this paper has a good imaging effect and solves the practical problems encountered by UAV-borne synthetic aperture radar.
  • Verification 2 Change the platform speed to 10m/s and the platform flight height to 350m in the conditions in Table 1, which are the test conditions for the actual flight test.
  • the test platform is a KWT-65 multi-rotor UAV, carrying a Ku-band small Frequency-modulated continuous wave synthetic aperture radar, using this platform to fly about 100 sorties, each flight distance is about 1km, and collected multiple sets of echo data in a certain area on the ground.
  • FIG. 9 is the imaging effect of the traditional imaging method
  • FIG. 9 (b) is the imaging effect of the method of this embodiment
  • FIG. 9 (c) is the optical model corresponding to the imaging area.
  • FIG. 9 shows that the focusing effect of the method of this embodiment is better than that of the traditional imaging method.
  • the target presents an "I" shape.
  • the imaging result of the method of the present invention has a clearer outline of the word "I” than the traditional method. Taking the left side of the "I"-shaped building in the image as an example, when the radar illuminates the glass, it will penetrate the glass and enter the interior of the room, and the echo cannot be Return to the radar, so the imaging result will appear black, and the edge of the window will have obvious reflections due to the existence of right angles and metals, forming bright spots.
  • the floor interval can be clearly observed through the black areas between the bright spots.
  • the bright spots are all blurred into a ball, and the glass between the floors cannot be clearly observed.
  • the imaging effect of the method of the invention is also clearer for the outline of the garden around the "I" shaped building.
  • the overall signal-to-noise ratio of the image is also better than traditional methods.
  • this method has a better focusing effect on 90% of the echo data imaging, while the traditional method has a better focusing effect on about 35% of the images.
  • the imaging time of this method is about 320s, and the imaging time of the traditional imaging method is about 5200s, and the speed is increased by about 16 times.
  • Verification 3 Process the test data of a certain test under the conditions of verification 2, and get the image while estimating the trajectory of the platform.
  • the comparison between the trajectory estimated by the method of the present invention and the trajectory collected by the inertial navigation equipment is shown in Fig. 10 .
  • (a) of Figure 10 is the comparison effect of the trajectory estimated by the method of the present invention and the trajectory collected by the inertial navigation equipment from top to bottom in the azimuth dimension, the distance dimension, and the dimension respectively, wherein the blue curve represents the measurement result of the inertial navigation equipment, The red curve is the estimated result of the method of the present invention. It can be seen that the two basically coincide in the azimuth.
  • the high sampling rate of the multi-rotor platform ensures a small error in the azimuth.
  • the method mainly includes: 1) acquiring target echoes according to the UAV-borne SAR system; 2) estimating the maneuvering platform according to the echo signals 3) Segment the echo signal according to the motion state of the platform; 4) Compensate the echo signal of each segment according to the motion state of the platform; 5) Perform two-dimensional Fourier transform on the echo signal after compensation of each segment Obtain the two-dimensional spectrum, and use the series inversion method to decompose the two-dimensional spectrum to obtain each section of the phase filter; 6) multiply each section of the two-dimensional spectrum by the corresponding phase filter of the section, and then perform two-dimensional inverse Fourier on the two-dimensional spectrum Leaf transformation to obtain each section of image; 7) Perform geometric correction on each section of image, and then splicing each section of image to obtain full aperture imaging result.
  • the method of estimating the platform motion parameters based on the phase history of the echo signal and segmenting the echo signal is used to accurately compensate the motion error of the platform, so as to achieve good imaging focusing effect and high imaging success rate, and improve the performance of multi-rotor UAVs.
  • the efficiency of platform data acquisition can effectively perform high-resolution imaging on the multi-rotor UAV platform synthetic aperture radar system.
  • this method can calculate the three-dimensional coordinates of the platform trajectory while imaging to achieve the effect of platform positioning, which can be applied to UAV navigation and provide the possibility for future intelligent and integrated detection systems.
  • Synthetic aperture radar imaging system Compared with the traditional airborne synthetic aperture radar imaging algorithm, it takes into account the influence of squint caused by the change of attitude angle of the UAV platform, compensates the coupling phase of distance and azimuth, improves the image focusing effect, and is suitable for applications without antenna servo. Synthetic aperture radar imaging system;
  • phase filter multiplication method is used instead of the interpolation method, which improves the imaging speed of a single imaging in each segment;
  • the image is obtained by using parallel imaging of each segment and then splicing into a complete image, which improves the imaging speed of the complete image.
  • the experimental verification shows that the Segmented Aperture Imaging (SAI) method proposed in this embodiment has a high imaging resolution and a fast calculation speed. At the same time, it can estimate the platform trajectory and reduce the The requirements of the algorithm on the hardware equipment are clarified, and it is shown that the invention can be effectively applied to the imaging system of the small maneuverable platform synthetic aperture radar.
  • SAI Segmented Aperture Imaging
  • This embodiment deduces the corresponding relationship between the platform motion parameters and the echo phase in detail, and for the second-order expansion of the two-dimensional frequency spectrum under squint, motion compensation and imaging can be performed without relying on the inertial navigation equipment.
  • the point spread function curves of each section of the method of this embodiment, the complete image and the traditional imaging algorithm were compared; the actual measurement verification compared the imaging result graphs of the method of the present invention and the traditional algorithm, proving that this embodiment can effectively detect the multi-rotor unmanned
  • the aircraft platform synthetic aperture radar system performs high-resolution imaging.

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Physics & Mathematics (AREA)
  • Electromagnetism (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • General Physics & Mathematics (AREA)
  • Signal Processing (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法，该方法包括：1）根据无人机载合成孔径雷达系统获取目标回波；2）根据回波信号估计机动平台的运动状态；3）根据平台运动状态对回波信号分段；4）根据平台运动状态对各段回波信号运动补偿；5）对各段补偿后的回波信号进行二维傅里叶变换得到二维频谱，采用级数反演法分解二维频谱得到各段相位滤波器；6）对各段二维频谱乘该段对应的相位滤波器，然后对二维频谱进行二维逆傅里叶变换得到各段图像；7）对各段图像进行几何校正，然后拼接各段图像，得到全孔径成像结果。8）对各段平台轨迹进行拼接，得到平台完整轨迹坐标。

Description

一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法 技术领域
本发明属于信号处理技术领域，具体涉及一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法。
背景技术
随着遥感技术广泛应用于地图测绘、地物表面形变检测、天气观测等民事应用，遥感设备的小型化已经成为必然趋势。自上世纪90年代以来，无人机平台由于其小型便携、灵活性高、操作简单等特点逐渐成为合成孔径雷达的重要平台，众多学者从系统设计和成像算法两个方面对无人机平台合成孔径雷达系统进行了大量研究。近年来随着无人机导航与控制技术的发展，各种型号的多旋翼无人机大量进入市场，但是还未被广泛使用于遥感领域，主要是因为传统机载成像算法对多旋翼无人机载合成孔径雷达的成像效果不稳定，受制于其携带的惯导设备性能，因此研究不依赖于惯导设备的多旋翼无人机载合成孔径雷达成像算法具有十分重要的意义。
目前国内外的研究机构采用的机载成像算法主要依赖于高精度惯导系统获取的数据进行运动补偿并成像，并且因为大型平台所携带的伺服系统可以实时补偿平台姿态角变化带来的影响，而无需考虑平台角度变化，因此无法解决多旋翼无人机载合成孔径雷达系统存在的三个问题：1)多旋翼平台的轨迹通常非常不稳定，因为它容易受到诸如风等环境因素的影响，可能会导致严重的抖动和急转弯；2)多旋翼独有的飞行原理会带来由于旋翼振动和机身倾斜操作引起的高频误差和入射角误差；3)高机动轨迹与低成本惯导系统相结合，导致用于运动补偿的运动状态和姿态角数据相当差。并且因为无人机载重的限制无法携带高精度惯导设备，因而无人机导航系统的定位精度也会有所下降，若飞机在飞行过程之中根据成像结果反演出飞机飞行轨迹则可为其导航系统提供依据，进行更智能化的航迹规划。因此研究不依赖于惯导设备的多旋翼无人机载合成孔径雷达成像及定位算法具有十分重要的意义。
发明内容
本发明是针对上述现有机载成像算法的缺点和不足而进行的，目的在于提供一种多旋翼无人机载合成孔径雷达分段孔径成像(SAI)及定位方法，适用于无惯导设备或低精度惯导设备的合成孔径雷达系统，提高雷达成像的聚焦效果与成像效率。
本发明提供了一种多旋翼无人机载合成孔径雷达分段孔径成像方法，用于根据合成孔径雷达的原始回波信号进行分段孔径成像，具有这样的特征，包括以下步骤：步骤1，对原始回波信号s(t，η)进行距离脉冲压缩得到距离脉压信号s RC(t，η)，其中t为距离快时间、η为方位慢时间，根据距离脉压信号s RC(t，η)中强散射点的相位历程
Figure PCTCN2021115882-appb-000001
计算机动平台的估计速度
Figure PCTCN2021115882-appb-000002
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000003
步骤2，根据估计速度
Figure PCTCN2021115882-appb-000004
的方向将距离脉压信号s RC(t，η)分段得到N个段及每个段对应的段脉压信号s RC，i(t，η)，i＝1...N；步骤3，根据估计速度
Figure PCTCN2021115882-appb-000005
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000006
计算每个段对应的段脉压信号s RC，i(t，η)的相位补偿量
Figure PCTCN2021115882-appb-000007
对段脉压信号s RC，i(t，η)乘以运动误差补偿滤波器
Figure PCTCN2021115882-appb-000008
得到N个已补偿信号为s MC，i(t，η)，虚数j＝
Figure PCTCN2021115882-appb-000009
步骤4，对已补偿信号s MC，i(t，η)进行二维傅里叶变换得到N个二维频谱s MC，i(f，f d)，采用级数反演法对二维频谱s MC，i(f，f d)进行分解，构建方位压缩滤波器H AC，i，其中，f为与距离快时间t对应的频率，f d为与方位慢时间η对应的多普勒频率；步骤5，对二维频谱s MC，i(f，f d)乘方位压缩滤波器H AC，i，然后进行二维逆傅里叶变换得到N个成像结果s IMG，i(t，η)；步骤6，依次对相邻的段对应的成像结果s IMG，i(t，η)重叠区域中的强聚焦点所在距离维包络对齐与相参积累，对不重叠区域拼接，得到最终成像结果S all。
在本发明提供的多旋翼无人机载合成孔径雷达分段孔径成像方法中，还可以具有这样的特征：步骤1中，对距离脉压信号s RC(t，η)中强散射点的相位历程
Figure PCTCN2021115882-appb-000010
进行二阶拟合得到该强散射点的相位历程
Figure PCTCN2021115882-appb-000011
根据二阶项系数β和一阶项系数α计算机动平台的估计速度
Figure PCTCN2021115882-appb-000012
以及估计波束中心斜视角
Figure PCTCN2021115882-appb-000013
其中，t为距离快时间，η为方位慢时间，o(η)为高阶相位误差，
Figure PCTCN2021115882-appb-000014
为常数相位项，λ为系统发射信号波长，
Figure PCTCN2021115882-appb-000015
为系统参考距离的估计值。
在本发明提供的多旋翼无人机载合成孔径雷达分段孔径成像方法中，还可以具有这样的特征：其中，步骤2中，根据估计速度
Figure PCTCN2021115882-appb-000016
的方向将距离脉压信号s RC(t，η)依次划分为速度方向一致的N个段，然后判断每个段的长度是否小于一个合成孔径长度，若小于，则将该段向两边扩展为一个合成孔径长度，最终得到分段后的N个段脉压信号s RC，i(t，η)，i＝1...N。
在本发明提供的多旋翼无人机载合成孔径雷达分段孔径成像方法中，还可以具有这样的特征：其中，步骤3中，相位补偿量
Figure PCTCN2021115882-appb-000017
式中，R 0为各段
Figure PCTCN2021115882-appb-000018
的均值，θ 0为N个估计波束中心斜视角
Figure PCTCN2021115882-appb-000019
的均值。
在本发明提供的多旋翼无人机载合成孔径雷达分段孔径成像方法中，还可以具有这样的特征：其中，步骤4中，方位压缩滤波器
Figure PCTCN2021115882-appb-000020
式中，f为与距离快时间t对应的频率，f d为与方位慢时间η对应的多普勒频率，f c为系统发射信号的载频，c为光速。
在本发明提供的多旋翼无人机载合成孔径雷达分段孔径成像方法中，还可以具有这样的特征：其中，步骤6中，先对相邻所述段对应的成像结果s IMG，i(t，η)进行几何校正得到已校正成像结果
Figure PCTCN2021115882-appb-000021
然后将已校正成像结果
Figure PCTCN2021115882-appb-000022
旋转
Figure PCTCN2021115882-appb-000023
度，得到斜距垂直于机动平台轨迹的已校正成像结果
Figure PCTCN2021115882-appb-000024
然后依次对相邻的段对应的已校正成像结果
Figure PCTCN2021115882-appb-000025
重叠区域中的强聚焦点所在距离维包络对齐与相参积累，对不重叠区域拼接，得到最终成像结果S all。
在本发明提供的多旋翼无人机载合成孔径雷达分段孔径成像方法中，还可以具有这样的特征：其中，步骤6对成像结果s IMG，i进行几何校正按如下子步骤进行：步骤6-1，对成像结果s IMG，i(t，η)进行方位向傅里叶变换，得到距离多普勒域图s IMG，i(t，f d)；步骤6-2，根据傅里叶变换特性和目标空间几何结构，构建校正图像倾斜的倾斜校正滤波器的表达式：
Figure PCTCN2021115882-appb-000026
Figure PCTCN2021115882-appb-000027
式中，l为典型建筑物的距离向尺度；步骤6-3，将距离多普勒域图s IMG，i(t，f d)乘以倾斜校正滤波器H GC-1，得到已倾斜校正频域图
Figure PCTCN2021115882-appb-000028
步骤6-4，对已倾斜校正频域图
Figure PCTCN2021115882-appb-000029
进行方位向逆傅里叶变换，得到已倾斜校正时域图
Figure PCTCN2021115882-appb-000030
步骤6-5，根据目标空间几何结构，得到拉伸/压缩因子表达式：
Figure PCTCN2021115882-appb-000031
步骤6-6，将拉伸/压缩因子表达式代入已倾斜校正时域图
Figure PCTCN2021115882-appb-000032
之中，得到已形变校正时域图
Figure PCTCN2021115882-appb-000033
步骤6-7，对已形变校正时域图
Figure PCTCN2021115882-appb-000034
进行距离向傅里叶变换，得到已形变校正频域图
Figure PCTCN2021115882-appb-000035
步骤6-8，根据傅里叶变换特性和目标空间几何结构，构建校正图像平移的二位置校正滤波器的表达式：
Figure PCTCN2021115882-appb-000036
步骤6-9，将已形变 校正频域图
Figure PCTCN2021115882-appb-000037
乘以位置校正滤波器H GC-3，得到已几何校正频域图
Figure PCTCN2021115882-appb-000038
步骤6-10，对已几何校正频域图
Figure PCTCN2021115882-appb-000039
进行距离向逆傅里叶变换，得到各段完成几何校正操作的已校正成像结果
Figure PCTCN2021115882-appb-000040
在本发明提供的多旋翼无人机载合成孔径雷达分段孔径成像方法中，还可以具有这样的特征：其中，步骤1中，估计速度
Figure PCTCN2021115882-appb-000041
及估计波束中心斜视角
Figure PCTCN2021115882-appb-000042
按如下步骤计算：步骤1-1，根据多普勒调频斜率K a和多普勒中心f dc的定义，得到K a与f dc的计算式：
Figure PCTCN2021115882-appb-000043
Figure PCTCN2021115882-appb-000044
式中，
Figure PCTCN2021115882-appb-000045
表示(·)对方位慢时间η求导；步骤1-2，将原始回波信号所构成的空间称为信号空间，对距离脉压信号s RC(t，η)中强散射点的相位历程
Figure PCTCN2021115882-appb-000046
进行二阶拟合得到该强散射点在信号空间中的表达式：
Figure PCTCN2021115882-appb-000047
式中，β为二阶项系数，α为一阶项系数，
Figure PCTCN2021115882-appb-000048
为常数相位项，o(η)为高阶相位误差；步骤1-3，将相位历程
Figure PCTCN2021115882-appb-000049
的表达式分别代入步骤1-1中K a的计算式与f dc的计算式之中，得到信号空间中K a与f dc的表达式：
Figure PCTCN2021115882-appb-000050
Figure PCTCN2021115882-appb-000051
步骤1-4，将目标和机动平台实际位置所构成的空间称为目标空间，根据空间几何关系，得到目标空间中的相位历程
Figure PCTCN2021115882-appb-000052
的表达式：
Figure PCTCN2021115882-appb-000053
式中，v表示机动平台运动速度，θ表示机动平台运动引起的波束中心斜视角，R 0表示目标和机动平台的最近距离，λ为系统发射信号波长；步骤1-5，将目标空间中的相位历程
Figure PCTCN2021115882-appb-000054
的表达式分别代入步骤1-1中K a的计算式与f dc的计算式之中，得到目标空间中K a与f dc的表达式：
Figure PCTCN2021115882-appb-000055
步骤1-6，分别对比信号空间中K a的表达式与目标空间中K a表达式，信号空间中f dc的表达式与目标空间中f dc表达式的右侧，得到估计速度、估计波束中心斜视角、距离估计值：
Figure PCTCN2021115882-appb-000056
Figure PCTCN2021115882-appb-000057
在本发明提供的多旋翼无人机载合成孔径雷达分段孔径成像方法中，还可以具有这样的特征：其中，步骤4中，采用级数反演法对二维频谱s MC，i(f，f d)分解得到方位压缩滤波器表达式H AC，i的过程按如下子步骤进行：步骤4-2-1，根据驻相点法求得段对应的二维频谱s MC，i(f，f d)的表达式：
Figure PCTCN2021115882-appb-000058
式中，f c表示系统发射信号的载频，c表示光速；步骤4-2-2，采用级数反演法对步骤4-2-1中的二维频谱s MC，i(f，f d)的表达式进行分解，得到消除f和f d的耦合项的二维频谱表达式：
Figure PCTCN2021115882-appb-000059
Figure PCTCN2021115882-appb-000060
步骤4-2-3，根据步骤4-2-2中的消除f和f d的耦合项的二维频谱表达式构建理想相位滤波器表达式：
Figure PCTCN2021115882-appb-000061
Figure PCTCN2021115882-appb-000062
步骤4-2-4，将估计速度
Figure PCTCN2021115882-appb-000063
估计波束中心斜视角
Figure PCTCN2021115882-appb-000064
距离估计值
Figure PCTCN2021115882-appb-000065
代入理想相位滤波器表达式之中，得到方位压缩滤波器表达式：
Figure PCTCN2021115882-appb-000066
Figure PCTCN2021115882-appb-000067
本发明还提供了一种多旋翼无人机载合成孔径雷达分段孔径成像方法，用于根据无人机载合成孔径雷达的原始回波信号，计算无人机飞行轨迹其特征在于，具有这样的特征，包括以下步骤：步骤S1，对原始回波信号s(t，η)进行距离脉冲压缩得到距离脉压信号s RC(t，η)，根据所述距离脉压信号s RC(t，η)中强散射点的相位历程
Figure PCTCN2021115882-appb-000068
计算机动平台的估计速度
Figure PCTCN2021115882-appb-000069
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000070
步骤S2，根据估计速度
Figure PCTCN2021115882-appb-000071
的方向将距离脉压信号s RC(t，η)分段得到N个段及每个段对应的段脉压信号s RC，i(t，η)，i＝1...N，并且根据估计速度
Figure PCTCN2021115882-appb-000072
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000073
计算第i个段对应的平台轨迹坐标
Figure PCTCN2021115882-appb-000074
为各段方位向长度，
Figure PCTCN2021115882-appb-000075
Figure PCTCN2021115882-appb-000076
式中，θ in表示无人机载合成孔径雷达波束与地面 法向的夹角；步骤S3，在第i个和第i-1个段中的相邻区域，提取三个强散射点，第i个段中强散射点的坐标分别为
Figure PCTCN2021115882-appb-000077
第i-1个段中强散射点的坐标分别为
Figure PCTCN2021115882-appb-000078
步骤S4，根据强散射点的坐标
Figure PCTCN2021115882-appb-000079
和
Figure PCTCN2021115882-appb-000080
计算第i个和第i-1个段的旋转矩阵γ，
Figure PCTCN2021115882-appb-000081
步骤S5，将第i个所述段的所述平台轨迹坐标
Figure PCTCN2021115882-appb-000082
旋转至以第i-1个所述段的所述平台轨迹的坐标系下，使相邻段的平台轨迹坐标对齐，即：
Figure PCTCN2021115882-appb-000083
步骤S6，对第i个和第i-1个段重叠区域的平台轨迹坐标进行相参积累，对不重叠区域的平台轨迹坐标进行拼接，得到第i个和第i-1个段的拼接后轨迹坐标
Figure PCTCN2021115882-appb-000084
Figure PCTCN2021115882-appb-000085
步骤S7，重复步骤S3-S6，直至得到所有相邻段的拼接后轨迹坐标[P x，P y，P z]，从而得到平台的最终轨迹坐标
Figure PCTCN2021115882-appb-000086
发明的作用与效果
根据本发明所涉及的多旋翼无人机载合成孔径雷达分段孔径成像及定位方法，采用基于回波信号相位历程估计平台运动参数和对回波信号分段的方法对平台运动误差精确补偿，达到成像聚焦效果好和成像成功率高的效果，提高了多旋翼无人机平台数据采集的效率，可以有效地对多旋翼无人机平台合成孔径雷达系统进行高分辨成像，并且在成像时可以对平台轨迹进行有效估计，得到平台相对于测量区域的坐标信息，达到对平台定位的效果。
本发明与现有技术相比具有如下优点：
与传统的机载合成孔径雷达成像算法相比，本发明考虑到了平台运动与信号相位之间的关系，可适用于无惯导设备或者低精度惯导设备的合成孔径雷达成像系统进行成像；
与传统的机载合成孔径雷达成像算法相比，本发明考虑到了无人机平台姿态角变化带来的斜视影响，补偿了距离和方位的耦合相位，提高了图像聚焦效果，可适用于无天线伺服的合成孔径雷达成像系统；
与传统的机载合成孔径雷达成像算法相比，本发明采用相位滤波器相乘的方法，替代插值方法，提高了各段内单次成像的成像速度；
与传统的机载合成孔径雷达成像算法相比，本发明采用各段并行成像然后拼接为完整图像的办法获取图像，提高了完整图像的成像速度。
与传统的机载合成孔径雷达成像算法相比，本发明可以在成像同时估计平台的位置，达到了对平台定位的效果，便于建立自动化、智能化、一体化的无人机探测系统。
附图说明
图1本发明的实施例中的流程图；
图2本发明的实施例中的机载合成孔径雷达系统目标空间和信号空间示意图；
图3本发明的实施例中的几何校正流程图；
图4本发明的实施例中的轨迹估计示意图；
图5本发明的实施例中的实测轨迹与姿态角变化曲线图；
图6本发明的实施例中的平台轨迹与点阵目标空间位置关系示意图；
图7本发明的实施例中的各段成像结果、全孔径成像结果与传统成像方法成像结果对比图；
图8本发明的实施例中的本发明与传统成像方法点扩散函数对比结果；
图9本发明的实施例中的传统成像方法成像结果、本方法成像结果、光学模型对比图；
图10本发明的实施例中的对平台估计轨迹曲线与惯导测量轨迹曲线对比及误差曲线图。
具体实施方式
为了使本发明实现的技术手段、创作特征、达成目的与功效易于明白了解，以下结合实施例及附图对本发明一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法作具体阐述。
在本发明中，无人机起飞后，无人机载合成孔径雷达通过发射天线发射载频为f c的线性调频信号，发射信号经目标散射后由雷达通过接收天线接收，接收到的信号为原始回波信号s(t，η)，其中t为距离快时间、η为方位慢时间。合成孔径长度L s＝R×θ BW，其中，R为参考距离，θ BW为方位波束宽度。
本发明提供的多旋翼无人机载合成孔径雷达分段孔径成像方法，用于根据无人机载合成孔径雷达的原始回波信号进行分段孔径成像，包括以下步骤：
步骤1，对原始回波信号s(t，η)进行距离脉冲压缩得到距离脉压信号s RC(t，η)，其中t为距离快时间、η为方位慢时间，根据距离脉压信号s RC(t，η)中强散射点的相位历程
Figure PCTCN2021115882-appb-000087
计算机动平台的估计速度
Figure PCTCN2021115882-appb-000088
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000089
步骤2，根据估计速度
Figure PCTCN2021115882-appb-000090
的方向将距离脉压信号s RC(t，η)分段得到N个段及每个段对应的段脉压信号s RC，i(t，η)，i＝1...N；
步骤3，根据估计速度
Figure PCTCN2021115882-appb-000091
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000092
计算每个段对应的段脉压信号s RC，i(t，η)的相位补偿量
Figure PCTCN2021115882-appb-000093
对段脉压信号s RC，i(t，η)乘以运动误差补偿滤波器
Figure PCTCN2021115882-appb-000094
得到N个已补偿信号为s MC，i(t，η)，虚数
Figure PCTCN2021115882-appb-000095
步骤4，对已补偿信号s MC，i(t，η)进行二维傅里叶变换得到N个二维频谱s MC，i(f，f d)，采用级数反演法对二维频谱s MC，i(f，f d)进行分解，构建方位压缩滤波器H AC，i，其中，f为与距离快时间t对应的频率，f d为与方位慢时间η对应的多普勒频率；
步骤5，对二维频谱s MC，i(f，f d)乘方位压缩滤波器H AC，i，然后进行二维逆傅里叶变换得到N个成像结果s IMG，i(t，η)；
步骤6，依次对相邻的段对应的成像结果s IMG，i(t，η)重叠区域中的强聚焦点所在距离维包络对齐与相参积累，对不重叠区域拼接，得到最终成像结果S all。
其中，步骤1中，对距离脉压信号s RC(t，η)中强散射点的相位历程
Figure PCTCN2021115882-appb-000096
进行二阶拟合得到该强散射点的相位历程
Figure PCTCN2021115882-appb-000097
根据二阶项系数β和一阶项系数α计算机动平台的估计速度
Figure PCTCN2021115882-appb-000098
以及估计波束中心斜视角
Figure PCTCN2021115882-appb-000099
其中，t为距离快时间，η为方位慢时间，o(η)为高阶相位误差，
Figure PCTCN2021115882-appb-000100
为常数相位项，λ为系统发射信号波长，
Figure PCTCN2021115882-appb-000101
为系统参考距离的估计值。
步骤2中，根据估计速度
Figure PCTCN2021115882-appb-000102
的方向将距离脉压信号s RC(t，η)依次划分为速度方向一致的N个段，然后判断每个段的长度是否小于一个合成孔径长度，若小于，则将该段向两边扩展为一个合成孔径长度，最终得到分段后的N个段脉压信号s RC，i(t，η)，i＝1...N。
步骤3中，相位补偿量
Figure PCTCN2021115882-appb-000103
式中，R 0为各段
Figure PCTCN2021115882-appb-000104
的均值，θ 0为N个估计波束中心斜视角
Figure PCTCN2021115882-appb-000105
的均值。
步骤4中，方位压缩滤波器
Figure PCTCN2021115882-appb-000106
式中，f为与距离快时间t对应的频率，f d为与方位慢时间η对应的多普勒频率，f c为系统发射信号的载频，c为光速。
步骤6中，先对相邻所述段对应的成像结果s IMG，i(t，η)进行几何校正得到已校正成像结果
Figure PCTCN2021115882-appb-000107
然后将已校正成像结果
Figure PCTCN2021115882-appb-000108
旋转
Figure PCTCN2021115882-appb-000109
度，得到斜距垂直于机动平台轨迹的已校正成像结果
Figure PCTCN2021115882-appb-000110
然后依次对相邻的段对应的已校正成像结果
Figure PCTCN2021115882-appb-000111
重叠区域中的强聚焦点所在距离维包络对齐与相参积累，对不重叠区域拼接，得到最终成像结果S all。
本发明还提供了一种多旋翼无人机载合成孔径雷达分段孔径成像方法，用于根据无人机载合成孔径雷达的原始回波信号，计算无人机飞行轨迹其特征在于，具有这样的特征，包括以下步骤：步骤S1，对原始回波信号s(t，η)进行距离脉冲压缩得到距离脉压信号s RC(t，η)，根据所述距离脉压信号s RC(t，η)中强散射点的相位历程
Figure PCTCN2021115882-appb-000112
计算机动平台的估计速度
Figure PCTCN2021115882-appb-000113
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000114
步骤S2，根据估计速度
Figure PCTCN2021115882-appb-000115
的方向将距离脉压信号s RC(t，η)分段得到N个段及每个段对应的段脉压信号s RC，i(t，η)，i＝1...N，并且根据估计速度
Figure PCTCN2021115882-appb-000116
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000117
计算第i个段对应的平台轨迹坐标
Figure PCTCN2021115882-appb-000118
为各段方位向长度，
Figure PCTCN2021115882-appb-000119
Figure PCTCN2021115882-appb-000120
Figure PCTCN2021115882-appb-000121
式中，θ in表示无人机载合成孔径雷达波束与地面法向的夹角；
步骤S3，在第i个和第i-1个段中的相邻区域，提取三个强散射点，第i个段中强散射点的坐标分别为
Figure PCTCN2021115882-appb-000122
第i-1个段中强散射点的坐标分别为
Figure PCTCN2021115882-appb-000123
步骤S4，根据强散射点的坐标
Figure PCTCN2021115882-appb-000124
和
Figure PCTCN2021115882-appb-000125
计算第i个和第i-1个段的旋转矩阵γ，
Figure PCTCN2021115882-appb-000126
步骤S5，将第i个所述段的所述平台轨迹坐标
Figure PCTCN2021115882-appb-000127
旋转至以第i-1个所述段的所述平台轨迹的坐标系下，使相邻段的平台轨迹坐标对齐，即：
Figure PCTCN2021115882-appb-000128
步骤S6，对第i个和第i-1个段重叠区域的平台轨迹坐标进行相参积累，对不重叠区域的平台轨迹坐标进行拼接，得到第i个和第i-1个段的拼接后轨迹坐标[P x，P y，P z]，
Figure PCTCN2021115882-appb-000129
步骤S7，重复步骤S3-S6，直至得到所有相邻段的拼接后轨迹坐标[P x，P y，P z]，从而得到平台的最终轨迹坐标
Figure PCTCN2021115882-appb-000130
<实施例>
在本实施例中，机动平台指多旋翼无人机，机载合成孔径雷达为Ku波段的线性调频连续波雷达，目标指需要进行合成孔径成像的地面区域。
如图1所示，多旋翼无人机载合成孔径雷达分段孔径成像方法包括如下步骤：
步骤1，机动平台起飞后，机载合成孔径雷达通过发射天线发射载频为f c的线性调频信号，发射信号经目标散射后由雷达通过接收天线接收，接收到的信号为原始回波信号s(t，η)，其中t为距离快时间、η为方位慢时间。对原始回波信号进行距离脉冲压缩得到脉冲压缩信号s RC(t，η)，根据脉冲压缩信号s RC(t，η)中强散射点的相位历程
Figure PCTCN2021115882-appb-000131
计算机动平台的估计速度
Figure PCTCN2021115882-appb-000132
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000133
具体按照如下子步骤进行：
步骤1-1，根据多普勒调频斜率K a和多普勒中心f dc的定义，得到K a与f dc的计算式：
Figure PCTCN2021115882-appb-000134
Figure PCTCN2021115882-appb-000135
式中，
Figure PCTCN2021115882-appb-000136
表示(·)对方位慢时间η求导，η为方位慢时间，
Figure PCTCN2021115882-appb-000137
为某个强散射点的相位历程；
步骤1-2，参照图2，将原始回波信号所构成的空间称为信号空间，对信号空间之中的原始回波信号进行距离脉冲压缩，得到脉冲压缩信号s RC(t，η)，其中，t为距离快时间，对脉冲压缩信号s RC(t，η)中强散射点的相位历程
Figure PCTCN2021115882-appb-000138
进行二阶拟合得到该强散射点在信号空间中的表达式：
Figure PCTCN2021115882-appb-000139
式中，β为方位慢时间二阶项系数，α为方位慢时间的一阶项系数，o(η)为高阶相位误差，
Figure PCTCN2021115882-appb-000140
为常数相位项；
步骤1-3，将上述式<3>分别代入式<1>与式<2>之中，得到信号空间中K a与f dc的表达式：
Figure PCTCN2021115882-appb-000141
Figure PCTCN2021115882-appb-000142
步骤1-4，将目标和机动平台实际位置所构成的空间称为目标空间，根据目标空间之中的几何关系与余弦定理，得到目标和平台之间距离R的表达式：
Figure PCTCN2021115882-appb-000143
式中，v表示机动平台运动速度，θ表示机动平台运动引起的波束中心斜视角，R 0表示目标和机动平台的最近距离；
将上述式<6>在η＝0时泰勒展开，保留至二阶项，得到R的表达式：
Figure PCTCN2021115882-appb-000144
根据相位计算公式及定义，得到目标空间中相位历程
Figure PCTCN2021115882-appb-000145
的表达式：
Figure PCTCN2021115882-appb-000146
式中，λ为系统发射信号波长；
步骤1-5，将上述式<8>分别代入式<1>与式<2>之中，得到目标空间中K a与f dc的表达式：
Figure PCTCN2021115882-appb-000147
Figure PCTCN2021115882-appb-000148
步骤1-6，分别对比上述式<4>与式<9>，式<5>与式<10>的等式右侧，得到估计速度、估计波束中心斜视角、距离估计值的表达式分别为：
Figure PCTCN2021115882-appb-000149
Figure PCTCN2021115882-appb-000150
Figure PCTCN2021115882-appb-000151
步骤2，根据估计速度
Figure PCTCN2021115882-appb-000152
的方向将脉冲压缩信号s RC(t，η)分段得到N个段及每个段对应的段压缩信号s RC，i(t，η)，i＝1...N，具体如下：
根据估计速度
Figure PCTCN2021115882-appb-000153
的正负方向初步将脉冲压缩信号s RC(t，η)依次划分为速度方向一致的段，然后检查每段长度是否小于一个合成孔径长度L s＝R×θ BW，若小于，则将该段向两边扩展为一个合成孔径长度，否则不进行处理。最终，脉冲压缩信号s RC(t，η)被划分为N个段，每段对应的回波信号表示为s RC，i(t，η)，i＝1...N。
步骤3，根据估计速度
Figure PCTCN2021115882-appb-000154
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000155
计算每个段对应的段压缩信号s RC，i(t，η)的相位补偿量
Figure PCTCN2021115882-appb-000156
对段压缩信号s RC，i(t，η)乘以运动误差补偿滤波器
Figure PCTCN2021115882-appb-000157
得到N个已补偿信号为s MC，i(t，η)，虚数
Figure PCTCN2021115882-appb-000158
实现对每个段压缩信号s RC，i(t，η)进行运动补偿。具体按照如下子步骤进行：
步骤3-1，计算各段平台运动状态变化带来的斜距变化δR对相位的影响，表达式为：
Figure PCTCN2021115882-appb-000159
式中，v 0为
Figure PCTCN2021115882-appb-000160
的均值，θ 0为
Figure PCTCN2021115882-appb-000161
的均值；
步骤3-2，构建各段对应的段压缩信号s RC，i(t，η)的相位补偿量的表达式：
Figure PCTCN2021115882-appb-000162
步骤3-3，构建运动误差补偿滤波器表达式：
Figure PCTCN2021115882-appb-000163
其中，j表示虚数，
Figure PCTCN2021115882-appb-000164
步骤3-4，对段压缩信号s RC，i(t，η)乘以如上述式<14>所示运动误差补偿滤波器，得到已补偿信号s MC，i(t，η)，完成对每个段压缩信号的运动补偿。
步骤4，对已补偿信号s MC，i(t，η)进行二维傅里叶变换得到N个二维频谱s MC，i(f，f d)，采用级数反演法对二维频谱s MC，i(f，f d)进行分解，构建方位压缩滤波器H AC，i，其中，f为与距离快时间t对应的频率，f d为与方位慢时间η对应的多普勒频率。具体按照如下步骤进行：
步骤4-1，对已补偿信号s MC，i(t，η)进行二维傅里叶变换得到各段对应的二维频谱s MC，i(f，f d)，其中，f为与距离快时间t对应的频率，f d为与方位慢时间η对应的多普勒频率；
步骤4-2-1，根据驻相点法求得各段二维频谱s MC，i(f，f d)的表达式：
Figure PCTCN2021115882-appb-000165
其中，f c表示系统发射信号的载频，c表示光速；
步骤4-2-2，采用级数反演法对上述式<17>分解，得到消除f和f d的耦合项的二维频谱表达式：
Figure PCTCN2021115882-appb-000166
步骤4-2-3，根据上述式<15>构建理想相位滤波器表达式：
Figure PCTCN2021115882-appb-000167
步骤4-2-4，将估计速度
Figure PCTCN2021115882-appb-000168
估计波束中心斜视角
Figure PCTCN2021115882-appb-000169
距离估计值
Figure PCTCN2021115882-appb-000170
代入上述式<19>之中，得到相位滤波器表达式：
Figure PCTCN2021115882-appb-000171
步骤5，对二维频谱s MC，i(f，f d)乘方位压缩滤波器H AC，i，然后进行二维逆傅里叶变换得到N个成像结果s IMG，i(t，η)。
步骤6，参照图3，依次对相邻的段对应的成像结果s IMG，i(t，η)重叠区域中的强聚焦点所在距离维包络对齐与相参积累，对不重叠区域拼接，得到最终成像结果S all。具体按照如下子步骤进行：
步骤6-1，对成像结果s IMG，i(t，η)进行方位向傅里叶变换，得到距离多普勒域图s IMG，i(t，f d)；
步骤6-2，根据傅里叶变换特性和目标空间几何结构，构建校正图像倾斜的倾斜校正滤波器的表达式：
Figure PCTCN2021115882-appb-000172
式中，l为典型建筑物的距离向尺度；
步骤6-3，将距离多普勒域图s IMG，i(t，f d)乘以倾斜校正滤波器H GC-1，得到已倾斜校正频域图
Figure PCTCN2021115882-appb-000173
步骤6-4，对已倾斜校正频域图
Figure PCTCN2021115882-appb-000174
进行方位向逆傅里叶变换，得到已倾斜校正时域图
Figure PCTCN2021115882-appb-000175
步骤6-5，根据目标空间几何结构，得到拉伸/压缩因子表达式：
Figure PCTCN2021115882-appb-000176
步骤6-6，将上述式<22>代入已倾斜校正时域图
Figure PCTCN2021115882-appb-000177
之中，得到已形变校正时域图
Figure PCTCN2021115882-appb-000178
步骤6-7，对已形变校正时域图
Figure PCTCN2021115882-appb-000179
进行距离向傅里叶变换，得到已形变校正频域图
Figure PCTCN2021115882-appb-000180
步骤6-8，根据傅里叶变换特性和目标空间几何结构，构建校正图像平移的二位置校正滤波器的表达式：
Figure PCTCN2021115882-appb-000181
步骤6-9，将已形变校正频域图
Figure PCTCN2021115882-appb-000182
乘以上述式<23>中的位置校正滤波器H GC-3，得到已几何校正频域图
Figure PCTCN2021115882-appb-000183
步骤6-10，对已几何校正频域图
Figure PCTCN2021115882-appb-000184
进行距离向逆傅里叶变换，得到各段完成几何校正操作的已几何校正时域图
Figure PCTCN2021115882-appb-000185
步骤6-11，将已几何校正时域图
Figure PCTCN2021115882-appb-000186
沿逆时针旋转
Figure PCTCN2021115882-appb-000187
度，得到斜距垂直于机动平台轨迹的待拼接时域图
Figure PCTCN2021115882-appb-000188
步骤6-12，依次对相邻的待拼接时域图
Figure PCTCN2021115882-appb-000189
重叠区域中的强聚焦点所在距离维包络对齐；
步骤6-13，对重叠区域进行相参积累，对不重叠区域依次连接，完成子孔径拼接，得到全孔径成像结果S all。
旋翼无人机载合成孔径雷达分段孔径成像及定位方法的定位部分包括如下步骤：
步骤S1，对原始回波信号s(t，η)进行距离脉冲压缩得到脉冲压缩信号s RC(t，η)，根据所述脉冲压缩信号s RC(t，η)中强散射点的相位历程
Figure PCTCN2021115882-appb-000190
计算机动平台的估计速度
Figure PCTCN2021115882-appb-000191
和估计波束中心斜视角
Figure PCTCN2021115882-appb-000192
本步骤的与成像方法部分的步骤1完全相同，在此不再赘述。
步骤S2，根据估计速度
Figure PCTCN2021115882-appb-000193
的方向将脉冲压缩信号s RC(t，η)分段得到N个段及每个段对应的段压缩信号s RC，i(t，η)，i＝1...N。并且根据估计波束中心斜视角
Figure PCTCN2021115882-appb-000194
估计第i个段对应的平台轨迹坐标
Figure PCTCN2021115882-appb-000195
其中，M为各段方位向长度。
Figure PCTCN2021115882-appb-000196
Figure PCTCN2021115882-appb-000197
Figure PCTCN2021115882-appb-000198
其中，θ in表示无人机载合成孔径雷达波束与地面法向的夹角。
步骤S3，在第i个和第i-1个段中的相邻区域，提取三个强散射点，第i个段中强散射点的坐标分别为
Figure PCTCN2021115882-appb-000199
第i-1个段中强散射点的坐标分别为
Figure PCTCN2021115882-appb-000200
步骤S4，根据强散射点的坐标
Figure PCTCN2021115882-appb-000201
和
Figure PCTCN2021115882-appb-000202
计算第i个和第i-1个段的旋转矩阵
Figure PCTCN2021115882-appb-000203
步骤S5，将第i个所述段的所述平台轨迹坐标
Figure PCTCN2021115882-appb-000204
旋转至以第i-1个所述段的所述平台轨迹的坐标系下，使相邻段的平台轨迹坐标对齐
Figure PCTCN2021115882-appb-000205
步骤S6，对第i个和第i-1个段重叠区域的平台轨迹坐标进行相参积累，对不重叠区域的平台轨迹坐标进行拼接，得到第i和第i-1段的拼接后轨迹坐标[P x，P y，P z]。
Figure PCTCN2021115882-appb-000206
步骤S7，重复N-1次步骤S3至S6，直至得到所有相邻段的拼接后轨迹坐标[P x，P y，P z]，从而得到平台的最终轨迹坐标
Figure PCTCN2021115882-appb-000207
下面结合具体的对比验证结果进一步说明本方法的优点。
1.仿真条件：
如表1所列：
表1仿真参数
Figure PCTCN2021115882-appb-000208
2.仿真与试验
验证1：在表1条件下，仿真实测轨迹得到的点阵目标原始回波信号s r(t，η)，实测轨迹来自于某次试验惯导采集的实测数据，其轨迹如图5(a)、(b)所示。其中，图5(a)表示惯导记录的机动平台飞行时的姿态角随机动平台移动的变化曲线，从上至下三幅图分别为横滚角，俯仰角，偏航角相对于平台移动的变化曲线，图中横坐标表示机动平台飞过的路径从-200m飞行至200m。图5(b)表示惯导记录的机动平台运动轨迹在三个维度的曲线，从上至下三幅图分别为方位维度、距离维度、高度维度的变化曲线。
点阵目标为距离、方位间隔分别为25m的7×7个点目标，点阵目标与平台轨迹的空间位置关系如图6所示。图6所示坐标系X、Y、Z分别对应于方位、距离、高度维度，图中R ref表示系统在当前入射角θ in下的参考距离，H表示机动平台飞行的平均高度，θ a表示系统方位维的波束宽度，θ bw表示系统距离维的波束宽度，S表示待测区域的宽度，v表示机动平台的运动速度，L s表示一个合成孔径长度。图右侧的点阵表示仿真目标的位置关系，点阵中心位于图片左侧R ref与Y轴的交点。点阵的长宽分别均匀排列了7个点，点间距为25m，其余各参数取值参考表1。
对s r(t，η)距离脉冲压缩后得到s RD(t，η)，对s RD(t，η)分别采用SAI算法、传统成像方法进行成像对比，结果如图7、图8所示。其中，该路径被本文提出的SAI算法分割为三段(段 1、段2、段3)，图7中的(a)、(b)、(c)分别为本实施例中段1、段2、段3几何校正后的成像结果，由于每个段只能对部分区域成像，因此各段成像结果只包含了目标区域的一部分，另外，在方位维边缘的地方，可以看到点目标的轮廓有扩散的情况。图7中的(d)为本实施例中的各段拼接后的全孔径成像结果，图7中的(e)为传统成像方法的成像结果。对比图7(d)与(e)可以看出聚焦点轮廓清晰度不同，并且明暗程度有差异，其中以方框框出的几个点对比最为明显。图8为本实施例与传统成像方法点扩散函数对比结果。图8之中上图为段1、段2、段3第三行点目标的点扩散函数曲线对比图，下图为本文提出的SAI方法与传统成像方法第三行点目标的点扩散函数曲线对比图。对段1、段2、段3、SAI方法、传统成像方法的点扩散函数曲线每个尖峰的-3dB位置处主瓣展宽进行分析，可以得到分辨率相对于理论值扩展的幅度，根据图中数据统计结果如表2A)所示，主瓣展宽系数越小，图像聚焦效果越好；对比这五条曲线每个峰值的主瓣和副瓣的幅度高低，统计结果如表2B)所示，主副瓣比越大表示图像信噪比越大，即图像质量越好；对比这五条曲线的峰值幅度损失可以让我们观察到不完整孔径对图像的影响，虽然各段峰值幅度损失在边缘处较大，但是SAI最终全孔径图像由于各段的互补关系使其峰值幅度损失大大降低；衡量成像结果的某一行点的点扩散函数曲线指标以后，通过图像熵判断完整图像所包含的信息量，图像聚焦效果越高所包含的信息量越多，图像熵就越小，SAI方法最终全孔径图像是熵最小的。从主瓣展宽系统、主副瓣比、峰值幅度损失、图像熵这四个方面对比，本文提出的SAI方法相比于传统成像方法具有更高的分辨率、更大的信噪比、并且峰值幅度损失和原始方法相当，所包含的信息更多，因此本文提出的SAI具备很好的成像效果，解决了无人机载合成孔径雷达所遇到的实际问题。
表2图像质量指标对比
A)主瓣展宽系数对比
Figure PCTCN2021115882-appb-000209
B)主副瓣比(dB)
Figure PCTCN2021115882-appb-000210
C)峰值幅度损失(dB)
Figure PCTCN2021115882-appb-000211
Figure PCTCN2021115882-appb-000212
D)图像熵
Figure PCTCN2021115882-appb-000213
验证2：将表1条件之中平台速度改为10m/s，平台飞行高度改为350m，即为实测飞行试验的试验条件，试验平台为KWT-65型号多旋翼无人机，携带Ku波段小型调频连续波合成孔径雷达，使用该平台飞行约100架次，每架次飞行距离约1km，采集了地面某区域的多组回波数据。采用传统方法和本发明方法对多组回波数据成像后，聚焦效果对比如图9所示。其中，图9的(a)为传统成像方法的成像效果；图9的(b)为本实施例的方法的成像效果；图9的(c)为对该成像区域对应的光学模型。
图9表明了本实施例的方法在聚焦效果上优于传统成像方法的聚焦效果。根据待测目标的光学图像可以看出目标呈现“工”字形。采用本发明方法的成像结果“工”字轮廓要比采用传统方法更为清晰，以图像“工”型楼的左边为例，雷达照射玻璃时，会穿透玻璃进入房间内部，无法将回波返回至雷达，因此成像结果会呈现出黑色，而窗户边缘因为存在直角和金属会有明显的反射，形成亮点。本发明方法所对应的成像结果，可以通过亮点之间黑色区域清晰观察得到楼层间隔，传统方法所对应的成像结果中亮斑都模糊成一团，无法清晰观察到楼层之间的玻璃。“工”型楼左上方存在一条马路，处采用本发明方法可以清晰看到马路边缘，传统方法成像得到的马路边缘轮廓不甚清晰，还有多余的纵向条纹。并且“工”型楼周围的花园轮廓本发明方法成像效果也更为清晰。图像整体的信噪比也优于传统方法。此外，经过多次/多地飞行试验验证本方法对90％的回波数据成像有较好的聚焦效果，而传统方法聚焦效果较好的图像约为35％。对于同样尺度的图像(12500×8192)在相同计算机平台，本方法成像时间约为320s，传统成像方法成像时间约为5200s，速度约提升了16倍。
验证3：对验证2条件下某次实测试验数据进行处理，得到图像的同时，可以估计得到平台的轨迹。采用本发明方法估计轨迹和用惯导设备采集的轨迹对比如图10所示。其中，图10的(a)为本发明方法估计轨迹和惯导设备采集轨迹从上至下分别是在方位维度、距离维度、维度的对比效果，其中蓝色曲线表示惯导设备的测量结果，红色曲线为本发明方法的估计结果，可以看出在方位上两者基本重合，使用多旋翼平台的高采样率保证了方位上的小误差，在距离维两者略有偏差，蓝色曲线在红色曲线的边缘处，尤其起伏较大的三个尖峰的地方，蓝色曲线会有偏差，在高度维可以看到蓝色曲线和红色曲线偏差会更大一些，由于本发明方法基于二维图像估计，因此部分误差会积累于第三个维度即高度维之上；图9的(b)为本实施例的方法的估计误差曲线，该曲线显示在方位维、距离维、高度维的偏差分别约为厘米级、厘米级、分米级，达到和惯导接近的水平。
实施例的作用与效果
根据本实施例提供的多旋翼无人机载合成孔径雷达分段孔径成像方法，方法主要包括：1)根据无人机载合成孔径雷达系统获取目标回波；2)根据回波信号估计机动平台的运动状 态；3)根据平台运动状态对回波信号分段；4)根据平台运动状态对各段回波信号运动补偿；5)对各段补偿后的回波信号进行二维傅里叶变换得到二维频谱，采用级数反演法分解二维频谱得到各段相位滤波器；6)对各段二维频谱乘该段对应的相位滤波器，然后对二维频谱进行二维逆傅里叶变换得到各段图像；7)对各段图像进行几何校正，然后拼接各段图像，得到全孔径成像结果。8)对各段平台轨迹进行拼接，得到平台完整轨迹坐标。本实施例采用基于回波信号相位历程估计平台运动参数和对回波信号分段的方法对平台运动误差精确补偿，达到成像聚焦效果好和成像成功率高的效果，提高了多旋翼无人机平台数据采集的效率，可以有效地对多旋翼无人机平台合成孔径雷达系统进行高分辨成像。同时，本方法可以在成像同时计算平台轨迹的三维坐标，达到对平台定位的效果，可应用于无人机导航，为以后智能化、一体化探测系提供可能。
与现有技术相比具有如下优点：
与传统的机载合成孔径雷达成像算法相比，考虑到了平台运动与信号相位之间的关系，可适用于无惯导设备或者低精度惯导设备的合成孔径雷达成像系统进行成像；
与传统的机载合成孔径雷达成像算法相比，考虑到了无人机平台姿态角变化带来的斜视影响，补偿了距离和方位的耦合相位，提高了图像聚焦效果，可适用于无天线伺服的合成孔径雷达成像系统；
与传统的机载合成孔径雷达成像算法相比，采用相位滤波器相乘的方法，替代插值方法，提高了各段内单次成像的成像速度；
与传统的机载合成孔径雷达成像算法相比，采用各段并行成像然后拼接为完整图像的办法获取图像，提高了完整图像的成像速度。
另外，通过实验验证表明：本实施例提出的多旋翼无人机载合成孔径雷达分段孔径成像(SAI)方法成像分辨率很高，并且运算速度快，同时，可以对平台轨迹进行估计，降低了算法对硬件设备的要求，表明本发明能有效应用于小型机动平台合成孔径雷达成像系统。
本实施例详细推导了平台运动参数与回波相位之间的对应关系，并对斜视下二维频谱的二阶展开式，可以不依赖于惯导设备进行运动补偿和成像。仿真时对比了本实施例方法各段、完整图像和传统成像算法的点扩散函数曲线；实测验证对比了本发明方法与传统算法的成像结果图，证明本实施例可以有效地对多旋翼无人机平台合成孔径雷达系统进行高分辨成像。
上述实施方式为本发明的优选案例，并不用来限制本发明的保护范围。

Claims (10)

  1. 一种多旋翼无人机载合成孔径雷达分段孔径成像方法，用于根据无人机载合成孔径雷达的原始回波信号进行分段孔径成像，其特征在于，包括以下步骤：
    步骤1，对所述原始回波信号s(t，η)进行距离脉冲压缩得到距离脉压信号s RC(t，η)，其中t为距离快时间、η为方位慢时间，根据所述距离脉压s RC(t，η)中强散射点的相位历程
    Figure PCTCN2021115882-appb-100001
    计算机动平台的估计速度
    Figure PCTCN2021115882-appb-100002
    和估计波束中心斜视角
    Figure PCTCN2021115882-appb-100003
    步骤2，根据所述估计速度
    Figure PCTCN2021115882-appb-100004
    的方向将所述距离脉压信号s RC(t，η)分段得到N个段及每个所述段对应的段脉压信号s RC，i(t，η)，i＝1...N；
    步骤3，根据所述估计速度
    Figure PCTCN2021115882-appb-100005
    和所述估计波束中心斜视角
    Figure PCTCN2021115882-appb-100006
    计算每个所述段对应的所述段脉压信号s RC，i(t，η)的相位补偿量
    Figure PCTCN2021115882-appb-100007
    对所述段脉压信号s RC，i(t，η)乘以运动误差补偿滤波器
    Figure PCTCN2021115882-appb-100008
    得到N个已补偿信号为s MC，i(t，η)，虚数
    Figure PCTCN2021115882-appb-100009
    步骤4，对所述已补偿信号s MC，i(t，η)进行二维傅里叶变换得到N个二维频谱s MC，i(f，f d)，采用级数反演法对所述二维频谱s MC，i(f，f d)进行分解，构建方位压缩滤波器H AC，i，其中，f为与距离快时间t对应的频率，f d为与方位慢时间η对应的多普勒频率；
    步骤5，对所述二维频谱s MC，i(f，f d)乘所述方位压缩滤波器H AC，i，然后进行二维逆傅里叶变换得到N个成像结果s IMG，i(t，η)；
    步骤6，依次对相邻的所述段对应的所述成像结果s IMG，i(t，η)重叠区域中的强聚焦点所在距离维包络对齐与相参积累，对不重叠区域拼接，得到最终成像结果S all。
  2. 根据权利要求1所述的多旋翼无人机载合成孔径雷达分段孔径成像方法，其特征在于：
    其中，所述步骤1中，对所述距离脉压信号s RC(t，η)中所述强散射点的所述相位历程
    Figure PCTCN2021115882-appb-100010
    进行二阶拟合得到该强散射点的相位历程
    Figure PCTCN2021115882-appb-100011
    根据二阶项系数β和一阶项系数α计算所述机动平台的所述估计速度
    Figure PCTCN2021115882-appb-100012
    以及所述估计波束中心斜视角
    Figure PCTCN2021115882-appb-100013
    其中，t为距离快时间，η为方位慢时间，o(η)为高阶相位误差，
    Figure PCTCN2021115882-appb-100014
    为常数相位项，λ为系统发射信号波长，
    Figure PCTCN2021115882-appb-100015
    为系统参考距离的估计值。
  3. 根据权利要求1所述的多旋翼无人机载合成孔径雷达分段孔径成像方法，其特征在于：
    其中，所述步骤2中，根据所述估计速度
    Figure PCTCN2021115882-appb-100016
    的方向将所述距离脉压信号s RC(t，η)依次划分为速度方向一致的N个所述段，然后判断每个所述段的长度是否小于一个合成孔径长度，若小于，则将该段向两边扩展为一个所述合成孔径长度，最终得到分段后的N个所述段脉压信号s RC，i(t，η)，i＝1...N。
  4. 根据权利要求1所述的多旋翼无人机载合成孔径雷达分段孔径成像方法，其特征在于：
    其中，所述步骤3中，所述相位补偿量
    Figure PCTCN2021115882-appb-100017
    式中，R 0为各段
    Figure PCTCN2021115882-appb-100018
    的均值，θ 0为N个所述估计波束中心斜视角
    Figure PCTCN2021115882-appb-100019
    的均值。
  5. 根据权利要求1所述的多旋翼无人机载合成孔径雷达分段孔径成像方法，其特征在于：
    其中，所述步骤4中，所述方位压缩滤波器
    Figure PCTCN2021115882-appb-100020
    式中，f为与距离快时间t对应的频率，f d为与方位慢时间η对应的多普勒频率，f c为所述系统发射信号的载频，c为光速。
  6. 根据权利要求1所述的多旋翼无人机载合成孔径雷达分段孔径成像方法，其特征在于：
    其中，所述步骤6中，先对相邻所述段对应的所述成像结果s IMG，i(t，η)进行几何校正得到已校正成像结果
    Figure PCTCN2021115882-appb-100021
    然后将所述已校正成像结果
    Figure PCTCN2021115882-appb-100022
    旋转
    Figure PCTCN2021115882-appb-100023
    度，得到斜距垂直于所述机动平台轨迹的所述已校正成像结果
    Figure PCTCN2021115882-appb-100024
    然后依次对相邻的所述段对应的所述已校正成像结果
    Figure PCTCN2021115882-appb-100025
    重叠区域中的所述强聚焦点所在距离维包络对齐与相参积累，对不重叠区域拼接，得到所述最终成像结果S all。
  7. 根据权利要求6所述的多旋翼无人机载合成孔径雷达分段孔径成像方法，其特征在于：
    其中，所述步骤6对所述成像结果s IMG，i进行几何校正按如下子步骤进行：
    步骤6-1，对所述成像结果s IMG，i(t，η)进行方位向傅里叶变换，得到距离多普勒域图s IMG，i(t，f d)；
    步骤6-2，根据傅里叶变换特性和目标空间几何结构，构建校正图像倾斜的倾斜校正滤波器的表达式：
    Figure PCTCN2021115882-appb-100026
    式中，l为典型建筑物的距离向尺度；
    步骤6-3，将所述距离多普勒域图s IMG，i(t，f d)乘以所述倾斜校正滤波器H GC-1，得到已倾斜校正频域图
    Figure PCTCN2021115882-appb-100027
    步骤6-4，对所述已倾斜校正频域图
    Figure PCTCN2021115882-appb-100028
    进行方位向逆傅里叶变换，得到已倾斜校正时域图
    Figure PCTCN2021115882-appb-100029
    步骤6-5，根据目标空间几何结构，得到拉伸/压缩因子表达式：
    Figure PCTCN2021115882-appb-100030
    步骤6-6，将所述拉伸/压缩因子表达式代入所述已倾斜校正时域图
    Figure PCTCN2021115882-appb-100031
    之中，得到已形变校正时域图
    Figure PCTCN2021115882-appb-100032
    步骤6-7，对所述已形变校正时域图
    Figure PCTCN2021115882-appb-100033
    进行距离向傅里叶变换，得到已形变校正频域图
    Figure PCTCN2021115882-appb-100034
    步骤6-8，根据傅里叶变换特性和目标空间几何结构，构建校正图像平移的二位置校正滤波器的表达式：
    Figure PCTCN2021115882-appb-100035
    步骤6-9，将所述已形变校正频域图
    Figure PCTCN2021115882-appb-100036
    乘以所述位置校正滤波器H GC-3，得到已几何校正频域图
    Figure PCTCN2021115882-appb-100037
    步骤6-10，对所述已几何校正频域图
    Figure PCTCN2021115882-appb-100038
    进行距离向逆傅里叶变换得到各段完成几何校正操作的已几何校正时域图
    Figure PCTCN2021115882-appb-100039
    步骤6-11，将所述已几何校正时域图
    Figure PCTCN2021115882-appb-100040
    沿逆时针旋转
    Figure PCTCN2021115882-appb-100041
    度，得到斜距垂直于机动平台轨迹的待拼接时域图
    Figure PCTCN2021115882-appb-100042
    步骤6-12，依次对相邻的所述待拼接时域图
    Figure PCTCN2021115882-appb-100043
    重叠区域中的强聚焦点所在距离维包络对齐；
    步骤6-13，对重叠区域进行相参积累，对不重叠区域依次连接，完成子孔径拼接，得到所述全孔径成像结果S all。
  8. 根据权利要求1所述的多旋翼无人机载合成孔径雷达分段孔径成像方法，其特征在于：
    其中，所述步骤1中，所述估计速度
    Figure PCTCN2021115882-appb-100044
    及所述估计波束中心斜视角
    Figure PCTCN2021115882-appb-100045
    按如下步骤计算：
    步骤1-1，根据多普勒调频斜率K a和多普勒中心f dc的定义，得到K a与f dc的计算式：
    Figure PCTCN2021115882-appb-100046
    Figure PCTCN2021115882-appb-100047
    式中，
    Figure PCTCN2021115882-appb-100048
    表示(·)对方位慢时间η求导；
    步骤1-2，将所述原始回波信号所构成的空间称为信号空间，对所述距离脉压信号s RC(t，η)中所述强散射点的所述相位历程
    Figure PCTCN2021115882-appb-100049
    进行二阶拟合得到该强散射点在所述信号空间中的表达式：
    Figure PCTCN2021115882-appb-100050
    式中，β为二阶项系数，α为一阶项系数，
    Figure PCTCN2021115882-appb-100051
    为常数相位项，o(η)为高阶相位误差；
    步骤1-3，将所述相位历程
    Figure PCTCN2021115882-appb-100052
    的表达式分别代入步骤1-1中K a的计算式与f dc的计算式之中，得到所述信号空间中K a与f dc的表达式：
    Figure PCTCN2021115882-appb-100053
    Figure PCTCN2021115882-appb-100054
    步骤1-4，将目标和所述机动平台实际位置所构成的空间称为目标空间，根据空间几何关系，得到所述目标空间中的所述相位历程
    Figure PCTCN2021115882-appb-100055
    的表达式：
    Figure PCTCN2021115882-appb-100056
    式中，v表示所述机动平台运动速度，θ表示所述机动平台运动引起的波束中心斜视角，R 0表示所述目标和所述机动平台的最近距离，λ为系统发射信号波长；
    步骤1-5，将所述目标空间中的所述相位历程
    Figure PCTCN2021115882-appb-100057
    的表达式分别代入步骤1-1中K a的计算式与f dc的计算式之中，得到所述目标空间中K a与f dc的表达式：
    Figure PCTCN2021115882-appb-100058
    Figure PCTCN2021115882-appb-100059
    步骤1-6，分别对比所述信号空间中K a的表达式与所述目标空间中K a表达式，所述信号空间中f dc的表达式与所述目标空间中f dc表达式的右侧，得到所述估计速度、所述估计波束中心斜视角、所述距离估计值：
    Figure PCTCN2021115882-appb-100060
    Figure PCTCN2021115882-appb-100061
    Figure PCTCN2021115882-appb-100062
  9. 根据权利要求1所述的多旋翼无人机载合成孔径雷达分段孔径成像方法，其特征在于：
    其中，所述步骤4中，采用级数反演法对所述二维频谱s MC，i(f，f d)分解得到所述方位压缩滤波器表达式H AC，i的过程按如下子步骤进行：
    步骤4-2-1，根据驻相点法求得所述段对应的所述二维频谱s MC，i(f，f d)的表达式：
    Figure PCTCN2021115882-appb-100063
    式中，f c表示系统发射信号的载频，c表示光速；
    步骤4-2-2，采用级数反演法对步骤4-2-1中的所述二维频谱s MC，i(f，f d)的表达式进行分解，得到消除f和f d的耦合项的二维频谱表达式：
    Figure PCTCN2021115882-appb-100064
    步骤4-2-3，根据步骤4-2-2中的所述消除f和f d的耦合项的二维频谱表达式构建理想相位滤波器表达式：
    Figure PCTCN2021115882-appb-100065
    步骤4-2-4，将所述估计速度
    Figure PCTCN2021115882-appb-100066
    所述估计波束中心斜视角
    Figure PCTCN2021115882-appb-100067
    所述距离估计值
    Figure PCTCN2021115882-appb-100068
    代入所述理想相位滤波器表达式之中，得到所述方位压缩滤波器表达式：
    Figure PCTCN2021115882-appb-100069
  10. 一种多旋翼无人机载合成孔径雷达分段孔径定位方法，用于根据无人机载合成孔径雷达的原始回波信号，计算无人机飞行轨迹其特征在于，包括以下步骤：
    步骤S1，对所述原始回波信号s(t，η)进行距离脉冲压缩得到距离脉压信号s RC(t，η)，根据所述距离脉压信号s RC(t，η)中强散射点的相位历程
    Figure PCTCN2021115882-appb-100070
    计算机动平台的估计速度
    Figure PCTCN2021115882-appb-100071
    和估计波束中心斜视角
    Figure PCTCN2021115882-appb-100072
    步骤S2，根据所述估计速度
    Figure PCTCN2021115882-appb-100073
    的方向将所述距离脉压信号s RC(t，η)分段得到N个段及每个所述段对应的段脉压信号s RC，i(t，η)，i＝1...N，并且根据所述估计速度
    Figure PCTCN2021115882-appb-100074
    和所述估计波束中心斜视角
    Figure PCTCN2021115882-appb-100075
    计算第i个所述段对应的平台轨迹坐标
    Figure PCTCN2021115882-appb-100076
    M为各所述段方位向长度，
    Figure PCTCN2021115882-appb-100077
    Figure PCTCN2021115882-appb-100078
    Figure PCTCN2021115882-appb-100079
    式中，θ in表示所述无人机载合成孔径雷达波束与地面法向的夹角；
    步骤S3，在第i个和第i-1个所述段中的相邻区域，提取三个强散射点，第i个所述段中所述强散射点的坐标分别为
    Figure PCTCN2021115882-appb-100080
    第i-1个所述段中所述强散射点的坐标分别为
    Figure PCTCN2021115882-appb-100081
    步骤S4，根据所述强散射点的坐标
    Figure PCTCN2021115882-appb-100082
    和
    Figure PCTCN2021115882-appb-100083
    计算第i个和第i-1个所述段的旋转矩阵γ，
    Figure PCTCN2021115882-appb-100084
    步骤S5，将第i个所述段的所述平台轨迹坐标
    Figure PCTCN2021115882-appb-100085
    旋转至以第i-1个所述段的所述平台轨迹的坐标系下，使相邻所述段的所述平台轨迹坐标对齐，即：
    Figure PCTCN2021115882-appb-100086
    步骤S6，对第i个和第i-1个所述段重叠区域的所述平台轨迹坐标进行相参积累，对不重叠区域的所述平台轨迹坐标进行拼接，得到第i个和第i-1个所述段的拼接后轨迹坐标[P x，P y，P z]，
    Figure PCTCN2021115882-appb-100087
    步骤S7，重复步骤S3至S6，直至得到所有相邻所述段的拼接后轨迹坐标[P x，P y，P z]，从而得到平台的最终轨迹坐标
    Figure PCTCN2021115882-appb-100088
PCT/CN2021/115882 2021-08-13 2021-09-01 一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法 WO2023015623A1 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
US18/576,667 US20240319364A1 (en) 2021-08-13 2021-09-01 Segmented aperture imaging and positioning method of multi-rotor unmanned aerial vehicle-borne synthetic aperture radar

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN202110929813.9A CN115704892A (zh) 2021-08-13 2021-08-13 一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法
CN202110929813.9 2021-08-13

Publications (1)

Publication Number Publication Date
WO2023015623A1 true WO2023015623A1 (zh) 2023-02-16

Family

ID=85181155

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2021/115882 WO2023015623A1 (zh) 2021-08-13 2021-09-01 一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法

Country Status (3)

Country Link
US (1) US20240319364A1 (zh)
CN (1) CN115704892A (zh)
WO (1) WO2023015623A1 (zh)

Cited By (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116482686A (zh) * 2023-06-21 2023-07-25 中国科学院空天信息创新研究院 一种基于方位向自适应分块的高分辨率isar成像方法
CN116540203A (zh) * 2023-07-04 2023-08-04 西安电子科技大学 基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法
CN116626629A (zh) * 2023-07-21 2023-08-22 中国科学院空天信息创新研究院 一种星载斜视sar系统模糊综合性能评估方法
CN116879895A (zh) * 2023-07-05 2023-10-13 南京航空航天大学 一种基于改进ωk的条带SAR斜视成像方法
CN117148351A (zh) * 2023-10-31 2023-12-01 西安电子科技大学 基于卫星sar图像的弹载sar图像成像方法及装置
CN117289277A (zh) * 2023-11-27 2023-12-26 中山大学 一种基于子带分割合成的多频雷达三维成像方法及系统
CN117310633A (zh) * 2023-11-28 2023-12-29 中国科学院空天信息创新研究院 一种多维度sar系统几何定标方法
CN117347971A (zh) * 2023-09-05 2024-01-05 武汉大学 一种激光雷达建筑物边界估测方法及计算机可读介质
CN117572401A (zh) * 2023-11-24 2024-02-20 北京理工大学 一种基于测量子孔径多普勒相位差的无源定位方法
CN117784133A (zh) * 2024-01-17 2024-03-29 中国科学院空天信息创新研究院 一种合成孔径无源定位迭代运动补偿方法和系统
CN117872373A (zh) * 2024-03-11 2024-04-12 中国科学院空天信息创新研究院 一种面向轻小型无人机的实时处理微型sar系统
CN118068332A (zh) * 2024-04-25 2024-05-24 中国石油大学(华东) 适合于调频连续波的合成孔径雷达距离多普勒成像方法
CN118112566A (zh) * 2024-04-25 2024-05-31 中国石油大学(华东) 一种无人机sar成像方法
CN118131238A (zh) * 2024-05-07 2024-06-04 中国人民解放军战略支援部队航天工程大学 基于pfa改进的isar成像算法的图像横向定标方法
CN118350201A (zh) * 2024-03-11 2024-07-16 中国科学院声学研究所 一种合成孔径声呐回波仿真方法及系统

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2012093257A (ja) * 2010-10-27 2012-05-17 Mitsubishi Electric Corp レーダ画像処理装置
CN104020471A (zh) * 2014-05-21 2014-09-03 北京无线电测量研究所 一种基于分块处理的sar实时成像方法及系统
CN108344982A (zh) * 2018-02-07 2018-07-31 成都电科智达科技有限公司 基于长时间相参积累的小型无人机目标雷达检测方法
CN109655831A (zh) * 2018-12-27 2019-04-19 中国科学院电子学研究所 一种无人机sar成像方法及装置
CN111190181A (zh) * 2020-01-09 2020-05-22 西安电子科技大学 颠簸平台无人机载sar实时成像处理方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2012093257A (ja) * 2010-10-27 2012-05-17 Mitsubishi Electric Corp レーダ画像処理装置
CN104020471A (zh) * 2014-05-21 2014-09-03 北京无线电测量研究所 一种基于分块处理的sar实时成像方法及系统
CN108344982A (zh) * 2018-02-07 2018-07-31 成都电科智达科技有限公司 基于长时间相参积累的小型无人机目标雷达检测方法
CN109655831A (zh) * 2018-12-27 2019-04-19 中国科学院电子学研究所 一种无人机sar成像方法及装置
CN111190181A (zh) * 2020-01-09 2020-05-22 西安电子科技大学 颠簸平台无人机载sar实时成像处理方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
WANG JIA, YANG LEI; GONG MING: "An Imaging Algorithm for High Squint Missile-borne SAR Based on Block Processing", GUIDANCE & FUZE, vol. 34, no. 3, 30 September 2013 (2013-09-30), XP093034319, ISSN: 1671-0576 *

Cited By (23)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116482686A (zh) * 2023-06-21 2023-07-25 中国科学院空天信息创新研究院 一种基于方位向自适应分块的高分辨率isar成像方法
CN116482686B (zh) * 2023-06-21 2023-08-15 中国科学院空天信息创新研究院 一种基于方位向自适应分块的高分辨率isar成像方法
CN116540203A (zh) * 2023-07-04 2023-08-04 西安电子科技大学 基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法
CN116540203B (zh) * 2023-07-04 2023-09-22 西安电子科技大学 基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法
CN116879895A (zh) * 2023-07-05 2023-10-13 南京航空航天大学 一种基于改进ωk的条带SAR斜视成像方法
CN116626629A (zh) * 2023-07-21 2023-08-22 中国科学院空天信息创新研究院 一种星载斜视sar系统模糊综合性能评估方法
CN116626629B (zh) * 2023-07-21 2023-09-22 中国科学院空天信息创新研究院 一种星载斜视sar系统模糊综合性能评估方法
CN117347971A (zh) * 2023-09-05 2024-01-05 武汉大学 一种激光雷达建筑物边界估测方法及计算机可读介质
CN117347971B (zh) * 2023-09-05 2024-06-11 武汉大学 一种激光雷达建筑物边界估测方法及计算机可读介质
CN117148351A (zh) * 2023-10-31 2023-12-01 西安电子科技大学 基于卫星sar图像的弹载sar图像成像方法及装置
CN117148351B (zh) * 2023-10-31 2024-02-06 西安电子科技大学 基于卫星sar图像的弹载sar图像成像方法及装置
CN117572401A (zh) * 2023-11-24 2024-02-20 北京理工大学 一种基于测量子孔径多普勒相位差的无源定位方法
CN117289277A (zh) * 2023-11-27 2023-12-26 中山大学 一种基于子带分割合成的多频雷达三维成像方法及系统
CN117289277B (zh) * 2023-11-27 2024-01-30 中山大学 一种基于子带分割合成的多频雷达三维成像方法及系统
CN117310633B (zh) * 2023-11-28 2024-02-23 中国科学院空天信息创新研究院 一种多维度sar系统几何定标方法
CN117310633A (zh) * 2023-11-28 2023-12-29 中国科学院空天信息创新研究院 一种多维度sar系统几何定标方法
CN117784133A (zh) * 2024-01-17 2024-03-29 中国科学院空天信息创新研究院 一种合成孔径无源定位迭代运动补偿方法和系统
CN117872373A (zh) * 2024-03-11 2024-04-12 中国科学院空天信息创新研究院 一种面向轻小型无人机的实时处理微型sar系统
CN117872373B (zh) * 2024-03-11 2024-05-07 中国科学院空天信息创新研究院 一种面向轻小型无人机的实时处理微型sar系统
CN118350201A (zh) * 2024-03-11 2024-07-16 中国科学院声学研究所 一种合成孔径声呐回波仿真方法及系统
CN118068332A (zh) * 2024-04-25 2024-05-24 中国石油大学(华东) 适合于调频连续波的合成孔径雷达距离多普勒成像方法
CN118112566A (zh) * 2024-04-25 2024-05-31 中国石油大学(华东) 一种无人机sar成像方法
CN118131238A (zh) * 2024-05-07 2024-06-04 中国人民解放军战略支援部队航天工程大学 基于pfa改进的isar成像算法的图像横向定标方法

Also Published As

Publication number Publication date
CN115704892A (zh) 2023-02-17
US20240319364A1 (en) 2024-09-26

Similar Documents

Publication Publication Date Title
WO2023015623A1 (zh) 一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法
CN109031219B (zh) 基于相位测距的宽带雷达弹道目标微动几何参数估计方法
CN102654576B (zh) 基于sar图像和dem数据的图像配准方法
CN111142105A (zh) 复杂运动目标isar成像方法
CN103293520B (zh) 基于等效距离方程的sar地面运动目标成像方法
CN108710111B (zh) 一种机载双基前视sar方位相位的二维空变校正方法
CN107991659B (zh) 基于字典学习的米波雷达低仰角目标测高方法
CN104597447B (zh) 一种子孔径SAR大斜视改进Omega‑K成像方法
CN109444882B (zh) 基于变斜视椭圆波束同步模型的双站sar成像方法
CN102645651A (zh) 一种sar层析超分辨成像方法
CN102866393B (zh) 一种基于pos与dem数据的sar多普勒参数估计方法
CN103278820A (zh) 临近空间慢速平台sar动目标检测方法及成像方法
CN102226841A (zh) 基于高阶多项式距离方程的同步轨道sar成像方法
CN102778681A (zh) 固定发射站双基地前视合成孔径雷达成像方法
CN106291548B (zh) 基于惯导信息和回波数据的Ka调频连续波SAR运动补偿方法
CN108427115A (zh) 合成孔径雷达对运动目标参数的快速估计方法
CN110187337B (zh) 一种基于ls和neu-ecef时空配准的高机动目标跟踪方法及系统
CN112859018B (zh) 一种基于图像几何校正的视频sar成像方法
CN110133682A (zh) 星载全方位sar自适应目标三维重建方法
CN104407349B (zh) 一站固定式双站低频超宽带sar的频域成像方法
CN102736073B (zh) 一种通用模式下星载sar距离向模糊度的计算方法
CN113671492A (zh) 一种面向机动平台前视成像的samp重构方法
CN110161503A (zh) 基于三维等距圆模型的近距宽域sar高分辨成像方法
CN109946682B (zh) 基于ICESat/GLAS的GF3数据基线估计方法
CN115685200A (zh) 一种高精度大前斜视sar成像运动补偿与几何校正方法

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 21953233

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE