Background
Synthetic Aperture Radars (SAR) can realize all-day, all-weather and large-area ground observation and high-resolution imaging, and have important application value and potential in the fields of military, economy, environment and the like. Common imaging algorithms for synthetic aperture radar are the Range-Doppler (RD) and the transform-Chirp Scaling (CS) algorithms. The CS algorithm starts from the echo directly, accurately deduces the expression of the echo signal in the range-Doppler domain, and avoids the defect that the distance migration correction interpolation calculation amount is increased greatly under the condition that the synthetic aperture time is longer in the RD algorithm.
As shown in fig. 1, the transform chirp scaling algorithm includes: inputting original data, transposing distance and azimuth data, transforming azimuth Fourier, transforming distance Fourier, inverse transforming distance Fourier, transposing azimuth distance data, inverse transforming azimuth Fourier, and outputting complex images. The whole flow of the algorithm is as follows: the Fourier transform is firstly carried out in the azimuth direction, the original data is transformed to a range-Doppler domain, and the range-displacement curves of all range gates are corrected to a certain preselected reference curve by carrying out complex multiplication with a Chirp Scaling (CS) factor, so that the range-displacement curves of all range gates are the same. Then, distance Fourier transform is carried out, the data are transformed to a two-dimensional frequency domain, complex multiplication of distance compensation factors is carried out on the data, and distance displacement correction and distance compression are completed; and performing inverse Fourier transform on the range to obtain a range Doppler domain, performing complex multiplication on the range Doppler domain and the range Doppler domain, performing azimuth compression, and performing inverse Fourier transform on the azimuth compressed range Doppler domain to obtain an SAR image.
It can be seen from the above flow: a CS compensation factor needs to be added between the azimuth fourier transform and the range fourier transform, a distance compensation factor needs to be added between the range fourier transform and the range fourier inverse transform, and an azimuth compensation factor needs to be added between the azimuth distance data transpose and the azimuth fourier inverse transform.
The calculation of the above three factors has already been known as a calculation formula, which is as follows: CS compensation factor: phi1(τ,f,Rref)=exp{-jπBr0(f,Rref)Cs(f)[τ-τref(f)]2}(1)
Distance compensation factor: <math><mrow>
<msub>
<mi>Φ</mi>
<mn>2</mn>
</msub>
<mrow>
<mo>(</mo>
<msub>
<mi>f</mi>
<mi>c</mi>
</msub>
<mo>,</mo>
<mi>f</mi>
<mo>,</mo>
<msub>
<mi>R</mi>
<mi>ref</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>exp</mi>
<mrow>
<mo>{</mo>
<mi>j</mi>
<mfrac>
<mrow>
<mi>π</mi>
<msubsup>
<mi>f</mi>
<mi>c</mi>
<mn>2</mn>
</msubsup>
</mrow>
<mrow>
<msub>
<mi>B</mi>
<mrow>
<mi>r</mi>
<mn>0</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>,</mo>
<msub>
<mi>R</mi>
<mi>ref</mi>
</msub>
<mo>)</mo>
</mrow>
<mrow>
<mo>[</mo>
<mn>1</mn>
<mo>+</mo>
<msub>
<mi>C</mi>
<mi>s</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>]</mo>
</mrow>
</mrow>
</mfrac>
<mo>}</mo>
</mrow>
<mo>·</mo>
<mi>exp</mi>
<mrow>
<mo>{</mo>
<mi>j</mi>
<mfrac>
<mrow>
<mn>4</mn>
<mi>π</mi>
</mrow>
<mi>c</mi>
</mfrac>
<msub>
<mi>f</mi>
<mi>τ</mi>
</msub>
<msub>
<mi>R</mi>
<mi>ref</mi>
</msub>
<msub>
<mi>C</mi>
<mi>s</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>}</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mo>)</mo>
</mrow>
</mrow></math>
azimuth compensation factor: <math><mrow>
<msub>
<mi>Φ</mi>
<mn>3</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>τ</mi>
<mo>,</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>exp</mi>
<mrow>
<mo>{</mo>
<mo>-</mo>
<mi>j</mi>
<mfrac>
<mrow>
<mn>4</mn>
<mi>πR</mi>
</mrow>
<mi>λ</mi>
</mfrac>
<mrow>
<mo>[</mo>
<mn>1</mn>
<mo>-</mo>
<mi>sin</mi>
<mi>φ</mi>
<mo>·</mo>
<mi>D</mi>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>]</mo>
</mrow>
<mo>+</mo>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mi>j</mi>
<mfrac>
<mrow>
<mn>2</mn>
<mi>π</mi>
<mi>Rf</mi>
<mi>cos</mi>
<mi>φ</mi>
</mrow>
<mi>v</mi>
</mfrac>
<mo>)</mo>
</mrow>
<mo>+</mo>
<mi>jΘ</mi>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>}</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>3</mn>
<mo>)</mo>
</mrow>
</mrow></math>
wherein,
<math><mrow>
<msub>
<mi>C</mi>
<mi>s</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mrow>
<mi>sin</mi>
<mi>φ</mi>
</mrow>
<msqrt>
<mn>1</mn>
<mo>-</mo>
<msup>
<mrow>
<mo>(</mo>
<mfrac>
<mi>λf</mi>
<mrow>
<mn>2</mn>
<mi>v</mi>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</msqrt>
</mfrac>
<mo>-</mo>
<mn>1</mn>
</mrow></math>
<math><mrow>
<msub>
<mi>B</mi>
<mrow>
<mi>r</mi>
<mn>0</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>;</mo>
<msub>
<mi>R</mi>
<mi>ref</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mi>b</mi>
<mrow>
<mn>1</mn>
<mo>+</mo>
<mi>b</mi>
<msub>
<mi>R</mi>
<mi>ref</mi>
</msub>
<mi>sin</mi>
<mi>φ</mi>
<mfrac>
<mrow>
<mn>2</mn>
<mi>λ</mi>
</mrow>
<msup>
<mi>c</mi>
<mn>2</mn>
</msup>
</mfrac>
<mfrac>
<msup>
<mrow>
<mo>(</mo>
<mfrac>
<mi>λf</mi>
<mrow>
<mn>2</mn>
<mi>v</mi>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<msup>
<mrow>
<mo>[</mo>
<mn>1</mn>
<mo>-</mo>
<msup>
<mrow>
<mo>(</mo>
<mfrac>
<mi>λf</mi>
<mrow>
<mn>2</mn>
<mi>v</mi>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>]</mo>
</mrow>
<mfrac>
<mn>3</mn>
<mn>2</mn>
</mfrac>
</msup>
</mfrac>
</mrow>
</mfrac>
</mrow></math>
<math><mrow>
<msub>
<mi>τ</mi>
<mi>ref</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mn>2</mn>
<mi>c</mi>
</mfrac>
<msub>
<mi>R</mi>
<mi>ref</mi>
</msub>
<mrow>
<mo>[</mo>
<mn>1</mn>
<mo>+</mo>
<msub>
<mi>C</mi>
<mi>s</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>]</mo>
</mrow>
</mrow></math>
Rreffor reference distances, the image center point is usually used;
phi is an included angle between the radar wave velocity center and the platform motion direction of the radar;
b is the transmission pulse modulation frequency;
c is the speed of light;
λ is the radar emission wavelength;
v is the radar ground speed;
the calculation formula of the variables is as follows:
<math><mrow>
<mi>v</mi>
<mo>=</mo>
<msqrt>
<mfrac>
<mrow>
<mi>λ</mi>
<msub>
<mi>R</mi>
<mi>ref</mi>
</msub>
<msub>
<mi>f</mi>
<mi>dr</mi>
</msub>
</mrow>
<mn>2</mn>
</mfrac>
<mo>+</mo>
<msup>
<mrow>
<mo>(</mo>
<mfrac>
<mrow>
<mi>λ</mi>
<msub>
<mi>f</mi>
<mi>dc</mi>
</msub>
<mi></mi>
</mrow>
<mn>2</mn>
</mfrac>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</msqrt>
</mrow></math>
<math><mrow>
<mi>φ</mi>
<mo>=</mo>
<mi>arccos</mi>
<mrow>
<mo>[</mo>
<mo>-</mo>
<mfrac>
<msub>
<mi>λf</mi>
<mi>dc</mi>
</msub>
<mrow>
<mn>2</mn>
<mi>v</mi>
</mrow>
</mfrac>
<mo>]</mo>
</mrow>
</mrow></math>
<math><mrow>
<mi>D</mi>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msup>
<mrow>
<mo>[</mo>
<mn>1</mn>
<mo>-</mo>
<msup>
<mrow>
<mo>(</mo>
<mfrac>
<mi>λf</mi>
<mrow>
<mn>2</mn>
<mi>v</mi>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>]</mo>
</mrow>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
</msup>
</mrow></math>
<math><mrow>
<mi>Θ</mi>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mrow>
<mn>4</mn>
<mi>π</mi>
</mrow>
<msup>
<mi>c</mi>
<mn>2</mn>
</msup>
</mfrac>
<msub>
<mi>B</mi>
<mrow>
<mi>r</mi>
<mn>0</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>,</mo>
<msub>
<mi>R</mi>
<mi>ref</mi>
</msub>
<mo>)</mo>
</mrow>
<mrow>
<mo>[</mo>
<mn>1</mn>
<mo>+</mo>
<msub>
<mi>C</mi>
<mi>s</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<mo>]</mo>
</mrow>
<msub>
<mi>C</mi>
<mi>s</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>f</mi>
<mo>)</mo>
</mrow>
<msup>
<mrow>
<mo>[</mo>
<mi>R</mi>
<mfrac>
<mrow>
<mi>sin</mi>
<mi>φ</mi>
</mrow>
<mrow>
<mi>sin</mi>
<msub>
<mi>φ</mi>
<mi>ref</mi>
</msub>
</mrow>
</mfrac>
<mo>-</mo>
<msub>
<mi>R</mi>
<mi>ref</mi>
</msub>
<mo>]</mo>
</mrow>
<mn>2</mn>
</msup>
</mrow></math>
in the above basic variables, f is an azimuth variable, and τ, fτAnd R is a distance variable.
In the expression of three factors, reference distance RrefThe frequency b of the transmitted pulse modulation, the speed of light c and the wavelength lambda of the radar transmitted are known or can be obtained directly by an instrument; and the ground speed v of the radar and the included angle phi between the center of the radar beam and the motion direction of the radar platform need to be calculated. In the formula for v, fdrAnd fdcRespectively representing the Doppler modulation frequency and the Doppler center frequency, f, τ, f in other formulaeτRespectively represent azimuth Doppler frequency, range time and range frequency in a two-dimensional frequency domain. From the basic variables mentioned above, C can be obtaineds(f)、Br0(f) And the like. For further information on the calculation of the three factors, reference is made to document 1: xanthate, "research on imaging processing technology of high-resolution spaceborne synthetic aperture radar", doctor paper of Beijing university of aerospace, 6 months 1999.
In the existing CS algorithm for synthetic aperture radar imaging, when the above calculation formula is used to calculate three factors, namely a Chirp Scaling (CS) factor, a distance compensation factor, and an orientation compensation factor, the calculation formula itself is complex, so the amount of calculated data is huge, the requirements on the storage resources, network resources, and calculation resources of the system are high, and the real-time performance of the calculation is not easily guaranteed. Therefore, the existing factor calculation method is usually realized by adopting a software programming method, and the method has complex calculation and low speed and is difficult to meet the requirement of real-time calculation. If the variables in the three factor formulas can be combined and analyzed based on the principle of numerical analysis, and the variables and the combination thereof are simplified and processed by adopting a numerical approximation calculation method within a certain precision range, the calculation efficiency can be improved, the hardware realization becomes possible, and the real-time performance required during the factor calculation is ensured. In addition, in the prior art, different software and hardware are required to realize the calculation of the related factors aiming at different factors, so that the calculation cost is high, the operation is complex, and no method or device can be simultaneously suitable for the calculation of the three factors at present. Therefore, it is necessary to develop a special factor reusable calculation method and apparatus suitable for calculating three factors with fast calculation speed.
Disclosure of Invention
The invention aims to unify three factor calculation processes of CS compensation, distance compensation and azimuth compensation in a transformation linear frequency modulation scale imaging algorithm, provide a reusable factor calculation method and a reusable factor calculation device, and reduce the design complexity of a factor calculation system.
The invention also aims to simplify the factor calculation method and realize the real-time calculation of the factors in the transformation linear frequency modulation scale imaging algorithm.
In order to achieve the above object, the present invention provides a reusable computing device for transforming a chirp-scale imaging algorithm factor, comprising: the system comprises an input interface module 1, a pre-calculation module 2, a line correlation quantity calculation module 3, a point correlation quantity calculation module 4, an output interface module 5 and a control module 6; the input interface module 1 is used for inputting initial variable values required by factor calculation; the pre-calculation module 2 is used for pre-calculating parameters to be used in the factor calculation; the line correlation quantity calculating module 3 is used for calculating the weak real-time line correlation quantity; the point correlation quantity calculation module 4 is used for realizing the calculation of strong real-time point correlation quantity; the output interface module 5 is used for outputting the result of the factor calculation; the control module 6 is used to control the workflow of the modules and synchronization and communication issues between the modules.
In the above technical solution, the function of the pre-calculation module 2 is implemented by a software programming method.
In the above technical solution, the point correlation quantity calculation module 4 adopts a four-stage pipelined floating-point adder.
A method for performing reusable calculation by using a reusable calculation device for transforming a linear frequency modulation scale imaging algorithm factor comprises the following steps:
1) simplifying the calculation formulas of the three factors to obtain a unified calculation model of factor calculation: a. the1×R1+A2×R2;
2) For one factor of the three factors, the original calculation formula before the factor is simplified is utilized to calculate the model A1×R1+A2×R2Component A of1Solving a plurality of sampling points, and utilizing the sampling points to realize the numerical approximation of a certain polynomial function to the components to obtain the form and the relevant parameters of the polynomial function;
3) replacing the original calculation formula of the component with a polynomial function, and evaluating the polynomial function to obtain the value of the component;
4) determining the other components R in the calculation model in a similar manner to steps 2) and 3)1、A2、R2A value of (d);
5) performing linear operation and data rounding processing on the obtained values of the components to obtain the value of a calculation model, namely obtaining the value of the factor;
6) the values of the other factors are evaluated in the same manner.
In the step 1), A in the model is calculated
1、A
2For weak real-time computation of quantities, R
1、R
2Calculating the quantity for strong real time; for the CS compensation factor, A
1Is (tau-tau)
ref)
2,R
1Is B
r0C
S,A
2And R
2Are all 0; for the distance compensation factor, A
1Is 1/[ B ]
r0(1+C
S)],R
1Is f
τ 2And A is
2Is R
refC
S(f)/c,R
2Is f
τ(ii) a To the orientation compensation factor, A
1Is 4B
r0C
s(1+C
s)/c
2,R
1Is composed of
A
2Is 4
R
2Is R.
In step 2), when the polynomial function is used for numerical approximation, the piecewise polynomial function is preferably used for numerical approximation.
In step 5), the data de-rounding process can be applied to A1×R1And A2×R2Respectively to or from A1×R1+A2×R2And (4) carrying out data trimming processing on the result.
The invention has the advantages that:
1. the invention realizes the multiplexing calculation of the conversion linear frequency modulation scale imaging algorithm factor, so that the same device can be used for the calculation of different factors, and the integration level of the factor calculation is improved.
2. The invention simplifies the design of a factor calculation system in the transformation linear frequency modulation scale imaging algorithm, and is beneficial to the realization of hardware design and software and hardware cooperation.
3. In the process of factor calculation, a polynomial approximation method or a segmented polynomial approximation method is adopted for complex variables, so that the complexity of factor calculation is reduced, and the efficiency of factor calculation is improved.
Detailed Description
The invention is described in further detail below with reference to the figures and the detailed description.
Taking the C band as an example, as shown in fig. 2, the operation steps of the reusable factor calculating method include:
in step 10, three types of basic variables in the calculation formula of the CS compensation factor, the distance compensation factor and the orientation compensation factor are separated: frame correlation quantity, line correlation quantity, point correlation quantity. The frame correlation quantity is a variable which needs to be updated in each frame and is a non-real-time calculation quantity, the non-real-time calculation quantity is a calculation quantity which only needs to be calculated for a fixed number of times when each frame of image data is processed, the value of the calculation quantity is kept as a constant in the data processing of each frame of image, the calculation complexity is a constant and is irrelevant to the size of the image; the line correlation quantity is a variable which needs to be updated when each line in the image is calculated, and is a weak real-time calculation quantity, the weak real-time calculation quantity is that any line of the image to be processed only needs to be calculated for a fixed number of times, and the calculation complexity is in direct proportion to the number of lines of the processed image; the point correlation quantity is a variable which needs to be updated for each pixel point in the image, and is a strong real-time calculation quantity, the strong real-time calculation quantity is that for any pixel of the image data to be processed, the value needs to be calculated independently, and the calculation complexity is in direct proportion to the number of points of the processed image.
The calculation direction is often expressed by a distance direction and an azimuth direction in the imaging process, the calculation direction of the factors is different in different factor calculation formulas, the calculation direction is the azimuth direction in the CS compensation factor, and the calculation direction is the distance direction in the distance compensation factor and the azimuth compensation factor. In the factor calculation, a variable perpendicular to the calculation direction is a weak real-time variable, and a variable in the same direction as the factor calculation direction is a strong real-time variable. Therefore, the azimuth direction variable in the CS compensation factor is a strong real-time point related variable, and the distance direction variable is a weak real-time line related variable; in the distance compensation factor and the azimuth compensation factor, the distance direction variable is a strong real-time point related variable, and the azimuth direction variable is a weak real-time line related variable.
Taking CS compensation factors as an example, the calculation direction is the azimuth direction, the azimuth direction frequency f is the azimuth direction variable, and the range direction time τ is the range direction variable, so f and its combination are calculated as strong real-time quantities, and τ and its combination are calculated as weak real-time quantities. And τref(f) The time is referred to in the azimuth direction, and the change range of the value in the whole azimuth direction is very small, so that the time can be used as a non-real-time change quantity and is approximate to the time in the same frameThe image data is kept constant.
In step 20, the basic variables in the computational expressions of the three factors are numerically analyzed, and some of the variables are simplified as appropriate. Within one frame of image data by counting the Doppler center frequency fdcDoppler frequency modulation fdrThe numerical analysis of the radar ground speed v, the included angle phi between the radar beam center and the radar platform motion direction and other variables shows that the slope of the variables is very low, and the influence on the factor calculation precision is not large after the variables are simplified into constants within a certain range, so that the variables can be regarded as constants.
In step 30, some variables in the factor calculation expression with the same calculation direction are combined to generate new variables, and the new variables are simplified by adopting a numerical approximation method for calculation. As in the CS compensation factor, the azimuth direction variable Br0And an azimuth direction variable CSBoth variables are variables related to f, and the calculation directions of the two variables are the same, so the two variables can be combined into a new azimuth variable. If for the new variable B synthesizedr0CSWhen the existing calculation formula is used for calculation, the existing calculation formula is too complex, and the calculation speed and the calculation resources are difficult to meet the actual requirements, so that the new variable B to be synthesized is requiredr0CSThe calculation formula of (2) is simplified. From the numerical analysis of the last step, it can be known that in the C-band SAR imaging transformation chirp scaling algorithm, the f dynamic range is smaller, so that for Br0CSThe least square method can be adopted to carry out polynomial (usually, binomial is adopted, and numerical analysis shows that the linear function is selected for C-band data to obtain quite good approximation effect) function numerical approximation, and the polynomial function used in numerical approximation is Br0CSNew calculation formula (2). Compared with the original calculation formula, the new calculation formula is greatly simplified, and the accuracy of the calculation result is equivalent.
In step 40, the calculation formula is simplified according to the calculation formula of the three factors to obtain a unified calculation model of factor calculationForm A1×R1+A2×R2. Wherein A is1、A2For weak real-time computation of quantities, R1、R2The amount is calculated in strong real time. For three factors, A1、A2、R1、R2The meanings of the expressions are different, and A can be obtained by simplifying the calculation formula of the three factors1、A2、R1、R2The respective expressions: for the CS compensation factor, A1Is (tau-tau)ref)2,R1Is Br0CSIn which τ isrefCan be regarded as a constant within a certain precision range, and A2And R2Are all 0; for the distance compensation factor, A1Is 1/[ B ]r0(1+CS)],R1Is fτ 2And A is2Is RrefCS(f)/c,R2Is fτ(ii) a To the orientation compensation factor, A1Is 4Br0Cs(1+Cs)/c2,R1Is composed ofAnd A is2Is 4R2Is R.
In step 50, to reduce the complexity of the calculation, model A is calculated1×R1+A2×R2Of a single component, e.g. A1、A2、R1、R2Respectively realizing the approximation of each component by adopting a polynomial function based on a least square method, and finally performing linear operation on each component to obtain a calculation model A1×R1+A2×R2To obtain the values of the three factors.
By R in the CS compensation factor1=Br0CSFor example, the implementation method comprises the following steps:
a) to R1By throwingPerforming approximate calculation on the object-line function, and setting the approximate formula as R1(f)=a1f2+b1f+c1The approximation is a parabolic function with the doppler frequency f as a variable.
b) According to CS compensation factor Br0CSIs calculated by
<math><mrow>
<msub>
<mi>B</mi>
<mrow>
<mi>r</mi>
<mn>0</mn>
</mrow>
</msub>
<msub>
<mi>C</mi>
<mi>s</mi>
</msub>
<mo>=</mo>
<mrow>
<mo>(</mo>
<mfrac>
<mrow>
<mi>sin</mi>
<mi>φ</mi>
</mrow>
<msqrt>
<mn>1</mn>
<mo>-</mo>
<msup>
<mrow>
<mo>(</mo>
<mfrac>
<mi>λf</mi>
<mrow>
<mn>2</mn>
<mi>v</mi>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
</msqrt>
</mfrac>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<mo>·</mo>
<mrow>
<mo>(</mo>
<mfrac>
<mi>b</mi>
<mrow>
<mn>1</mn>
<mo>+</mo>
<mi>b</mi>
<msub>
<mi>R</mi>
<mi>ref</mi>
</msub>
<mi>sin</mi>
<mi>φ</mi>
<mfrac>
<mrow>
<mn>2</mn>
<mi>λ</mi>
</mrow>
<msup>
<mi>c</mi>
<mn>2</mn>
</msup>
</mfrac>
<mfrac>
<msup>
<mrow>
<mo>(</mo>
<mfrac>
<mi>λf</mi>
<mrow>
<mn>2</mn>
<mi>v</mi>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<msup>
<mrow>
<mo>[</mo>
<mn>1</mn>
<mo>-</mo>
<msup>
<mrow>
<mo>(</mo>
<mfrac>
<mi>λf</mi>
<mrow>
<mn>2</mn>
<mi>v</mi>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mn>2</mn>
</msup>
<mo>]</mo>
</mrow>
<mfrac>
<mn>3</mn>
<mn>2</mn>
</mfrac>
</msup>
</mfrac>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
</mrow></math>
To obtain R1More than 3 sampling points, substituting values of the sampling points with respect to R1In the approximation of (a), thereby obtaining a1、b1、c1The numerical value of (c). A is to1、b1、c1Substitution of the value of (into R)1Can be used to conveniently obtain R in real-time calculation1The value of (c).
c) Solving for A by a similar method1、A2、R2And finally according to the calculation model A1×R1+A2×R2The value of the orientation compensation factor is found.
The orientation compensation factor and the distance compensation factor can be evaluated in the same manner.
In the method of the present invention, a polynomial function is used in steps 20 and 50 to achieve numerical approximations of complex variables in the computational formula and various components in the computational model. In calculating these polynomial functions, an adder may be used for implementation. But in the process of using the adder to calculateThe single step error accumulation of the polynomial function is amplified. In order to eliminate the accumulated error possibly brought by the adder, a piecewise function is adopted to carry out numerical approximation. For convenience of description, taking a one-dimensional linear function y as ex + t as an example, when performing hardware calculation on the function, one multiplier and one adder are required; if the function uses the recursion formula y in the calculationi+1-yi=e(xi+1-xi) When processing the SAR image, the difference dx between the independent variables of two adjacent points is equal to (x)i+1-xi) Is constant, independent of the position i of the point, so when calculating yi+1=yiOnly one adder is needed for + dy (where dy is e × dx), which simplifies hardware implementation and reduces system complexity compared to a method that requires both a multiplier and an adder. However, the adder inevitably has an error in each accumulation of the linear function, and the adder takes the result of the previous accumulation as an input value for the next accumulation, so that the addition error is cumulatively amplified as the number of accumulations increases. To eliminate such accumulated error, a piecewise linear function is used to the linear function yi=eixi+ti(where i represents different segments) by approximating the variables using a plurality of fold lines, the slope e of the linear function in each segment of the linear functioniWith a slight difference. The accumulation length of the adder is limited to the quotient between the point number in the calculation direction and the segment number, and meanwhile, the error of the variable approximation by adopting a plurality of broken lines is further reduced compared with the error of the variable approximation by adopting a single straight line. Similar to the one-dimensional linear function y ═ ex + t, when a polynomial function is used for numerical approximation, a piecewise polynomial function may also be used to achieve numerical approximation of variables or components, so as to improve the calculation accuracy. In addition, when the polynomial function is segmented, the dimensionality of the segmented function can be reduced, and the numerical approximation by adopting a lower function can achieve the effect equivalent to that of a high-dimensional function without segmentation.
While a one-dimensional linear function is illustrated in the present embodiment, it should be understood by those skilled in the art that the precision can be improved when a high-dimensional function is used to approximate a variable.
The method of the invention can adopt a uniform calculation model to realize the calculation of the three factors. Since the approximation of each component is realized by using a polynomial function based on the least square method in the factor calculation process, for each factor, the polynomial coefficients obtained by the least square method only vary with the frame-related parameters and are independent of the position in the frame, and the polynomial coefficients, such as a1、b1、c1Also considered as frame correlation quantity. For the calculation of the frame correlation quantity, the calculation is very flexible, and the data quantity needing to be calculated is very small, so that the calculation can be carried out by adopting a software pre-calculation method. In the factor calculation process, real-time quantities are calculated, and because the real-time requirements on the real-time quantities are high, the calculation of the real-time quantities is realized by adopting a hardware or software and hardware combined method. In summary, the reusable computing device of the factor is cooperatively realized by software and hardware layering.
As shown in fig. 3, the reusable computing device for transforming the chirp-scale imaging algorithm factors comprises: the device comprises an input interface module 1, a pre-calculation module 2, a line correlation quantity calculation module 3, a point correlation quantity calculation module 4, an output interface module 5 and a control module 6.
The function of the input interface module 1 is to input the initial variable values required for factor calculation, including: reference distance RrefFrequency b of transmitted pulse, speed of light c, radar transmission wavelength lambda, Doppler center frequency fdcDoppler frequency modulation fdr. The values of other variables required for the factor calculation can be obtained by the above-mentioned variable calculation.
The function of the pre-calculation module 2 is to pre-calculate parameters to be used in subsequent steps by a programmed method, including: the calculation is used for realizing the calculation of each component A in the calculation model1、A2、R1、R2Polynomial functions used for approximating numerical values, e.g. R in CS compensation factors1(f)=a1f2+b1f+c1(ii) a Correlation coefficient of polynomial function, e.g. a1、b1、c1(ii) a Frame-related parameters required in performing line-related or point-related parameter calculations, e.g. R in CS compensation factor1The initial value of f of (a), and the like. As known from the factor calculation method, the correlation coefficient of the polynomial function is obtained by calculating a certain number of sampling points for a variable using an existing complex formula and substituting the values of the sampling points into the polynomial function. The above operations have high computational complexity, take long time, and are difficult to implement in hardware due to the complex formula. However, most of these parameters are non-real-time frame-related quantities, and the requirement on real-time performance is low, so that the calculation of these parameters can be realized in software by a programming method.
The function of the line correlation quantity calculation module 3 is to realize the calculation of the line correlation quantity in weak real time, and the line correlation quantity comprises: factor unified computing model A1×R1+A2×R2A in (A)1And A2For different factors, A1And A2The meaning indicated is different, but the calculation of different factors can be carried out with the same module. The line correlation amount calculation module 3 calculates the required line correlation amount by using the frame correlation amount obtained in the pre-calculation module 2. Because the line correlation quantity is weak real-time quantity and the requirement on the real-time performance is not very high, the line correlation quantity calculating module 3 can adopt a method of combining software and hardware, and realize the calculation of the line correlation quantity by using less hardware resources and a method of controlling software and state. The hardware resources only need one multiplier, one adder, and some registers and state machines at minimum. According to the concrete condition of hardware resource and polynomial function for numerical value approximation, different software is edited for controlling hardware to implement calculation of line correlation quantity. Compared with a method of using hardware only, the method of combining software and hardware can achieve the purpose of replacing space and hardware resources with time.
The function of the point correlation quantity calculation module 4 is to realize the calculation of the point correlation quantity in strong real time, and the point correlation quantity comprises a factor unified calculation moduleForm A1×R1+A2×R2R in (1)1And R2And a multiply-add operation on the four components.
For different factors, R1And R2The meaning indicated is different, but the calculation of different factors can be carried out with the same module. The point correlation quantity calculating module 4 inputs the frame correlation quantity obtained by the pre-calculating module 2 and the line correlation quantity obtained by the line correlation quantity calculating module 3, and calculates the point correlation quantity by using hardware. Because the point correlation quantity is strong in real time, the calculation of the point correlation quantity is mainly realized by adopting a hardware method. At point-to-point correlation quantity R1And R2The method of floating-point calculation is adopted in the calculation, and the real-time calculation requires that one result is output per beat, so in the embodiment, the floating-point adder can adopt a four-stage pipelined floating-point adder as shown in fig. 4. Four stages in the four-stage pipelined floating-point adder mean that four beats are required in the floating-point adder to complete one addition operation, and the pipelined means that although four beats are required to complete one addition operation, one result can be output per beat, and the time overlapping is completed by using the expansion of space like a pipelined in industrial generation. For example, when the first addition operation is at the second stage of the four-stage pipeline, the second addition operation is at the first stage of the four-stage pipeline. Four beats in the four-stage pipelined floating-point adder respectively complete complement pair-order, complement addition/subtraction, complement normalization and rounding processing. The four-stage pipeline floating-point adder has two inputs, one is the calculation result output by the adder, and the other is the accumulated pipeline characteristic parameter. The accumulated pipeline characteristic parameter is a frame-dependent parameter that includes an initial value of the pipeline and a gain value of the pipeline, the gain value of the pipeline being related to a slope of the polynomial function. The present adder is a four-stage pipeline, and therefore requires four initial values, which are assumed to be represented by s0, s1, s2, and s3, respectively. The input sequence is s0, s1, s2, s3, s0+4, s1+4 and s2+4 … respectively, and the output sequence of the adder is s0+4, s1+4, s2+4, s3+4 and s0+4+4 … respectively. It is clear that neither s0, s1, s2, s3 can be obtained from the output of the pipeline, only by means of pre-predictionThe computation is obtained, i.e. the pipeline is initialized with the result of the pre-computation, and once the pipeline is running, its output can be continuously run as input. In the foregoing method, it is mentioned that the error is smaller when the piecewise function is used to implement the numerical approximation of the variable than when the simple function is used to implement the numerical approximation, and the gain value of the pipeline is related to the slope of the function, so that when the piecewise function is used to implement the numerical approximation of the variable, the gain value of the pipeline is also changed accordingly. In the process of pipeline calculation, under the control of the control state machine, the gain value of the pipeline can be reset so as to improve the precision of the factor calculation result.
Unified calculation model A for solving factors1×R1+A2×R2R in (1)1And R2Then, A is obtained by a multiplier1×R1And A2×R2Then, the obtained data is subjected to rounding processing, and finally, the rounded data is added by an adder to obtain a final value, so that the value of the correlation factor is obtained. The data rounding processing means: since the factor is calculated as an angle value, if the angle exceeds 360 °, the value is meaningless, and the trimming process is to remove the portion of the angle exceeding 360 ° and control the angle value within 360 °. In the invention, the calculation adopts a floating point format, so that the data is only required to be shifted according to the exponent of the floating point data during the rounding processing, and then the integer part of the data is discarded. In the prior art, the rounding-off process is performed after all computations are completed, and in the method of the present invention, the rounding-off process is completed before the last addition operation, which can be used for adding A to A1×R1And A2×R2The original floating-point adder for addition processing is simplified into a fixed-point adder.
The output interface module 5 is used for outputting the factor calculation result obtained in the previous step. The output interface module 5 includes a corresponding output interface circuit.
The control module 6 is used for controlling the workflow of the module and the synchronization and communication problem among the modules, and under the action of the control signal of the control module 6, the system can realize the calculation of different factors. Meanwhile, the control module is responsible for controlling and monitoring the state of the whole system and for the synchronization and communication problems between the system and other systems.
The factor calculation system has universality on the CS compensation factor, the distance compensation factor and the orientation compensation factor, and different factors can be calculated by inputting different initial values into the factor calculation system.
Partial imaging results realized by using the C-band SAR imaging transformation chirp scaling algorithm based on the present invention can be seen in fig. 5 and 6, where the imaging area in fig. 5 is the beijing palace and the imaging area in fig. 6 is the beijing palace.
For SAR imaging conversion chirp scaling algorithms of other wave bands, only slight differences exist on the internal structure of specific hardware components, and the calculation method and the system structure have no essential differences. The method of the present invention is therefore widely adaptable.