CN108983169B - Meter wave radar terrain correction method based on digital elevation model - Google Patents

Meter wave radar terrain correction method based on digital elevation model Download PDF

Info

Publication number
CN108983169B
CN108983169B CN201810781175.9A CN201810781175A CN108983169B CN 108983169 B CN108983169 B CN 108983169B CN 201810781175 A CN201810781175 A CN 201810781175A CN 108983169 B CN108983169 B CN 108983169B
Authority
CN
China
Prior art keywords
meter
wave radar
data
time
flying target
Prior art date
Legal status (The legal status 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 status listed.)
Active
Application number
CN201810781175.9A
Other languages
Chinese (zh)
Other versions
CN108983169A (en
Inventor
宋万杰
张萌
金莉
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Xidian University
Original Assignee
Xidian University
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 Xidian University filed Critical Xidian University
Priority to CN201810781175.9A priority Critical patent/CN108983169B/en
Publication of CN108983169A publication Critical patent/CN108983169A/en
Application granted granted Critical
Publication of CN108983169B publication Critical patent/CN108983169B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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
    • G01S7/00Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/40Means for monitoring or calibrating

Abstract

The invention discloses a terrain correction method for a meter-wave radar based on a digital elevation model, which belongs to the technical field of meter-wave radar signal processing and mainly comprises the following steps: establishing a spherical coordinate system, wherein a meter-wave radar exists in the spherical coordinate system, and a flying target exists in a detection range of the meter-wave radar; acquiring digital elevation model data, determining ADS-B data of a flying target, and then determining an important reflection point in a spherical coordinate system; acquiring the terminal data of the meter-wave radar, and obtaining the correlation completion result of the ADS-B data of the flight target and the terminal data of the meter-wave radar according to the ADS-B data of the flight target and the determined important reflection point, thereby obtaining the height curve of the flight target; obtaining a height measurement error curve of the meter wave radar according to the flight target height curve; and correcting the height measurement error curve of the meter-wave radar to obtain a final error correction result, and then using the final error correction result as a terrain correction result of the meter-wave radar based on the digital elevation model.

Description

Meter wave radar terrain correction method based on digital elevation model
Technical Field
The invention belongs to the technical field of meter-wave radar signal processing, and particularly relates to a terrain correction method for a meter-wave radar based on a digital elevation model, which is suitable for the research of radar height measurement and angle measurement refinement.
Background
The meter-wave radar is also called as a very high frequency radar, belongs to a long-wave radar, and therefore has the advantages of long detection distance, stealth reflection, interference resistance and the like; meanwhile, due to the characteristic that the wave beam is wide, on one hand, the direct wave and the reflected wave are in the same wave beam and cannot be distinguished; on the other hand, the problem of beam landing exists, so that lobes are split, and target elevation angle estimation is deviated when the elevation angle is low; aiming at the multipath effect problem mainly existing in the meter wave radar at present, the problem is mainly solved from two aspects of algorithm improvement and algorithm correction through a large amount of analysis and research; specifically, on one hand, an angle measurement and height measurement algorithm is improved and optimized so as to obtain a more high-precision angle measurement result, on the other hand, the influence of topographic relief on angle measurement and height measurement is researched, specific influence factors are analyzed, a topographic compensation method is provided, and an angle measurement error is corrected.
According to the early angle and height measurement algorithm of the radar, although a plurality of targets can be measured, the formed lobe is wide, so that the angle resolution is low, and when two targets are in the same beam width, the two targets cannot be distinguished, so that the angle measurement accuracy is reduced; in order to solve the problem of insufficient estimation precision of a low elevation angle, two types of super-resolution algorithms, namely a subspace-type algorithm and a maximum likelihood algorithm, are developed, a Multiple Signal Classification algorithm (MUSIC) is a commonly used subspace-type algorithm, the MUSIC algorithm requires that a Signal subspace and a noise subspace have an orthogonal relationship, namely the Signal subspace and the noise subspace must be incoherent, however, a direct wave Signal received by a radar antenna at a low elevation angle is approximately the same as the Doppler frequency of a multipath Signal, namely, the two signals have a coherent relationship, and further the angle measurement precision of the MUSIC algorithm is reduced; the realization process of the Maximum Likelihood algorithm (ML) requires multidimensional search of a space domain, so that the ML algorithm has large computation load; although the two algorithms have high angular resolution, the two algorithms cannot meet the requirement of radar for processing signals in real time due to the problems of large calculation amount and low calculation efficiency.
The learners provide a Synthetic Steering Vector (SVML) algorithm through long-term analysis and research, the algorithm adopts an equivalent spherical model, the earth curvature and the information of the terrain near the radar are taken into consideration, and the elevation angle of a target can be accurately estimated; in addition, the scholars have proposed an Alternative Projection (AP) algorithm which can accurately calculate the elevation angle of an object without depending on the terrain, but has problems that the amount of calculation is large and the convergence cannot be determined.
In conclusion, it can be seen that although the problem of multipath effect can be solved by the existing partial super-resolution algorithm, the target elevation angle can be accurately measured only by matching with a proper terrain; however, the ADS-B (Automatic Dependent Surveillance-Broadcasting) technology adopted at present is difficult to perform terrain correction without civil aviation and high mobility.
Disclosure of Invention
Aiming at the defects of the prior art, namely, the invention aims to provide a terrain correction method of a meter-wave radar based on a digital elevation model, which is convenient to operate, can form complementation with an ADS-B technology and can still correct the terrain without civil aviation and high mobility; in engineering application, the height and angle measurement of the meter-wave radar is more accurate by the terrain correction method of the meter-wave radar based on the digital elevation model, and the measurement accuracy is improved.
In order to achieve the technical purpose, the invention is realized by adopting the following technical scheme.
A terrain correction method for a meter-wave radar based on a digital elevation model comprises the following steps:
step 1, establishing a spherical coordinate system, wherein a meter-wave radar exists in the spherical coordinate system, and a flying target exists in a detection range of the meter-wave radar; acquiring digital elevation model data, determining ADS-B data of a flying target, and then determining an important reflection point in a spherical coordinate system;
step 2, acquiring the terminal data of the meter-wave radar, and obtaining the correlation completion result of the ADS-B data of the flying target and the terminal data of the meter-wave radar according to the ADS-B data of the flying target and the determined important reflection point;
step 3, obtaining a flight target height curve according to the correlation completion result of the ADS-B data and the terminal data;
step 4, obtaining a height measurement error curve of the meter wave radar according to the flight target height curve;
and 5, correcting the height measurement error curve of the meter-wave radar to obtain a final error correction result, wherein the final error correction result is a terrain correction result of the meter-wave radar based on the digital elevation model.
Compared with the prior art, the invention has the advantages that:
the method can form complementation with an ADS-B method, can directly correct the elevation error generated when the radar measures the target through the DEM, can correct the terrain under the conditions of no civil aviation and high maneuverability, greatly simplifies the correction process and reduces the correction cost, and has better target estimation effect than that without correction.
Drawings
The present invention will be described in further detail with reference to the accompanying drawings and specific embodiments.
FIG. 1 is a flow chart of a terrain correction method for a meter-wave radar based on a digital elevation model according to the present invention;
FIG. 2 is a schematic diagram of ADS-B coordinate transformation;
FIG. 3 is a schematic diagram illustrating formation of ground reflected waves;
FIG. 4 is a flow chart of the specific modification steps in the present invention;
FIG. 5(a) is a comparison graph of a target height and a target true height obtained by a first set of measured data through SVLM algorithm simulation;
FIG. 5(b) is a comparison graph of the target height and the target true height obtained by the second set of measured data through SVLM algorithm simulation;
FIG. 6(a) is a schematic view of SVLM altimetry error and corresponding elevation for a first set of measured data;
FIG. 6(b) is a schematic view of SVLM altimetry error and corresponding elevation for a second set of measured data;
FIG. 7(a) is a plot of the terrain corrected height measurement error for the SVLM height measurement error for the first set of measured data;
fig. 7(b) is a graph of the height measurement error of the SVLM height measurement error of the second set of measured data after terrain correction.
Detailed Description
Step 1, format conversion is carried out on ADS-B data and the ADS-B data and DEM data are stored in the same array.
Determining a meter-wave radar, wherein a flying target exists in a detection range of the meter-wave radar, and acquiring ADS-B data of the flying target, wherein the ADS-B data of the flying target is obtained after the flying process of the flying target is continuously monitored by using the meter-wave radar through a broadcast type automatic correlation monitoring technology, and the ADS-B data of the flying target comprises a flight number Num, time, longitude and latitude (omega, upsilon) and altitude real _ H of the flying target; firstly, coordinate conversion is carried out on ADS-B data of a flying target obtained through an ADS-B technology, and distance and azimuth information of the flying target when a meter-wave radar is used as a reference point are obtained.
The coordinate transformation firstly needs to establish a coordinate system as shown in fig. 2, and establishes a spherical coordinate system oxy with the geocentric O as the origin, wherein the line passing through the geocentric O and the east longitude 0 degree is taken as the x axis of the abscissa, the east direction is taken as the positive direction of the x axis, and the line passing through the geocentric O and the east longitude O is taken as the positive direction of the x axisThe connecting line of the north pole C is a y-axis, and the north direction is the positive direction of the y-axis; the actual position of the meter-wave radar in the spherical coordinate system oxy is point A, the position of the flying target is point B, and the meter-wave radar and the flying target are different two points on the spherical surface of the spherical coordinate system oxy, and the longitude and the latitude of the position of the meter-wave radar are (omega) AA ) The corresponding longitude and latitude of the position of the flying target is (omega) BB ),υ A Denotes the latitude, upsilon, of the meter wave radar B Representing the latitude of the flying target, and the value range is-90 degrees to 90 degrees; omega A Representing longitude, omega, of meter-wave radar B Representing the longitude of the flying target, and the value range is-180 degrees to 180 degrees; north latitude is positive, south latitude is negative, east longitude is positive, west longitude is negative, and the units are degrees.
Recording the spherical distance between the meter-wave radar and the flying target as L AB The distance is the distance of a flying target relative to the meter-wave radar measured by the meter-wave radar, and the radius of the earth is R; the spherical included angle between the spherical connecting line of the meter-wave radar and the flying target and the spherical connecting line of the meter-wave radar and the north pole C is recorded as the azimuth angle A of the flying target relative to the meter-wave radar, namely the azimuth angle of the flying target relative to the meter-wave radar, and the value range is [0 DEG, 360 DEG](ii) a The spherical angle between the spherical connecting line of the flying target and the meter-wave radar and the spherical connecting line of the north pole C is recorded as the azimuth angle B of the meter-wave radar, the spherical angle between the spherical connecting line of the north pole C and the meter-wave radar and the spherical connecting line of the north pole C and the flying target is recorded as the azimuth angle C of the north pole C, the geocentric angle of the spherical connecting line of the north pole C and the flying target is alpha, the geocentric angle of the spherical connecting line of the meter-wave radar and the north pole C is beta, and the geocentric angle of the spherical connecting line of the meter-wave radar and the flying target is phi.
(1a) And (4) calculation of azimuth angles:
firstly, knowing the longitude and latitude of the positions of a meter-wave radar and a flying target, and solving the azimuth angle of the flying target relative to the meter-wave radar; according to the trigonometric cosine formula:
cos(φ)=cosα·cosβ+sinα·sinβ·cosθ (1-1)
wherein, sin tableShowing a sine function and cos showing a cosine function; theta represents the dihedral angle of the plane AOC and the plane BOC in the spherical coordinate system; difference ω between longitudes of flying target by meter-wave radar BA And (6) calculating to obtain.
The known data is then substituted into the equation:
cos(φ)=cos(90-υ B )·cos(90-υ A )+sin(90-υ B )·sin(90-υ A )·cos(ω BA ) (1-2)
and further solving the sine value of phi as follows:
Figure BDA0001732603180000041
according to the spherical sine formula
Figure BDA0001732603180000042
Solving the azimuth angle A of the flight target relative to the meter wave radar as follows:
Figure BDA0001732603180000043
wherein arcsin represents an arcsine function.
Finally, considering the difference of the positions of the flying target relative to the meter-wave radar, carrying out different processing on the calculation result to obtain a final azimuth angle; the method comprises the following steps of taking a meter-wave radar as a coordinate origin, determining an azimuth angle according to the position of a flying target relative to the meter-wave radar, and representing the azimuth angle of the meter-wave radar as A _ fangwei:
if the flying target is in the first quadrant of the spherical coordinate system, namely the direction of the flying target is the northeast direction of the meter-wave radar, the azimuth angle A _ fangwei of the meter-wave radar is & lt A.
If the flying target is in the second quadrant of the spherical coordinate system, namely the direction of the flying target is in the northwest direction of the meter-wave radar, the azimuth angle A _ fangwei of the meter-wave radar is (360 degrees +/-A).
If the flying target is in the third quadrant or the fourth quadrant of the spherical coordinate system, namely the direction of the flying target is the southwest or southeast direction of the meter-wave radar, the azimuth angle A _ fangwei of the meter-wave radar is (180 degrees-A).
(1b) Calculating the distance:
calculating the spherical distance between the meter-wave radar and the flying target through the longitude and latitude of two points; the specific process is that the inverse cosine is obtained according to the formula (1-2) in the formula (1a), phi is the earth center included angle of the spherical connecting line of the meter wave radar and the flying target, and the unit is degree; then converts it into radian phi rad Finally multiplying by the earth radius R to obtain the spherical distance L between the meter-wave radar and the flying target AB The calculation formula is as follows:
Figure BDA0001732603180000051
(1c) and acquiring digital elevation model DEM data.
Acquiring Digital Elevation Model (DEM) data, wherein the DEM data is USGS (U.S. geological Survey) DEM data, the DEM data comprises longitude and latitude data and elevation data of all positions on the ground, each longitude and latitude corresponds to one elevation data, and the longitude and latitude and the elevation data in the DEM data are in one-to-one correspondence; the elevation data corresponding to the longitude and latitude can be found through the known longitude and latitude.
Because only the ground elevations corresponding to the longitude and latitude where the flying target is located are not taken into consideration in the correction process, the correction error is still large at this time, so that the following correction process takes important reflection points into consideration, wherein the important reflection points are main ground reflection points when the flying target reflected signals are received by the meter-wave radar through ground secondary reflection.
Firstly, establishing a model of an important reflection point, as shown in fig. 3, using a point A to represent a meter-wave radar, using a point A' to represent the ground position of the meter-wave radar, using a point B to represent the position of a flying target, using a point M to represent the ground surface, using a certain ground position between the meter-wave radar and the flying target as an important reflection point D, determining the position of the important reflection point D according to a formula (1-7), using elevation information used for correcting height measurement errors as elevation data of the important reflection point D, so that the longitude and latitude of the important reflection point are required to be obtained, and then using the longitude and latitude as an index to find the elevation data of the important reflection point D in DEM data.
Firstly, the distance and azimuth angle of the important reflection point D relative to the meter-wave radar are required to be calculated, the azimuth angle is the same as the azimuth angle of the flying target relative to the meter-wave radar, the distance L between the important reflection point D and the point A 'is calculated in (1a), and then the distance L between the important reflection point D and the point A' is calculated A′D The method comprises the following steps:
Figure BDA0001732603180000061
in the formula, lambda represents the signal wavelength transmitted by the meter-wave radar, and the distance between the important reflection point D and the meter-wave radar is only related to the height of the antenna frame and the signal wavelength transmitted by the meter-wave radar; l is A′D Namely the distance between the important reflection point D and the point A', and the azimuth of the important reflection point D is the same as the azimuth angle A _ fangwei of the meter-wave radar; at the moment, the longitude and latitude where the meter wave radar is located are known to be (omega) AA ) And the distance L of the significant reflection point D from the point A A′D And the azimuth angle A _ fangwei and R of the meter-wave radar represent the radius of the earth, and the longitude and latitude of the important reflection point D are calculated through the following formula, and the principle is the same as the calculation of the azimuth angle in the step (1 a).
Firstly, calculating an included angle D between a connecting line of the meter wave radar and the geocenter O and a connecting line of the important reflection point D and the geocenter O:
Figure BDA0001732603180000062
then, solving an included angle a between a connecting line of the north pole C and the geocenter O of the earth and a connecting line of the important reflection point D and the geocenter O:
a=arccos(cos(90-υ A )·cos(d)+sin(90-υ A )·sin(d)·cos(A_fangwei)) (1-9)
wherein arccos represents an inverse cosine function.
And finally, solving a spherical angle xi of the meter wave radar, the north pole C and the important reflection point D at the north pole C:
Figure BDA0001732603180000063
obtaining the longitude and latitude (upsilon) of the important reflection point D DD ) Comprises the following steps:
υ D =υ A +ξ (1-11)
ω D =90-a (1-12)
wherein, ω is D Longitude, upsilon, representing the significant reflection point D D Representing the latitude of the important reflection point D; and after the longitude and latitude of the important reflection point D are obtained, the elevation information corresponding to the important reflection point D is searched in the digital elevation model data through the longitude and latitude.
In summary, the detailed data obtained by the ADS-B comprises a flight number Num of the flight target, time, longitude and latitude (omega, upsilon) and altitude real _ H, and the azimuth A _ fangwei and distance information L of the flight target relative to the meter wave radar at a certain moment are obtained by converting the data obtained by the ADS-B technology through (1a) and (1B) AB And (1c) taking the important reflection point into account to obtain the longitude and latitude (upsilon) of the important reflection point D DD ) And searching the digital elevation model DEM data by taking the digital elevation model DEM data as an index to obtain a digital elevation, namely the elevation data delta h of the important reflection point D.
Processing and converting ADS-B data of N flying targets at different moments according to (1a), (1B) and (1c), and arranging flying target information into L according to the sequence of distance, azimuth angle, altitude (namely the real height of the flying target), time and digital elevation AB (i),A_fangwei(i),real_H(i),time(i),△h(i)]For use in step 2; wherein L is AB (i) Representing a spherical distance curve between the meter-wave radar and the flying target at the time (i), A _ fangwei (i) representing an azimuth angle curve of the meter-wave radar at the time (i), real _ H (i) representing a real altitude curve of the flying target at the time (i), and delta h (i) representing the time(i) An elevation data curve of the time important reflection point D, wherein i is 1,2,3, N, N represents the total number of different times of ADS-B data of the selected flight target.
And 2, matching the ADS-B data with the terminal data.
After format conversion is carried out on the ADS-B data of the flying target in the step 1, the flying target information is arranged and placed into [ L ] in the ADS-B data of the flying target according to the sequence of distance, azimuth, altitude and time AB (i),A_fangwei(i),real_H(i),time(i),△h(i)]N, N represents the total number of different times of the ADS-B data of the selected flight targets.
2.1, acquiring the terminal data of the meter-wave radar, wherein the terminal data of the meter-wave radar selects M time ' (1) for the meter-wave radar terminal, and the M time ' (M) respectively measures the spherical distance between the meter-wave radar and the flying target and the azimuth angle of the meter-wave radar to respectively obtain the spherical distance L between the meter-wave radar and the flying target measured at the time of the time ' (1) AB The spherical distance L between the meter-wave radar and the flying target measured from the moment (1) to the moment (M) AB ' (M) and the azimuth angle a _ fangwei ' (1) of the meter-wave radar measured at time ' (1) to the azimuth angle a _ fangwei ' (M) of the meter-wave radar measured at time ' (M).
And placing the meter-wave radar terminal data according to the same sequence and format to obtain:
[L AB ′(j),A_fangwei′(j),time′(j)],L AB ' (j) represents the spherical distance between the meter-wave radar and the flying target at the time of time ' (j), A _ fangwei ' (j) represents the azimuth angle curve of the meter-wave radar, j is 1,2, 3.
2.2 in order to synchronize the ADS-B data of the flight target and the terminal data of the meter-wave radar, the ADS-B data and the terminal data of the meter-wave radar need to be correlated; the specific correlation method is a nearest neighbor method, a proper screening threshold is set, proper echo data are screened according to the screening threshold, and certain probability of a receiving terminal is required to be ensured when the screening threshold is selected.
Specifically, the azimuth angle and the distance are selectedThe value is a measured value, and the measured value of the flying target detected by the meter-wave radar at the time' (j) is recorded as Z j ,Z j =[L AB ′(j),A_fangwei′(j)](j 1, 2.. said., M), and the measured value of the ADS-B data of the flight target at the time (j) is recorded as P j ,P j =[L AB (j),A_fangwei(j)](j=1,2,...,M),M≤N。
Then calculate Z j And P j Is the Euclidean distance S j
Figure BDA0001732603180000081
Z j (k) Measurement value Z representing flying target detected by meter-wave radar at time' (j) j The k-th element of (1), P j (k) Measurement value P of ADS-B data representing flight target at time' (j) j The kth element, k ═ 1, 2.
2.3 treatment of S j Comparing with a screening threshold gamma (taking an empirical value of 10), if S j <Gamma, then Z j And P j Correlating, and measuring the value Z of the flying target detected by the meter-wave radar at the time' (j) j And obtaining associated data of the jth ' effective measurement value, wherein the associated data of the jth ' effective measurement value is meter wave radar terminal data corresponding to the jth ' effective measurement value and ADS-B data corresponding to the jth ' effective measurement value, and the meter wave radar terminal data corresponding to the jth ' effective measurement value is [ L ] AB ′(j'),A_fangwei′(j'),time′(j')]The ADS-B data corresponding to the jth effective measurement value is [ L ] AB (j'),A_fangwei(j'),real_H(j'),time(j'),△h(j')]The initial value of j 'is 1 and then the value of j' is incremented by 1.
If S j ≥γ,Z j And P j Can not be correlated, the measured value Z of the flying target detected by the meter-wave radar at the time (j) j Corresponding meter wave radar terminal data L AB (j '), A _ fangwei (j '), real _ H (j '), time (j '), Δ H (j ') are discarded, the correlation step is completed.
2.4, repeatedly executing 2.2 and 2.3 to obtain the correlation data from the 1 st effective measurement value to the mth effective measurement value, and recording the correlation data as the correlation completion result of the ADS-B data of the flight target and the terminal data of the meter wave radar; wherein, j' is 1,2,3, and M is less than or equal to M.
And 3, calculating the correlation completion result of the ADS-B data and the terminal data by adopting a terrain-sensitive SVML algorithm to obtain a time ' (j) flying target height curve SVLM _ H (j '), wherein j ' is 1,2, 3.
And 4, the altitude in the ADS-B data is the real height of the flying target, and for the flying target measured by the meter-wave radar, the result obtained by subtracting the real height curve real _ H (j ') of the flying target at the time (j) from the real height curve SVLM _ H (j') of the flying target at the time (j) calculated in the step 3 is the height Error curve Error (Error j '), j' is 1,2, 3.
And 5, performing terrain correction according to the relation between an elevation data curve delta h (j '), j' 1,2,3, and the j 'of the meter wave radar height measurement Error curve Error (j'), j '1, 2,3, and j', m, and the elevation data curve delta h (j ') of the important reflection point D at the time of time' (j).
According to the correction method based on the DEM, an elevation data curve delta h (j '), j' is 1,2,3, a. Therefore, according to the similarity of DEM data and the change characteristics of the error curve, the whole correction process is segmented, points with larger curve changes are divided, a proper change threshold value is selected, when the curve changes and exceeds the change threshold value, the search is stopped, and the section with smaller fluctuation is corrected; and after the correction is completed, the backward search is continued until all searches are completed.
Specific correction process referring to fig. 4, fig. 4 is a flowchart of a specific correction step in the present invention, which includes the following sub-steps:
(5a) first, a maximum value Δ h _ max ═ max [ [ Δ h (j ') ], (j' ═ 1,2,... multidot.m), ] and a minimum value Δ h _ min [ [ Δ h (j ') ], (j'. 1,2,. multidot.m) ] of the elevation data and an average value of all the elevation data are calculated
Figure BDA0001732603180000091
(5b) An appropriate threshold value of 2000 is selected according to actual conditions and simulation experience.
Determining two different time moments from a time ' (j) time moment metric wave radar height measurement Error curve Error (j '), j ' is 1,2,3 1 And search point x 2 ,x 1 =1,2,3,...,m,x 2 =1,2,3,...,m,x 1 Is the 1 st time, x, of the m times 2 Is the 2 nd time of the m times.
(5c) Computing a search point x 1 And search point x 2 Is equal to abs (Error (x) 2 )-Error(x 1 ))。
(5d) Searching point x 1 And search point x 2 Is compared with a selected threshold value, if d is<value, then let x 2 Adds 1 to the value of (1), returns (5 c);
if d is more than or equal to value, beginning to perform k correction, and firstly calculating a search point x after the k correction 1 To search point x 2 Corresponding to the average value delta h _ mean _ k of the DEM data, the calculation expression is as follows:
Figure BDA0001732603180000092
where Δ h (j ") represents the elevation data of the significant reflection point D at the time j".
Then, the digital elevation information of the target and the previously calculated delta h _ max, delta h _ min, delta h _ mean are processed on the xth-th position 1 Elevation data Δ h (x) of the temporally significant reflection point D 1 ) And x 2 Elevation data Δ h (x) of the temporally significant reflection point D 2 ) Comparing, selecting the value closest to delta h _ mean _ k, and recording as delta h _ suit k Then for the search point x 1 To search point x 2 The data in the segment adopts delta h _ suit k The value is corrected for the k time to obtain a search point x 1 To search point x 2 The k-th correction result of the k-th data
Figure BDA0001732603180000101
Figure BDA0001732603180000102
And k is 1, and a correction result is stored, so that correction is performed mainly based on the characteristic that elevation data of the terrain within a certain range has small difference and smooth change, and then (5e) is performed.
(5e) Judgment of x 2 Whether or not it is equal to m; if not, let x 1 =x 2 +1 and adding 1 to the value of k, returning to (5 c); otherwise, the search is complete, now x is obtained 1 To x 2 The 1 st correction result of the 1 st data
Figure BDA0001732603180000103
To x 1 To x 2 The k-th correction result of the k-th data
Figure BDA0001732603180000104
Then will be
Figure BDA0001732603180000105
Splicing in sequence to obtain the final error correction result, i.e. [ result 1 ,result 2 ,...,result m ],k≤m。
Wherein the content of the first and second substances,
Figure BDA0001732603180000106
in (1)
Figure BDA0001732603180000107
Corresponding result 1
Figure BDA0001732603180000108
In (1)
Figure BDA0001732603180000109
Corresponding result m And finishing the whole operation process, wherein the final error correction result is a terrain correction result of the meter wave radar based on the digital elevation model.
The effect of the present invention will be further described with reference to the simulation diagram.
Simulating conditions and processes.
The simulation tool is Matlab 2014a, the simulation data are measured data, and the measured data are simulated by an SVML algorithm to obtain the flying target height; when the ADS-B data is converted, format conversion is carried out by taking the longitude and latitude where the millimeter wave radar is located as the origin of coordinates and the earth equivalent radius of 6378137; and then, carrying out difference on the experimental result obtained by the SVLM and the real target height obtained by the ADS-B data to obtain a height measurement error, correcting the height measurement error according to DEM data segments, wherein the correction method is that the segment difference is carried out on the height measurement error and the DEM data, selecting a proper change threshold value in the simulation process, stopping backward search when the curve fluctuation exceeds the change threshold value, and correcting the segment with smaller fluctuation. And after the correction is completed, continuing to search backwards until all searches are completed.
(II) simulation results and analysis
Fig. 5(a) is a comparison graph of a target height and a target true height obtained by a first group of measured data through simulation of an SVLM algorithm, and fig. 5(b) is a comparison graph of a target height and a target true height obtained by a second group of measured data through simulation of an SVLM algorithm, and it can be seen from the measurement results of the two groups of measured data that the error of the target height measured by the SVLM algorithm is large and certain correction is required; FIG. 6(a) is a comparison graph of height measurement errors of a first set of measured data SVLM and corresponding elevations, FIG. 6(b) is a comparison graph of height measurement errors of a second set of measured data SVLM and corresponding elevations, it can be seen from FIG. 6(a) that the change processes of DEM elevation data and SVML algorithm height measurement errors have certain similarities, when a large change occurs in an elevation curve, the height measurement errors are obviously increased, and FIG. 6(b) has the same phenomenon, that is, the terrain elevation change is a main factor affecting the height measurement errors; fig. 7(a) is a curve obtained by terrain correction of SVLM height measurement errors of a first set of measured data, fig. 7(b) is a curve obtained by terrain correction of SVLM height measurement errors of a second set of measured data, and it can be seen from fig. 7(a) that the SVLM height measurement errors obtained after terrain correction are significantly smaller than the errors directly measured by using the SVML algorithm, the maximum error is reduced by half from nearly 6 km before correction, and the corresponding heights of all target points are corrected, and fig. 7(b) also has the same characteristics, thus proving the effectiveness of the terrain correction of the SVML height measurement errors.
The above description is only for the specific embodiments of the present invention, but the scope of the present invention is not limited thereto, and any person skilled in the art can easily conceive of the changes or substitutions within the technical scope of the present invention, and all the changes or substitutions should be covered within the scope of the present invention. Therefore, the protection scope of the present invention shall be subject to the protection scope of the claims.

Claims (1)

1. A meter-wave radar terrain correction method based on a digital elevation model is characterized by comprising the following steps:
step 1, establishing a spherical coordinate system, wherein a meter-wave radar exists in the spherical coordinate system, and a flying target exists in a detection range of the meter-wave radar; acquiring digital elevation model data, determining ADS-B data of a flying target, and then determining an important reflection point in a spherical coordinate system;
step 2, acquiring the terminal data of the meter-wave radar, and obtaining the correlation completion result of the ADS-B data of the flying target and the terminal data of the meter-wave radar according to the ADS-B data of the flying target and the determined important reflection point;
step 3, obtaining a flight target height curve according to the correlation completion result of the ADS-B data and the terminal data;
step 4, obtaining a height measurement error curve of the meter wave radar according to the flight target height curve;
step 5, correcting the height measurement error curve of the meter-wave radar to obtain a final error correction result, wherein the final error correction result is a terrain correction result of the meter-wave radar based on a digital elevation model;
in step 1, the spherical coordinate system specifically includes:
a spherical coordinate system oxy which takes a connecting line passing through the geocenter and the east longitude 0 degrees as an x axis of a horizontal coordinate, takes the east direction as the positive direction of the x axis, takes a connecting line passing through the geocenter and the north pole as a y axis and takes the north direction as the positive direction of the y axis;
the digital elevation model data is USGSDEM data, the digital elevation model DEM data comprises longitude and latitude data and elevation data of all positions on the ground, each longitude and latitude corresponds to one elevation data, and the longitude and latitude and the elevation data in the digital elevation model data correspond to one another;
the ADS-B data of the flying target is specifically data obtained after the flying process of the flying target is continuously monitored by using a meter-wave radar through a broadcast type automatic correlation monitoring technology, and comprises a flight number Num, time, longitude and latitude (omega, upsilon) and a real height real _ H of the flying target;
the important reflection point is determined by the following process:
first, the distance between the important reflection point and the point A' is recorded as L A′D
Figure FDA0003623660010000011
The point A' represents the ground position of the meter-wave radar, the lambda represents the signal wavelength transmitted by the meter-wave radar, and the h represents the height of the meter-wave radar to the ground;
then, calculating the longitude and latitude (upsilon) of the important reflection point DD ) Comprises the following steps:
υ D =υ A
ω D =90-a
wherein upsilon is A Which represents the latitude of the meter-wave radar,
Figure FDA0003623660010000012
r represents the earth radius, A _ fangwei represents the azimuth angle of the meter-wave radar, ω D Longitude, upsilon, representing the significant reflection point D D Denotes the latitude of the significant reflection point, sin denotes positiveChord function, arcsin represents an arcsine function; a-arccos (cos (90-upsilon) A )·cos(d)+sin(90-υ A ) Sin (d) cos (A _ fangwei)), cos representing a cosine function, arccos representing an inverse cosine function, upsilon A Representing the latitude of the meter-wave radar;
through the longitude and latitude (upsilon) of the important reflection point DD ) Searching the digital elevation model data to obtain digital elevation, namely the elevation data delta h of the important reflection point;
obtaining the longitude and latitude of the important reflection point and the distance L of the important reflection point relative to the point A A′D Determining the important reflection points in the spherical coordinate system after the longitude and latitude of the important reflection points and the elevation data delta h of the important reflection points;
the A _ fangwei represents the azimuth angle of the meter-wave radar, and the determination process is as follows:
setting the meter-wave radar and the flying target as two different points on the spherical surface of the spherical coordinate system oxy, and recording the spherical distance between the meter-wave radar and the flying target as L AB
Figure FDA0003623660010000021
If the flying target is in a first quadrant of a spherical coordinate system, namely the direction of the flying target is the northeast direction of the meter-wave radar, the azimuth angle A _ fangwei of the meter-wave radar is less than A;
if the flying target is in the second quadrant of the spherical coordinate system, namely the direction of the flying target is the northwest direction of the meter-wave radar, the azimuth angle A _ fangwei of the meter-wave radar is (360 degrees +/-A);
if the flying target is in the third quadrant or the fourth quadrant of the spherical coordinate system, namely the direction of the flying target is the southwest or southeast direction of the meter-wave radar, the azimuth angle A _ fangwei of the meter-wave radar is (180 degrees-A);
wherein R represents the earth radius, and angle A represents the azimuth angle of the flying target relative to the meter-wave radar,
Figure FDA0003623660010000022
ω A meter wave radarLongitude, ω of B Representing the longitude, upsilon, of the flight object B Representing the latitude of the flight target; phi represents the included angle of the geocentric of the spherical connecting line of the meter-wave radar and the flying target, and satisfies the following conditions: cos (phi) ═ cos (90-upsilon) B )·cos(90-υ A )+sin(90-υ B )·sin(90-υ A )·cos(ω BA ),ω A Representing the longitude of the meter-wave radar;
the substep of step 2 is:
2.1 determining ADS-B data arrangement of N flying targets at different moments;
[L AB (i),A_fangwei(i),real_H(i),time(i),Δh(i)]for use in step 2; wherein L is AB (i) Representing a spherical distance curve between a time meter-wave radar and a flying target, A _ fangwei (i) representing an azimuth angle curve of the time meter-wave radar, real _ H (i) representing a real height curve of the time meter-wave radar, delta h (i) representing an elevation data curve of an important reflection point D at the time (i), wherein i is 1,2,3,.
The method comprises the steps of obtaining meter-wave radar terminal data, wherein the meter-wave radar terminal data are obtained by selecting M time ' (1) for a meter-wave radar terminal, and measuring the spherical distance between a meter-wave radar and a flying target and the azimuth angle of the meter-wave radar respectively through the time ' (M) to obtain the spherical distance L between the meter-wave radar and the flying target measured at the time of the time ' (1) respectively AB The spherical distance L between the meter-wave radar and the flying target measured from the moment (1) to the moment (M) AB ' (M) and azimuth angle a _ fangwei ' (1) of the meter-wave radar measured at time ' (1) to azimuth angle a _ fangwei ' (M) of the meter-wave radar measured at time ' (M);
then, the meter-wave radar terminal data are placed to obtain:
[L AB ′(j),A_fangwei′(j),time′(j)],L AB ' (j) represents the spherical distance between the meter-wave radar and the flying target at the time of time ' (j), A _ fangwei ' (j) represents the azimuth angle curve of the meter-wave radar, j is 1,2,3, and M represents the total number of different times selected for acquiring the meter-wave radar terminal data;
2.2 recording the measured value of the flying target detected by the meter-wave radar at the time (j) as Z j
Z j =[L AB ′(j),A_fangwei′(j)](j 1, 2.. said., M), and the measured value of the ADS-B data of the flight target at the time (j) is recorded as P j ,P j =[L AB (j),A_fangwei(j)](j=1,2,...,M),M≤N;
Then calculate Z j And P j Is the Euclidean distance S j
Figure FDA0003623660010000031
Z j (k) Measurement value Z representing flying target detected by meter-wave radar at time' (j) j The k-th element of (1), P j (k) Measurement value P of ADS-B data representing flight target at time' (j) j The kth element, k ═ 1, 2;
2.3 treatment of S j Comparing with a set screening threshold gamma if S j < gamma, then Z j And P j Correlating, and measuring the value Z of the flying target detected by the meter-wave radar at the time' (j) j And obtaining associated data of the jth ' effective measurement value, wherein the associated data of the jth ' effective measurement value is meter wave radar terminal data corresponding to the jth ' effective measurement value and ADS-B data corresponding to the jth ' effective measurement value, and the meter wave radar terminal data corresponding to the jth ' effective measurement value is [ L ] AB ′(j'),A_fangwei′(j'),time′(j')]The ADS-B data corresponding to the jth effective measurement value is [ L ] AB (j'),A_fangwei(j'),real_H(j'),time(j'),Δh(j')]The initial value of j 'is 1, then the value of j' is added with 1;
if S j ≥γ,Z j And P j Can not be correlated, the measured value Z of the flying target detected by the meter-wave radar at the time (j) j Corresponding meter wave radar terminal data L AB (j '), A _ fangwei (j '), real _ H (j '), time (j '), Δ H (j ') is discarded;
2.4, repeatedly executing 2.2 and 2.3 to obtain the correlation data from the 1 st effective measurement value to the mth effective measurement value, and recording the correlation data as the correlation completion result of the ADS-B data of the flight target and the terminal data of the meter wave radar; wherein j' is 1,2,3, a, M, M is less than or equal to M;
in step 3, the flying target altitude curve, specifically the flying target altitude curve at time' (j), is obtained by the following process:
calculating the association completion result of the ADS-B data and the terminal data by adopting a terrain-sensitive SVML algorithm to obtain a time ' (j) flying target height curve SVLM _ H (j '), wherein j ' is 1,2,3,. and m;
in step 4, the height measurement error curve of the meter-wave radar, specifically the height measurement error curve of the meter-wave radar at the time of time' (j), is obtained by the following steps:
a result obtained by subtracting the time ' (j) time flying target height curve SVLM _ H (j '), j ' is 1,2,3, and is the time (j) time meter-wave radar height measurement Error curve Error (j '), j ' is 1,2,3, and.
The substep of step 5 is:
(5a) first, a maximum elevation data value Δ h _ max ═ max [ Δ h (j ') ], (j' 1,2,.. multidot.m), a minimum elevation data value Δ h _ min ═ min [ Δ h (j ') ], (j' 1,2,. multidot.m), and an average value of all elevation data are calculated
Figure FDA0003623660010000041
(5b) Setting a threshold value;
determining two different time moments from a time ' (j) time moment metric wave radar height measurement Error curve Error (j '), j ' is 1,2,3 1 And search point x 2 ,x 1 =1,2,3,...,m,x 2 =1,2,3,...,m,x 1 Is the 1 st time, x, of the m times 2 Is the 2 nd time among the m times;
(5c) computing a search point x 1 And search point x 2 Is equal to abs (Error (x) 2 )-Error(x 1 ));
(5d) Searching point x 1 And search point x 2 Is compared with a threshold value, if d < value, let x 2 Adds 1 to the value of (c) and returns to (5 c);
if d is more than or equal to value, beginning to perform k correction, and firstly calculating a search point x after the k correction 1 To search point x 2 Corresponding to the average value delta h _ mean _ k of the DEM data, the calculation expression is as follows:
Figure FDA0003623660010000051
wherein, Δ h (j ") represents the elevation data of the important reflection point D at the j' th moment;
then, the digital elevation information of the delta h _ mean _ k, the delta h _ max, the delta h _ min, the delta h _ mean and the target is processed in the xth 1 Altitude data Δ h (x) of the time-critical reflection point D 1 ) And x 2 Elevation data delta h (x) of the time critical reflection point D 2 ) Comparing, and selecting the value closest to delta h _ mean _ k, and marking as delta h _ suit k Then for the search point x 1 To search point x 2 Data in this segment is represented by Δ h _ suit k The value is corrected for the k time to obtain a search point x 1 To search point x 2 The k-th correction result of the k-th data
Figure FDA0003623660010000052
Figure FDA0003623660010000053
The initial value of k is 1, and then (5e) is executed;
(5e) judgment of x 2 Whether or not it is equal to m; if not, let x 1 =x 2 +1 and adding 1 to the value of k, returning to (5 c); otherwise, the search is complete, now x is obtained 1 To x 2 The 1 st correction result of the 1 st data
Figure FDA0003623660010000054
To x 1 To x 2 The k-th correction result of the k-th data
Figure FDA0003623660010000055
Then will be
Figure FDA0003623660010000056
Splicing in sequence to obtain the final error correction result, i.e. [ result 1 ,result 2 ,...,result m ],k≤m;
Wherein the content of the first and second substances,
Figure FDA0003623660010000057
in (1)
Figure FDA0003623660010000058
Corresponding result 1
Figure FDA0003623660010000059
In (1)
Figure FDA00036236600100000510
Corresponding result m
CN201810781175.9A 2018-07-17 2018-07-17 Meter wave radar terrain correction method based on digital elevation model Active CN108983169B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810781175.9A CN108983169B (en) 2018-07-17 2018-07-17 Meter wave radar terrain correction method based on digital elevation model

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810781175.9A CN108983169B (en) 2018-07-17 2018-07-17 Meter wave radar terrain correction method based on digital elevation model

Publications (2)

Publication Number Publication Date
CN108983169A CN108983169A (en) 2018-12-11
CN108983169B true CN108983169B (en) 2022-08-02

Family

ID=64549750

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810781175.9A Active CN108983169B (en) 2018-07-17 2018-07-17 Meter wave radar terrain correction method based on digital elevation model

Country Status (1)

Country Link
CN (1) CN108983169B (en)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110007336A (en) * 2019-04-04 2019-07-12 首都师范大学 A method of earthquake is monitored based on Law of DEM Data
CN110633447B (en) * 2019-08-30 2020-09-04 中国人民解放军军事科学院国防科技创新研究院 Spherical distance fixed-point calculation method based on FPGA and calculation device thereof
CN112098953A (en) * 2020-09-21 2020-12-18 中国人民解放军63921部队 Rapid iteration method and device for calculating atmospheric refraction correction
CN113189562B (en) * 2021-07-02 2021-09-07 成都众享天地网络科技有限公司 Terrain detection algorithm based on elevation

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101482616A (en) * 2008-08-13 2009-07-15 中国科学院电子学研究所 Topographic survey method
CN105466391A (en) * 2015-11-18 2016-04-06 中国能源建设集团江苏省电力设计院有限公司 A tower base section generating method by utilization of a digital elevation model and field-data correction
CN105954746A (en) * 2016-04-29 2016-09-21 西安电子科技大学 Landform correction meter wave radar height measurement method based on broadcast automatic mutual supervisory signals
CN206223983U (en) * 2016-11-22 2017-06-06 无锡孚嘉航海科技有限公司 A kind of radar object height estimating system
CN106990401A (en) * 2017-05-24 2017-07-28 武汉大学 Based on the class vertical error modification method of Full wave shape airborne laser radar data two

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101482616A (en) * 2008-08-13 2009-07-15 中国科学院电子学研究所 Topographic survey method
CN105466391A (en) * 2015-11-18 2016-04-06 中国能源建设集团江苏省电力设计院有限公司 A tower base section generating method by utilization of a digital elevation model and field-data correction
CN105954746A (en) * 2016-04-29 2016-09-21 西安电子科技大学 Landform correction meter wave radar height measurement method based on broadcast automatic mutual supervisory signals
CN206223983U (en) * 2016-11-22 2017-06-06 无锡孚嘉航海科技有限公司 A kind of radar object height estimating system
CN106990401A (en) * 2017-05-24 2017-07-28 武汉大学 Based on the class vertical error modification method of Full wave shape airborne laser radar data two

Also Published As

Publication number Publication date
CN108983169A (en) 2018-12-11

Similar Documents

Publication Publication Date Title
CN108983169B (en) Meter wave radar terrain correction method based on digital elevation model
CN110426690B (en) Automatic calibration method for airborne weather radar beam pointing
CN108614258B (en) Underwater positioning method based on single underwater sound beacon distance measurement
CN103926572B (en) A kind of clutter suppression method of battle array radar self adaptation subspace, airborne anon-normal side
CN108872932B (en) Beyond-visual-range target direct positioning result deviation rectifying method based on neural network
CN110646782B (en) Satellite-borne laser on-orbit pointing calibration method based on waveform matching
CN102914776A (en) Multichannel SAR (synthetic aperture radar) mobile object localization method on the basis of fuzzy-c-mean algorithm
CN108845325A (en) Towed linear-array sonar submatrix error misfits estimation method
CN110703259B (en) Underwater acoustic array channel phase consistency calibration method based on moving sound source
CN109507635A (en) Utilize the array amplitude phase error evaluation method of two unknown orientation auxiliary sources
CN109856605A (en) A kind of while formation of the digital multiple beam quadratic fit curve is directed toward modification method
CN110488283B (en) Error correction method for multi-channel HRWS-SAR channel
CN111199280B (en) Multi-station target source geographic coordinate estimation method combining signal complex envelope and carrier phase information in presence of short wave channel model error
CN109738902B (en) High-precision autonomous acoustic navigation method for underwater high-speed target based on synchronous beacon mode
CN108710111A (en) A kind of two-dimentional space-variant bearing calibration of airborne biradical Forward-looking SAR orientation phase
CN109856616B (en) Method for correcting error of radar positioning relative system
CN107861096A (en) Least square direction-finding method based on voice signal reaching time-difference
CN110132281A (en) A kind of autonomous acoustic navigation method of underwater high-speed target with high precision based on inquiry answer-mode
CN108398659A (en) A kind of Wave arrival direction estimating method that pencil of matrix is combined with rooting MUSIC
CN104515974A (en) Processing method of microwave landing airborne equipment angle and ranging data
CN110208741B (en) Beyond-visual-range single target direct positioning method based on multi-circle array phase measurement
Yu In-situ calibration of transceiver alignment for a high-precision USBL system
CN117146830A (en) Self-adaptive multi-beacon dead reckoning and long-baseline tightly-combined navigation method
CN115826004B (en) Three-star cooperative direct positioning method based on two-dimensional angle and time difference combination
CN102062851A (en) Direction finding method based on improved L array star-carrying broadband multipurpose

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant