Closed hfattahi closed 1 year ago
This seems good in my testing so far. There seems to be an very noticeable difference between the MLE-only version and this new version on L-band data (higher long-term coherence) with larger stack size (needs more multilook samples to have a positive definite covariance matrix).
In this example, the left side of the image is a rocky mountainous area with the highest coherence, and the bottom right is city. The MLE version is dropping most of the good pixels in those areas, while this PR looks much better.
This seems good in my testing so far. There seems to be an very noticeable difference between the MLE-only version and this new version on L-band data (higher long-term coherence) with larger stack size (needs more multilook samples to have a positive definite covariance matrix).
In this example, the left side of the image is a rocky mountainous area with the highest coherence, and the bottom right is city. The MLE version is dropping most of the good pixels in those areas, while this PR look much better.
![]()
Thanks @scottstanie for trying this PR.
This PR adds a new workflow to estimate wrapped phase series. This module named "phase-link" tries to estimate the wrapped phase series for a given neighborhood based on EMI estimation. If the coherence matrix is positive definite and computing inverse is doable, the EMI estimation proceeds. If the coherence is not positive definite and computing its inverse is not doable, instead of giving up with any estimate, here we switch the gear to EVD (Eigenvalue and Vector Decomposition) of Covariance matrix. This significantly improves the estimation and gives wrapped phase estimates for many pixels which EMI alone may miss. In practice the positive definiteness check for the coherence matrix can fail even for cases with very high coherence. I experienced many examples where the coherence matrix is not positive definite and therefore the EMI estimation is not doable while falling back to evd leads to reliable phase estimate with very high temporal coherence of more than 0.9.
Here is an example of a coherence matrix which EMI fails to give an estimate because the coherence is not positive definite. Obviously the coherence is very high and we should provide an estimate.
The results of the new module (EMI + EVD) compared with the EMI estimation for one example dataset can be seen here where the first plot shows EMI only and the second plot shows EMI + EVD from phase-link module:
Note that the light blue pixels are all those pixels that no estimation is provided for by EMI. The example is for one ministack. This improvement is even more critical for workflows which would run on many mini-stacks. If randomly bunch of pixels are unnecessarily masked in each mini-stack estimation, the final results will have a lot more pixels excluded.
In this PR I have not bothered to regularize the coherence matrix when it is not positive definite. Down the line we need to evaluate how efficient the regularization is and what is the gain compared to simply falling back to EVD when coherence is not positive definite. Nevertheless we know that regularization is not always doable and we will still need to fall back to EVD in many cases.
Currently a lot of code is duplicated from evd. I had to make a choice with code duplication or modifying evd. I didn't want to make evd over complicated. Moreover I have relaxed the checks on the coherence and covariance matrices. The only checks that I have kept in phase-link are checking the positive-definiteness of the coherence matrix and check if Lapack can estimate the EVD of Covariance matrix. The rest of the checks seem un-necessary to me and were leading to skipping the estimation for many highly coherent neighborhoods. I think it is ok for short term to have the code duplication. Over long term it would be nice to have EMI and EVD core modules which can be called from phase-link.cpp or evd.cpp. Such changes coupled with handling rasters outside the modules will be a future clean up. Anyways for now I chose the code duplication for further testing the new workflow against the existing estimators in evd.cpp.