-
Notifications
You must be signed in to change notification settings - Fork 1
/
Potential_changes.py
143 lines (122 loc) · 5.18 KB
/
Potential_changes.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
#!/usr/bin/env python
# -*- coding: utf-8 -*-
# @Time : 2018/7/27 14:07
# @Author : Deyu.Tian
# @Site : ChangGuang Satellite
# @File : Potential_changes.py
# @Software: PyCharm Community Edition
import Image
import Image_process
import numpy as np
import sys
import os
import Config
import Visulazation as visul
def getPCC(diff_map):
pass
def Yuzhi_band(bandpath):
Image_process.draw_hist(bandpath)
pass
def normalize_arr(inarr):
"""
normalize array to [0, 1]
:param inarr:
:return:
"""
return (inarr-np.max(inarr))/(np.max(inarr)-np.min(inarr))
pass
def diff_map_matrix(imsatpath, imuavpath, w=0):
"""
use matrix to faster computation
:param imsatpath:卫星图像的路径
:param imuavpath:无人机图像的路径
:param w:搜索邻域块的大小
:return:两幅图像的差值图
"""
if w == 0 or w % 2 == 0:
print("use this function like below: \n"
"diff_map(imsatpath, imuavpath, w)\n need param w of search box,"
" its must be odd like 1, 3, 5, 7, 9,.etc.")
sys.exit(1)
imsatarr, imuavarr = np.moveaxis(Image.img2array(imsatpath), 0, 2), np.moveaxis(Image.img2array(imuavpath), 0, 2)
print(imsatarr.shape, imuavarr.shape)
visul.visul_arr_rgb(imsatarr)
visul.visul_arr_rgb(imuavarr)
imsat_desc, imuav_desc = descriptor(imsatarr), descriptor(imuavarr)
print("features of IMUAV and IMSAT have generated successed, "
"now compute difference maps with matrix ways!")
diffarr = np.zeros((imuavarr.shape[0], imuavarr.shape[1]), dtype='f')
for i in range((w-1)/2, imuav_desc.shape[0]-(w-1)/2):
for j in range((w-1)/2, imuav_desc.shape[1]-(w-1)/2):
xmin, xmax = j-(w-1)/2, j+(w-1)/2+1
ymin, ymax = i-(w-1)/2, i+(w-1)/2+1
v = imsat_desc[ymin:ymax, xmin:xmax]
v_2d = v.reshape((v.shape[0] * v.shape[1]), v.shape[2])
u = imuav_desc[i, j].reshape(1, imuav_desc[i, j].shape[0])
diffarr[i, j] = np.min(np.linalg.norm(v_2d - u, axis=1))
#print(v.shape, v_2d.shape, u.shape, diffarr[i, j])
#print(imuav_desc[i, j], u, diffarr[i, j])
np.save(os.path.join(Config.data, "diffmap.npy"), diffarr)
imggt = Image.read_tif_metadata(imsatpath)
Image.array2rasterUTM(os.path.join(Config.data, "diffmap.tif"), imggt, diffarr)
return diffarr
def diff_map(imsatpath, imuavpath, w=0):
"""
input imuav and imsat, return different map
:param imsatpath:
:param imuavpath:
:param w:search box in imsat, w must be odd
:return:
"""
if w == 0 or w % 2 == 0:
print("use this function like below: \n"
"diff_map(imsatpath, imuavpath, w)\n need param w of search box,"
" its must be odd like 1, 3, 5, 7, 9,.etc.")
sys.exit(1)
imsatarr, imuavarr = np.moveaxis(Image.img2array(imsatpath), 0, 2), np.moveaxis(Image.img2array(imuavpath), 0, 2)
print(imsatarr.shape, imuavarr.shape)
imsat_desc, imuav_desc = descriptor(imsatarr), descriptor(imuavarr)
print("descriptors of IMUAV and IMSAT have generated successed, now compute difference maps!")
diffarr = np.zeros((imuavarr.shape[0], imuavarr.shape[1]), dtype='f')
for i in range(imuav_desc.shape[0]):
for j in range(imuav_desc.shape[1]):
xmin, xmax = int(j-(w-1)/2), int(j+(w-1)/2+1)
ymin, ymax = int(i-(w-1)/2), int(i+(w-1)/2+1)
distsset = []
for x in range(xmin, xmax):
for y in range(ymin, ymax):
if imuav_desc.shape[1] > x > 0 and 0 < y < imuav_desc.shape[0]:
distsset.append(np.linalg.norm(imsat_desc[x, y] - imuav_desc[i, j])) #norm1 of vector/matrix
diffarr[i, j] = min(distsset)
np.save(os.path.join(Config.data, "diffmap.npy"), diffarr)
return diffarr
def descriptor(rgbarr):
"""
input rgbarr, return relative 36d features set
:param rgbarr:
:return:
"""
desc = np.zeros((rgbarr.shape[0], rgbarr.shape[1], 36), dtype='i')
print(desc.shape)
# visul.visul_arr_rgb(rgbarr[:, :, 0])
desc[:, :, 0:9] = Image_process.band_9_neibour_layers_ignoreedge(rgbarr[:, :, 0])
desc[:, :, 9:18] = Image_process.band_9_neibour_layers_ignoreedge(rgbarr[:, :, 1])
desc[:, :, 18:27] = Image_process.band_9_neibour_layers_ignoreedge(rgbarr[:, :, 2])
luminate = rgbarr[:, :, 0] * 0.299 + rgbarr[:, :, 1] * 0.587 + rgbarr[:, :, 2] * 0.114
desc[:, :, 27:36] = getIG(np.uint8(luminate)) #if( sdepth == CV_16S && ddepth == CV_32F)
print("max of IG is ", np.max(desc[:, :, 27:36]))
return desc
def getIG(lumimgarr):
"""
get IG(stacked layers of pixel-and-its 9-neibours' gradients )
:param lumimgarr:array of Y
:return:stacked band layers of gradients
"""
gradient = Image_process.sobel(lumimgarr, kernel_size=7)
return Image_process.band_9_neibour_layers_ignoreedge(gradient)
if __name__ == '__main__':
diff_map_matrix(Config.ImSat, Config.ImUAV, 21)
# imsatarr, imuavarr = np.moveaxis(Image.img2array(Config.ImSat), 0, 2), \
# np.moveaxis(Image.img2array(Config.ImUAV), 0, 2)
# print(imsatarr.shape, imuavarr.shape)
# descriptor(imsatarr)