forked from insarlab/MintPy
-
Notifications
You must be signed in to change notification settings - Fork 0
/
variance.py
97 lines (79 loc) · 3.19 KB
/
variance.py
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
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
"""Utilities to calculate structural functions."""
############################################################
# Program is part of MintPy #
# Copyright (c) 2013, Zhang Yunjun, Heresh Fattahi #
# Author: Zhang Yunjun, 2017 #
############################################################
# Recommend usage:
# from mintpy.simulation import variance as var
import numpy as np
import pyproj
from mintpy.utils import ptime
def sample_data(lat, lon, mask=None, num_sample=500):
''''''
## Flatten input data
lat = lat.flatten()
lon = lon.flatten()
## Check number of samples and number of pixels
num_pixel = len(lat)
if num_sample > num_pixel:
print('Number of samples > number of pixels, fix number of samples to number of pixels.')
num_sample = num_pixel
# Check input mask
if mask is None:
mask = np.ones((num_pixel), dtype=np.bool_)
mask = mask.flatten()
# Random select samples
idx = np.arange(num_pixel)[mask]
rng = np.random.default_rng()
idx_sample = rng.choice(idx, size=int(num_sample))
lat_sample = lat[idx_sample]
lon_sample = lon[idx_sample]
return idx_sample, lat_sample, lon_sample
def get_distance(lat, lon, i):
'''Return the distance of all points in lat/lon from its ith point'''
lat1 = lat[i]*np.ones(lat.shape)
lon1 = lon[i]*np.ones(lon.shape)
g = pyproj.Geod(ellps='WGS84')
dist = g.inv(lon1, lat1, lon, lat)[2]
return dist
def structure_function(data, lat, lon, step=5e3, min_pair_num=100e3, print_msg=True):
num_sample = len(data)
distance = np.zeros(num_sample**2)
variance = np.zeros(num_sample**2)
if print_msg:
prog_bar = ptime.progressBar(maxValue=num_sample)
for i in range(num_sample):
distance[i*num_sample:(i+1)*num_sample] = get_distance(lat, lon, i)
variance[i*num_sample:(i+1)*num_sample] = np.square(data - data[i])
if print_msg:
prog_bar.update(i+1, every=10)
if print_msg:
prog_bar.close()
(bin_dist,
bin_struct_func,
bin_struct_func_std) = bin_variance(distance, variance,
step=step,
min_pair_num=min_pair_num,
print_msg=print_msg)
return bin_dist, bin_struct_func, bin_struct_func_std
def bin_variance(distance, variance, step=5e3, min_pair_num=100e3, print_msg=True):
x_steps = np.arange(0,np.max(distance),step)
num_step = len(x_steps)
var = np.zeros(x_steps.shape)
var_std = np.zeros(var.shape)
p_num = np.zeros(x_steps.shape)
if print_msg:
prog_bar = ptime.progressBar(maxValue=num_step)
for i in range(num_step):
x = x_steps[i]
idx = (distance > max(0, x-step/2.)) * (distance < x+step/2.)
p_num[i] = np.sum(idx)
var[i] = np.mean(variance[idx])
var_std[i] = np.std(variance[idx])
if print_msg:
prog_bar.update(i+1, every=10)
if print_msg:
prog_bar.close()
max_step_idx = int(max(np.argwhere(p_num > min_pair_num)))
return x_steps[0:max_step_idx], var[0:max_step_idx], var_std[0:max_step_idx]