forked from GengPeiLin/gmtsar2stamps
-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathmt_extract_cands_gmtsar
executable file
·71 lines (60 loc) · 2.22 KB
/
mt_extract_cands_gmtsar
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
#!/bin/bash
#
# extract data for candidiate pixels
#
threshold="0.4"
region="0/2000/0/10000"
reg="0/2000/0/10000"
gmt grdmath scatter.grd $threshold LT 0 NAN scatter.grd MUL = cands0.grd
gmt grdsample -T cands0.grd -R$region -Gcands.grd
cp cands.grd cands_old.grd
gmt grdedit cands.grd -R$reg
input="cands.grd"
width=`gmt grdinfo -C $input | awk '{print $10}'`
echo $width > width.txt # range
length=`gmt grdinfo -C $input | awk '{print $11}'`
echo $length > len.txt # azimuth
# Get the dispersion index for PS candidates
gmt grd2xyz cands.grd -s > cands.dat
gmt gmtconvert cands.dat -o2 > pscands.1.da
awk '{printf("%d %d %d\n", NR,$2,$1)}' cands.dat > pscands.1.ij
# Retrieve lon/lat for PS candidates
gmt grd2xyz cands_old.grd -s > cands_old.dat
proj_ra2ll_ascii.csh trans.dat cands_old.dat cands.ll.dat
rm ralt.grd raln.grd
gmt gmtconvert cands.ll.dat -o0,1 -bos > pscands.1.ll
# Retrieve hgt for PS candidates
gmt grdtrack cands.ll.dat -Gdemfilt.grd -Z -nb -bos > pscands.1.hgt
# Retrieve look angles
gmt grdsample cands_old.grd -I50+/50+ -Gtmp.grd
gmt grd2xyz tmp.grd | awk '{print $1,$2,1}' > tmp.xy
proj_ra2ll_ascii.csh trans.dat tmp.xy tmp.ll
#rm ralt.grd raln.grd
awk '{printf("%f %f %d\n%f %f %d\n",$1,$2,0,$1,$2,1000)}' tmp.ll > tmp.llt
ln -s ../topo/master.PRM .
SAT_look master.PRM < tmp.llt > tmp.lltn
gmt gmtmath -C5 -o5 tmp.lltn ACOS 180 MUL PI DIV = look_angle.1.in
rm tmp.xy tmp.ll tmp.llt tmp.lltn
# Retrieve perpendicular baseline for each interferogram
cat list | while read name prm date
do
cp ../raw/*PRM .
SAT_baseline master.PRM $prm.PRM >> $prm.PRM
bperp $prm.PRM tmp.grd bperp.grd
gmt grd2xyz bperp.grd | awk '{print $3}' > bperp0
value=`grep earth_radius $prm.PRM | awk '{print $3+1000}'`
update_PRM.csh $prm.PRM earth_radius $value
bperp $prm.PRM tmp.grd bperp.grd
gmt grd2xyz bperp.grd | awk '{print $3}' > bperp1000
paste bperp0 bperp1000 > tmp
awk '{printf("%f\n%f\n", $1,$2)}' tmp > bperp_$date.1.in
done
rm tmp bperp0 bperp1000 bperp.grd tmp.grd
# Retrieve phase for PS candidates
rm pscands.1.ph
cat list | while read name prm date
do
gmt grdtrack cands_old.dat -Gre_$name.grd -Z -bos >> pscands.1.ph
gmt grdtrack cands_old.dat -Gim_$name.grd -Z -bos >> pscands.1.ph
done
rm cands_old.grd