CN103454636A - Differential interferometric phase estimation method based on multi-pixel covariance matrixes - Google Patents

Differential interferometric phase estimation method based on multi-pixel covariance matrixes Download PDF

Info

Publication number
CN103454636A
CN103454636A CN2013104049900A CN201310404990A CN103454636A CN 103454636 A CN103454636 A CN 103454636A CN 2013104049900 A CN2013104049900 A CN 2013104049900A CN 201310404990 A CN201310404990 A CN 201310404990A CN 103454636 A CN103454636 A CN 103454636A
Authority
CN
China
Prior art keywords
deformation
interferometric phase
differential interferometry
phase
interferometric
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.)
Granted
Application number
CN2013104049900A
Other languages
Chinese (zh)
Other versions
CN103454636B (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 CN201310404990.0A priority Critical patent/CN103454636B/en
Publication of CN103454636A publication Critical patent/CN103454636A/en
Application granted granted Critical
Publication of CN103454636B publication Critical patent/CN103454636B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Radar Systems Or Details Thereof (AREA)

Abstract

The invention discloses a differential interferometric phase estimation method based on multi-pixel covariance matrixes. The differential interferometric phase estimation method mainly solves the problem that an existing differential interferometric phase generating algorithm is not stable in registration errors. The differential interferometric phase estimation method includes the first step of inputting images and parameters, the second step of registering interferometric phase diagrams, the third step of obtaining the interferometric phase diagrams through secondary interference, the fourth step of obtaining image windows, the fifth step of constructing repeated interferometric phase joint data vectors, the sixth step of estimating the covariance matrixes, the seventh step of estimating correlation coefficient matrixes, the eighth step of carrying out matrix characteristic decomposition, the ninth step of constructing cost functions, and the tenth step of estimating differential interferometric phases. The method has the advantages of restoring pixel information in a self-adaptation mode under the condition that the interferometric phase diagrams have the registration errors, accurately estimating the differential interferometric phases, and being capable of effectively reducing the influence of the registration errors on the differential interferometric phases.

Description

Differential interferometry phase estimation method based on many pixels covariance matrix
Technical field
The invention belongs to communication technical field, further relate to the differential interferometry phase estimation method based on many pixels covariance matrix in the radar exploration technique field.The present invention can the condition very poor for interferometric phase image registration accuracy in differential interferometry synthetic-aperture radar (Differential Interferometric Synthetic Aperture Radar, D-InSAR) under, sane estimation difference interferometric phase.
Background technology
D-InSAR is the very fast remote sensing technology of development in recent years, there is round-the-clock, round-the-clock, wide area high resolution detection Ground Deformation, with respect to existing GPS (Global Position System, GPS) and the measurement of the level means, huge advantage and wide application prospect are arranged, be widely used in the fields such as city settlement monitoring, earthquake deformation measurement, volcanic explosion.
Image registration accuracy has determined the order of accuarcy of the differential interferometry phase place that D-InSAR obtains, and the differential interferometry phase place has directly been reacted the deformation quantity on ground, and the registration accuracy of interferometric phase image is to affect the key factor that the differential interferometry phase place generates.Existing differential interferometry phase estimation method is mainly based on the conventional InSAR treatment scheme, conventional InSAR processing requirements registration accuracy reaches in 1/10 pixel, and need complicated Interpolation Process, when registration error surpasses 1/10 pixel, will obtain inaccurate differential interferometry phase place.
M.Massonnet, the people such as R.M.Goldstein and H.A.Zebker is at article " Mapping small elevation changes over large areas:differential interferometry " (Journal of Geophysical Research, 1989, vol.94:9183-9191) the middle method of estimation that proposes a kind of differential interferometry phase place generated based on existing InSAR interferometric phase.The method is first carried out smart registration operation to SAR image and interferometric phase image, then the interferometric phase image obtained is done to second order difference, thereby finally resulting filtering is as a result obtained to the differential interferometry phase place.The weak point that the method exists is: when interferometric phase image has larger registration error to exist, can't obtain result accurately, and need to carry out filtering when obtaining the differential interferometry phase place, registration and filtering are carried out step by step, registration error can be lost differential interferometry phase estimation precision, and the phase filtering step can not be recovered the loss of significance that registration error is brought, its process errors is constantly accumulation.
Liao Mingsheng, Lu Lijun, Wang Yan and Li Deren are in article " research of the InSAR technology for detection earth's surface miniature deformation of analyzing based on point target " (" urban geology ", 2006,1 (2): proposed a kind of method of estimating the difference interferometric phase based on the coherent point target analysis 38-41).The weak point that the method exists is: only have interferometric phase image before and after deformation all during smart registration, just can obtain differential interferometry phase place accurately, and, for the interferometric phase image before and after the deformation of thick registration, the noise that the differential interferometry phase diagram that the method obtains contains is larger.
The patented claim that He'nan University proposes " the earth's surface three-dimensional deformation monitoring method based on InSAR and gps data the fusion " (date of application: on 02 05th, 2010, application number: 201010106794, publication number: CN101770027A) in a kind of differential interferometry phase place preparation method merged in conjunction with gps data is disclosed.The method is revised the SAR satellite orbital error by the coordinate transformation relation of setting up GPS and SAR data, then obtains geography information by the interferometric phase image of unwrapping accuracy registration, finally merges InSAR and the GPS deformation data obtains the differential interferometry phase place.The weak point that the method exists is: while with gps data, revising the SAR satellite orbital error, corresponding pixel is registration fully, when merging InSAR and GPS deformation data, and the necessary registration fully of pixel that same landform is corresponding.
Summary of the invention
While the present invention is directed to the generation differential interferometry phase place of above-mentioned prior art existence, the defect of the necessary accuracy registration of the interferometric phase image before and after deformation, proposed a kind of differential interferometry phase estimation method based on many pixels covariance matrix.The present invention builds covariance matrix by the information consolidation of a plurality of pixels, again estimated covariance matrix feature decomposition is obtained to signal subspace, finally utilize the method for signal subspace fitting to estimate the difference interferometric phase, guaranteed under interferometric phase image before and after the deformation condition in registration accuracy very poor, the sane estimation of differential interferometry phase place completes the differential interferometry phase filtering simultaneously.
For achieving the above object, key step of the present invention is as follows:
(1) input picture and parameter:
Before target area deformation 1a) respectively the differential interferometry synthetic-aperture radar obtained and the level land interferometric phase image that goes after deformation be input to system;
1b) systematic parameter of differential interferometry synthetic-aperture radar is input to system;
(2) interferometric phase image registration:
2a) deduct rear data line with the every data line in the interferometric phase image of level land that goes before deformation respectively, obtain removing level land row differential interferometry phase diagram before deformation;
2b) deduct a rear column data with each column data in the row differential interferometry phase diagram of level land that goes before deformation respectively, obtain the front interferometric phase gradient map of deformation;
2c) deduct rear data line with the every data line in the interferometric phase image of level land that goes after deformation respectively, obtain removing level land row differential interferometry phase diagram after deformation;
2d) respectively with after deformation go each column data in the row differential interferometry phase diagram of level land to deduct after a column data, obtain the interferometric phase gradient map after deformation;
The interferometric phase gradient map of 2e) take before deformation is reference picture, and the interferometric phase gradient map after deformation is carried out to registration process, before obtaining the deformation of registration and the side-play amount when interferometric phase gradient map after deformation and registration;
2f) by all pixels in the interferometric phase image of level land of going after deformation, a side-play amount during registration of translation, obtain before the deformation of thick registration interferometric phase image after interferometric phase image and deformation;
(3) obtain the interferometric phase image that secondary is interfered:
3a) by the data of all pixels in interferometric phase image before the deformation of thick registration, be multiplied by scale factor, obtain the front mark interferometric phase image that becomes of deformation of thick registration;
3b) by the data of all pixels in interferometric phase image after the deformation of thick registration, deduct the front data that become corresponding pixel points in the mark interferometric phase image of deformation of thick registration, obtain the interferometric phase image that secondary is interfered;
(4) obtain image window:
In the interferometric phase image of interfering at secondary, centered by pixel to be estimated, regular length is the piece radius, obtains a square window;
(5) build multiple interferometric phase associating data vector:
5a) choose arbitrarily a pixel in square window, the data of neighbor pixel around selected pixel and its are formed a line, obtain the phase vectors of selected pixel;
5b) by each element in the phase vectors of selected pixel, ask after multiplying each other with imaginary unit and take the index that natural Exponents e is the end, obtain the index vector of selected pixel;
5c) before each element of the index vector of selected pixel, insert positive integer 1, obtain the multiple interferometric phase associating data vector of pixel to be estimated;
5d) travel through all pixels in square window, judge whether to obtain the multiple interferometric phase associating data vector of all pixels, if so, perform step (6), otherwise, execution step 5a);
(6) adopt following formula, the estimate covariance matrix:
C = 1 M 2 Σ i = 1 M Σ j = 1 M l ( i , j ) × l ( i , j ) H
Wherein, C means estimated covariance matrix, and l (i, j) means the multiple interferometric phase associating data vector of the capable j column element of i in square window, i=1,2 ..., M, j=1,2,, M, M means the radius of square window, H means to do the conjugate transpose operation;
(7) estimate the coefficient of coherence matrix:
All elements value in covariance matrix is taken absolute value, obtain the coefficient of coherence matrix of pixel to be estimated;
(8) the coefficient of coherence matrix is carried out to the feature decomposition operation, obtain eigenwert and the proper vector of coefficient of coherence matrix;
(9) adopt following formula, the structure cost function:
Wherein, J means constructed cost function,
Figure BDA0000378885740000042
mean the differential interferometry phase place
Figure BDA0000378885740000043
corresponding differential interferometry phase vectors, N means the number of large eigenwert in the eigenwert of coefficient of coherence matrix, β plarge eigenwert characteristic of correspondence vector in the proper vector of expression coefficient of coherence matrix, p=1,2 ..., N, β qlarge eigenwert characteristic of correspondence vector in the proper vector of expression coefficient of coherence matrix, q=1,2 ..., N, H means to do the conjugate transpose operation,
Figure BDA0000378885740000046
mean that Hadamard Hadamard is long-pending;
(10) estimate the difference interferometric phase:
(π, π] scope in, search differential interferometry phase place
Figure BDA0000378885740000044
, the maximal value of searching cost function, by the corresponding differential interferometry phase place of cost function maximal value
Figure BDA0000378885740000045
differential interferometry phase place as pixel to be estimated.
The present invention has the following advantages compared with prior art:
First, the present invention utilizes the cost function of the proper vector structure of covariance matrix to estimate the difference interferometric phase, overcome estimating in difference interferometric phase method based on the coherent point target analysis of prior art, in the thick registration situation of interferometric phase image, the inaccurate defect of differential interferometry phase place obtained, make the present invention there is good robustness to the interferometric phase image registration error, have in the situation that there is very large registration error in interferometric phase image, still can obtain the advantage of differential interferometry phase place accurately.
Second, the present invention utilizes many pixels joint mapping covariance matrix, overcome in the differential interferometry phase estimation method generated based on the InSAR interferometric phase of prior art, the defect that the differential interferometry phase loss that registration error causes can't recover by phase filtering, make the present invention can adaptive recovery differential interferometry phase information, have the phase information utilization abundant, the differential interferometry phase estimation is advantage accurately.
The accompanying drawing explanation
Fig. 1 is process flow diagram of the present invention;
Fig. 2 removes the level land interferometric phase image for what the InSAR interference technique emulation that adopts prior art generated;
Fig. 3 is for adopting the present invention and prior art to process the differential interferometry phase diagram of the Generation of simulating data of interferometric phase image accuracy registration;
Fig. 4 is for adopting the present invention and prior art to process the differential interferometry phase diagram of the Generation of simulating data that the interferometric phase image registration error is 1 pixel;
Fig. 5 is for adopting the present invention and prior art to process the differential interferometry phase diagram that measured data generates.
Embodiment
Below in conjunction with accompanying drawing, the present invention will be further described.
With reference to accompanying drawing 1, the specific embodiment of the invention step is as follows:
Step 1, input picture and parameter.
Before the target area deformation that the differential interferometry synthetic-aperture radar is obtained with deformation after remove the level land interferometric phase image, the number of times of the imaging to target area and the base length of differential interferometry synthetic-aperture radar are input to system.The base length that the level land interferometric phase image is corresponding of going after deformation is less than scarcely than the base length that the level land interferometric phase image is corresponding of going before deformation.
Step 2, the interferometric phase image registration.
The registration of interferometric phase image, be that the gradient map of two width interferometric phase images is carried out to registration, then according to the side-play amount of gradient map, interferometric phase image carried out to translation.In the present invention, the registration accuracy of needs two a width interferometric phase image reaches a pixel and gets final product, and has greatly alleviated the difficulty of image registration.
The concrete steps of interferometric phase image registration are as follows:
The first step, deduct rear data line with the every data line in the interferometric phase image of level land that goes before deformation, obtains removing level land row differential interferometry phase diagram before deformation;
Second step, deduct a rear column data with each column data in the row differential interferometry phase diagram of level land that goes before deformation, obtains the front interferometric phase gradient map Q of deformation;
The 3rd step, deduct rear data line with the every data line in the interferometric phase image of level land that goes after deformation, obtains removing level land row differential interferometry phase diagram after deformation;
The 4th step, with after deformation go each column data in the row differential interferometry phase diagram of level land to deduct after a column data, obtain the interferometric phase gradient map S after deformation;
The 5th step, utilize following formula, obtains the front interferometric phase gradient map of deformation and the cross-correlation matrix of the interferometric phase gradient map after deformation:
L=IFFT2(FFT2(Q)×conj(FFT2(S)))
Wherein, L means interferometric phase gradient map before deformation and the cross-correlation matrix of the interferometric phase gradient map after deformation, IFFT2 () means to do two-dimentional inverse Fourier transform operation, FFT2 () means to do the two-dimensional Fourier transform operation, Q means the interferometric phase gradient map data matrix before deformation, S means the interferometric phase gradient map data matrix after deformation, and conj () means to get conjugate operation;
The 6th step, maximizing from all elements of cross-correlation matrix L, obtain coordinate figure corresponding to maximal value element in cross-correlation matrix L;
The 7th step, with coordinate figure corresponding to maximal value element in cross-correlation matrix L, deduct coordinate figure corresponding to central element of cross-correlation matrix L, the side-play amount while obtaining registration;
The 8th step, the side-play amount of the interferometric phase gradient map during according to registration after deformation, do integral translation to the level land interferometric phase image that goes after deformation, obtained going before the deformation of thick registration and removed the level land interferometric phase image after level land interferometric phase image and deformation.
Step 3, obtain the interferometric phase image that secondary is interfered.
Interferometric phase image to two width registrations is done Dual-interferometry, the interferometric phase image of take after deformation is reference, interferometric phase image before deformation is adjusted on identical with it ratio, and scale factor and radar are relevant to the base length of the imaging number of times of target area and radar.
The concrete steps of the interferometric phase image that the acquisition secondary is interfered are as follows:
The first step, when the difference interference synthetic aperture radar equals 2 to the imaging number of times of target area, scale factor equals 1;
Second step, when the difference interference synthetic aperture radar is greater than 2 to the imaging number of times of target area, according to the following formula, calculate scale factor:
k = B A
Wherein, k means scale factor, and B means the base length of going the differential interferometry synthetic-aperture radar that the level land interferometric phase image is corresponding after deformation, and A means the base length of going the differential interferometry synthetic-aperture radar that the level land interferometric phase image is corresponding before deformation;
The 3rd step, by the data of removing all pixels in the interferometric phase image of level land before deformation, with scale factor, k multiplies each other, and obtains the front change mark of deformation and removes the level land interferometric phase image;
The 4th step, by the data of removing all pixels in the interferometric phase image of level land after deformation, deduct the data that the front change mark of deformation goes corresponding pixel points in the interferometric phase image of level land, obtains the interferometric phase image that secondary is interfered.
Step 4, obtain image window.
In the interferometric phase image of interfering at secondary, centered by pixel to be estimated, 3~7 pixels of regular length of take are the piece radius, obtain a square window.The square window that to have chosen size in example of the present invention be 9 * 9.
Step 5, build multiple interferometric phase associating data vector.
What while building multiple interferometric phase associating data vector, use is current pixel point and pixel on every side thereof, has utilized fully the information of current pixel point, thus make while generating the differential interferometry phase place can adaptive reply current pixel point information.
The concrete steps that build multiple interferometric phase associating data vector are as follows:
5a) choose arbitrarily a pixel in square window, the data of 1~8 neighbor pixel around selected pixel and its are formed a line, obtain the phase vectors of selected pixel, in the invention of this example, 8 neighbor pixels around having adopted build jointly;
5b) by each element in the phase vectors of selected pixel, ask after multiplying each other with imaginary unit and take the index that natural Exponents e is the end, obtain the index vector of selected pixel;
5c) before each element of the index vector of selected pixel, insert positive integer 1, obtain the multiple interferometric phase associating data vector of pixel to be estimated;
5d) all pixels in the traversal square window, judge whether the multiple interferometric phase associating data vector of all pixels obtains, if all obtain, performs step 6, otherwise, execution step 5a).
Step 6, the estimate covariance matrix.
Generally, can't directly obtain covariance matrix, so utilize the space sample of covariance matrix on average to replace statistical average.Adopt following formula estimate covariance matrix in example of the present invention:
C = 1 M 2 Σ i = 1 M Σ j = 1 M l ( i , j ) × l ( i , j ) H
Wherein, C means estimated covariance matrix, and l (i, j) means the multiple interferometric phase associating data vector of the capable j column element of i in square window, i=1,2 ..., M, j=1,2,, M, M means the radius of square window, M=4 in the invention of this example, H means to do the conjugate transpose operation.
Step 7, estimate the coefficient of coherence matrix.
In interferometric phase image, noise power is much smaller than signal power, and the loading of very little diagonal angle is very little on eigenwert and the proper vector impact of matrix, so utilize covariance matrix to estimate the coefficient of coherence matrix.All elements value in covariance matrix C is taken absolute value, obtain the coefficient of coherence matrix R of pixel to be estimated.
Step 8, to the coefficient of coherence matrix, R carries out the feature decomposition operation, obtains eigenwert and the proper vector of coefficient of coherence matrix R.
Step 9, estimate the difference interferometric phase.
In the invention of this example, by the method for signal subspace fitting, estimate the difference interferometric phase, first by the base vector of signal subspace, construct cost function J.
Adopt following formula, the structure cost function:
Figure BDA0000378885740000081
Wherein, J means constructed cost function, and N means the number of large eigenwert in the eigenwert of coefficient of coherence matrix, N=9 in the invention of this example, β plarge eigenwert characteristic of correspondence vector in the proper vector of expression coefficient of coherence matrix, p=1,2 ..., 9, β qlarge eigenwert characteristic of correspondence vector in the proper vector of expression coefficient of coherence matrix, q=1,2 ..., 9, H means to do the conjugate transpose operation,
Figure BDA0000378885740000089
mean that Hadamard Hadamard is long-pending,
Figure BDA0000378885740000082
mean the differential interferometry phase place
Figure BDA0000378885740000083
corresponding differential interferometry phase vectors, the concrete steps that build the differential interferometry phase vectors are as follows:
The first step, by the differential interferometry phase place of pixel to be estimated
Figure BDA0000378885740000084
, ask after multiplying each other with the j of imaginary unit and take the index that natural Exponents e is the end, obtain the exponential quantity of pixel to be estimated
Figure BDA0000378885740000085
;
Second step, the exponential quantity by positive integer 1 with pixel to be estimated order forms a line, and obtains the single differential interferometry phase vectors of pixel to be estimated;
The 3rd step, choose 9 single differential interferometry phase vectors, and all elements in 9 single differential interferometry vectors is formed a line, and obtains the differential interferometry phase vectors of pixel to be estimated
Figure BDA0000378885740000087
.
Step 10, estimate the difference interferometric phase.
(π, π] scope in, search differential interferometry phase place
Figure BDA0000378885740000088
, the maximal value of searching cost function J, by the corresponding differential interferometry phase place of cost function J maximal value
Figure BDA0000378885740000091
differential interferometry phase place as pixel to be estimated.
Below in conjunction with emulation experiment, effect of the present invention is described further.
1, simulated conditions:
With reference to accompanying drawing 1, by the differential interferometry synthetic-aperture radar, areal is carried out to simulation imaging and obtain three width SAR images, after adopting the InSAR of prior art to interfere to process, before obtaining respectively deformation and after deformation, remove level land interferometric phase image, the input picture using it as system.
Fig. 2 removes the level land interferometric phase image for what the InSAR interference technique emulation that adopts prior art generated, wherein, remove the level land interferometric phase image before the deformation that Fig. 2 (a) generates for the InSAR interference technique emulation that adopts prior art, after the deformation that Fig. 2 (b) generates for the InSAR interference technique emulation emulation that adopts prior art, remove the level land interferometric phase image.Because maximum landform shape becomes 0.1m, so go the level land interferometric phase image approximate identical in Fig. 2 (a) and Fig. 2 (b), simulation parameter arranges as follows:
Systematic parameter Parameter value
Wavelength 0.03m
Antenna length 2m
Pulse repetition rate 500Hz
Platform speed 80m/s
The radar downwards angle of visibility 45°
Pulse width 5us
The platform flying height 6000m
Base length 0.8m
Range resolution 1m
The largest deformation amount 0.1m
2, analysis of simulation result:
Fig. 3 is for adopting the present invention and prior art to process the differential interferometry phase diagram of the Generation of simulating data of interferometric phase image accuracy registration, wherein, Fig. 3 (a) is for adopting prior art to process the differential interferometry phase diagram of the Generation of simulating data of interferometric phase image accuracy registration, abscissa axis mean the distance to, unit is pixel, axis of ordinates mean orientation to, unit is pixel, in the image in left side, the gray-scale value of pixel means the size of differential interferometry phase place, the size of differential interferometry phase fluctuation in figure, reflected the number of differential interferometry phase noise, the image strip on right side means the corresponding relation of gray-scale value in differential interferometry phase value and left-side images, the size of the numeric representation differential interferometry phase place of image right, Fig. 3 (b) is for adopting the present invention to process the differential interferometry phase diagram of the Generation of simulating data of interferometric phase image accuracy registration.
Contrast accompanying drawing 3(a) and accompanying drawing 3(b), large than Fig. 3 (b) of the fluctuation of differential interferometry phase place corresponding to deformation region in Fig. 3 (a), the differential interferometry phase noise produced in Fig. 3 (a) than Fig. 3 (b) greatly.Thereby can obtain, when the interferometric phase image accuracy registration, the differential interferometry phase diagram quality that the present invention obtains, higher than prior art, has self-adaptive recovery differential interferometry phase information, and the differential interferometry phase estimation is advantage accurately.
Fig. 4 adopts the present invention and prior art to process the differential interferometry phase diagram of the Generation of simulating data that the interferometric phase image registration error is 1 pixel, wherein, Fig. 4 (a) is for adopting prior art to process the differential interferometry phase diagram of the Generation of simulating data that the interferometric phase image registration error is 1 pixel, and Fig. 4 (b) is for adopting the present invention to process the differential interferometry phase diagram of the Generation of simulating data that the interferometric phase image registration error is 1 pixel.
Contrast accompanying drawing 4(a) and accompanying drawing 4(b), in Fig. 4 (a), all by noise, formed, do not comprise any useful information, and Fig. 4 (b) still can obtain differential interferometry phase diagram clearly, and differential interferometry phase place contrast of fringes is more clear.Thereby can obtain, when between interferometric phase image, having the error of 1 pixel, the present invention can obtain the good differential interferometry phase diagram of quality, significantly is better than prior art.Contrast accompanying drawing 3(b) and accompanying drawing 4(b), in the differential interferometry phase diagram of Fig. 4 (b), the differential interferometry phase fluctuation, only less times greater than Fig. 3 (b), can obtain thus, and the present invention has good robustness to the interferometric phase image registration error.
Below in conjunction with the measured data result, the present invention is described further.
Process measured data used, the ENVISAT satellite SAR imaging results to centralItaly L ' Aquila area on February 23rd, 2009 and May 4 respectively of European Space Agency, and the digital elevation diagram data that publishes of US National Aeronautics and Space Administration (NASA).
Fig. 5 is for adopting the present invention and prior art to process the differential interferometry phase diagram that measured data generates, wherein, Fig. 5 (a) is for adopting prior art to process the differential interferometry phase diagram that measured data generates, and Fig. 5 (b) is for adopting the present invention to process the differential interferometry phase diagram that measured data generates.
Contrast accompanying drawing 5(a) and accompanying drawing 5(b), in the differential interferometry phase diagram of Fig. 5 (a) only to see fuzzy differential interferometry phase place striped, and contain a lot of phase noises in figure, differential interferometry phase place striped is clearly arranged in the interferometric phase image of Fig. 5 (b), the differential interferometry phase place is comparatively level and smooth, only contain a small amount of interferometric phase noise, recovered the partial information of differential interferometry phase place.Can obtain thus, the present invention is in the situation that process the unknown of interferometric phase image registration error, still can obtain the differential interferometry phase diagram than prior art better quality, completed the operation of differential interferometry phase filtering simultaneously, illustrated that the present invention can adaptive recovery differential interferometry phase information, further verified the robustness of the present invention to registration error, and practicality of the present invention and validity have been described simultaneously.

Claims (8)

1. the differential interferometry phase estimation method based on many pixels covariance matrix, comprise the steps:
(1) input picture and parameter:
Before target area deformation 1a) respectively the differential interferometry synthetic-aperture radar obtained and the level land interferometric phase image that goes after deformation be input to system;
1b) systematic parameter of differential interferometry synthetic-aperture radar is input to system;
(2) interferometric phase image registration:
2a) deduct rear data line with the every data line in the interferometric phase image of level land that goes before deformation respectively, obtain removing level land row differential interferometry phase diagram before deformation;
2b) deduct a rear column data with each column data in the row differential interferometry phase diagram of level land that goes before deformation respectively, obtain the front interferometric phase gradient map of deformation;
2c) deduct rear data line with the every data line in the interferometric phase image of level land that goes after deformation respectively, obtain removing level land row differential interferometry phase diagram after deformation;
2d) respectively with after deformation go each column data in the row differential interferometry phase diagram of level land to deduct after a column data, obtain the interferometric phase gradient map after deformation;
The interferometric phase gradient map of 2e) take before deformation is reference picture, and the interferometric phase gradient map after deformation is carried out to registration process, before obtaining the deformation of registration and the side-play amount when interferometric phase gradient map after deformation and registration;
2f) by all pixels in the interferometric phase image of level land of going after deformation, a side-play amount during registration of translation, obtain before the deformation of thick registration interferometric phase image after interferometric phase image and deformation;
(3) obtain the interferometric phase image that secondary is interfered:
3a) by the data of all pixels in interferometric phase image before the deformation of thick registration, be multiplied by scale factor, obtain the front mark interferometric phase image that becomes of deformation of thick registration;
3b) by the data of all pixels in interferometric phase image after the deformation of thick registration, deduct the front data that become corresponding pixel points in the mark interferometric phase image of deformation of thick registration, obtain the interferometric phase image that secondary is interfered;
(4) obtain image window:
In the interferometric phase image of interfering at secondary, centered by pixel to be estimated, regular length is the piece radius, obtains a square window;
(5) build multiple interferometric phase associating data vector:
5a) choose arbitrarily a pixel in square window, the data of neighbor pixel around selected pixel and its are formed a line, obtain the phase vectors of selected pixel;
5b) by each element in the phase vectors of selected pixel, ask after multiplying each other with imaginary unit and take the index that natural Exponents e is the end, obtain the index vector of selected pixel;
5c) before each element of the index vector of selected pixel, insert positive integer 1, obtain the multiple interferometric phase associating data vector of pixel to be estimated;
5d) travel through all pixels in square window, judge whether to obtain the multiple interferometric phase associating data vector of all pixels, if so, perform step (6), otherwise, execution step 5a);
(6) adopt following formula, the estimate covariance matrix:
C = 1 M 2 Σ i = 1 M Σ j = 1 M l ( i , j ) × l ( i , j ) H
Wherein, C means estimated covariance matrix, and l (i, j) means the multiple interferometric phase associating data vector of the capable j column element of i in square window, i=1,2 ..., M, j=1,2,, M, M means the radius of square window, H means to do the conjugate transpose operation;
(7) estimate the coefficient of coherence matrix:
All elements value in covariance matrix is taken absolute value, obtain the coefficient of coherence matrix of pixel to be estimated;
(8) the coefficient of coherence matrix is carried out to the feature decomposition operation, obtain eigenwert and the proper vector of coefficient of coherence matrix;
(9) adopt following formula, the structure cost function:
Figure FDA0000378885730000022
Wherein, J means constructed cost function, mean the differential interferometry phase place
Figure FDA0000378885730000024
corresponding differential interferometry phase vectors, N means the number of large eigenwert in the eigenwert of coefficient of coherence matrix, β plarge eigenwert characteristic of correspondence vector in the proper vector of expression coefficient of coherence matrix, p=1,2 ..., N, β qlarge eigenwert characteristic of correspondence vector in the proper vector of expression coefficient of coherence matrix, q=1,2 ..., N, H means to do the conjugate transpose operation,
Figure FDA0000378885730000025
mean that Hadamard Hadamard is long-pending;
(10) estimate the difference interferometric phase:
(π, π] scope in, search differential interferometry phase place find the maximal value of cost function, by the corresponding differential interferometry phase place of cost function maximal value
Figure FDA0000378885730000032
differential interferometry phase place as pixel to be estimated.
2. the differential interferometry phase estimation method based on many pixels covariance matrix according to claim 1, it is characterized in that: the systematic parameter of the differential interferometry synthetic-aperture radar step 1b) comprises imaging number of times and base length; Described imaging number of times refers to the number of times of differential interferometry synthetic-aperture radar to the target area imaging, and base length refers to respectively with before deformation and the length of removing the differential interferometry synthetic-aperture radar baseline that the level land interferometric phase image is corresponding after deformation.
3. the differential interferometry phase estimation method based on many pixels covariance matrix according to claim 1, it is characterized in that: the concrete steps of the gradient map of interferometric phase step 2e) registration are as follows:
The first step, utilize following formula, obtains the front interferometric phase gradient map of deformation and the cross-correlation matrix of the interferometric phase gradient map after deformation:
L=IFFT2(FFT2(R)×conj(FFT2(S)))
Wherein, L means interferometric phase gradient map before deformation and the cross-correlation matrix of the interferometric phase gradient map after deformation, IFFT2 () means to do two-dimentional inverse Fourier transform operation, FFT2 () means to do the two-dimensional Fourier transform operation, R means the interferometric phase gradient map data matrix before deformation, S means the interferometric phase gradient map data matrix after deformation, and conj () means to get conjugate operation;
Second step, maximizing from all elements of cross-correlation matrix L, obtain coordinate figure corresponding to maximal value element in cross-correlation matrix L;
The 3rd step, with coordinate figure corresponding to maximal value element in cross-correlation matrix L, deduct coordinate figure corresponding to central element of cross-correlation matrix L, the side-play amount while obtaining registration;
The 4th step, by the interferometric phase gradient map after deformation, take the side-play amount of pixel during as a registration of unit translation, complete before deformation with deformation after the registration of interferometric phase gradient map.
4. the differential interferometry phase estimation method based on many pixels covariance matrix according to claim 1, it is characterized in that: the computing method of scale factor step 3a) are as follows:
The first step, when the difference interference synthetic aperture radar equals 2 to the imaging number of times of target area, scale factor equals 1;
Second step, when the difference interference synthetic aperture radar is greater than 2 to the imaging number of times of target area, according to the following formula, calculate scale factor:
k = B A
Wherein, k means scale factor, and B means the base length of going the differential interferometry synthetic-aperture radar that the level land interferometric phase image is corresponding after deformation, and A means the base length of going the differential interferometry synthetic-aperture radar that the level land interferometric phase image is corresponding before deformation.
5. the differential interferometry phase estimation method based on many pixels covariance matrix according to claim 1, it is characterized in that: regular length described in step (4) is 3~7 pixels.
6. the differential interferometry phase estimation method based on many pixels covariance matrix according to claim 1 is characterized in that: the number that is used for building the neighbor pixel of multiple interferometric phase associating data vector step 5a) is 1~8.
7. the differential interferometry phase estimation method based on many pixels covariance matrix according to claim 1, it is characterized in that: the preparation method of the differential interferometry phase vectors that differential interferometry phase place described in step (9) is corresponding is as follows:
The first step, by the differential interferometry phase place of pixel to be estimated, ask after multiplying each other with imaginary unit and take the index that natural Exponents e is the end, obtains the exponential quantity of pixel to be estimated;
Second step, sequentially form a line positive integer 1 and the exponential quantity of pixel to be estimated, obtains the single differential interferometry phase vectors of pixel to be estimated;
The 3rd step, choose the identical single differential interferometry phase vectors of number of using pixel when building multiple interferometric phase associating data vector, all elements in selected single differential interferometry vector is formed a line, obtain the differential interferometry phase vectors of pixel to be estimated.
8. the differential interferometry phase estimation method based on many pixels covariance matrix according to claim 1 is characterized in that: the number of large eigenwert described in step (9) is with to build the number that multiple interferometric phase combines data vector pixel used identical.
CN201310404990.0A 2013-09-08 2013-09-08 Differential interferometric phase estimation method based on multi-pixel covariance matrixes Active CN103454636B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201310404990.0A CN103454636B (en) 2013-09-08 2013-09-08 Differential interferometric phase estimation method based on multi-pixel covariance matrixes

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201310404990.0A CN103454636B (en) 2013-09-08 2013-09-08 Differential interferometric phase estimation method based on multi-pixel covariance matrixes

Publications (2)

Publication Number Publication Date
CN103454636A true CN103454636A (en) 2013-12-18
CN103454636B CN103454636B (en) 2015-05-20

Family

ID=49737211

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201310404990.0A Active CN103454636B (en) 2013-09-08 2013-09-08 Differential interferometric phase estimation method based on multi-pixel covariance matrixes

Country Status (1)

Country Link
CN (1) CN103454636B (en)

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105469368A (en) * 2015-11-30 2016-04-06 中国人民解放军国防科学技术大学 Interferometric phase filtering method
CN105974411A (en) * 2016-04-28 2016-09-28 电子科技大学 DINSAR (differential interferometric synthetic aperture radar)-based high-voltage transmission tower top tilt displacement monitoring method
CN103996175B (en) * 2014-05-13 2017-02-15 西安电子科技大学 Forest or urban area high-resolution interference phase filtering method
CN108303735A (en) * 2018-01-30 2018-07-20 单新建 The earthquake deformation acquisition methods of subband interferometry based on optimized parameter setting
CN108627832A (en) * 2018-05-11 2018-10-09 电子科技大学 A method of passway for transmitting electricity Ground Deformation is extracted based on multiple timings SAR image
CN108663678A (en) * 2018-01-29 2018-10-16 西北农林科技大学 More baseline InSAR phase unwrapping algorithms based on mixed integer optimization model
CN109541596A (en) * 2018-11-28 2019-03-29 中国电子科技集团公司电子科学研究院 InSAR image processing method and device based on deep learning algorithm
CN111595293A (en) * 2020-05-29 2020-08-28 山东交通学院 Surface deformation monitoring method and system based on multi-source monitoring data fusion
CN118549934A (en) * 2024-07-30 2024-08-27 云南大学 Phase estimation optimization method, device and equipment

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101551455A (en) * 2009-05-13 2009-10-07 西安电子科技大学 3D terrain imaging system of interferometric synthetic aperture radar and elevation mapping method thereof
CN101770027A (en) * 2010-02-05 2010-07-07 河海大学 Ground surface three-dimensional deformation monitoring method based on InSAR and GPS data fusion

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101551455A (en) * 2009-05-13 2009-10-07 西安电子科技大学 3D terrain imaging system of interferometric synthetic aperture radar and elevation mapping method thereof
CN101770027A (en) * 2010-02-05 2010-07-07 河海大学 Ground surface three-dimensional deformation monitoring method based on InSAR and GPS data fusion

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
ZHENFANG LI等: "Image autocoregistration and InSAR interferogram estimation using joint subspace projection", 《IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING》 *
廖明生等: "INSAR干涉条纹图的复数空间自适应滤波", 《遥感学报》 *

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103996175B (en) * 2014-05-13 2017-02-15 西安电子科技大学 Forest or urban area high-resolution interference phase filtering method
CN105469368A (en) * 2015-11-30 2016-04-06 中国人民解放军国防科学技术大学 Interferometric phase filtering method
CN105974411A (en) * 2016-04-28 2016-09-28 电子科技大学 DINSAR (differential interferometric synthetic aperture radar)-based high-voltage transmission tower top tilt displacement monitoring method
CN108663678A (en) * 2018-01-29 2018-10-16 西北农林科技大学 More baseline InSAR phase unwrapping algorithms based on mixed integer optimization model
CN108303735A (en) * 2018-01-30 2018-07-20 单新建 The earthquake deformation acquisition methods of subband interferometry based on optimized parameter setting
CN108627832A (en) * 2018-05-11 2018-10-09 电子科技大学 A method of passway for transmitting electricity Ground Deformation is extracted based on multiple timings SAR image
CN109541596A (en) * 2018-11-28 2019-03-29 中国电子科技集团公司电子科学研究院 InSAR image processing method and device based on deep learning algorithm
CN111595293A (en) * 2020-05-29 2020-08-28 山东交通学院 Surface deformation monitoring method and system based on multi-source monitoring data fusion
CN118549934A (en) * 2024-07-30 2024-08-27 云南大学 Phase estimation optimization method, device and equipment

Also Published As

Publication number Publication date
CN103454636B (en) 2015-05-20

Similar Documents

Publication Publication Date Title
CN103454636B (en) Differential interferometric phase estimation method based on multi-pixel covariance matrixes
CN110058237B (en) InSAR point cloud fusion and three-dimensional deformation monitoring method for high-resolution SAR image
CN104123464B (en) Method for inversion of ground feature high elevation and number of land subsidence through high resolution InSAR timing sequence analysis
Goel et al. A distributed scatterer interferometry approach for precision monitoring of known surface deformation phenomena
JP6106285B2 (en) Ion distortion correction method and apparatus for satellite radar interference degree
CN103885059B (en) A kind of multi-baseline interference synthetic aperture radar three-dimensional rebuilding method
CN103439708B (en) Polarized InSAR interferogram estimation method based on generalized scattering vector
CN106772342A (en) A kind of Timing Difference radar interference method suitable for big gradient surface subsidence monitoring
CN109212522B (en) High-precision DEM inversion method and device based on double-base satellite-borne SAR
Ren et al. Calculating vertical deformation using a single InSAR pair based on singular value decomposition in mining areas
CN109031301A (en) Alpine terrain deformation extracting method based on PSInSAR technology
CN104007439B (en) Interferential circular SAR elevation estimation processing method
CN109239710B (en) Method and device for acquiring radar elevation information and computer-readable storage medium
CN102680972A (en) Method and device for monitoring surface deformation and data processing equipment
CN111650579B (en) InSAR mining area three-dimensional deformation estimation method and device for rock migration parameter adaptive acquisition and medium
CN104482877A (en) Motion compensation method and system in three-dimensional imaging of dynamic object
CN115201825B (en) Atmospheric delay correction method in InSAR (interferometric synthetic aperture radar) inter-seismic deformation monitoring
CN103018741A (en) Interferometric synthetic aperture radar (InSAR) imaging and flat ground removing integral method based on back projection
CN110109105A (en) A method of the InSAR technical monitoring Ground Deformation based on timing
CN103308914B (en) One-station fixed bistatic interference synthetic aperture radar (SAR) processing method
CN103018740A (en) InSAR (Interferometric Synthetic Aperture Radar) imaging method based on curved surface projection
CN114858238A (en) Underground water storage and transformation monitoring and early warning system
CN112051572A (en) Method for monitoring three-dimensional surface deformation by fusing multi-source SAR data
CN103616682A (en) Multi baseline InSAR (Interferometric Synthetic Aperture Radar) processing method based on curved surface projection
CN107144213A (en) The big magnitude three-D sequential deformation method of estimation in mining area and device based on SAR intensity images

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant