get_sqear_prec_by_geo_coord_polygon.py
9.68 KB
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
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
import urllib.request
import json
import os
import datetime
from imageio import imread
import cv2
import numpy as np
png_path = r"d:\PYTHON\tests\test1\png"
png_result_path = r"d:\PYTHON\tests\test1\cutted"
lenght_lan = 79497.3714882
lenght_lon = 111162.6
base_pixel_len = 152.8740566
str_limiter = '||'
d1plan_g = base_pixel_len / lenght_lan
d1plon_g = base_pixel_len / lenght_lon
x286lon = 38.72287
x286lan = 44.417242
base_lon = x286lon - (128 * d1plon_g)
min_lon = x286lon - (128 * d1plon_g)
max_lon = x286lon + (128 * d1plon_g)
base_lan = x286lan + (128 * d1plan_g)
min_lan = x286lan - (128 * d1plan_g)
max_lan = x286lan + (128 * d1plan_g)
# Генерировать измерители
def GenerateMeasurements():
result = []
l_dlan = d1plan_g * 48
l_dlon = d1plon_g * 48
'''
result.append(GetMeas(1, 'X286|sp_dm', x286lan+l_dlan, x286lon-l_dlon))
result.append(GetMeas(2, 'X286|sp_d0', x286lan+l_dlan, x286lon-0))
result.append(GetMeas(3, 'X286|sp_dp', x286lan+l_dlan, x286lon+l_dlon))
result.append(GetMeas(4, 'X286|s0_dm', x286lan+0, x286lon-l_dlon))
result.append(GetMeas(5, 'X286|s0_d0', x286lan+0, x286lon-0))
result.append(GetMeas(6, 'X286|s0_dp', x286lan+0, x286lon+l_dlon))
result.append(GetMeas(7, 'X286|sm_dm', x286lan-l_dlan, x286lon-l_dlon))
result.append(GetMeas(8, 'X286|sm_d0', x286lan-l_dlan, x286lon-0))
result.append(GetMeas(9, 'X286|sm_dp', x286lan-l_dlan, x286lon+l_dlon))
'''
result.append(GetMeas(0, 'TEST_00', x286lan, x286lon))
result.append(GetMeas(1, 'TEST_01', 44.3251277777778, 38.6788694444444))
result.append(GetMeas(2, 'TEST_02', 44.3727555555556, 38.7005))
result.append(GetMeas(3, 'TEST_03', 44.4186305555556, 38.6651361111111))
result.append(GetMeas(4, 'TEST_04', 44.4545888888889, 38.6328638888889))
result.append(GetMeas(5, 'TEST_05', 44.4508222222222, 38.6833333333333))
result.append(GetMeas(6, 'TEST_06', 44.4329277777778, 38.7327722222222))
result.append(GetMeas(7, 'TEST_07', 44.4106166666667, 38.7633277777778))
result.append(GetMeas(8, 'TEST_08', 44.3846138888889, 38.7533722222222))
result.append(GetMeas(9, 'TEST_09', 44.3416583333333, 38.735175))
return result
# получить измеритель (имя, номер, х-коорд картинки и у-коорд картинки) из широты и долготы
def GetMeas(_number, _name, _lan, _lon):
if _lan <= min_lan or _lan >= max_lan:
print("MeasUnit named {} have not current lan ({})".format(_name, _lan))
return None
if _lon <= min_lon or _lon >= max_lon:
print("MeasUnit named {} have not current lon ({})".format(_name, _lon))
return None
lan_delta = _lan - base_lan
lon_delta = _lon - base_lon
y_delta = lan_delta * lenght_lan
x_delta = lon_delta * lenght_lon
my = y_delta // base_pixel_len
mx = x_delta // base_pixel_len
result = dict(number = _number, name = _name, x = abs(mx), y = abs(my))
return result
# Получить осадки по значению dbZ
def GetPrecByDBZ_2(_dbz):
result = 0
if _dbz == 0 or _dbz > 127:
return result
step1 = (_dbz - 32) / 10
step2 = 10 ** step1
step3 = step2 / 200
step4 = step3 ** 0.625
result = step4
return result
# Сгенерировать полигон по набору геокоординат
def GeneratePolygons(_meas_list, _number, _name):
# Максимальные значения координат Х и У для полигона
max_x = 0
max_y = 0
min_x = 256
min_y = 256
#Временный полигон
tmp_pol = []
# Перебрать измерители последовательно и сформировать полигон как набор пар (х, у) плюс данные о границах полигона (прямоугольных)
for imeas in _meas_list:
mx =imeas["x"]
my = imeas["y"]
# Ищем прямоугольные границы полигона
if mx > max_x:
max_x = mx
if my > max_y:
max_y = my
if mx < min_x:
min_x = mx
if my < min_y:
min_y = my
tmp_pol.append((mx, my))
result = dict(number = _number, name = _name, max_x = max_x, max_y = max_y, min_x = min_x, min_y = min_y, pol = tmp_pol)
result["pol"] = tmp_pol
return result
# Получить "обрезанную" матрицу по полигону
def GetCuttedMatrix(_image, _mpolygon):
matrix = _mpolygon["pol"]
max_x = int(_mpolygon["max_x"])
max_y = int(_mpolygon["max_y"])
min_x = int(_mpolygon["min_x"])
min_y = int(_mpolygon["min_y"])
mask = np.zeros(_image.shape, dtype=np.uint8)
roi_corners = np.array([matrix], dtype=np.int32)
channel_count = _image.shape[2] # i.e. 3 or 4 depending on your image
ignore_mask_color = (255,)*channel_count
cv2.fillPoly(mask, roi_corners, ignore_mask_color)
masked_image = cv2.bitwise_and(_image, mask)
masked_image = masked_image[min_y:max_y, min_x:max_x]
return masked_image
# Получаем все файлы из директории с сохранёнными изображениями
png_files = os.listdir(png_path)
print(png_files)
print("Lan border = {} - {}".format(min_lan, max_lan))
print("Lon border = {} - {}".format(min_lon, max_lon))
# Сформировать массив "матрица картинки-время" (время из названия файла)
store = []
for png_f in png_files:
dt_value = png_f.replace('dj_286m_', '').replace('.png', '')
dt_value = int(dt_value)
dt_value = datetime.datetime.fromtimestamp(dt_value)
date_diff = datetime.timedelta(hours=2)
dt_value = dt_value + date_diff
file_data = imread(png_path + '\\' + png_f)
st_unit = dict(
datetime = dt_value,
matrix = file_data
)
store.append(st_unit)
# Генерируем тестовые измерители
prec_meas = GenerateMeasurements()
# Формируем данные
main_data = ""
################################################################
# TEST - генерация обрезанных файлов
main_pol = GeneratePolygons([prec_meas[0], prec_meas[3], prec_meas[4], prec_meas[5], prec_meas[6]], 1, "Verh")
main_pol = GeneratePolygons([prec_meas[0], prec_meas[6], prec_meas[7], prec_meas[8], prec_meas[9], prec_meas[1], prec_meas[2]], 2, "Nige_Polkovnichego")
pol_store = []
pol_store.append(main_pol)
cutted_store = []
'''
for m_image in store:
image = m_image["matrix"]
masked_image = GetCuttedMatrix(image, main_pol)
path_to_save = "d:\\PYTHON\\tests\\test1\\cutted\\{name}_{dt}.png".format(name = main_pol["name"], dt = m_image["datetime"].timestamp())
#cv2.imwrite(path_to_save, masked_image)
cutted_store.append(masked_image)
'''
################################################################
# Последниё интервал расчёта
last_datetime = store[0]["datetime"]
# Дельта (служебная) между интервалами расчёта
dt = datetime.timedelta(minutes=10)
# Пробегаемся по всем картинкам в хранилище
for prec_data in store:
# Формируем локальные переменные для цикла
lmatrix = prec_data["matrix"]
ldatetime = prec_data["datetime"]
date_diff = ldatetime - last_datetime
local_data = ""
# Фиксим (пустой строкой) промежутки между измерениями, если таковые имеются
while date_diff.seconds >= 1200:
local_data = ""
last_datetime = last_datetime + dt
local_data += last_datetime.__str__() + str_limiter
main_data += local_data + "\n"
date_diff = ldatetime - last_datetime
# Формируем для каждого полигона осадки, согласно его прямоугольным координатам относительно картинки
local_data = ""
local_data += ldatetime.__str__() + str_limiter
for pol in pol_store:
masked_image = GetCuttedMatrix(lmatrix, pol)
sum_prec = 0
for yy in masked_image:
for xx in yy:
lm_value = xx[0]
m_prec = GetPrecByDBZ_2(lm_value)
sum_prec += m_prec
local_data += sum_prec.__str__() + str_limiter
# Вываливаем данные
main_data += local_data + "\n"
last_datetime = ldatetime
# Заменяем временные разделители и прочие символы для правильного отобраежния в Эксель
main_data = main_data.replace("||", "\t")
main_data = main_data.replace(".", ",")
# Пишем результат в текстовый файл
f = open('result.txt', 'w')
f.write(main_data)
f.close()
"""
1. Для каждого измерителя ищем разницу в широте и долготе лот базовых координат
1.1 Проверяем на выход за предельные значения
2. По здначениям дельты и длинам градусов ищем разницу в метрах
3. По разнице в метрах и базовому размеру пикселя ищем нужный квадрат как количество целых делений дельты в метрах на базовую длину градуса
"""