forked from Narayana-Rao/PolSAR-tools
-
Notifications
You must be signed in to change notification settings - Fork 0
/
mod_dop_cp.py
150 lines (109 loc) · 4.92 KB
/
mod_dop_cp.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
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
from qgis.PyQt.QtCore import *
from qgis.PyQt.QtWidgets import *
from qgis.PyQt.QtGui import *
from qgis.PyQt import *
from qgis.core import *
import requests
import numpy as np
import multiprocessing
# Initialize Qt resources from file resources.py
from .resources import *
# Import the code for the dialog
from .SAR_Tools_dialog import MRSLabDialog
import os.path
from osgeo import gdal
import time
import os.path
##############################################################################################
class dop_cp(QtCore.QObject):
'''DOP CP '''
def __init__(self,iFolder,C2,ws,tau):
QtCore.QObject.__init__(self)
self.iFolder = iFolder
self.C2 = C2
self.ws=ws
self.tau=tau
self.killed = False
# self.mainObj = MRSLab()
def conv2d(self,a, f):
filt = np.zeros(a.shape)
wspad = int(f.shape[0]/2)
s = f.shape + tuple(np.subtract(a.shape, f.shape) + 1)
strd = np.lib.stride_tricks.as_strided
subM = strd(a, shape = s, strides = a.strides * 2)
filt_data = np.einsum('ij,ijkl->kl', f, subM)
filt[wspad:wspad+filt_data.shape[0],wspad:wspad+filt_data.shape[1]] = filt_data
return filt
def run(self):
finish_cond = 0
try:
def dopcp_fn(C2_stack,ws):
if self.tau==0:
chi_in = -45.0
else:
chi_in = 45.0
kernel = np.ones((ws,ws),np.float32)/(ws*ws)
c11_T1 = C2_stack[:,:,0]
c12_T1 = C2_stack[:,:,1]
c21_T1 = C2_stack[:,:,2]
c22_T1 = C2_stack[:,:,3]
c11_T1r = self.conv2d(np.real(c11_T1),kernel)
c11_T1i = self.conv2d(np.imag(c11_T1),kernel)
c11s = c11_T1r+1j*c11_T1i
c12_T1r = self.conv2d(np.real(c12_T1),kernel)
c12_T1i = self.conv2d(np.imag(c12_T1),kernel)
c12s = c12_T1r+1j*c12_T1i
self.pBar.emit(25)
c21_T1r = self.conv2d(np.real(c21_T1),kernel)
c21_T1i = self.conv2d(np.imag(c21_T1),kernel)
c21s = c21_T1r+1j*c21_T1i
c22_T1r = self.conv2d(np.real(c22_T1),kernel)
c22_T1i = self.conv2d(np.imag(c22_T1),kernel)
c22s = c22_T1r+1j*c22_T1i
self.pBar.emit(50)
# c2_det = (c11s*c22s-c12s*c21s)
# c2_trace = c11s+c22s
# t2_span = t11s*t22s
# m1 = np.real(np.sqrt(1.0-(4.0*c2_det/np.power(c2_trace,2))))
# Stokes Parameter
s0 = c11s + c22s;
s1 = c11s - c22s;
s2 = (c12s + c21s);
if (chi_in >= 0):
s3 = (1j*(c12s - c21s)); # The sign is according to RC or LC sign !!
if (chi_in < 0):
s3 = -(1j*(c12s - c21s)); # The sign is according to RC or LC sign !!
dop= np.sqrt(np.power(s1,2) + np.power(s2,2) + np.power(s3,2))/(s0);
self.pBar.emit(90)
self.progress.emit('->> Write files to disk...')
"""Write files to disk"""
infile = self.iFolder+'/C11.bin'
ofiledop = self.iFolder+'/DOP_CP.bin'
write_bin(ofiledop,dop,infile)
self.pBar.emit(100)
self.progress.emit('->> Finished DOP calculation!!')
def write_bin(file,wdata,refData):
ds = gdal.Open(refData)
[cols, rows] = wdata.shape
driver = gdal.GetDriverByName("ENVI")
outdata = driver.Create(file, rows, cols, 1, gdal.GDT_Float32)
outdata.SetGeoTransform(ds.GetGeoTransform())##sets same geotransform as input
outdata.SetProjection(ds.GetProjection())##sets same projection as input
outdata.SetDescription(file)
outdata.GetRasterBand(1).WriteArray(wdata)
# outdata.GetRasterBand(1).SetNoDataValue(np.NaN)##if you want these values transparent
outdata.FlushCache() ##saves to disk!!
# self.dop_fp(self.T3)
dopcp_fn(self.C2,self.ws)
finish_cond = 1
except Exception as e:
# forward the exception upstream
self.error.emit(e, traceback.format_exc())
self.finished.emit(finish_cond)
# def kill(self):
# self.killed = True
"""***************************************"""
finished = QtCore.pyqtSignal(object)
error = QtCore.pyqtSignal(Exception, str)
progress = QtCore.pyqtSignal(str)
pBar = QtCore.pyqtSignal(int)