---
lang: zh-tw
tags: timeseries, Stamps, isce,tutorial, tomosar
title: ISCE-stamps/MTI-TomoSAR
author: Jun-Yan Chen
time: 20210711
---

# ISCE - Stamps/MT-InSAR/PSDS InSAR

:::warning
- 閱讀前請先建立影像堆，可參考前篇: [TOPS-stack for timeseries](/NfHR48A1Q8u2LiceJH84OQ)
- 建議安裝isce2.3 (master/slave)版本，才可以正常在StaMPS中運作，如果想要在stamps中運行ISCE2.4 (reference/secondary) 版本，需要將所有資料夾與檔案改為master/slave，可能會需要修改stamps的matlab程式碼。
- 執行前需要安裝好StaMPS，安裝請參考：[StaMPS/MTI](/hxddWhmLRuS5M6MGoi9W8Q)
- 下文主要參考[isce2stamps manual](https://usermanual.wiki/Document/Manual.467280380/help)，並經過測試後修正流程
- 以下僅介紹Stamps基本執行與需要使用到ISCE的步驟，相關詳細Stamps指令與繪圖、輸出，請參考 [StaMPS/MTI](/hxddWhmLRuS5M6MGoi9W8Q)
:::

### Step 0: check ISCE result and processing concept

在開始之前，須注意mereged/SLC與merged/baselines 資料夾中的檔案數目是否相同，若否（通常baselines較SLC多)可以先重新處理沒有產生的檔案，或是移除多出來的baselines/yyyymmdd檔案。

如果沒有執行，會導致後續的嚴重錯誤！！！！！！！！！！！！！！！！
e.g. stamps step-1 load data error
e.g. ps_sb_merge error

- isce的可能失敗原因:沒有提取出足夠的burst，進而無法正確對位
- 可以用下列指令來確定哪個檔案不存在！
``for i in baselines/*; do if [ ! -d SLC/`basename $i` ]; then echo $i ; fi ; done``

如果已經執行完make_single_reference_isce處理完，此時可以直接修改檔案 bperp.1.in 與 day.1.in (此方法僅限PS，小基線請將上述的檔案刪除後重跑)

# update day file

for i in $(cat slcs.list)
do
awk '{print $1 }' $i/baseline 
done > day.1.in


# update bperp file
for i in $(cat slcs.list)
do
awk '{print $2 }' $i/baseline 
done > bperp.1.in

#### workflow
左邊是標準PS流程，右邊則是處理MT-InSAR之流程(青色代表在matlab中執行)。

### Step 1: crop stack to AOI (optional)

使用 contrib/timeseries/prepStackToStaMPS/run_SLCcropStack.csh 指令檔來裁切。
將其複製到目標資料夾內(任意位置即可)，並編輯行：
> NOTE
> 若研究區域在山區如果建議要設定較大的裁切區域，因為crop_rdr.py是假設雷達座標與地理座標為線性變換，因此誤差可能會相當大，或是修改crop_rdr.py之演算法。

set overwrite = 0
# 改為isce的merged資料夾
set data_path = /u/k-data/dbekaert/HMA_nepal/Sentinel1/track_019/processing_1/merged

#改為任意新資料夾
set proc_dir = /u/k-data/dbekaert/HMA_nepal/Sentinel1/track_019/processing_1/crop_testing

# check if the procesing dir already exists
if (! -e $proc_dir) then
   mkdir $proc_dir
endif

# getting the crop extend
cd $data_path/geom_reference

# 設定裁切範圍
crop_rdr.py -b '27.86 28.2 85.1 85.4' > $proc_dir/crop_log.txt

# -b: SNWE 裁切範圍 
# stripmap裁切可能需要指定：
# -lon lon.rdr -lat lat.rdr 

#### NO changes required below ######

# getting the files to crop
cd $proc_dir

# 若使用stripmap 影像，可能要手動刪除.full的檔案後綴
ls -1 $data_path/geom_reference/*.full  > $proc_dir/geomFiles2crop.txt
ls -1 $data_path/SLC/2*/2*.slc.full > $proc_dir/slcFiles2crop.txt
ls -1 $data_path/baselines/2*/2*.full.vrt >  $proc_dir/baselineFiles2crop.txt

修改結束後，執行 `run_SLCcropStack.csh`
產生crop_testing 資料夾，進入後便可以進入下一步。

---

### Step2: Prepare ISCE Stack (required)

建立input_file，包含內容如下:

source_data slc_stack
slc_stack_path $PROJECT/merged/SLC 
slc_stack_reference yyyymmdd
slc_stack_geom_path $PROJECT/merged/geom_reference
slc_stack_baseline_path $PROJECT/merged/baselines
range_looks 9
azimuth_looks 3
aspect_ratio 3
lambda 0.0555
slc_suffix .full
geom_suffix .full

[注意]
- 執行前修改日期、路徑
- 不同的sensor，可能會需要改變lambda (波長)
- 使用stripmapStack所產製的影像，不會有suffix，可空白。

並執行指令:
- `make_single_master_stack_isce` 或 `make_single_reference_stack_isce`
- 結束後會產生INSAR_yyyymmdd的資料夾
- 如果想要執行single-master PSInSAR 可以進入step 3，若想要進行MT-InSAR則進入step 2-2。

### Step 2-2(optional): Small baseline prepare
:::info
本步驟是在選取SDFT pixel，可參考:
Hooper, A. (2008). A multi-temporal InSAR method incorporating both persistent scatterer and small baseline approaches. Geophysical Research Letters, 35(16). https://doi.org/10.1029/2008GL034654
:::
1. 進入INSAR_yyyymmdd 資料夾，執行:
`mt_extract_info_isce`
2. 進入matlab執行，執行後會繪圖顯示基線的連結狀況，調整到符合預期即可。

%%%% matlab %%%%%
ps_load_info;
sb_find(rho, Ddiff_max, Bdiff_max);
% rho 為 baseline coherence (Default 0.5)
% Ddiff_max: 最大時間基線
% Bdiff_max: 最大空間基線,對Sentinel-1影像來說可以到非常小。

3. 在Bash中執行:

`make_small_baselines_isce`
- 結束後會產生SMALL_BASELINES 資料夾
- 進入SMALL_BASELINES，後執行下方步驟(同single-master PSInSAR)

### Step3: PS Candidate Selection (Required)
- 進入INSAR_yyyymmdd(小基線進入SMALL_BASELINES) 資料夾後執行：

mt_prep_isce threshold [rg_patch az_patch rg_overlap az_overalp]
# threshold: 振幅閾值(必要), PS: 0.4, SBAS: 0.6
# rg_patch: 處理時是否要在range方向分割
# az_patch: 處理時是否要在azimuth方向分割
# rg_overlap
# az_overlap

* 當干涉影像數目超過1024張時，會產生`core dump (segmentation fault)`，可以透過調高可以系統開檔數量來解決(通常是小基線才會發生):
`ulimit -n 2000`

### Step4-1: run stamps in matlab
*詳細步驟說明可以參考[StaMPS/MTI](/hxddWhmLRuS5M6MGoi9W8Q)*
- 可以在一般PS在InSAR_xxxxxxxx,小基線在SMALL_BSELINE中，開啟matlab ，執行:

stamps %run all steps
stamps(start,end) %run from step start to end (number)

- 背景執行(bash):
`$ nohup matlab -nodesktop -nosplash -r 'stamps; exit;' >nohup.out 2>nohup.err &`

### Step4-2: Merge PS-SBAS (optional)

- 如果想要合併PS點與SDFP點，可以在各自執行完Step 5之後，後，在INSAR_yyyymmdd中執行matlab指令:`ps_sb_merge`
- 執行結束後產生MERGED資料夾，並進入該資料夾中完成剩下的步驟
`stamps(6,8)`

*可能的失敗原因:一開始SLC與baselines的數量沒有檢查所導致的錯誤...

### Step5: Plot data
*請參考[StaMPS/MTI](/hxddWhmLRuS5M6MGoi9W8Q)*

# TomoSAR Processing
:::info Update: 2022/7/2
Code can be download from this github page: [DinhHoTongMinh/TomoSAR](https://github.com/DinhHoTongMinh/TomoSAR)

這個程式碼**只能**用來處理**小區域**的變形(如山崩、山區等等)，如果影像範圍太大，在matlab中會發生記憶體不足的問題

> 未來需要有熱誠的人把整個程式移植到python，使用H5PY格式來改寫讀檔才可能有效處理，與盡可能的平行化

:::spoiler Reference
Dinh Ho Tong Minh and Yen Nhi Ngo. "Compressed SAR Interferometry in the Big Data Era". Remote Sensing. 2022, 14, 390. [https://www.mdpi.com/2072-4292/14/2/390](https://www.mdpi.com/2072-4292/14/2/390)
:::

### Step 1: 裁切isce stack
*tomosar 一定要進行裁切才能夠執行*

### Step 2: make_single_master_stack_isce
同上一小節內所述，建立input_file，執行指令。

### Step 3: TomoSAR 參數設定

進入TomoSAR的資料夾內，修改matalb指令檔Parameter_input.m

% These follow parameters can be good for most areas.
Alpha = 0.05;  % significance level. 0.05 for 5% significance.
BroNumthre = 20; % less than 20 is likely PS
Cohthre = 0.25; % threshold to select DS in which its phase variance is mostly less than 20 degree.
Cohthre_slc_filt = 0.05; % less than 0.05 is mostly water

%%%%%%%%%%%%%%%%%%%%%%%%% For ComSAR analysis  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PSDSInSAR is heavy on memory use due to full covariance estimation. A rough approximation for
% PSDSInSAR RAM requirement is 1.5*Nslc*Nslc*Nline*Nwidth/2.7e8 (GB)
% ComSAR is much friendly Big Data processing. A rough approximation for
% ComSAR RAM requirement is 0.3*Nslc*Nslc*Nline*Nwidth/2.7e8 (GB)
% i.e., 200 images of 500x2000 size, 220 GB is for PSDS, but for ComSAR it requires only 45 GB.
ComSAR_flag = true; % true is for ComSAR, false is for PSDSInSAR
miniStackSize = 5; % 5 (or 10) can help to reduce up to 80% (or 90%) computation.
Unified_flag = true; % true is for full time series ComSAR, false is just for compressed version
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

InSAR_processor = 'isce'; % snap or isce
...
InSAR_path='/home/junyan1998/R103/Taichung/PSTC01/crop_merged/INSAR_20210320'
reference_date='20210320'

將processor改為isce，設定InSAR_path與reference_date後即可，其他參數可後續在做調整。

- 完成後可直接在matlab中執行：
`PSDS_main`
- 或背景執行：
`$ nohup matlab -nodesktop -nosplash -r 'PSDS_main; exit;' >nohup.out 2>nohup.err &`
- ComSAR 主要有四個步驟: 分別為SHP_select, phase linking, phase filtering, despeckle

> **Problem issue**
> 在PSDS_main執行時，他會想要寫出YYYYMMDD/secondary.slc.csar以及reference/reference.slc.casr，但未知雲因導致漏掉檔案(原因不明，也可能是自m)，但是會導致在後續的跑不動
> Also, code不知為何會想要寫入reference date的資料夾，但是理論上不應該存在

### Step 4: mt_prep_isce_comsar

與mt_prep_isce 相同